std::accumulate() 仅是复数 std::vector 的实部

Posted 2023-02-21

【中文标题】std::accumulate() 仅是复数 std::vector 的实部

【英文标题】：std::accumulate() only the real part of a complex std::vector

【发布时间】：2016-07-20 18:31:30

【问题描述】：

我曾经取对称（Hermitian）矩阵的总和（该矩阵在std::vector 中），这是一个巨大的浪费，因为虚部总是会加到零（我说的是巨大的矩阵n>1000 的边长，所以我认为这会有所不同。

所以现在我只添加矩阵的上三角部分。但我想进一步优化这一点，避免添加复杂的部分，因为我不需要它。

我目前使用的是：

std::real(std::accumulate(myMatrix.begin(), myMatrix.end(), std::complex<Real>(0,0)));

这会将myMatrix 的所有元素添加到std::complex<Real>(0,0)，得到我需要的总和。

但这会添加我的向量的实部和虚部，这是一种浪费！如何编写仅添加此矩阵的实部的最优化版本？

---

更新：

虽然我接受了有效的答案，但我发现它比对矩阵的实部和虚部求和要慢。对于边长为 128 的矩阵，它慢了 5%-10%。这令人惊讶。任何其他更快的建议都将受到高度赞赏。

如果您需要更多信息，请询问。

【问题讨论】：

如果你-1这个帖子请解释一下。我将不胜感激。

您是否考虑过基于 std::thread 的并行或 cuda？另外，对于一个对称矩阵，你只需要访问 n(n-1)/2 个元素，而不是 n*n 个元素来进行累加。

@FengWang 实际上我正在使用 OpenMP 进行并行化，并且我已经访问了 n(n-1)/2 个元素 :)

对于大型数组的简单聚合，性能瓶颈确实在于内存带宽。 Real 是多少字节？如果Real 定义为 4 字节浮点数，则 64 位机器加载 8 个对齐字节的速度可能与其前 4 个字节一样快。您可以通过将 myMatrix 分成实部和虚部来测试这一点，两者都由std::vector<float> 表示，然后只聚合实部，但如果您不必在不同的表示之间进行大量转换，这将是有益的。跨度>

【参考方案1】：

Real real_sum = std::accumulate(
    myMatrix.cbegin(), myMatrix.cend(), Real,
    [](Real const acc, std::complex<Real> const& c)  return acc + std::real(c); 
);

【讨论】：

顺便说一句，我刚刚发现这在计算上比较慢:(

@TheQuantumPhysicist ：试试[](Real& acc, std::complex<Real> const c) -> Real& return acc += c.real(); ——它的表现有什么不同吗？

仍然……这比实部和虚部之和慢 5%-10%……我真的很惊讶。

@TheQuantumPhysicist 可能复杂的加法版本得到了 SIMD

【参考方案2】：

std::accumulate 有两个重载，其中一个采用运算符：

template< class InputIt, class T, class BinaryOperation >
T accumulate( InputIt first, InputIt last, T init,
              BinaryOperation op );

只需提供您自己的而不是默认为+：

std::accumulate(myMatrix.begin(), myMatrix.end(), Real0,
    [](Real const& sum, std::complex<Real> const& next)
        return sum + std::real(next);
    );

---

或者，你可以做一些有趣的事情，比如使用boost::transform_iterator：

auto make_real = [](std::complex<Real> const& c) 
    return std::real(c);
;

std::accumulate(
    boost::make_transform_iterator(myMatrix.begin(), make_real),
    boost::make_transform_iterator(myMatrix.end(), make_real),
    Real0);

或者range-v3:

accumulate(myMatrix,
    Real0,
    std::plus<>,
    [](std::complex<Real> const& c)  return c.real(); 
);

如果real 没有重载，上述两种方法都会更好，您可以在boost 示例中提供std::real<Real>，在第二个示例中提供&std::complex<Real>::real。

【参考方案3】：

绝对有必要使用std::accumulate 吗？

    Real acc = 0;
    for(auto c : myMatrix)
       acc += real(c);

不要误会我的意思，我赞成在适合的情况下使用标准算法，但在可读性方面似乎很难超越这个循环。

这与我的 g++-4.8.4 安装附带的实现相比：

  template<typename _InputIterator, typename _Tp>
    inline _Tp
    accumulate(_InputIterator __first, _InputIterator __last, _Tp __init)
    
      // concept requirements
      __glibcxx_function_requires(_InputIteratorConcept<_InputIterator>)
      __glibcxx_requires_valid_range(__first, __last);

      for (; __first != __last; ++__first)
    __init = __init + *__first;
      return __init;
    

所以你可以看到他们几乎在做同样的事情。

【讨论】：

std 算法针对处理器进行了优化。所以是的，如果你想要最好的结果，你应该尽可能地使用标准算法。

通过使用auto c 而不是auto& c，您正在制作副本。这是你的意图吗？

@kfsone 我相信编译器会做出避免复制的优化。