Eigencode für die Matrixmultiplikation läuft langsamer als Schleifenmultiplikation mit std::vector

Question

Ich lerne sowohl C++ als auch maschinelles Lernen, also habe ich mich entschieden, die Eigenbibliothek für die Matrizenmultiplikation zu verwenden. Ich habe ein Perzeptron darauf trainiert, eine Ziffer aus der MNIST-Datenbank zu erkennen. Für die Trainingsphase setze ich die Anzahl der Trainingszyklen (bzw. Epochen) auf T = 100.

Die „Trainingsmatrix“ ist eine 10000 x 785-Matrix. Das nullte Element jeder Zeile enthält das 'Label', das die Ziffer identifiziert, auf die die Eingabedaten (die verbleibenden 784 Elemente der Zeile) abgebildet werden.

Es gibt auch einen 784 x 1 „Gewichts“-Vektor, der die Gewichte für jedes der 784 Merkmale enthält. Der Gewichtungsvektor würde mit jedem Eingangsvektor (einer Zeile der Trainingsmatrix ohne das nullte Element) multipliziert und bei jeder Iteration aktualisiert, und dies würde T-mal für jede der 10000 Eingaben geschehen.

Ich habe das folgende Programm geschrieben (das die Essenz dessen einfängt, was ich tue), in dem ich den "Vanille" -Ansatz der Multiplikation der Zeilen einer Matrix mit dem Gewichtsvektor (unter Verwendung von std::vector und Schleifen) mit dem verglich, was ich fühlte das Beste, was ich mit einem Eigen-Ansatz tun konnte. Es ist nicht wirklich eine Multiplikation einer Matrix mit einem Vektor, ich schneide tatsächlich die Zeile der Trainingsmatrix und multipliziere sie mit dem Gewichtsvektor.

Die Zeitdauer für die Trainingsperiode betrug beim std::vector-Ansatz 160,662 ms und beim Eigen-Verfahren meist über 10.000 ms.

Ich kompiliere das Programm mit dem folgenden Befehl:

clang++ -Wall -Wextra -pedantic -O3 -march=native -Xpreprocessor -fopenmp permute.cc -o perm -std=c++17

Ich verwende ein "Mitte" 2012 MacBook Pro mit macOS Catalina und 2,5 GHz Dual Core i5.

#include <iostream>
#include <algorithm>
#include <random>
#include <Eigen/Dense>
#include <ctime>
#include <chrono>
using namespace Eigen;

int main() {
    Matrix<uint8_t, Dynamic, Dynamic> m = Matrix<uint8_t, Dynamic, Dynamic>::Random(10000, 785);
    Matrix<double, 784, 1> weights_m = Matrix<double, 784, 1>::Random(784, 1);
    Matrix<uint8_t, 10000, 1> y_m, t_m;

    std::minstd_rand rng;
    rng.seed(time(NULL));
    std::uniform_int_distribution<> dist(0,1); //random integers between 0 and 1
    for (int i = 0; i < y_m.rows(); i++) {
        y_m(i) = dist(rng);
        t_m(i) = dist(rng);
    }

    int T = 100;
    int err;
    double eta;
    eta = 0.25; //learning rate
    Matrix<double, 1, 1> sum_wx_m;

    auto start1 = std::chrono::steady_clock::now(); //start of Eigen Matrix loop

    for (int iter = 0; iter < T; iter++) {
        for (int i = 0; i < m.rows(); i++) {
            sum_wx_m = m.block(i, 1, 1, 784).cast<double>() * weights_m;
        
            //some code to update y_m(i) based on the value of sum_wx_m which I left out
        
            err = y_m(i) - t_m(i);
            if (fabs(err) > 0) { //update the weights_m matrix if there's a difference between target and predicted
                weights_m = weights_m - eta * err * m.block(i, 1, 1, 784).transpose().cast<double>();
            } 
        }
    }

    auto end1 = std::chrono::steady_clock::now();
    auto diff1 = end1 - start1;
    std::cout << "Eigen matrix time is "<<std::chrono::duration <double, std::milli> (diff1).count() << " ms" << std::endl;

    //checking how std::vector form performs;

    std::vector<std::vector<uint8_t>> v(10000);
    std::vector<double> weights_v(784);
    std::vector<uint8_t> y_v(10000), t_v(10000);

    for (unsigned long i = 0; i < v.size(); i++) {
        for (int j = 0; j < m.cols(); j++) {
            v[i].push_back(m(i, j));
        }
    }

    for (unsigned long i = 0; i < weights_v.size(); i++) {
        weights_v[i] = weights_m(i);
    }

    for (unsigned long i = 0; i < y_v.size(); i++) {
        y_v[i] = dist(rng);
        t_v[i] = dist(rng);
    }

    double sum_wx_v;

    auto start2 = std::chrono::steady_clock::now(); //start of vector loop

    for (int iter = 0; iter < T; iter++) {
        for(unsigned long j = 0; j < v.size(); j++) {
            sum_wx_v = 0.0;
            for (unsigned long k = 1; k < v[0].size() ; k++) {
                sum_wx_v += weights_v[k - 1] * v[j][k];
            }
        
            //some code to update y_v[i] based on the value of sum_wx_v which I left out
        
            err = y_v[j] - t_v[j];
            if (fabs(err) > 0) {//update the weights_v matrix if there's a difference between target and predicted
                for (unsigned long k = 1; k < v[0].size(); k++) {
                    weights_v[k - 1] -= eta * err * v[j][k];
                }
            }
        }
    }

    auto end2 = std::chrono::steady_clock::now();
    auto diff2 = end2 - start2;
    std::cout << "std::vector time is "<<std::chrono::duration <double, std::milli> (diff2).count() << " ms" << std::endl;
}

Welche Änderungen sollte ich vornehmen, um bessere Laufzeiten zu erhalten?

Answer 1

Ist vielleicht nicht die beste Lösung, aber du kannst es versuchen:

• Da die Standarddatenreihenfolge von Eigen Column-Major ist, können Sie Ihre Trainingsmatrix 785 x 10000 haben, sodass jedes Trainingslabel/Datenpaar im Speicher zusammenhängend ist (ändern Sie auch die Zeile, in der sum_wx_m berechnet wird).
• Verwenden Sie eine Version von Blockoperationen mit fester Größe, dh Sie können m.block(i, 1, 1, 784) durch m.block<1,784>(i, 1) ersetzen (in umgekehrter Reihenfolge, wenn Sie Ihre Trainingsmatrix bereits gewechselt haben ). oder Sie können einfach den Datenteil Ihrer Trainingsmatrix abbilden und die .col()-Referenz verwenden [siehe Beispiel unten])

Hier ist Ihr Code, der basierend auf diesen Ideen modifiziert wurde:

#include <iostream>
#include <algorithm>
#include <random>
#include <Eigen/Dense>
#include <ctime>
#include <chrono>
using namespace Eigen;

int main() {
    Matrix<uint8_t, Dynamic, Dynamic> m = Matrix<uint8_t, Dynamic, Dynamic>::Random(785, 10000);
    Map<Matrix<uint8_t, Dynamic, Dynamic>> m_data(m.data() + 785, 784, 10000);

    Matrix<double, 784, 1> weights_m = Matrix<double, 784, 1>::Random(784, 1);
    Matrix<uint8_t, 10000, 1> y_m, t_m;

    std::minstd_rand rng;
    rng.seed(time(NULL));
    std::uniform_int_distribution<> dist(0,1); //random integers between 0 and 1
    for (int i = 0; i < y_m.rows(); i++) {
        y_m(i) = dist(rng);
        t_m(i) = dist(rng);
    }

    int T = 100;
    int err;
    double eta;
    eta = 0.25; //learning rate
     Matrix<double, 1, 1> sum_wx_m;

    auto start1 = std::chrono::steady_clock::now(); //start of Eigen Matrix loop

    for (int iter = 0; iter < T; iter++) {
        for (int i = 0; i < m.cols(); i++) {
            sum_wx_m = weights_m.transpose() * m_data.col(i).cast<double>();
        
            //some code to update y_m(i) based on the value of sum_wx_m which I left out
        
            err = y_m(i) - t_m(i);
            if (fabs(err) > 0) { //update the weights_m matrix if there's a difference between target and predicted
                weights_m = weights_m - eta * err * m_data.col(i).cast<double>();
            } 
        }
    }

    auto end1 = std::chrono::steady_clock::now();
    auto diff1 = end1 - start1;
    std::cout << "Eigen matrix time is "<<std::chrono::duration <double, std::milli> (diff1).count() << " ms" << std::endl;

    //checking how std::vector form performs;

    std::vector<std::vector<uint8_t>> v(10000);
    std::vector<double> weights_v(784);
    std::vector<uint8_t> y_v(10000), t_v(10000);

    for (unsigned long i = 0; i < v.size(); i++) {
        for (int j = 0; j < m.rows(); j++) {
            v[i].push_back(m(j, i));
        }
    }

    for (unsigned long i = 0; i < weights_v.size(); i++) {
        weights_v[i] = weights_m(i);
    }

    for (unsigned long i = 0; i < y_v.size(); i++) {
        y_v[i] = dist(rng);
        t_v[i] = dist(rng);
    }

    double sum_wx_v;

    auto start2 = std::chrono::steady_clock::now(); //start of vector loop

    for (int iter = 0; iter < T; iter++) {
        for(unsigned long j = 0; j < v.size(); j++) {
            sum_wx_v = 0.0;
            for (unsigned long k = 1; k < v[0].size() ; k++) {
                sum_wx_v += weights_v[k - 1] * v[j][k];
            }
        
            //some code to update y_v[i] based on the value of sum_wx_v which I left out
        
            err = y_v[j] - t_v[j];
            if (fabs(err) > 0) {//update the weights_v matrix if there's a difference between target and predicted
                for (unsigned long k = 1; k < v[0].size(); k++) {
                    weights_v[k - 1] -= eta * err * v[j][k];
                }
            }
        }
    }

    auto end2 = std::chrono::steady_clock::now();
    auto diff2 = end2 - start2;
    std::cout << "std::vector time is "<<std::chrono::duration <double, std::milli> (diff2).count() << " ms" << std::endl;
}

Ich habe diesen Code in meinem Ubuntu Desktop mit i7-9700K kompiliert:

g++ -Wall -Wextra -O3 -std=c++17
====================================
Eigen matrix time is 110.523 ms
std::vector time is 117.826 ms


g++ -Wall -Wextra -O3 -march=native -std=c++17
=============================================
Eigen matrix time is 66.3044 ms
std::vector time is 71.2296 ms

Answer 2

Nach Diskussionen mit den Benutzern J. Schultke und puhu habe ich die folgenden Änderungen in meinem Code vorgenommen:

1. Ich habe alle Aufrufe von m.block(i, 1, 1, 784) in m.block<1, 784>(i, 1) geändert , dies reduziert die für die Schleife der Eigenmatrix benötigte Zeit um ein Drittel . (erstmals vorgeschlagen von J. Schultke)
2. Ich habe meine m-Matrix als in RowMajor- Reihenfolge gespeichert deklariert. Dies liegt daran, dass die Eigenmatrizen standardmäßig in der Reihenfolge ColMajor (Spalte-Major) gespeichert werden. Dies würde dazu führen, dass jeder Eintrag in einer Reihe zusammenhängend gespeichert wird. Jetzt würden also die m.block() -Aufrufe, die ich verwende, um auf einen Teil einer Zeile in der m-Matrix zu verweisen, einfach den gesamten Speicher auf einmal abrufen und die Zeit der "Eigen-Matrix" auf unter die "std" reduzieren: :Vektorzeit. (empfohlen von puhu)

Die durchschnittlichen Laufzeiten sind jetzt

cpp:Pro$ ./perm
Eigen matrix time is 134.76 ms
std::vector time is 155.574 ms

und der geänderte Code lautet:

#include <iostream>
#include <algorithm>
#include <random>
#include <Eigen/Dense>
#include <chrono>
#include <ctime>
using namespace Eigen;
int main() {
    Matrix<uint8_t, Dynamic, Dynamic, RowMajor> m = Matrix<uint8_t, Dynamic, Dynamic, RowMajor>::Random(10000, 785);
    Matrix<double, 784, 1> weights_m = Matrix<double, 784, 1>::Random(784, 1);
    Matrix<uint8_t, 10000, 1> y_m, t_m;
    std::minstd_rand rng;
    rng.seed(time(NULL));
    std::uniform_int_distribution<> dist(0,1); //random integers between 0 and 1
    for (int i = 0; i < y_m.rows(); i++) {
        y_m(i) = dist(rng);
        t_m(i) = dist(rng);
    }

    int T = 100;
    int err;
    double eta;
    eta = 0.25; //learning rate
    Matrix<double, 1, 1> sum_wx_m;

    auto start1 = std::chrono::steady_clock::now(); //start of Eigen Matrix loop

    for (int iter = 0; iter < T; iter++) {
        for (int i = 0; i < m.rows(); i++) {
            auto b = m.block<1, 784>(i, 1).cast<double>();
            sum_wx_m = b * weights_m;
    
            //some code to update y_m(i) based on the value of sum_wx_m which I left out
    
            err = y_m(i) - t_m(i);
            if (fabs(err) > 0) { //update the weights_m matrix if there's a difference between target and predicted
                weights_m = weights_m - eta * err * b.transpose();
            } 
        }
    }

    auto end1 = std::chrono::steady_clock::now();
    auto diff1 = end1 - start1;
    std::cout << "Eigen matrix time is "<<std::chrono::duration <double, std::milli> (diff1).count() << " ms" << std::endl;

    //checking how std::vector form performs;

    std::vector<std::vector<uint8_t>> v(10000);
    std::vector<double> weights_v(784);
    std::vector<uint8_t> y_v(10000), t_v(10000);

    for (unsigned long i = 0; i < v.size(); i++) {
        for (int j = 0; j < m.cols(); j++) {
            v[i].push_back(m(i, j));
        }
    }

    for (unsigned long i = 0; i < weights_v.size(); i++) {
        weights_v[i] = weights_m(i);
    }

    for (unsigned long i = 0; i < y_v.size(); i++) {
        y_v[i] = dist(rng);
        t_v[i] = dist(rng);
    } 

    double sum_wx_v;

    auto start2 = std::chrono::steady_clock::now(); //start of vector loop

    for (int iter = 0; iter < T; iter++) {
        for(unsigned long j = 0; j < v.size(); j++) {
            sum_wx_v = 0.0;
            for (unsigned long k = 1; k < v[0].size() ; k++) {
                sum_wx_v += weights_v[k - 1] * v[j][k];
            }
    
            //some code to update y_v[i] based on the value of sum_wx_v which I left out
    
            err = y_v[j] - t_v[j];
            if (fabs(err) > 0) {//update the weights_v matrix if there's a difference between target and predicted
                for (unsigned long k = 1; k < v[0].size(); k++) {
                    weights_v[k - 1] -= eta * err * v[j][k];
                }
            }
        }
    }

    auto end2 = std::chrono::steady_clock::now();
    auto diff2 = end2 - start2;
    std::cout << "std::vector time is "<<std::chrono::duration <double, std::milli> (diff2).count() << " ms" << std::endl;
}