高性能矩阵乘法

矩阵乘法

xiafei@

2026/5/61并行算法优化研究相对于传统面向对象串行算法的4个挑战：同步：两个或者多个线程协调其行为的过程通信：与线程之间交换数据相关的带宽和延迟问题负载均衡：多个线程之间工作量分布的情况，给各个线程（执行核）分配均匀的工作可扩展性：衡量在性能更加强劲的系统上运行软件时能否有效利用更多线程的指标， 观察应用程序在更高级的平台上运行 4核到8核线性增长2026/5/62多线程（核）设计主要分解模式任务分解： 对程序根据其执行的功能进行分解的过程数据分解： 将应用程序根据各任务所处理的数据而非按任务的天然特性来进行分解数据流分解： 研究数据在诸任务之间如何流动，根据任务之间的数据流关系对问题 进行分解模式分解方式任务级并行模式任务分解DivideandConquer任务/数据分解几何分解模式数据分解流水线模式数据流分解波峰（wavefront）模式数据流分解2026/5/63多线程（核）设计主要分解模式任务分解： 对程序根据其执行的功能进行分解的过程数据分解： 将应用程序根据各任务所处理的数据而非按任务的天然特性来进行分解数据流分解： 研究数据在诸任务之间如何流动，根据任务之间的数据流关系对问题 进行分解

分解方式设计说明任务分解不同的程序行为采用不同的线程实现常用于GUI应用程序数据分解多个线程对不同的数据块执行相同的操作常用于音频、图像处理和科学计算应用程序数据流分解一个线程的输出作为另一个线程的输入尤其应注意尽量消除启动和排空延迟2026/5/64矩阵乘法算法探讨

在工程科学计算中，矩阵乘积是最基本的运算 典型的n阶稠密方阵乘积算法的时间复杂度是O(n3)。 目前对大型矩阵乘积运算的处理主要是采用分治思想，将矩阵分布在多个节点上，但每个结点上的小矩阵仍要立方级乘法次数。 基于分之思想的两种划分策略：条形划分和块状（棋盘）划分的6种常见分布式矩阵乘法并行算法。

2026/5/65基于不同划分策略的矩阵乘法算法探讨

1、条形（stripedpartitioning）划分的矩阵乘法并行算法

行条划分 列条划分两两组合：行列、行行、列列、列行2026/5/66基于不同划分策略的矩阵乘法算法探讨

2、块状划分(checkerboardpartitioning)的矩阵乘法并行算法

称为棋盘划分2026/5/67Cannon

DescriptionforimplementationofMPIprogramtocomputeMatrixMatrixMultiplicationusingblockcheckerboardpartitioning

andCannonAlgorithm

2026/5/68Cannon

Objective

Computingthematrix-matrixmultiplicationonSMPSystem.UseblockcheckerboardpartitioningofthematricesandCannon'sAlgorithm.

AssumptionSizeofthesquarematricesp=q2andthesizeofsquarematricesAandBisevenlydivisiblebyq.

Itisassumedthatthenumberofblocksareequaltothenumberofprocessors.2026/5/69Cannon

Cannon'salgorithmisbasedoncartesianvirtualtopologyAandBaresquarematricesofsizenandCbetheoutput

matrix.Thesematricesaredivedintoblocksorsubmatricestoperformmatrix-matrixoperationsinparallelnxnmatrixAcanberegardedasqxqarrayofblocksAi,j(0<=i<q,0<=j<q)suchthateachblockisan(n/q)x(n/q)submatrixWeusep

processorstoimplementtheblockversionofmatrixmultiplicationinparallelbychoosingqasasquarerootof

p

andcomputeadistinctblockCi,joneachprocessor.2026/5/610传统并行

2026/5/611传统并行

ThematricesAandBare

partitionedintop

blocks,A

i,j

andBi,j

(0<=i<q,0<=j<q)ofsize

(n/qxn/q)

oneachprocess.Theseblocksaremappedontoaqxq

logicalmeshofprocesses.TheprocessesarelabeledfromP0,0toPq-1,q-1.2026/5/612传统并行

ProcessPi,jinitiallystoreblockmatricesAi,jandBi,jandcomputesblockCi,jofresultmatrix.TocomputesubmatrixCi,j,weneedallsubmatrices,Ai,kand

Bk,j(0

<=k<q).Toacquirealltherequiredblocks,anall-to-allbroadcastofmatrixAi,j's

isperformedineachrowandsimilarlyineachcolumnofmatrixBi,j's.MPIcollectivecommunicationisusedtoperformthisoperations.2026/5/613传统并行

AfterPi,jacquires,Ai,0,Ai,1,Ai,2

,Ai,q-1andB0,j

,B1,j

,B2,j

,Bq-1,j

,itperformstheserialblockmatrixtomatrixmultiplicationandaccumulatesthepartialblock

matrixCi,jofmatrixC.ToobtaintheresultantproductmatrixC,processeswithrank0gathersalltheblockmatricesbyusingMPI_Gather

collectivecommunicationoperation.2026/5/614Cannon

pprocessorsarrangedinqxqsquaregridofprocessorsandtheinputmatrices.AandBaredistributedamongtheprocessesincheckerboardfashion.ItresultsinconstructingpblockmatricesofAandB.Itusesonlypoint-to-pointcommunication

forcircularlyshiftingblocksofmatrixAandmatrixBamongpprocesses.

2026/5/615Cannon-inital

Thealgorithmproceedsinq

stages.ThefirststepinthisalgorithmistoperforminitialalignmentoftheblockmatrixAandblockmatrixB.TheblocksofmatrixAarecircularlyshiftedtothe

i

positionstoleftintherowofthesquaregridofprocesses,wherei

istherow

numberoftheprocessinthemesh.TheblocksofmatrixBarecircularlyshiftedj

positionsupwards,wherejisthecolumn

numberoftheprocessintheprocessesmesh.2026/5/616Cannon-inital2026/5/617Cannon-runningThealgorithmperformsthefollowingstepsineachstage:

1.MultiplytheblockofmatrixAandmatrixBandaddtheresultantmatrixtogettheblockmatrixC,whichisinitiallysettozero.

2.CircularlyshifttheblocksofmatrixAtoleftintherowsoftheprocessesandtheblocksofmatrixBupwardsinthecolumnsofthesquaregridofprocessesinawraparoundmanner.

2026/5/618Cannon-running2026/5/619Cannon-running2026/5/620书中Cannon-bug

MPI_SendandMPI_Recvisnotusedforpoint-to-pointcommunicationbecauseifalltheprocessescallMPI_SendorMPI_Recv

indifferentorderthedeadlockedsituationmayarise.

Howtofix?指派一个缓冲区，使用MPI_Irecv/MPI_Isend非阻塞式通讯函数，MPI_wait.MPI_Sendrecv.

2026/5/621Cannon-bug

死锁的问题 问题来源于main_shift()这个函数中MPI函数的使用。在Cannon-mpi代码的main_shift()模块中，文献中算法使用的是MPI的阻塞通信函数:MPI_Send/MPI_Recv,这就使得Cannon算法在执行循环左移和循环上移时，矩阵规模超过共享buff的容量时出现循环等待的死锁状况。 在曙光4000集群系统上，该算法的发生死锁的矩阵下限规模是200×200的浮点型矩阵。2026/5/622Cannon-bug原始（阻塞式）的main_shift模块：

voidmain_shift() { … /*将分块b左移位*/

MPI_Send(a,dl2,MPI_FLOAT,get_index(my_row,my_col-1, sp),1, MPI_COMM_WORLD);

MPI_Recv(a,dl2,MPI_FLOAT,get_index(my_row, my_col+1,sp),1, MPI_COMM_WORLD,&status); /*将分块b上移位*/

MPI_Send(b,dl2,MPI_FLOAT,get_index(my_row-1,my_col, sp),1, MPI_COMM_WORLD);

MPI_Recv(b,dl2,MPI_FLOAT,get_index(my_row+1, my_col,sp),1, MPI_COMM_WORLD,&status);

}2026/5/623Cannon-bug改进（非阻塞式）的main_shift模块

…

c[i*dl+j]+=a[i*dl+k]*b[j*dl+k]; //改进了的Cannon按行存取

/*将分块a左移位*/

MPI_Isend(a,dl2,MPI_FLOAT,get_index(my_row,my_col-1,sp),1,MPI_COMM_WORLD,&myrequest_s);

MPI_Irecv(buf,dl2,MPI_FLOAT,get_index(my_row,my_col+1,sp),1,MPI_COMM_WORLD,&myrequest_r);

MPI_Wait(&myrequest_s,&status);

MPI_Wait(&myrequest_r,&status);

memcpy(a,buf,sizeof(float)*dl2); /*将分块b上移位*/

MPI_Isend(b,dl2,MPI_FLOAT,get_index(my_row-1,my_col,sp),1,MPI_COMM_WORLD,&myrequest_s);

MPI_Irecv(buf,dl2,MPI_FLOAT,get_index(my_row+1,my_col,sp),1,MPI_COMM_WORLD,&myrequest_r);

MPI_Wait(&myrequest_s,&status);

MPI_Wait(&myrequest_r,&status);

memcpy(b,buf,sizeof(float)*dl2);2026/5/624Cannon-bugMPI_Irecv仅仅初始化接受操作,在与之对应的MPI_Wait函数的调用返回之前，将不能访问bufferMPI_Irecv函数返回时，handle指向一个MPI_Request对象，它代表了一个已近初始化了的通信操作。这个函数并不返回一个指向MPI_Status对象的指针，因为实际的接受操作并未完成。MPI_Wait会一直阻塞，直至参数handle所关联的操作完成，对发送来说，此时就可以向缓冲区写入新的值。而对接收来 说，便可以从