c++ - ¿Por qué es tan lento el valarray?

Question

Disculpe por mi pregunta por el valarray de nuevo. Estoy tratando de usarlo, ya que se parece mucho a matlab mientras opera el vector y las matrices. La primera vez que hice una comprobación de rendimiento y descubrí que valarray no puede lograr el rendimiento declarado como en el lenguaje de programación c ++ del libro por stroustrup.

El programa de prueba en realidad hizo 5M multiplicación de dobles. Pensé que c = a * b sería al menos comparable a la multiplicación de elementos de tipo doble de bucle for, pero estoy totalmente equivocado. Probado en varias computadoras y vc6.0 y vs2008.

Por cierto, lo probé en matlab usando el siguiente código:

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

Y el resultado es de 46ms. Este tiempo no es de alta precisión, solo funciona como referencia.

El código es:

#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;
}

Los resultados en ejecución: (modo de lanzamiento con optimización de velocidad máxima)

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
------------------------------------------------------

Modo de depuración con la misma optimización:

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
------------------------------------------------------

Answer 1

Creo que la respuesta de Michael Burr es correcta. Y puede ser que pueda crear un tipo virtual como el tipo del valor de retorno del operador + y volver a cargar otro operador = para este tipo virtual como operador = (tipo virtual & v) {& valarray = & v; v = NULL;} (en términos generales) . Por supuesto, es difícil implementar la idea en valarray. Pero cuando creas una nueva clase, puedes probar esta idea. Y luego, la eficiencia para el operador + es casi la misma que la del operador + =

Answer 2

El objetivo principal de valarray es ser rápido en máquinas vectoriales, que las máquinas x86 simplemente no lo son. Una buena implementación en una máquina no vectorial debe ser capaz de igualar el rendimiento que se obtiene con algo como:
for (i=0; i < N; ++i) c1[i] = a1[i] * b1[i];

y uno malo, por supuesto, no lo hará. A menos que haya algo en el hardware para acelerar el procesamiento en paralelo, eso será lo más cercano a lo mejor que puede hacer.

Answer 3

Estoy compilando la versión x64, VS 2010. Cambié su código muy ligeramente:

    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" ;

Aquí puedes ver que usé * = en lugar de c = a * b . En bibliotecas matemáticas más modernas, se utilizan mecanismos de plantillas de expresión muy complejos que eliminan este problema. En este caso, en realidad obtuve resultados un poco más rápidos de valarray, aunque probablemente sea porque el contenido ya estaba en la memoria caché. La sobrecarga que está viendo es simplemente temporarios redundantes y nada intrínseco a valarray, específicamente, vería el mismo comportamiento con algo como std::string .

Answer 4

Finalmente conseguí esto a través del uso de la evaluación retrasada. El código puede ser feo ya que estoy empezando a aprender estos conceptos avanzados de c ++. Corrígeme si tienes mejor idea por favor. Muchas gracias por toda su ayuda. Aquí está el código:

#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 temp first, assigned to c and delete the temp
//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;
}

El resultado en ejecución en Visual Studio es:

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

Answer 5

Lo probé en un sistema Linux x86-64 (CPU Sandy Bridge):

gcc 4.5.0:

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

Intel ICC 12.0.2:

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

En ambos casos, acabo de usar -O3 y ninguna otra marca relacionada con la optimización.

Parece que el compilador de MS C ++ y / o la implementación de valarray apestan.

Aquí está el código del OP modificado para Linux:

#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;
}

Answer 6

Sospecho que la razón c = a*b es mucho más lenta que realizar las operaciones, un elemento a la vez es que

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

el operador debe asignar memoria para poner el resultado, luego lo devuelve por valor.

Incluso si se utiliza una "swaptimization" para realizar la copia, esa función todavía tiene la sobrecarga de

• valarray el nuevo bloque para el valarray resultante valarray
• inicializando el nuevo valarray (es posible que esto pueda optimizarse)
• poniendo los resultados en el nuevo valarray
• paginación en la memoria para el nuevo valarray medida que se inicializa o se establece con valores de resultados
• desasignar el antiguo valarray que es reemplazado por el resultado

Answer 7

hmm ... probé bombardeo y es igual que valarray ... y más bombardeo ++ [] operatpr es muy lento

 #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" ;
        }
    }