使用 pagefun
提高 GPU 上小矩阵问题的性能
此示例说明如何使用 pagefun
来提高对多维数组中排列的多个矩阵施加独立运算的性能。
多维数组是二维矩阵的扩展,使用额外的下标进行索引。例如,三维数组使用三个下标。前两个维度表示矩阵,第三个维度表示元素页面(有时称为切片)。有关详细信息,请参阅多维数组。
虽然 GPU 可以有效地将小型独立操作应用于大型矩阵,但当这些操作以串行方式应用时(例如,当在 for
循环中应用这些操作时),性能并不理想。为了避免串行处理,arrayfun
函数在 GPU 上并行对数组的每个元素应用标量运算。类似地,pagefun
函数将一个函数应用于多维 GPU 数组的每个页面。
pagefun
函数支持应用大多数元素函数和许多支持 GPU 数组输入的矩阵运算。MATLAB® 还提供了许多专用的页面函数,包括 pagemtimes
、pagemldivide
、pagemrdivide
、pagetranspose
、pagectranspose
、pageinv
、pagenorm
和 pagesvd
。根据任务的不同,这些函数可能会简化您的代码或提供比使用 pagefun
更好的性能。
在这个例子中,机器人正在导航一个已知地图,该地图包含大量机器人可以使用其传感器识别的特征。机器人通过测量这些特征的相对位置和方向并将它们与地图位置进行比较来确定地图上的定位。假设机器人没有完全迷失,它可以利用两者之间的任何差异来纠正其位置,例如通过使用卡尔曼滤波器。这个例子展示了一种计算相对于机器人的特征位置的有效方法。
设置地图
定义包含许多特征的房间的尺寸。
roomDimensions = [50 50 5];
本示例末尾提供了支持函数 randomTransforms
,并用随机值初始化 N
变换,提供一个结构作为输出。它使用 3×1 向量 T
和 3×3 旋转矩阵 R
来表示位置和方向。N
平移被打包成一个 3×N
矩阵,而旋转被打包成一个 3×3×N
数组。
使用 randomTransforms
函数设置 1000
特征的地图以及机器人的起始位置。
numFeatures = 1000; Map = randomTransforms(numFeatures,roomDimensions); Robot = randomTransforms(1,roomDimensions);
plotRobot
函数作为本示例的支持文件提供,绘制了房间的俯视图以及机器人和附近特征的近距离视图。机器人用一个带轮子的蓝色盒子表示,特征用红色圆圈表示,伴随的线条表示它们的方向。要使用此函数,请将示例作为 livescript 打开。
调用 plotRobot
函数。
plotRobot(Robot,Map)
定义方程
为了正确识别地图中的特征,机器人需要转换地图以将其传感器置于原点。然后,它可以通过比较所看到的内容与期望看到的内容来找到地图特征。
对于地图特征 ,我们可以通过转换其全局地图位置来找到它相对于机器人 的位置和方向 :
其中 和 是机器人的位置和方向, 和 表示地图数据。等效的 MATLAB 代码如下所示:
Rrel(:,:,i) = Rbot' * Rmap(:,:,i) Trel(:,i) = Rbot' * (Tmap(:,i) - Tbot)
使用 for
循环在 CPU 上执行矩阵变换
本示例末尾提供了支持函数 loopingTransform
,依次循环遍历所有变换,将每个特征变换到其相对于机器人的位置。注意 zeros
函数的 like
名称-值参量,它使函数返回与原型数组相同数据类型的零数组。例如,如果原型数组是 gpuArray
,那么 zeros
将返回 gpuArray
。这使得您可以在下一节中在 GPU 上使用相同的代码。
使用timeit
函数对计算进行计时。timeit
函数对 loopingTransform
的执行进行多次计时,并返回测量值的中值。由于 timeit
需要一个没有参量的函数,因此使用 @()
语法来创建正确形式的匿名函数。
cpuTime = timeit(@()loopingTransform(Robot,Map,numFeatures))
cpuTime = 0.0042
使用 for
循环在 GPU 上执行矩阵变换
要在 GPU 上运行相同的代码,只需将输入数据作为 gpuArray
传递给函数。gpuArray
代表存储在 GPU 内存中的数组。MATLAB 和其他工具箱中的许多函数都支持 gpuArray
对象,允许您以最少的代码更改在 GPU 上运行代码。有关详细信息,请参阅在 GPU 上运行 MATLAB 函数。
确保您想要的 GPU 可用且已被选择。
gpu = gpuDevice;
disp(gpu.Name + " GPU selected.")
NVIDIA RTX A5000 GPU selected.
创建包含机器人的位置和方向以及地图中的特征的 GPU 数组。
gMap.R = gpuArray(Map.R); gMap.T = gpuArray(Map.T); gRobot.R = gpuArray(Robot.R); gRobot.T = gpuArray(Robot.T);
使用gputimeit
函数对计算进行计时。对于包含 GPU 计算的代码,gputimeit
函数相当于 timeit
。它确保在记录时间之前所有 GPU 操作都已完成。
gpuTime = gputimeit(@()loopingTransform(gRobot,gMap,numFeatures))
gpuTime = 0.1905
使用 pagefun
在 GPU 上执行矩阵变换
GPU 版本非常慢,因为尽管所有计算都是独立的,但它们都是串行运行的。使用 pagefun
我们可以并行运行所有计算。
本示例末尾提供了支持函数 pagefunTransform
,它应用与 loopingTransform
函数相同的变换,但使用 pagefun
而不是 for
循环。第一个计算是旋转的计算。这涉及矩阵乘法,转换为函数 mtimes
(*
)。将其与要相乘的两组旋转一起传递给 pagefun
:
Rel.R = pagefun(@mtimes,Robot.R',Map.R);
Robot.R'
是一个 3×3 矩阵,而 Map.R
是一个 3×3×N 数组。pagefun
函数将地图中的每个独立矩阵与相同的机器人旋转进行匹配,并为我们提供所需的 3×3×N 输出。
平移计算也涉及矩阵乘法,但矩阵乘法的正常规则允许它在循环之外进行而不进行任何更改:
Rel.T = Robot.R' * (Map.T - Robot.T);
使用gputimeit
函数对计算进行计时。
gpuPagefunTime = gputimeit(@()pagefunTransform(gRobot,gMap))
gpuPagefunTime = 3.3066e-04
比较结果
绘制时序结果。
figure labels = categorical(["CPU Execution","GPU Execution","GPU Execution with \fontname{consolas}pagefun"]); bar(labels,[cpuTime,gpuTime,gpuPagefunTime]) ylabel("Execution Time (s)") set(gca,YScale="log")
计算使用 pagefun
的执行比 CPU 和简单的 GPU 执行快多少。
fprintf("Executing the transforms on the GPU using pagefun is %3.2f times faster than on the CPU.\n", ... cpuTime/gpuPagefunTime);
Executing the transforms on the GPU using pagefun is 12.65 times faster than on the CPU.
fprintf("Executing the transforms on the GPU using pagefun is %3.2f times faster than using for-loops on the GPU.\n", ... gpuTime/gpuPagefunTime);
Executing the transforms on the GPU using pagefun is 576.17 times faster than using for-loops on the GPU.
使用多个可能的机器人位置来定位丢失的机器人
如果机器人位于地图的未知部分,它可以使用全局搜索算法来定位自己。该算法通过执行上述计算并寻找机器人传感器看到的特征与机器人期望在该位置看到的特征之间的良好对应关系来测试许多可能的位置。
现在机器人有多个可能的位置以及多个特征。N 个特征和 M 个机器人需要 N*M
转换。为了区分“机器人空间”和“特征空间”,使用第 4 维进行旋转,使用第 3 维进行平移。这意味着机器人的旋转将是 3×3×1×M,平移将是 3×1×M。
使用十个随机机器人位置初始化搜索。好的搜索算法会使用拓扑或其他线索来更智能地进行搜索。
numRobots = 10; Robot = randomTransforms(numRobots,roomDimensions); Robot.R = reshape(Robot.R,3,3,1,[]); % Spread along the 4th dimension Robot.T = reshape(Robot.T,3,1,[]); % Spread along the 3rd dimension
本例末尾定义了一个支持函数 loopingTransform2
,它使用两个嵌套循环执行循环转换,以循环遍历机器人和特征。
使用 timeit
来计算时间。
cpuTime = timeit(@()loopingTransform2(Robot,Map,numFeatures,numRobots))
cpuTime = 0.0759
创建包含机器人旋转和平移的 GPU 数组。
gRobot.R = gpuArray(Robot.R); gRobot.T = gpuArray(Robot.T);
使用 gputimeit
对 GPU 上的计算进行计时。
gpuTime = gputimeit(@() loopingTransform2(gRobot,gMap,numFeatures,numRobots))
gpuTime = 2.1059
与以前一样,循环版本在 GPU 上的运行速度要慢得多,因为它不并行进行计算。
本示例末尾提供了一个支持函数 pagefunTransform2
,它使用两个 pagefun
调用而不是嵌套的 for
循环应用与 loopingTransform2
函数相同的转换。此函数需要将 transpose
运算符以及 mtimes
合并到对 pagefun
的调用中。该函数还将 squeeze
函数应用于转置的机器人方向,以将机器人的分布放入第 3 维,从而匹配平移。尽管如此,最终的代码还是更加紧凑。
pagefun
函数会适当扩展维度,因此当我们将 3×3×1×M 矩阵 Rt
与 3×3×N×1 矩阵 Map.R
相乘时,我们得到一个 3×3×N×M 矩阵。
使用 gputimeit
对 GPU 上的计算进行计时。
gpuPagefunTime = gputimeit(@()pagefunTransform2(gRobot,gMap))
gpuPagefunTime = 0.0014
比较结果
绘制时序结果。
labels = categorical(["CPU Execution","GPU Execution","GPU Execution with \fontname{consolas}pagefun"]); bar(labels,[cpuTime,gpuTime,gpuPagefunTime]) ylabel("Execution Time (s)") set(gca,YScale="log")
fprintf("Executing the transforms on the GPU using pagefun is %3.2f times faster than on the CPU.\n", ... cpuTime/gpuPagefunTime);
Executing the transforms on the GPU using pagefun is 55.85 times faster than on the CPU.
fprintf("Executing the transforms on the GPU using pagefun is %3.2f times faster than using nested for-loops on the GPU.\n", ... gpuTime/gpuPagefunTime);
Executing the transforms on the GPU using pagefun is 1549.19 times faster than using nested for-loops on the GPU.
结论
pagefun
函数支持许多二维运算,以及 arrayfun
支持的大多数标量运算。通过这些函数,您可以对涉及矩阵代数和数组操作的一系列计算进行向量化,从而无需循环并可获得巨大的性能提升。
无论您在循环中对 GPU 数据进行小规模计算,都应该考虑以这种方式转换为向量化实现。这也是一个利用 GPU 来提高性能的机会,而之前 GPU 并没有带来任何性能提升。
支持函数
随机变换函数
randomTransforms
函数在指定尺寸的房间中创建定义 N
随机变换的矩阵。每个变换包括一个随机平移 T
和一个随机旋转 R
。该函数可用于设置房间特征地图以及机器人的起始位置和方向。
function Tform = randomTransforms(N,roomDimensions) % Preallocate matrices. Tform.T = zeros(3,N); Tform.R = zeros(3,3,N); for i = 1:N % Create random translation. Tform.T(:,i) = rand(3,1) .* roomDimensions'; % Create random rotation by extracting an orthonormal % basis from a random 3-by-3 matrix. Tform.R(:,:,i) = orth(rand(3,3)); end end
循环变换函数
loopingTransform
函数通过依次循环变换将每个特征转换为其相对于机器人的位置。
function Rel = loopingTransform(Robot,Map,numFeatures) % Preallocate matrices. Rel.R = zeros(size(Map.R),like=Map.R); Rel.T = zeros(size(Map.T),like=Map.T); for i = 1:numFeatures % Find orientation of map feature relative to the robot. Rel.R(:,:,i) = Robot.R' * Map.R(:,:,i); % Find position of map feature relative to the robot. Rel.T(:,i) = Robot.R' * (Map.T(:,i) - Robot.T); end end
pagefun
变换函数
pagefunTransform
函数通过使用 pagefun
函数应用变换将每个特征变换到其相对于机器人的位置。
function Rel = pagefunTransform(Robot,Map) % Find orientation of map feature relative to the robot. Rel.R = pagefun(@mtimes,Robot.R', Map.R); % Apply translation. Rel.T = Robot.R' * (Map.T - Robot.T); end
嵌套循环转换函数
loopingTransform2
函数使用两个嵌套循环执行循环转换,循环遍历机器人以及特征。变换将每个特征映射到相对于每个机器人的位置。
function Rel = loopingTransform2(Robot,Map,numFeatures,numRobots) % Preallocate matrices. Rel.R = zeros(3,3,numFeatures,numRobots,like=Map.R); Rel.T = zeros(3,numFeatures,numRobots,like=Map.T); for i = 1:numFeatures for j = 1:numRobots % Find orientation of map feature relative to the robot. Rel.R(:,:,i,j) = Robot.R(:,:,1,j)' * Map.R(:,:,i); % Find position of map feature relative to the robot. Rel.T(:,i,j) = ... Robot.R(:,:,1,j)' * (Map.T(:,i) - Robot.T(:,1,j)); end end end
两次调用 pagefun
转换函数
pagefunTransform2
函数通过两次调用 pagefun
函数执行转换,将每个特征映射到相对于每个机器人的位置。
function Rel = pagefunTransform2(Robot,Map) % Find orientation of map feature relative to the robot. Rt = pagefun(@transpose,Robot.R); Rel.R = pagefun(@mtimes,Rt,Map.R); % Find position of map feature relative to the robot. Rel.T = pagefun(@mtimes,squeeze(Rt), ... (Map.T - Robot.T)); end
另请参阅
pagefun
| gpuArray
| arrayfun
| gputimeit