Matlab并行化计算及GPU计算教程
前置要求和设置
要求电脑CPU有超过2个核心,内存大于2G。
建议先调试好代码,再进行并行化计算。
查看并行化计算工具箱版本
>>> ver(parallel)
------------------------------------------------------------------------------------------------
MATLAB 版本: 9.12.0.1884302 (R2022a)
MATLAB 许可证编号: 40523737
操作系统: Microsoft Windows 11 专业版 Version 10.0 (Build 22000)
Java 版本: Java 1.8.0_202-b08 with Oracle Corporation Java HotSpot(TM) 64-Bit Server VM mixed mode
------------------------------------------------------------------------------------------------
Parallel Computing Toolbox 版本 7.6 (R2022a)
确保并行化工具箱正确配置
菜单栏 >> Parallel >> Manage Clusters >> local >> Validation >> 全部勾选,点击Validate,绿色对勾为正确配置
验证GPU是否可以使用
>>> gpuDevice
ans =
CUDADevice - 属性:
Name: 'NVIDIA GeForce RTX 3080'
Index: 1
ComputeCapability: '8.6' % 计算能力至少1.3以上,建议2.0以上
SupportsDouble: 1
DriverVersion: 11.6000
ToolkitVersion: 11.2000
MaxThreadsPerBlock: 1024
MaxShmemPerBlock: 49152
MaxThreadBlockSize: [1024 1024 64]
MaxGridSize: [2.1475e+09 65535 65535]
SIMDWidth: 32
TotalMemory: 1.2884e+10
AvailableMemory: 1.1487e+10
MultiprocessorCount: 70
ClockRateKHz: 1785000
ComputeMode: 'Default'
GPUOverlapsTransfers: 1
KernelExecutionTimeout: 1
CanMapHostMemory: 1
DeviceSupported: 1
DeviceAvailable: 1
DeviceSelected: 1
快速上手parfor
将for循环,改为parfor循环。并行化计算会自动进行。
parfor ...
statement
end
也可以使用左下角的符号,选择Start parallel pool,来手动控制并行化计算池的启动与关闭。
也可以使用命令行
parpool('local',N) % 启动并行化,local为预设的方案,N为线程数量
delete(gcp('nocreate')) % 关闭并行化
使用parfor的注意事项
- 每次循环的计算过程是独立的,当前循环不依赖历史循环的信息;
- 循环的顺序是没有限制的;
- parfor循环的变量必须是整数,如果是非整数,则需要在循环外部转换为整数。
% 循环变量不能是非整数
parfor x = 0:0.1:1
pass
end
% 循环变量,要修改为整数
xALL = 0:0.1:1;
parfor ii = 1:length(xAll)
x = xALL(ii)
pass
end
- parfor循环,不允许进行嵌套,但是可以通过函数的形式进行调用
% 不能进行循环嵌套
parfor x = 1:5
parfor y = 1:10
pass
end
end
% 利用函数调用,进行parfor嵌套
parfor x = 1:5
myprog(x)
end
function myprog(x)
parfor y = 1:10
pass
end
- 如果代码分析器,无法得出代码哪里出问题了,可以将parfor循环体,改为函数调用的形式
% 代码分析器,无法获取代码不能并行化的原因
function ex2_bad
data = rand(4,4);
means = zeros(1,4);
parfor I = 1:4
means(I) = z.mean;
end
disp(means)
% 尝试使用函数调用的方式,替换parfor循环体
function ex2_bad
data = rand(4,4);
means = zeros(1,4);
parfor I = 1:4
means(I) = computeMeans(data, I);
end
disp(means)
function means = computeMeans(data,I)
z.means = mean(data(:,I));
means = z.mean;
- 如果parfor加速效果不明显,可能是以下几个原因导致的:计算时间太短;内存或者文件读写速度限制;已经存在并行化计算的其他代码(可以使用任务管理器查看);个别迭代耗时远大于其他迭代。
批处理
交互式的并行化池,需要matlab客户端一直打开。
对于批处理,在使用集群计算时,可以关闭matlab客户端,处理其他工作。适合集群计算、耗时长的计算。
批处理的三个步骤:
- 提交任务
% 对于脚本的使用方法
job = batch("ascript")
% 对于函数的使用方法,(函数名称,返回输出的数量,输入参数)
job = batch(@ex_serial, 1, (50,4000));
% 使用集群计算,可以关闭matlab客户端,甚至关闭电脑
clust = parcluster('local');
% N = 4; % 计算机CPU核心数
N = clust.NumIdleWorkers % 可以调用所有可用核心
job = batch(clust, @ex_parallel, 1, (50,4000), 'pool', N-1)
% N-1,是因为需要一个线程负责任务的分发以及结果的收集,即存在一个线程对其他线程进行管理
- 查看任务是否完成
wait(job, 'finished') % 当以自动方式运行多个任务时,这个命令可以屏蔽其他任务,直到当前任务完成
get(job, 'state') % 获取当前任务的进度
% 也可以使用菜单栏的Job Monitor,以交互式的方法检查任务进度
- 获取任务结果
results = fetchOutputs(job) % 获取任务结果,以cell array的形式返回
outk = results{k} % 获取任务结果的第k个输出
delete(job) % 任务完成后,删除任务
对于不使用for循环的代码,可以参考parallel_template.m
扩展至集群和云
大型项目需要更多的计算能力和内存,此时可以考虑使用公司、组织、云端的集群进行计算,只需要对代码进行很小的修改即可实现。
local profile,并行化本地配置方案,可以使用单台计算机的多核心进行加速,也可以用于在使用其他并行化方案前,测试并行化代码。
-
添加新的并行化配置方案,home >> parallel >> manage cluster profiles >> add profile / discover clusters / import profile
-
matlab job scheduler (MJS) profile,可以利用安装了matlab distributed computing server(MDCS)上的计算资源,MDCS的使用方法这里不做叙述。
第三方的cluster schedulers,在现有的cluster scheduler中集成了MDCS,如果没有你想要的scheduler,选择generic即可。
-
检查当前并行化配置方案是否可用,home >> parallel >> manage cluster profiles,选择要测试的profile,点击validate,全部为绿色表示通过。
-
使用方法,简单的将’local’替换为你的并行化方案即可。
-
如果报错,无法找到文件,可将数据发送到集群。
% 对于并行化计算池
parpool(...,'AttachedFiles',{'C:\mydir'}) % 对于新的并行化计算池
addAttachedFiles(gcp, {'C:\mydir'}) % 更新现有的并行化计算池
% 对于批处理
batch(...,'AttachedFiles',{'C:\mydir'})
- 如果加速没有达到想要的效果,可能是文件传输占用了很多时间,要避免频繁的将数据复制到集群。
比parfor更高级的方法:spmd
spmd(single program multiple data),spmd代码块中的代码,在所有的线程运行,可以用于parfor循环之前加载数据。
spmd
% myfile.mat,在所有线程都可以访问
data = load('myfile.mat')
end
parfor I = 1:N
% 使用数据的代码块
end
也可以控制每个线程做什么工作。
% labindex返回worker的索引
% numlabs返回worker的总数
spmd
switch labindex
case 1
% 线程1的代码
case 2
% 线程2的代码
...
end
end
也可以控制线程之间的数据交互
% labSend(dataToSend,workerIdx)可以将需要发送的数据,发送至指定索引的线程
% labReceive(workerIdx)可以接受指定线程的数据
% labSend和labReceive需要成对使用
spmd
...
if labindex < numlabs
% 除了最后一个线程之外的线程,将数据发送至下一个线程
labSend(dataToSend, labindex+1);
else
% 最后一个线程数据发送至第一个线程
labSend(dataToSend, 1);
end
if labindex > 1
% 除了最后一个线程之外的线程,接受上一个线程的数据
dataReceived = labReceive(labindex-1);
else
% 第一个线程接受最后一个线程的数据
dataReceived = labReveice(numlabs);
end
...
end
% 所有线程的数据可以组合输出,通过索引可以访问不同线程的数据
dataOnClient = dataReceived(:);
分布式矩阵
可以划分到各个线程进行计算,仍可以被视为一个矩阵,这种矩阵称为分布式矩阵。
主要被用于集群计算,可以利用多个机器的内存。
- 创建分布式矩阵的方法
% 方法一
zeros(...,'distributed')
randn(...,'distributed')
% 如果所有线程的数据划分方式一样,可以直接使用distributed,但是只能按照列进行划分
a = reshape(1:16, 4, 4)
a = distributed(a)
% 方法二
% 从多个文件,或者是一个巨大到无法一次载入内存的文件,创建分布式矩阵
% 可以在每个线程载入数据,然后使用codistributed.build及codistributorld将每个线程的矩阵结合起来
% 需要指定分布式矩阵的size,(可选)指定各个线程的数据
spmd
N = 1;
mat = repmat([1;2;3], 1, N) + (labindex-1)*3;
codistr = codistributorld(codistributorld.unsetDimension,codistributorld.unsetPartition,[3,numlabs*N]);
distmat = codistributed.build(mat, codistr)
disp(getLocalPart(distmat))
end
% 使用gather可以将分布式矩阵,转为数值化矩阵
distsum = sum(distmat);
numericsum = gather(distsum);
GPU计算
-
将数据迁移至GPU显存,或者在GPU显存创建数据
% 使用gpuArray将数据迁移至GPU a = reshape(1:16, 4, 4); a = gpuArray(a); % 直接在GPU创建数据,这种方法更高效 z = zeros(4, 4, 'gpuArray') % 如需查看哪些操作支持gpuArray >> methods gpuArray
-
使用GPU代码处理数据
-
将数据迁移回CPU
% 使用CPU对矩阵进行滤波操作
A = magic(5000);
f = ones(1,20)/20;
tic;
B = filter(f, 1, A);
tCPU = toc;
disp(['Total time on CPU: ' num2str(tCPU)])
% 使用GPU对矩阵进行滤波操作
tic;
AonGPU = gpuArray(A);
BonGPU = filter(f, 1, AonGPU);
BonCPU = gather(BonGPU);
wait(gpuDevice)
tGPU = toc;
disp(['Total time on GPU: ' num2str(tGPU)])
% 只查看计算耗时,不考虑数据传输耗时
tic;
BonGPU = filter(f, 1, AonGPU);
wait(gpuDevice)
tCompGPU = toc;
disp(['Computation time on GPU: ' num2str(tCompGPU)])
% 使用arrayfun或者pagefun进一步加速
M = 300; % 行数
K = 500; % 矩阵内部维度
N = 100; % 列数
P = 200; % 页数
% CPU模式
tic
A = rand(M,K);
B = rand(K,N,P);
C = zeros(M,N,P);
for I = 1:P
C(:,:,I) = A*B(:,:,I);
end
t = toc;
disp(['CPU time: ' num2str(t)])
% GPU模式,gpuArray
tic
A = rand(M,K,'gpuArray');
B = rand(K,N,P,'gpuArray');
C = zeros(M,N,P,'gpuArray');
for I = 1:P
C(:,:,I) = A*B(:,:,I);
end
wait(gpuDevice)
t = toc;
disp(['GPU using gpuAarrys time: ' num2str(t)])
% GPU模式,pagefun
tic
A = rand(M,K,'gpuArray');
B = rand(K,N,P,'gpuArray');
C = pagefun(@mtimes,A,B);
wait(gpuDevice)
t = toc;
disp(['GPU using pagefun time: ' num2str(t)])