2022/05/13(金)Ubuntu 22.04のBLAS (dgemm) をベンチマーク
BLAS, LAPACKとは何か、個々のBLASの特徴などは、20.04のときの記事を見て下さい。前回と同様、Reference BLAS, ATLAS, OpenBLAS, BLIS, MKLの5つのパッケージをaptでインストールし、dgemmを呼び出して計算時間を計測するプログラムを使って、update-alternativesの機能でlibblas.soを切り替えながらベンチマークを取り、計算時間を比較します。前回と大雑把な傾向はほとんど変わっていないのですが、MKLのパッケージが怪しげな挙動を示した(ような気がした)ので備忘録として書いています。
各BLASのインストール方法
- Reference BLAS
apt install libblas-dev
- ATLAS
apt install libatlas-base-dev
- OpenBLAS
apt install libopenblas-base apt install libopenblas-dev
- BLIS
apt install libblis-dev
- MKL
apt install intel-mkl最後のMKLですが、aptの実行中に
Use libmkl_rt.so as the default alternative to BLAS/LAPACK? [yes/no]と聞かれます。これには「no」と答えました。「yes」と答えると、「-lblas」を付けてコンパイルしたものはlibblas.soがリンクされずにMKL(libmkl_rt.so)がリンクされてしまってupdate-alternativesでlibblas.soの実体を切り替える作戦が不可能になってしまいました。何か自分が勘違いしてるのかも知れません。
これらのライブラリは、
update-alternatives --config libblas.so-x86_64-linux-gnuとすると
There are 5 choices for the alternative libblas.so-x86_64-linux-gnu (providing /usr/lib/x86_64-linux-gnu/libblas.so). Selection Path Priority Status ------------------------------------------------------------ 0 /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so 100 auto mode 1 /usr/lib/x86_64-linux-gnu/atlas/libblas.so 35 manual mode 2 /usr/lib/x86_64-linux-gnu/blas/libblas.so 10 manual mode * 3 /usr/lib/x86_64-linux-gnu/blis-openmp/libblas.so 80 manual mode 4 /usr/lib/x86_64-linux-gnu/libmkl_rt.so 1 manual mode 5 /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so 100 manual mode Press <enter> to keep the current choice[*], or type selection number:と表示され、libblas.soの実体をsymbolic linkで切り替えることができます。libblas.so.3の方も同様に
update-alternatives --config libblas.so.3-x86_64-linux-gnu
There are 5 choices for the alternative libblas.so.3-x86_64-linux-gnu (providing /usr/lib/x86_64-linux-gnu/libblas.so.3). Selection Path Priority Status ------------------------------------------------------------ 0 /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 100 auto mode 1 /usr/lib/x86_64-linux-gnu/atlas/libblas.so.3 35 manual mode 2 /usr/lib/x86_64-linux-gnu/blas/libblas.so.3 10 manual mode * 3 /usr/lib/x86_64-linux-gnu/blis-openmp/libblas.so.3 80 manual mode 4 /usr/lib/x86_64-linux-gnu/libmkl_rt.so 1 manual mode 5 /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 100 manual mode Press <enter> to keep the current choice[*], or type selection number:と切り替えることができ、前回と同様、両方を同じものに切り替えて実験を行いました。
ベンチマーク
dgemmの実行時間計測に使ったプログラムは、以下のようなものです。前回の20.04のときと同じです。
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
// prototype declaration
void dgemm_(char *transA, char *transB, int *m, int *n, int *k, double *alpha, double *A, int *ldA, double *B, int *ldB, double *beta, double *C, int *ldC);
int main(int argc, char **argv)
{
int i, j;
int m, n, k;
int size;
double *a, *b, *c;
double alpha, beta;
int lda, ldb, ldc;
struct timespec ts1, ts2;
size = atoi(argv[1]);
m = size;
n = size;
k = size;
a = (double *)malloc(sizeof(double) * m * k); // m x k matrix
b = (double *)malloc(sizeof(double) * k * n); // k x n matrix
c = (double *)malloc(sizeof(double) * m * n); // m x n matrix
for (i=0; i<m; i++) {
for (j=0; j<k; j++) {
a[i + m * j] = rand() / (1.0 + RAND_MAX);
}
}
for (i=0; i<k; i++) {
for (j=0; j<n; j++) {
b[i + k * j] = rand() / (1.0 + RAND_MAX);
}
}
for (i=0; i<m; i++) {
for (j=0; j<n; j++) {
c[i + m * j] = 0;
}
}
alpha = 1.;
beta = 0.;
lda = m;
ldb = k;
ldc = m;
// dgemm_(TransA, TransB, M, N, K, alpha, A, LDA, B, LDB, beta, C, LDC)
// C = alpha * A * B + beta * C
// A=M*K, B=K*N, N=M*N
// Trans: "N"/"T"/"C"
// LDA = number of row of A
clock_gettime(CLOCK_REALTIME, &ts1);
dgemm_("N", "N", &m, &n, &k, &alpha, a, &lda, b, &ldb, &beta, c, &ldc);
clock_gettime(CLOCK_REALTIME, &ts2);
printf("%g\n", (ts2.tv_sec - ts1.tv_sec) + (ts2.tv_nsec - ts1.tv_nsec) / 1e9);
free(a);
free(b);
free(c);
return 0;
}
単なる比較のための、単純な三重forループによる行列積のプログラムは、以下の通り。これも前回と同じ。
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
double **alloc_matrix(int n, int m)
{
double **a;
int i;
a = (double **)malloc(sizeof(double *) * n);
a[0] = (double *)malloc(sizeof(double) * n * m);
for (i=1; i<n; i++) {
a[i] = a[0] + i * m;
}
return a;
}
void free_matrix(double **a)
{
free(a[0]);
free(a);
}
// a: n x m, b: m x s, c: n x s
void m_m_mul(double **a, double **b, double **c, int n, int m, int s)
{
int i, j, k;
int i1, j1, k1;
for (i=0; i<n; i++) {
for (j=0; j<s; j++) {
c[i][j] = 0.0;
}
}
for (i=0; i<n; i++) {
for (j=0; j<s; j++) {
for (k=0; k<m; k++) {
c[i][j] += a[i][k] * b[k][j];
}
}
}
}
int main(int argc, char **argv)
{
int i, j, n;
double **a;
double **b;
double **c;
struct timespec ts1, ts2;
n = atoi(argv[1]);
a = alloc_matrix(n,n);
b = alloc_matrix(n,n);
c = alloc_matrix(n,n);
for (i=0; i<n; i++) {
for (j=0; j<n; j++) {
a[i][j] = rand()/(1.0 + RAND_MAX);
b[i][j] = rand()/(1.0 + RAND_MAX);
}
}
clock_gettime(CLOCK_REALTIME, &ts1);
m_m_mul(a, b, c, n, n, n);
clock_gettime(CLOCK_REALTIME, &ts2);
printf("%g\n", (ts2.tv_sec - ts1.tv_sec) + (ts2.tv_nsec - ts1.tv_nsec) / 1e9);
free_matrix(a);
free_matrix(b);
free_matrix(c);
return 0;
}
これらは、cc -O3 dgemm.c -lblas cc -O3 forloop.cのようにコンパイルします。ただし、dgemm.cの方のコンパイルをするときに、update-alternativesの選択をMKL以外にする必要がありました。MKLの状態でコンパイルすると、libblas.soがリンクされずにlibmkl_rt.soがリンクされてしまって、update-alternativesによる事後切り替えが不可能なbinaryができてしまいました。これもまた自分の勘違いかも。
さて、ベンチマークの結果は以下のようになりました。使ったCPUはcore i7 1195G7 (ノートPC) で、前回と異なるので絶対値には意味がないかも。

前回と似たような結果ですね。最速はOpenBLAS, BLIS, MKL。小サイズでMKLが妙に遅いのも同じ。