为什么std::valarray那么慢（2）

原创作品，欢迎批评，转载请保留作者信息。

上一篇文章我发现valarray比普通的循环慢了很多，原因在哪里呢？

我们看看valarray是如何实现c=a*b的。

#define _VALOP(TYPE, LENGTH, RHS) \
valarray<TYPE> _Ans(LENGTH); \
for (size_t _I = 0; _I < _Ans.size(); ++_I) \
_Ans[_I] = RHS; \
return (_Ans)

template<class _Ty> inline
valarray<_Ty> operator*(const valarray<_Ty>& _L,
const valarray<_Ty>& _R)
{_VALOP(_Ty, _L.size(), _L[_I] * _R[_I]); }

用比较容易懂得C语言表示也就是：

valarray<T> operator*(const valarray<T>& a,const valarray<T>& b)

{

  valarray<T> ans(a.size());

  for(int i=0;i<ans.size();i++) ans[i]=a[i]*b[i];

  return ans;

}

也就是这个需要构造一个临时对象，并返回。

当遇到c=这个操作时，我们要调用operator=这个成员函数：

void operator=(const valarray<_Ty>& _R) const
{_SLOP(= _R[_I]); }

这个操作符把ans的每一个元素赋值给a, 然后销毁临时对象ans.

这个临时对象的生成与销毁就是对性能影响最大的原因。

找到了这个原因，是不是就没有办法了呢？办法就是Stroustrup《C++ Programming Language》中提到的Lazy Evaluation。其思想就是当执行a*b时，并不进行实际的乘法操作，而只是保存了a和b的引用。当遇到赋值运算符时在进行实际的运算。这样我们就要创造一个新的类来完成这个存储乘法操作信息。

假设我们新加的类叫VecMul，这个类要保存对a和b的引用，遇到a*b时就构造VecMul(a,b)这样一个对象。这样问题就变成了：

c=VecMul(a,b) 左边是valarray类型， 右边是VecMul类型，这样是无法通过的，而且valarray的operator=已经被定义好了。怎么办呢？

左边必须定义为一个新的类，我们从而可以重新定义operator=,假设这个新类是Vector。这个新类主要的功能是完成Vector=VecMul这个操作。

这个新的 Vector我们可以继承std::valarray,这样我们还可以使用valarray的所有特性，而不用再去进行内存管理。为了方便以及也可以存取Vector里面的元素，我还定义了[]这个运算符。

下面是我为改进效率所作的新程序，新程序可以保证valarray的效率可以与C循环媲美。感兴趣的可以一试，注意要在release模式下运行才可以全速运行。

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


using namespace std ;

class hptime
{
LARGE_INTEGER sys_freq;
public:
hptime(){QueryPerformanceFrequency(&sys_freq);}   


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


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


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()
{
hptime t0;
    enum { N = 5*1024*1024 };
    Vector 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=t0.gettime();
        for(  i=0 ; i<N ; ++i ) c1[i] = a1[i] * b1[i] ;
        dtime=t0.gettime()-dtime;
        cout << "double operator* " << dtime << " ms\n" ;


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


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


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

运行时间（release mode using vc6.0）：

double operator* 40.0126 ms
valarray operator* 40.5328 ms
valarray[i] operator* 40.0425 ms
------------------------------------------------------
double operator* 40.3596 ms
valarray operator* 41.2441 ms
valarray[i] operator* 41.1438 ms
------------------------------------------------------
double operator* 40.2104 ms
valarray operator* 40.3188 ms
valarray[i] operator* 39.8274 ms
------------------------------------------------------
double operator* 40.316 ms
valarray operator* 40.7231 ms
valarray[i] operator* 40.4965 ms
------------------------------------------------------
double operator* 40.3205 ms
valarray operator* 42.9521 ms
valarray[i] operator* 40.3596 ms
------------------------------------------------------
double operator* 40.451 ms
valarray operator* 40.77 ms
valarray[i] operator* 40.1738 ms
------------------------------------------------------
double operator* 39.8601 ms
valarray operator* 40.6535 ms
valarray[i] operator* 40.6479 ms
------------------------------------------------------
double operator* 39.688 ms
valarray operator* 40.7351 ms
valarray[i] operator* 39.7916 ms
------------------------------------------------------
Press any key to continue

注意，我们为什么不能在Vector里面同样使用Valarray的引用呢？如果这样，用户就必须同时使用valarray, vector, VecMul三个类，在执行赋值时把valarray转化为Vector类。这样我们就不能保持我们习惯的数学表达式了。

这个方法可以推广，这就是所谓的表达式模板，我将在另一篇文章中作详细的描述。