Warum ist Valarray so langsam?

Question 1

Ich versuche, Valarray zu verwenden, da es MATLAB beim Betrieb von Vektoren und Matrizen sehr ähnlich ist. Ich habe zuerst eine Leistungsprüfung durchgeführt und festgestellt, dass valarray die im Buch deklarierte Leistung nicht erreichen kann Programmiersprache C++ von Stroustrup.

Das Testprogramm führte tatsächlich eine 5-Millionen-Multiplikation von Doubles durch. Ich dachte, dass c = a*b zumindest vergleichbar wäre mit dem for Loop Double Typ Element Multiplikation, aber ich bin völlig falsch. Ich habe es auf mehreren Computern und Microsoft Visual C++ 6.0 und Visual Studio 2008 versucht.

Übrigens habe ich auf MATLAB mit folgendem Code getestet:

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

Und das Ergebnis ist 46 ms. Diese Zeit ist nicht hochpräzise; es funktioniert nur als Referenz.

Der Code lautet:

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

Die Laufergebnisse: (Freigabemodus mit maximaler Geschwindigkeitsoptimierung)

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

Debugging-Modus mit gleicher Optimierung:

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

Question 2

Ich habe es gerade auf einem Linux x86-64-System (Sandy Bridge-CPU) ausprobiert:

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

In beiden Fällen habe ich gerade verwendet -O3 und keine anderen optimierungsbezogenen Flags.

Es sieht so aus, als ob der MS C++-Compiler und/oder die Valarray-Implementierung scheiße sind.

Hier ist der für Linux modifizierte OP-Code:

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

Question 3

Ich vermute, dass der Grund c = a*b ist so viel langsamer als das Ausführen der Operationen Element für Element

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

Der Operator muss Speicher zuweisen, in den das Ergebnis eingefügt werden soll, und gibt es dann als Wert zurück.

Selbst wenn eine “Auslagerungsoptimierung” verwendet wird, um die Kopie durchzuführen, hat diese Funktion immer noch den Overhead von

Zuweisen des neuen Blocks für das Ergebnis valarray

Initialisieren des Neuen valarray (es ist möglich, dass dies wegoptimiert wird)
die Ergebnisse in das Neue übertragen valarray
im Gedächtnis nach Neuem blättern valarray wie es initialisiert oder mit Ergebniswerten gesetzt wird

Altes auflösen valarray die durch das Ergebnis ersetzt wird

Question 4

Der springende Punkt bei Valarray ist es, auf Vektormaschinen schnell zu sein, was x86-Maschinen einfach nicht sind.

Eine gute Implementierung auf einer Nicht-Vektor-Maschine sollte in der Lage sein, die Leistung zu erreichen, die Sie mit etwas wie erhalten

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

und ein schlechter natürlich nicht. Wenn die Hardware nicht etwas enthält, um die Parallelverarbeitung zu beschleunigen, wird das dem Besten, was Sie tun können, ziemlich nahe kommen.

Question 5

Ich habe es endlich durch die verzögerte Auswertung geschafft. Der Code kann hässlich sein, da ich gerade erst anfange, diese fortgeschrittenen C++-Konzepte zu lernen.

Hier ist der Code:

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

Das laufende Ergebnis in Visual Studio ist:

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

Question 6

Ich kompiliere in Version x64, Visual Studio 2010. Ich habe Ihren Code sehr leicht geändert:

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

Hier können Sie sehen, dass ich *= anstelle von verwendet habe c = a * b. In moderneren mathematischen Bibliotheken werden sehr komplexe Ausdrucksvorlagenmechanismen verwendet, die dieses Problem beseitigen. In diesem Fall habe ich tatsächlich etwas schnellere Ergebnisse von valarray erhalten, obwohl das wahrscheinlich nur daran liegt, dass sich der Inhalt bereits in einem Cache befand. Der Overhead, den Sie sehen, sind einfach redundante Temporäre und nichts, was valarray eigen ist, insbesondere – Sie würden das gleiche Verhalten mit etwas wie sehen std::string.

Question 7

Ich denke, die Antwort von Michael Burr ist richtig. Und vielleicht können Sie einen virtuellen Typ als Typ für den Rückgabewert des Operators erstellen +und laden Sie eine andere neu operator= für diesen virtuellen Typ wie operator=(virtual type& v){&valarray=&v;v=NULL;} (grob gesprochen).

Natürlich ist es schwierig, die Idee auf Valarray umzusetzen. Aber wenn Sie eine neue Klasse erstellen, können Sie diese Idee ausprobieren. Und dann die Effizienz für operator+ ist fast das gleiche wie operator+=.

Question 8

Hmm..ich habe es getestet Blitz++ und es ist dasselbe wie valarray … Und außerdem der Blitz ++ [] Bediener ist sehr langsam.

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