下面的代码执行操作上gpuArrays相同的操作a和b在两种不同的方式.第一部分计算(a'*(a*b)')',第二部分计算a*b*a.然后验证结果是相同的.
%function test
clear
rng('default');rng(1);
a=sprand(3000,3000,0.1);
b=rand(3000,3000);
a=gpuArray(a);
b=gpuArray(b);
tic;
c1=gather(transpose(transpose(a)*transpose(a*b)));
disp(['time for (a''*(a*b)'')'': ' , num2str(toc),'s'])
clearvars -except c1
rng('default');
rng(1)
a=sprand(3000,3000,0.1);
b=rand(3000,3000);
a=gpuArray(a);
b=gpuArray(b);
tic;
c2=gather(a*b*a);
disp(['time for a*b*a: ' , num2str(toc),'s'])
disp(['error = ',num2str(max(max(abs(c1-c2))))])
%end
但是,计算(a'*(a*b)')'速度大约是计算速度的4倍a*b*a.以下是R2018a上Nvidia K20上面脚本的输出(我尝试过不同的版本和不同的GPU,具有相似的行为).
>> test
time for (a'*(a*b)')': 0.43234s
time for a*b*a: 1.7175s
error = 2.0009e-11
甚至更奇怪的是,如果上述脚本的第一行和最后一行是未注释的(它变成一个函数),则两个取较长的时间量(〜1.7S代替〜0.4秒).以下是此案例的输出:
>> test
time for (a'*(a*b)')': 1.717s
time for a*b*a: 1.7153s
error = 1.0914e-11
我想知道是什么导致了这种行为,以及如何在matlab函数内而不是在脚本内部的较短时间内(即~0.4s而不是~1.7s)执行 …
我注意到,将数据传输到最近的高端GPU比将其收集回CPU更快.以下是使用由旧版Nvidia K20和最近使用PCIE的Nvidia P100运行的mathworks技术支持提供给我的基准测试功能的结果:
Using a Tesla P100-PCIE-12GB GPU.
Achieved peak send speed of 11.042 GB/s
Achieved peak gather speed of 4.20609 GB/s
Using a Tesla K20m GPU.
Achieved peak send speed of 2.5269 GB/s
Achieved peak gather speed of 2.52399 GB/s
我已经在下面附上了基准功能以供参考.P100不对称的原因是什么？这个系统是依赖还是近期高端GPU的标准？可以提高聚集速度吗？
gpu = gpuDevice();
fprintf('Using a %s GPU.\n', gpu.Name)
sizeOfDouble = 8; % Each double-precision number needs 8 bytes of storage
sizes = power(2, 14:28);
sendTimes = inf(size(sizes));
gatherTimes = inf(size(sizes));
for ii=1:numel(sizes)
numElements = sizes(ii)/sizeOfDouble;
hostData = randi([0 …我很好奇为什么将稀疏矩阵乘以稠密矩阵所需的时间与相反的时间不同。算法有明显不同吗？
这是 matlab 2018a 中的一个示例：
a=sprand(M,M,0.01);
b=rand(M);
tic;ref1=a*b;t_axb=toc
tic;ref2=b*a;t_bxa=toc
这是使用 1 个线程的 Eigen 3 和 C++ 的示例：
//prepare acol=MxM ColMajor Eigen sparse matrix with 0.01 density
...
Map<Matrix<double,M,M,ColMajor> > bcol (PR, M, M );
double tic,toc;
tic=getHighResolutionTime();
result=acol*bcol;
toc=getHighResolutionTime();
printf("\nacol*bcol time: %f seconds", (toc - tic));
tic=getHighResolutionTime();
result=bcol*acol;
toc=getHighResolutionTime();
printf("\nbcol*acol time: %f seconds\n", (toc - tic));
当 M=4000 时，结果为：
t_axb =
0.6877
t_bxa =
0.4803
acol*bcol time: 0.937590 seconds
bcol*acol time: 0.532622 seconds
当 M=10000 时，结果为
t_axb =
11.5649
t_bxa …在我的MEX文件的最后一行完成执行后,返回matlab命令行大约需要大约14秒.
从matlab定时MEX文件:
D=rand(14000)+rand(14000)*1i;
tic;
[A B C]=myMexFile(D);
toc
disp(datetime('now'));
输出是:
Elapsed time is 35.192704 seconds.
15-Sep-2018 16:51:35
使用以下最小工作示例从C内对MEX文件进行计时:
#include <mex.h>
#include <sys/time.h>
#include <time.h>
#include <cuComplex.h>
double getHighResolutionTime() {
struct timeval tod;
gettimeofday(&tod, NULL);
double time_seconds = (double) tod.tv_sec + ((double) tod.tv_usec / 1000000.0);
return time_seconds;
}
void double2cuDoubleComplex(cuDoubleComplex* p, double* pr, double* pi,int numElements){
for(int j=0;j<numElements;j++){
p[j].x=pr[j];
p[j].y=pi[j];
}
}
void cuDoubleComplex2double(cuDoubleComplex* p, double* pr, double* pi,int numElements){
for(int j=0;j<numElements;j++){
pr[j]= p[j].x;
pi[j]= p[j].y;
}
}
void mexFunction( int …该函数mkl_malloc类似于malloc但有一个额外的alignment参数。这是原型：
void* mkl_malloc (size_t alloc_size, int alignment);
我注意到具有不同值的不同性能alignment。除了反复试验之外，是否有一种规范的或记录在案的有条理的方法来决定 的最佳值alignment？即正在使用的处理器、正在调用的函数、正在执行的操作等。
这个问题广泛适用于任何使用 MKL 的人，所以我很惊讶它不在参考手册中。
更新：我已经尝试过mkl_sparse_spmm并且没有注意到将对齐设置为 2 的幂（最多 1024 字节）的性能有显着差异，之后性能往往会下降。我使用的是英特尔至强 E5-2683。
c=complex(a,b)在matlab中比做起来慢得多c=a+1i*b.
以下是Matlab 2018a
a=rand(15000);
b=rand(15000);
%
clear c;
tic; c=a+1i*b; toc
Elapsed time is 0.338525 seconds
%
clear c;
tic; c=complex(a,b); toc
Elapsed time is 2.542403 seconds.
是complex在任何情况下实际上有用吗？为什么这么慢？
以下是使用Eigen将密集数组g和G与matlab相乘的mex代码.当g稀疏时,我该怎么做？
#include <iostream>
#include <Eigen/Dense>
#include "mex.h"
using Eigen::MatrixXd;
using namespace Eigen;
/*gateway function*/
void mexFunction( int nlhs, mxArray *plhs[],
int nrhs, const mxArray *prhs[]) {
int nRows=(int)mxGetM(prhs[0]);
int nCols=nRows;
double* g=mxGetPr(prhs[0]);
double* Gr=mxGetPr(prhs[1]);
Map<MatrixXd> gmap (g, nRows, nCols );
Map<MatrixXd> Grmap (Gr, nRows, nCols );
plhs[0] = mxCreateDoubleMatrix(nRows, nCols, mxREAL);
Map<MatrixXd> resultmap (mxGetPr(plhs[0]), nRows, nCols);
resultmap = gmap*Grmap;
}
我目前正在将矩阵的实部和虚部分别从 Matlab 导入到 C++。然后，我还将实部和虚部分别映射到特征值。我也分别进行计算并绘制最终结果，如下图：
//import real and imaginary parts from matlab
mwSize M = mxGetM (prhs[1]);
mwSize N = mxGetN (prhs[1]);
double * PR = mxGetPr (prhs[1]);
double * PI = mxGetPi (prhs[1]);
//map real and imaginary parts to Eigen
Map<Matrix<double,Dynamic,Dynamic,ColMajor> > Br (PR, M, N );
Map<Matrix<double,Dynamic,Dynamic,ColMajor> > Bi (PI, M, N );
//map real and imaginary parts of result
plhs[0] = mxCreateDoubleMatrix(M, N, mxCOMPLEX);
Map<Matrix<double,Dynamic,Dynamic,ColMajor> > resultr (mxGetPr(plhs[0]), M, N);
Map<Matrix<double,Dynamic,Dynamic,ColMajor> > resulti (mxGetPi(plhs[0]), M, N);
//calculate real …