为什么 valarray 这么慢？答案

【问题标题】：Why is valarray so slow?为什么 valarray 这么慢？
【发布时间】：2011-10-14 14:45:24
【问题描述】：

我正在尝试使用 valarray，因为它在操作向量和矩阵时很像 MATLAB。我首先做了一些性能检查，发现 valarray 无法达到 Stroustrup 在 C++ programming language 书中声明的性能。

测试程序实际上做了 500 万次双精度数的乘法运算。我认为 c = a*b 至少可以与 for 循环双类型元素乘法相媲美，但我完全错了。我在几台计算机和 Microsoft Visual C++ 6.0 和 Visual Studio 2008 上进行了尝试。

顺便说一句，我在 MATLAB 上使用以下代码进行了测试：

len = 5*1024*1024;
a = rand(len, 1);
b = rand(len, 1);
c = zeros(len, 1);
tic;
c = a.*b;
toc;

结果是 46 毫秒。这个时间精度不高；仅供参考。

代码是：

#include <iostream>
#include <valarray>
#include <iostream>
#include "windows.h"

using namespace std;
SYSTEMTIME stime;
LARGE_INTEGER sys_freq;

double gettime_hp();

int main()
{
    enum { N = 5*1024*1024 };
    valarray<double> a(N), b(N), c(N);
    QueryPerformanceFrequency(&sys_freq);
    int i, j;
    for (j=0 ; j<8 ; ++j)
    {
        for (i=0 ; i<N ; ++i)
        {
            a[i] = rand();
            b[i] = rand();
        }

        double* a1 = &a[0], *b1 = &b[0], *c1 = &c[0];
        double dtime = gettime_hp();
        for (i=0 ; i<N ; ++i)
            c1[i] = a1[i] * b1[i];
        dtime = gettime_hp()-dtime;
        cout << "double operator* " << dtime << " ms\n";

        dtime = gettime_hp();
        c = a*b ;
        dtime = gettime_hp() - dtime;
        cout << "valarray operator* " << dtime << " ms\n";

        dtime = gettime_hp();
        for (i=0 ; i<N ; ++i)
            c[i] = a[i] * b[i];
        dtime = gettime_hp() - dtime;
        cout << "valarray[i] operator* " << dtime<< " ms\n";

        cout << "------------------------------------------------------\n";
    }
}

double gettime_hp()
{
    LARGE_INTEGER tick;
    extern LARGE_INTEGER sys_freq;
    QueryPerformanceCounter(&tick);
    return (double)tick.QuadPart * 1000.0 / sys_freq.QuadPart;
}

运行结果：（最大速度优化的发布模式）

double operator* 52.3019 ms
valarray operator* 128.338 ms
valarray[i] operator* 43.1801 ms
------------------------------------------------------
double operator* 43.4036 ms
valarray operator* 145.533 ms
valarray[i] operator* 44.9121 ms
------------------------------------------------------
double operator* 43.2619 ms
valarray operator* 158.681 ms
valarray[i] operator* 43.4871 ms
------------------------------------------------------
double operator* 42.7317 ms
valarray operator* 173.164 ms
valarray[i] operator* 80.1004 ms
------------------------------------------------------
double operator* 43.2236 ms
valarray operator* 158.004 ms
valarray[i] operator* 44.3813 ms
------------------------------------------------------

同样优化的调试模式：

double operator* 41.8123 ms
valarray operator* 201.484 ms
valarray[i] operator* 41.5452 ms
------------------------------------------------------
double operator* 40.2238 ms
valarray operator* 215.351 ms
valarray[i] operator* 40.2076 ms
------------------------------------------------------
double operator* 40.5859 ms
valarray operator* 232.007 ms
valarray[i] operator* 40.8803 ms
------------------------------------------------------
double operator* 40.9734 ms
valarray operator* 234.325 ms
valarray[i] operator* 40.9711 ms
------------------------------------------------------
double operator* 41.1977 ms
valarray operator* 234.409 ms
valarray[i] operator* 41.1429 ms
------------------------------------------------------
double operator* 39.7754 ms
valarray operator* 234.26 ms
valarray[i] operator* 39.6338 ms
------------------------------------------------------

【问题讨论】：

您是否运行了可执行文件？还是您在调试器中尝试过（通过 Visual Studio）？
您使用了哪些优化设置？
调试器和exe没有区别
VC6?真的吗？它已有 13 年历史，早于标准。
在带有-flto -O3 -march=native -std=c++0x 的 GCC 4.6.1 上，我在所有三种情况下的性能几乎相同，从第一名到第三名略有增加。

标签： c++ valarray

【解决方案1】：

我刚刚在 Linux x86-64 系统（Sandy Bridge CPU）上试了一下：

gcc 4.5.0：

double operator* 9.64185 ms
valarray operator* 9.36987 ms
valarray[i] operator* 9.35815 ms

英特尔 ICC 12.0.2：

double operator* 7.76757 ms
valarray operator* 9.60208 ms
valarray[i] operator* 7.51409 ms

在这两种情况下，我只使用了-O3，没有使用其他与优化相关的标志。

看起来 MS C++ 编译器和/或 valarray 实现很烂。

这是针对 Linux 修改的 OP 代码：

#include <iostream>
#include <valarray>
#include <iostream>
#include <ctime>

using namespace std ;

double gettime_hp();

int main()
{
    enum { N = 5*1024*1024 };
    valarray<double> a(N), b(N), c(N) ;
    int i,j;
    for(  j=0 ; j<8 ; ++j )
    {
        for(  i=0 ; i<N ; ++i )
        {
            a[i]=rand();
            b[i]=rand();
        }

        double* a1 = &a[0], *b1 = &b[0], *c1 = &c[0] ;
        double dtime=gettime_hp();
        for(  i=0 ; i<N ; ++i ) c1[i] = a1[i] * b1[i] ;
        dtime=gettime_hp()-dtime;
        cout << "double operator* " << dtime << " ms\n" ;

        dtime=gettime_hp();
        c = a*b ;
        dtime=gettime_hp()-dtime;
        cout << "valarray operator* " << dtime << " ms\n" ;

        dtime=gettime_hp();
        for(  i=0 ; i<N ; ++i ) c[i] = a[i] * b[i] ;
        dtime=gettime_hp()-dtime;
        cout << "valarray[i] operator* " << dtime<< " ms\n" ;

        cout << "------------------------------------------------------\n" ;
    }
}

double gettime_hp()
{
    struct timespec timestamp;

    clock_gettime(CLOCK_REALTIME, &timestamp);
    return timestamp.tv_sec * 1000.0 + timestamp.tv_nsec * 1.0e-6;
}

【讨论】：

很好 - 仅供参考，您能否添加构建中使用的选项（我今晚可能会玩这些东西......）
+1。我在 libc++ 实现上运行了这个基准测试。它没有 MS 慢，但没有 gcc 快（它与您为 ICC 报告的速度大致相同）。原来我在表达式模板引擎中缺少一个键赋值运算符。补充说。现在 libc++ 和 gcc 一样快。致 OP：感谢您的速度测试！（这个问题也是+1）:-)
谢谢 - 我已经添加了一个注释重新编译器开关（在这两种情况下都只是-O3）并且还附加了我用于这些测试的针对 Linux 修改的 OP 代码。

【解决方案2】：

我怀疑c = a*b 比一次执行一个元素的操作慢得多的原因是

template<class T> valarray<T> operator*
    (const valarray<T>&, const valarray<T>&);

操作员必须分配内存来存放结果，然后按值返回。

即使使用“swaptimization”来执行复制，该函数仍然具有

的开销

为产生的valarray分配新块
初始化新的valarray（这可能会被优化掉）
将结果放入新的valarray
在内存中为新的valarray 分页，因为它被初始化或使用结果值设置
释放被结果替换的旧 valarray

【讨论】：

我刚刚搜了一下，它实际上返回了一个引用：template inline valarray<_ty>& operator*=(valarray<_ty>& _L, const _Ty& _R) {_VALGOP2(*= _R ); }
您在上面的评论中发布的声明与问题中发布的代码中使用的声明有两个不同之处：1) operator*= 与使用 operator*() 后跟 @987654331 不同@ 和 2) 这是 *= 运算符的声明，它采用标量参数将 valarray 乘以
Michael，根据分析，如果需要性能，我们没有办法使用valarray。但是，根据本书，这个类是专门为提高性能而设计的。你能给我一些分数吗？还有其他方法可以让我将数组作为一个整体来处理，就像 c++ 中的 valarray 一样吗？谢谢
@Michael：您可能是对的，但请在单独的答案中查看我的 Linux benchamrks - 只要您使用体面的编译器和valarray 实现。
新闻快讯：使用 表达式模板 (en.wikipedia.org/wiki/Expression_templates) 允许但不要求使用 valarray 算术。使用表达式模板可以完全消除 OP 问题中的临时性问题，从而完全消除表达式c = a*b 的堆分配和释放。很明显 gcc 可以做到这一点（以及稍微更正的 libc++），而 MS C++ 没有。

【解决方案3】：

valarray 的全部意义在于在向量机上快速，而 x86 机器则不然。

在非向量机上的良好实现应该能够与您获得的性能相匹配

for (i=0; i < N; ++i) 
    c1[i] = a1[i] * b1[i];

当然不会。除非硬件中有东西可以加速并行处理，否则这将非常接近您可以做到的最佳水平。

【讨论】：

【解决方案4】：

我终于通过延迟评估得到了这个。由于我刚刚开始学习这些 C++ 高级概念，因此代码可能很难看。

代码如下：

#include <iostream>
#include <valarray>
#include <iostream>
#include "windows.h"

using namespace std;
SYSTEMTIME stime;
LARGE_INTEGER sys_freq;

double gettime_hp();

// To improve the c = a*b (it will generate a temporary first, assigned to 'c' and delete the temporary.
// Which causes the program really slow
// The solution is the expression template and let the compiler to decide when all the expression is known.


// Delayed evaluation
//typedef valarray<double> Vector;
class Vector;

class VecMul
{
    public:
        const Vector& va;
        const Vector& vb;
        //Vector& vc;
        VecMul(const Vector& v1, const Vector& v2): va(v1), vb(v2) {}
        operator Vector();
};

class Vector:public valarray<double>
{
    valarray<double> *p;

    public:
        explicit Vector(int n)
        {
            p = new valarray<double>(n);
        }
        Vector& operator = (const VecMul &m)
        {
            for(int i=0; i<m.va.size(); i++)
                (*p)[i] = (m.va)[i]*(m.vb)[i]; // Ambiguous
            return *this;
        }
        double& operator[](int i) const {return (*p)[i];} //const vector_type[i]
        int size()const {return (*p).size();}
};


inline VecMul operator*(const Vector& v1, const Vector& v2)
{
    return VecMul(v1, v2);
}


int main()
{
    enum {N = 5*1024*1024};
    Vector a(N), b(N), c(N);
    QueryPerformanceFrequency(&sys_freq);
    int i, j;
    for (j=0 ; j<8 ; ++j)
    {
        for (i=0 ; i<N ; ++i)
        {
            a[i] = rand();
            b[i] = rand();
        }

        double* a1 = &a[0], *b1 = &b[0], *c1 = &c[0];
        double dtime = gettime_hp();
        for (i=0 ; i<N ; ++i)
            c1[i] = a1[i] * b1[i];
        dtime = gettime_hp()-dtime;
        cout << "double operator* " << dtime << " ms\n";

        dtime = gettime_hp();
        c = a*b;
        dtime = gettime_hp()-dtime;
        cout << "valarray operator* " << dtime << " ms\n";

        dtime = gettime_hp();
        for (i=0 ; i<N ; ++i)
            c[i] = a[i] * b[i];
        dtime = gettime_hp() - dtime;
        cout << "valarray[i] operator* " << dtime << " ms\n";

        cout << "------------------------------------------------------\n";
    }
}

double gettime_hp()
{
    LARGE_INTEGER tick;
    extern LARGE_INTEGER sys_freq;
    QueryPerformanceCounter(&tick);
    return (double)tick.QuadPart*1000.0/sys_freq.QuadPart;
}

在Visual Studio上的运行结果是：

double operator* 41.2031 ms
valarray operator* 43.8407 ms
valarray[i] operator* 42.49 ms

【讨论】：

为什么你的class Vector 有两个 valarray 对象？一个作为基类，另一个分配在堆上？似乎您只使用堆上的那个，但这需要额外的分配。只需使用基类！

【解决方案5】：

我在 Visual Studio 2010 x64 版本中进行编译。我对您的代码稍作改动：

    double* a1 = &a[0], *b1 = &b[0], *c1 = &c[0];
    double dtime = gettime_hp();
    for (i=0 ; i<N ; ++i)
        a1[i] *= b1[i];
    dtime = gettime_hp() - dtime;
    cout << "double operator* " << dtime << " ms\n";

    dtime = gettime_hp();
    a *= b;
    dtime = gettime_hp() - dtime;
    cout << "valarray operator* " << dtime << " ms\n";

    dtime = gettime_hp();
    for (i=0 ; i<N ; ++i)
        a[i] *= b[i];
    dtime = gettime_hp() - dtime;
    cout << "valarray[i] operator* " << dtime<< " ms\n";

    cout << "------------------------------------------------------\n" ;

在这里您可以看到我使用 *= 而不是 c = a * b。在更现代的数学库中，使用了非常复杂的表达式模板机制来消除这个问题。在这种情况下，我实际上从 valarray 获得了稍微快一点的结果，尽管这可能只是因为内容已经在缓存中。您看到的开销只是多余的临时变量，而 valarray 没有任何内在的东西，特别是 - 您会看到与 std::string 类似的行为。

【讨论】：

我验证了你的结果。这种变化虽然不是微小的变化。许多复合表达式不能总是使用 *=, += /=
@shangping：在这种情况下，如果您为所需的每个临时变量分配一个新的结果数组，您会看到double 的速度与valarray 类似。跨度>
“在更现代的数学库中，使用了非常复杂的表达式模板机制来消除这个问题。” 在std::valarray 的高质量实现中也是如此。

【解决方案6】：

我认为 Michael Burr 的回答是正确的。也许你可以创建一个虚拟类型作为运算符+的返回值，然后为这个虚拟类型重新加载另一个operator=，比如operator=(virtual type& v){&valarray=&v;v=NULL;}（粗略地说）。

当然，在 valarray 上实现这个想法是很困难的。但是当你创建一个新类时，你可以试试这个想法。然后，operator+ 的效率几乎与operator+= 相同。

【讨论】：

这应该是对他的回答的评论，而不是新的回答。

【解决方案7】：

嗯..我测试了Blitz++，它和valarray一样......而且Blitz++ []操作符很慢。

#include <blitz/array.h>
#include <iostream>

#ifdef WIN32
#include "windows.h"
LARGE_INTEGER sys_freq;
#endif

#ifdef LINUX
<ctime>
#endif

using namespace std;
SYSTEMTIME stime;

__forceinline double gettime_hp();
double gettime_hp()
{
    #ifdef WIN32
        LARGE_INTEGER tick;
        extern LARGE_INTEGER sys_freq;
        QueryPerformanceCounter(&tick);
        return (double)tick.QuadPart * 1000.0 / sys_freq.QuadPart;
    #endif

    #ifdef LINUX
        struct timespec timestamp;

        clock_gettime(CLOCK_REALTIME, &timestamp);
        return timestamp.tv_sec * 1000.0 + timestamp.tv_nsec * 1.0e-6;
    #endif
}
BZ_USING_NAMESPACE(blitz)

int main()
{
    int N = 5*1024*1024;

    // Create three-dimensional arrays of double
    Array<double, 1> a(N), b(N), c(N);

    int i, j;

    #ifdef WIN32
        QueryPerformanceFrequency(&sys_freq);
    #endif

    for (j=0 ; j<8 ; ++j)
    {
        for (i=0 ; i<N ; ++i)
        {
            a[i] = rand();
            b[i] = rand();
        }

        double* a1 = a.data(), *b1 = b.data(), *c1 = c.data();
        double dtime = gettime_hp();
        for (i=0 ; i<N ; ++i)
            c1[i] = a1[i] * b1[i];
        dtime = gettime_hp() - dtime;
        cout << "double operator* " << dtime << " ms\n";

        dtime = gettime_hp();
        c = a*b;
        dtime = gettime_hp() - dtime;
        cout << "blitz operator* " << dtime << " ms\n";

        dtime = gettime_hp();
        for (i=0 ; i<N ; ++i)
            c[i] = a[i] * b[i];
        dtime = gettime_hp() - dtime;
        cout << "blitz[i] operator* " << dtime<< " ms\n";

        cout << "------------------------------------------------------\n";
    }
}

【讨论】：

...结果如何？