Kod własny dla mnożenia macierzy działa wolniej niż mnożenie zapętlone przy użyciu std :: vector

Aug 15 2020

Uczę się języka C ++, a także uczenia maszynowego, więc zdecydowałem się skorzystać z biblioteki Eigen do mnożenia macierzy. Uczyłem perceptron rozpoznawania cyfry z bazy danych MNIST. Dla fazy treningu ustawiłem liczbę cykli treningowych (lub epok) na T = 100.

„Macierz szkoleniowa” to macierz 10000 x 785. Zerowy element każdego wiersza zawiera „etykietę” identyfikującą cyfrę, na którą są odwzorowywane dane wejściowe (pozostałe 784 elementy wiersza).

Istnieje również wektor „wag” 784 x 1, który zawiera wagi dla każdej z 784 cech. Wektor wag byłby mnożony z każdym wektorem wejściowym (wiersz macierzy szkoleniowej z wyłączeniem elementu zerowego) i byłby aktualizowany w każdej iteracji, co miałoby miejsce T razy dla każdego z 10000 danych wejściowych.

Napisałem następujący program (który oddaje istotę tego, co robię), w którym porównałem "waniliowe" podejście mnożenia wierszy macierzy z wektorem wagi (używając std :: vector i pętli) do tego, co czułem najlepsze, co mogłem zrobić, stosując podejście własne. To nie jest tak naprawdę mnożenie macierzy z wektorem, tak naprawdę przecinam wiersz macierzy uczącej i mnożę go przez wektor wagi.

Czas trwania okresu treningowego dla podejścia std :: vector wynosił 160,662 ms, a dla metody Eigen wynosił zwykle ponad 10 000 ms.

Kompiluję program poleceniem:

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

Używam MacBooka Pro z połowy 2012 r. Z systemem macOS Catalina i dwurdzeniowym procesorem i5 2,5 GHz.

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

Jakie zmiany należy wprowadzić, aby uzyskać lepsze czasy działania?

Odpowiedzi

1 puhu Aug 16 2020 at 00:34

Może nie jest to najlepsze rozwiązanie, ale możesz spróbować:

  • Ponieważ domyślną kolejnością danych własnych jest kolumna główna, możesz pozwolić macierzy treningowej na 785x10000, tak aby każda etykieta szkoleniowa / para danych była ciągła w pamięci (zmień również wiersz, w którym obliczana jest suma_wx_m).
  • Użyj wersji stałej wielkości operacji blokowych, czyli, można zastąpić m.block (I, 1, 1, 784) z m.block <1,784> (i, 1) (w odwrotnej kolejności, jeśli już włączony swoją macierz szkoleniowy układ lub możesz po prostu zmapować część danych swojej macierzy treningowej i użyć odniesienia .col () [patrz poniższy przykład])

Oto twój kod zmodyfikowany na podstawie tych pomysłów:

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

Skompilowałem ten kod na moim Ubuntu Desktop z i7-9700K:

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
tf3 Aug 16 2020 at 10:22

Po rozmowach z użytkownikami J. Schultke i puhu dokonałem następujących zmian w moim kodzie:

  1. Zmieniłem wszystkie wywołania m.block (i, 1, 1, 784) na m.block <1, 784> (i, 1) , co zmniejsza czas potrzebny na pętlę macierzy Eigena o jedną trzecią . (pierwszy zasugerowany przez J.Schultke)
  2. Deklarowałem, że moja macierz m jest przechowywana w kolejności RowMajor . Dzieje się tak, ponieważ domyślnie macierze własne są przechowywane w kolejności ColMajor (kolumna główna). Spowodowałoby to, że każdy wpis w wierszu byłby przechowywany w sposób ciągły. Więc teraz wywołania m.block () , których używam w odniesieniu do wycinka wiersza w macierzy m, po prostu pobierałyby cały fragment pamięci naraz, skracając czas „macierzy własnej” poniżej „std: : wektor "czas. (sugerowane przez puhu)

Średnie czasy działania są teraz

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

a zmodyfikowany kod to:

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