涉及临时内存分配时避免 blas?
Avoid blas when involving temporary memory allocation?
我有一个程序可以重复计算矩阵乘积 x'Ay
。更好的做法是通过调用 MKL 的 blas 来计算这个,即 cblas_dgemv
和 cblas_ddot
,这需要将内存分配给一个临时向量,或者最好简单地取 x_i * a_ij * y_j
的总和?换句话说,MKL的blas理论上有没有增加什么价值?
我为我的笔记本电脑进行了基准测试。除了 g++_no_blas 的性能是其他测试的两倍之外(为什么?),每个测试几乎没有区别。 O2、O3和Ofast之间也没有区别。
- g++_blas_static 57ms
- g++_blas_dynamic 58ms
- g++_no_blas 100 毫秒
- icpc_blas_static 57ms
- icpc_blas_dynamic58ms
- icpc_no_blas58ms
util.h
#ifndef UTIL_H
#define UTIL_H
#include <random>
#include <memory>
#include <iostream>
struct rng
{
rng() : unif(0.0, 1.0)
{
}
std::default_random_engine re;
std::uniform_real_distribution<double> unif;
double rand_double()
{
return unif(re);
}
std::unique_ptr<double[]> generate_square_matrix(const unsigned N)
{
std::unique_ptr<double[]> p (new double[N * N]);
for (unsigned i = 0; i < N; ++i)
{
for (unsigned j = 0; j < N; ++j)
{
p.get()[i*N + j] = rand_double();
}
}
return p;
}
std::unique_ptr<double[]> generate_vector(const unsigned N)
{
std::unique_ptr<double[]> p (new double[N]);
for (unsigned i = 0; i < N; ++i)
{
p.get()[i] = rand_double();
}
return p;
}
};
#endif // UTIL_H
main.cpp
#include <iostream>
#include <iomanip>
#include <memory>
#include <chrono>
#include "util.h"
#include "mkl.h"
double vtmv_blas(double* x, double* A, double* y, const unsigned n)
{
double temp[n];
cblas_dgemv(CblasRowMajor, CblasNoTrans, n, n, 1.0, A, n, y, 1, 0.0, temp, 1);
return cblas_ddot(n, temp, 1, x, 1);
}
double vtmv_non_blas(double* x, double* A, double* y, const unsigned n)
{
double r = 0;
for (unsigned i = 0; i < n; ++i)
{
for (unsigned j = 0; j < n; ++j)
{
r += x[i] * A[i*n + j] * y[j];
}
}
return r;
}
int main()
{
std::cout << std::fixed;
std::cout << std::setprecision(2);
constexpr unsigned N = 10000;
rng r;
std::unique_ptr<double[]> A = r.generate_square_matrix(N);
std::unique_ptr<double[]> x = r.generate_vector(N);
std::unique_ptr<double[]> y = r.generate_vector(N);
auto start = std::chrono::system_clock::now();
const double prod = vtmv_blas(x.get(), A.get(), y.get(), N);
auto end = std::chrono::system_clock::now();
auto duration = std::chrono::duration_cast<std::chrono::milliseconds>(
end - start);
std::cout << "Result: " << prod << std::endl;
std::cout << "Time (ms): " << duration.count() << std::endl;
GCC no blas 很差,因为它不使用矢量化 SMID 指令,而其他的都使用。 icpc 将自动矢量化你的循环。
您没有显示矩阵大小,但通常 gemv 受内存限制。由于矩阵比临时向量大得多,因此消除它可能无法大幅提高性能。
我有一个程序可以重复计算矩阵乘积 x'Ay
。更好的做法是通过调用 MKL 的 blas 来计算这个,即 cblas_dgemv
和 cblas_ddot
,这需要将内存分配给一个临时向量,或者最好简单地取 x_i * a_ij * y_j
的总和?换句话说,MKL的blas理论上有没有增加什么价值?
我为我的笔记本电脑进行了基准测试。除了 g++_no_blas 的性能是其他测试的两倍之外(为什么?),每个测试几乎没有区别。 O2、O3和Ofast之间也没有区别。
- g++_blas_static 57ms
- g++_blas_dynamic 58ms
- g++_no_blas 100 毫秒
- icpc_blas_static 57ms
- icpc_blas_dynamic58ms
- icpc_no_blas58ms
util.h
#ifndef UTIL_H
#define UTIL_H
#include <random>
#include <memory>
#include <iostream>
struct rng
{
rng() : unif(0.0, 1.0)
{
}
std::default_random_engine re;
std::uniform_real_distribution<double> unif;
double rand_double()
{
return unif(re);
}
std::unique_ptr<double[]> generate_square_matrix(const unsigned N)
{
std::unique_ptr<double[]> p (new double[N * N]);
for (unsigned i = 0; i < N; ++i)
{
for (unsigned j = 0; j < N; ++j)
{
p.get()[i*N + j] = rand_double();
}
}
return p;
}
std::unique_ptr<double[]> generate_vector(const unsigned N)
{
std::unique_ptr<double[]> p (new double[N]);
for (unsigned i = 0; i < N; ++i)
{
p.get()[i] = rand_double();
}
return p;
}
};
#endif // UTIL_H
main.cpp
#include <iostream>
#include <iomanip>
#include <memory>
#include <chrono>
#include "util.h"
#include "mkl.h"
double vtmv_blas(double* x, double* A, double* y, const unsigned n)
{
double temp[n];
cblas_dgemv(CblasRowMajor, CblasNoTrans, n, n, 1.0, A, n, y, 1, 0.0, temp, 1);
return cblas_ddot(n, temp, 1, x, 1);
}
double vtmv_non_blas(double* x, double* A, double* y, const unsigned n)
{
double r = 0;
for (unsigned i = 0; i < n; ++i)
{
for (unsigned j = 0; j < n; ++j)
{
r += x[i] * A[i*n + j] * y[j];
}
}
return r;
}
int main()
{
std::cout << std::fixed;
std::cout << std::setprecision(2);
constexpr unsigned N = 10000;
rng r;
std::unique_ptr<double[]> A = r.generate_square_matrix(N);
std::unique_ptr<double[]> x = r.generate_vector(N);
std::unique_ptr<double[]> y = r.generate_vector(N);
auto start = std::chrono::system_clock::now();
const double prod = vtmv_blas(x.get(), A.get(), y.get(), N);
auto end = std::chrono::system_clock::now();
auto duration = std::chrono::duration_cast<std::chrono::milliseconds>(
end - start);
std::cout << "Result: " << prod << std::endl;
std::cout << "Time (ms): " << duration.count() << std::endl;
GCC no blas 很差,因为它不使用矢量化 SMID 指令,而其他的都使用。 icpc 将自动矢量化你的循环。
您没有显示矩阵大小,但通常 gemv 受内存限制。由于矩阵比临时向量大得多,因此消除它可能无法大幅提高性能。