Собственный код для матричного умножения работает медленнее, чем циклическое умножение с использованием std :: vector
Я изучаю C ++, а также машинное обучение, поэтому решил использовать библиотеку Eigen для умножения матриц. Я тренировал перцептрон распознавать цифры из базы данных MNIST. Для фазы обучения я установил количество тренировочных циклов (или эпох) равным T = 100.
«Матрица обучения» - это матрица размером 10000 x 785. Нулевой элемент каждой строки содержит «метку», идентифицирующую цифру, которой соответствуют входные данные (оставшиеся 784 элемента строки).
Существует также вектор «весов» 784 x 1, который содержит веса для каждой из 784 функций. Вектор весов будет умножен на каждый входной вектор (строка обучающей матрицы, исключая нулевой элемент) и будет обновляться каждую итерацию, и это будет происходить T раз для каждого из 10000 входов.
Я написал следующую программу (которая отражает суть того, что я делаю), в которой я сравнил "ванильный" подход умножения строк матрицы на вектор веса (с использованием std :: vector и циклов) с тем, что, по моему мнению, было лучшее, что я мог сделать с подходом Эйгена. На самом деле это не умножение матрицы на вектор, я на самом деле разрезаю строку обучающей матрицы и умножаю ее на вектор веса.
Продолжительность периода обучения для подхода std :: vector составляла 160,662 мс, а для метода Eigen обычно превышала 10 000 мс.
Я компилирую программу с помощью следующей команды:
clang++ -Wall -Wextra -pedantic -O3 -march=native -Xpreprocessor -fopenmp permute.cc -o perm -std=c++17
Я использую MacBook Pro середины 2012 года, работающий под управлением macOS Catalina и имеющий двухъядерный процессор i5 2,5 ГГц.
#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;
}
Какие изменения мне следует внести, чтобы улучшить время работы?
Ответы
Возможно, это не лучшее решение, но вы можете попробовать:
- Поскольку по умолчанию порядок данных для Eigen равен Column-Major, вы можете позволить вашей обучающей матрице быть 785x10000, чтобы каждая пара обучающей метки / данных была смежной в памяти (также измените строку, в которой вычисляется sum_wx_m).
- Используйте версию блочных операций фиксированного размера, т.е. вы можете заменить m.block (i, 1, 1, 784) на m.block <1,784> (i, 1) (в обратном порядке, если вы уже переключили свою обучающую матрицу макет или вы можете просто сопоставить часть данных вашей обучающей матрицы и использовать ссылку .col () [см. пример ниже])
Вот ваш код, измененный на основе этих идей:
#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;
}
Я скомпилировал этот код на своем рабочем столе Ubuntu с 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
После обсуждения с пользователями J. Schultke и puhu я внес в свой код следующие изменения:
- Я изменил все вызовы m.block (i, 1, 1, 784) на m.block <1, 784> (i, 1) , это на треть сокращает время, необходимое для цикла матрицы Eigen . (впервые предложено Дж. Шультке)
- Я объявил, что моя матрица m хранится в порядке RowMajor . Это связано с тем, что по умолчанию матрицы Eigen хранятся в порядке ColMajor (по столбцам). Это приведет к тому, что каждая запись в строке будет храниться непрерывно. Итак, теперь вызовы m.block () , которые я использую для ссылки на фрагмент строки в матрице m, будут просто извлекать сразу весь фрагмент памяти, уменьшая время «собственной матрицы» до уровня ниже «std: : vector "время. (предложено пуху)
Среднее время работы сейчас составляет
cpp:Pro$ ./perm
Eigen matrix time is 134.76 ms
std::vector time is 155.574 ms
и модифицированный код:
#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;
}