矩阵乘法的MPI并行计算

最新推荐文章于 2024-09-28 14:08:19 发布
最新推荐文章于 2024-09-28 14:08:19 发布
阅读量5k 收藏 26
点赞数 5
分类专栏： MPI 文章标签： MPI

转载自 http://www.cnblogs.com/fengfu-chris/p/4364142.html

1、问题描述

矩阵乘法问题描述如下：

　　给定矩阵A和B，其中A是m*p大小矩阵，B是p*n大小的矩阵。求C = A*B。

求解这个问题最简单的算法是遍历A的行和B的列，求得C的相应元素，时间复杂度O(mnp)，空间复杂度O(1)。

// 矩阵乘法的C++实现
for(int i=0; i<m; i++){
    for(int j=0; j<n; j++){
        float temp = 0.0;
        for(int k=0; k<p; k++){
            temp += A[i*p + k] * B[k*n + j];
        }
        C[i*n + j] = temp;
    }
}

2、最简单的并行方案

要改进上述算法为并行算法，需要了解到C++ MPI编程的特点：

　　a. 各个进程之间不能有依赖。这是因为各个进程可以以任意的时间顺序执行。

　　b. 数据是分布式存储的。也就是说，每个进程有自己独立的数据备份。

有了这两点认识后，一种最简单的并行方案就出来了：（假设开启np个进程）

　　(1). 首先将矩阵A和C按行分为np块；

　　(2). 进程号为 id 的进程读取A的第 id 个分块和B；

　　(3). 进程号为 id 的进程求解相应的C的第 id 个分块。

代码如下：

/* filename: matMultiplyWithMPI.cpp
 * parallel matrix multiplication with MPI
 * C(m,n) = A(m,p) * B(p,n)
 * input: three parameters - m, p, n
 * @copyright: fengfu-chris
 */
#include<iostream>
#include<mpi.h>
#include<math.h>
#include<stdlib.h>

void initMatrixWithRV(float *A, int rows, int cols);
void matMultiplyWithSingleThread(float *A, float *B, float *matResult, int m, int p, int n);

int main(int argc, char** argv)
{
    int m = atoi(argv[1]);
    int p = atoi(argv[2]);
    int n = atoi(argv[3]);

    float *A, *B, *C;
    float *bA, *bC;  

    int myrank, numprocs;

    MPI_Status status;

    MPI_Init(&argc, &argv);  // 并行开始
    MPI_Comm_size(MPI_COMM_WORLD, &numprocs); 
    MPI_Comm_rank(MPI_COMM_WORLD, &myrank); 

    int bm = m / numprocs;

    bA = new float[bm * p];
    B  = new float[p * n];
    bC = new float[bm * n];

    if(myrank == 0){
        A = new float[m * p];
        C = new float[m * n];

        initMatrixWithRV(A, m, p);
        initMatrixWithRV(B, p, n);
    }
    MPI_Barrier(MPI_COMM_WORLD);

　　/* step 1: 数据分配 */
    MPI_Scatter(A, bm * p, MPI_FLOAT, bA, bm *p, MPI_FLOAT, 0, MPI_COMM_WORLD);
    MPI_Bcast(B, p * n, MPI_FLOAT, 0, MPI_COMM_WORLD);

　　/* step 2: 并行计算C的各个分块 */
    matMultiplyWithSingleThread(bA, B, bC, bm, p, n);

    MPI_Barrier(MPI_COMM_WORLD);

　　/* step 3: 汇总结果 */
    MPI_Gather(bC, bm * n, MPI_FLOAT, C, bm * n, MPI_FLOAT, 0, MPI_COMM_WORLD);

　　/* step 3-1: 解决历史遗留问题（多余的分块） */
    int remainRowsStartId = bm * numprocs;
    if(myrank == 0 && remainRowsStartId < m){
        int remainRows = m - remainRowsStartId;
        matMultiplyWithSingleThread(A + remainRowsStartId * p, B, C + remainRowsStartId * n, remainRows, p, n);
    }

    delete[] bA;
    delete[] B;
    delete[] bC;

    if(myrank == 0){
        delete[] A;
        delete[] C;
    }

    MPI_Finalize(); // 并行结束

    return 0;
}

void initMatrixWithRV(float *A, int rows, int cols)
{
    srand((unsigned)time(NULL));
    for(int i = 0; i < rows*cols; i++){
        A[i] = (float)rand() / RAND_MAX;
    }
}
void matMultiplyWithSingleThread(float *A, float *B, float *matResult, int m, int p, int n)
{
    for(int i=0; i<m; i++){
        for(int j=0; j<n; j++){
            float temp = 0;
            for(int k=0; k<p; k++){
                temp += A[i*p+k] * B[k*n + j];
            }
            matResult[i*n+j] = temp;
        }
    }
}

编译：

$mpigxx matMultiplyWithMPI.cpp -o matMultiplyWithMPI

运行：

$mpirun -np 8 matMultiplyWithMPI 3000 2000 4000

这里假设m = 3000， p = 2000， n = 4000。另外，开启的进程数为8个。 np的个数可以大于CPU的个数。

一般来讲，只有当矩阵大小大于5000的量级时，开启几十上百个进程的威力才能凸显出来。尤其是当矩阵量级达到万维以上时，串行或是少数几个进程并行的矩阵乘法将变得特别耗时。

3、改进的并行方案：内存考虑

上面的并行方案有个很大的缺陷，那就是 B 的备份数和开启的进程数一致。这对于内存不是很充裕或矩阵很大的时候，会导致灾难！例如，假设 B 是10000*10000维的，用double类型存储大概占700M左右的内存。当开启的进程数达到128个时，单是 B 的备份占据的内存开销将达到 128 * 700 M = 90G。 这将耗掉巨大的内存！

有什么改进的方案呢？

必须了解MPI的第三个特点：

　　c. 进程之间可以很方便地通信，并且支持多种通信方案。

这样，就可以把 B 也同时分布式的存储到各个进程对应的内存中，然后利用进程之间的通信来轮换各个 B 的分块，从而达到减小内存开销的效果。当然，几乎和所有的程序一样，离不开时间与空间的trade-off。所以，这种方法虽然节省了内存，却要消耗大量的时间在进程之间的通信上。

下面给出改进的并行方案：

　　(1). 将A和C按行分为np块，将B按列分为np块（B可以按列存储）；

　　(2). 进程号为 id 的进程读取 A 和 B 的第id个分块；

　　(3). 循环np次：

　　　　<1>. 各个进程用各自的A、B分块求解C的分块；

　　　　<2>. 轮换B的分块（例如：id 号进程发送自己当前的B的分块到 id+1号进程）

代码如下：

/* filename: matMultiplyWithMPI_updated.cpp
 * parallel matrix multiplication with MPI: updated
 * C(m,n) = A(m,p) * B(p,n)
 * input: three parameters - m, p, n
 * @copyright: fengfu-chris
 */
#include<iostream>
#include<mpi.h>
#include<math.h>
#include<stdlib.h>

void initMatrixWithRV(float *A, int rows, int cols);
void copyMatrix(float *A, float *A_copy, int rows, int cols);
// A: m*p, B: p*n  ！！！ note that B is stored by column first
void matMultiplyWithTransposedB(float *A, float *B, float *matResult, int m, int n, int p);

int main(int argc, char** argv)
{
　 int m = atoi(argv[1]);
　 int n = atoi(argv[2]);
　 int p = atoi(argv[3]);

   float *A, *B, *C;
   float *bA, *bB_send, *bB_recv, *bC, *bC_send;
　 int myrank, numprocs;

    MPI_Status status;

    MPI_Init(&argc, &argv);  // 并行开始
    MPI_Comm_size(MPI_COMM_WORLD, &numprocs); 
    MPI_Comm_rank(MPI_COMM_WORLD, &myrank); 
 　　
　  int bm = m / numprocs;
    int bn = n / numprocs;

    bA = new float[bm * p];
    bB_send = new float[bn * p];
    bB_recv = new float[bn * p];
    bC = new float[bm * bn];
    bC_send = new float[bm * n];

    if(myrank == 0){
        A = new float[m * p];
        B = new float[n * p];
        C = new float[m * n];

        initMatrixWithRV(A, m, p);
        initMatrixWithRV(B, n, p);
    }

    MPI_Barrier(MPI_COMM_WORLD);
    MPI_Scatter(A, bm * p, MPI_FLOAT, bA, bm * p, MPI_FLOAT, 0, MPI_COMM_WORLD);
    MPI_Scatter(B, bn * p, MPI_FLOAT, bB_recv, bn * p, MPI_FLOAT, 0, MPI_COMM_WORLD);

    int sendTo = (myrank + 1) % numprocs;
    int recvFrom = (myrank - 1 + numprocs) % numprocs;

    int circle = 0;  
    do{
        matMultiplyWithTransposedB(bA, bB_recv, bC, bm, bn, p);
        int blocks_col = (myrank - circle + numprocs) % numprocs;
        for(int i=0; i<bm; i++){
            for(int j=0; j<bn; j++){
                bC_send[i*n + blocks_col*bn + j] = bC[i*bn + j];
            }
        }

        if(myrank % 2 == 0){
            copyMatrix(bB_recv, bB_send, bn, p);
            MPI_Ssend(bB_send, bn*p, MPI_FLOAT, sendTo, circle, MPI_COMM_WORLD);
            MPI_Recv(bB_recv, bn*p, MPI_FLOAT, recvFrom, circle, MPI_COMM_WORLD, &status);
        }else{
            MPI_Recv(bB_recv, bn*p, MPI_FLOAT, recvFrom, circle, MPI_COMM_WORLD, &status); 
            MPI_Ssend(bB_send, bn*p, MPI_FLOAT, sendTo, circle, MPI_COMM_WORLD);
            copyMatrix(bB_recv, bB_send, bn, p);    
        }

        circle++;
    }while(circle < numprocs);

　 MPI_Barrier(MPI_COMM_WORLD);
    MPI_Gather(bC_send, bm * n, MPI_FLOAT, C, bm * n, MPI_FLOAT, 0, MPI_COMM_WORLD);

    if(myrank == 0){
        int remainAStartId = bm * numprocs;
        int remainBStartId = bn * numprocs;

        for(int i=remainAStartId; i<m; i++){
            for(int j=0; j<n; j++){
                float temp=0;
                for(int k=0; k<p; k++){
                    temp += A[i*p + k] * B[j*p +k];
                }
                C[i*p + j] = temp;
            }
        }

        for(int i=0; i<remainAStartId; i++){
            for(int j=remainBStartId; j<n; j++){
                float temp = 0;
                for(int k=0; k<p; k++){
                    temp += A[i*p + k] * B[j*p +k];
                }
                C[i*p + j] = temp;
            }
        }
    }

    delete[] bA;
    delete[] bB_send;
    delete[] bB_recv;
    delete[] bC;
    delete[] bC_send;

    if(myrank == 0){
        delete[] A;
        delete[] B;
        delete[] C;
    }

    MPI_Finalize(); // 并行结束

    return 0;
}

void initMatrixWithRV(float *A, int rows, int cols)
{
    srand((unsigned)time(NULL));
    for(int i = 0; i < rows*cols; i++){
        A[i] = (float)rand() / RAND_MAX;
    }
}
void copyMatrix(float *A, float *A_copy, int rows, int cols)
{
    for(int i=0; i<rows*cols; i++){
        A_copy[i] = A[i];
    }
}

void matMultiplyWithTransposedB(float *A, float *B, float *matResult, int m, int p, int n)
{
    for(int i=0; i<m; i++){
        for(int j=0; j<n; j++){
            float temp = 0;
            for(int k=0; k<p; k++){
                temp += A[i*p+k] * B[j*p+k];
            }
            matResult[i*n+j] = temp;
        }
    }
}

这里最需要注意的地方就是B的轮换。 有两点需要注意：

(1) 防阻塞机制。这里采用奇偶原则：偶数号进程先发送，再接收；奇数号进程则相反。这样可以避免所有进程同时发送造成死锁的情况；

　(2) 数据备份。发送和接收的信息存储在不同的矩阵中，这样保证原来的信息不会被覆盖。

这种方法的优点是显而易见的。对于足够牛的服务器/计算机集群，开启成百上千个进程来并行完全不是问题。

__Sunny__
关注
专栏目录
矩阵乘法 mpi c语言,C语言基于MPI并行计算矩阵的乘法
weixin_35997546的博客
05-18 1470
浏览量：1,023(在苹果系统下，如果文章中的图片不能正常显示，请升级Safari浏览器到最新版本，或者使用Chrome、Firefox浏览器打开。)矩阵的计算是并行计算里的一个很重要的问题。矩阵是一种网格化的数据，是一组同类型数值的集合，矩阵的出现，使得代数系统更完善，对各种实际问题的求解产生了巨大的作用。但是其庞大的计算量往往令人生畏，稍微大一点的矩阵，计算就变得非常繁琐，不仅仅是人，而且串行...
并行编程简单的MPI编程hello world
07-06
基于vs2010环境编写的并行编程简单的MPI编程hello world
MATLABMPI：掌握MATLAB并行计算技术
最新发布
weixin_35750747的博客
09-28 1373
本文还有配套的精品资源，点击获取 简介：并行计算是加速复杂科学计算任务的关键技术。MATLABMPI工具箱，由麻省理工学院林肯实验室开发，结合了MPI标准与MATLAB易用性，使得用户可以在MATLAB环境下进行高效的并行编程。该工具箱支持利用多核处理器或集群资源，通过关键功能如并行启动管理、工作区共享、并行函数、MPI通信原语、数据分布和并行优化等，大幅提升计算效率。此外...
C++ MPICH
weixin_34292402的博客
03-09 98
假设一个C++的MPI程序在单机上能够跑。可是在多机上跑会报下面错误： Fatal error in MPI_Send: Unkown error class , error stack。 解决方法： 把machinefile里面的全部hostname换成相应的IP地址就能够了。 然后再又一次跑： mpirun -n 10 -f machinefile ./a.out ...
c语言 mpi程序,C语言使用MPI实现并行操作（MPI入门）
weixin_35815640的博客
05-18 1517
C语言使用MPI：#include "mpi.h"#include "stdio.h"#include "math.h"#include "unistd.h"#define n 10void printArr(int *arr){for(int i = 0; i < n; i++){for(int j = 0; j < n; j++){printf("%d ",arr[i * n + ...
MPI并行计算 矩阵乘法
05-31
MPI并行计算 矩阵乘法 如果是10000个线程的话可能会出问题
矩阵乘法MPI并行程序报告.pdf
06-06
【矩阵乘法MPI并行程序】是分布式计算中的一种常见应用，主要利用Message Passing Interface (MPI)库来实现。MPI是一种标准，允许程序员在多处理器系统或者计算机集群上编写并行程序，它提供了丰富的通信原语，如MPI...
矩阵乘法mpi实现 并行运算
05-27
矩阵乘法mpi实现 并行运算 在linux下成功运行 使用mpicc -o 编译 使用mpirun命令运行
矩阵乘法MPI并行程序报告.docx
07-12
【矩阵乘法MPI并行程序】是计算机科学领域中一种并行计算技术的应用，主要用于提高大规模矩阵乘法的计算效率。MPI（Message Passing Interface）是一种标准的并行编程接口，适用于分布式内存系统，如集群环境。在本...
MPI的环境搭建C++编程
06-24
MPI简介，环境的搭建，VC的配置，简单的程序示例，并行算法.....................
《并行程序设计导论》 mpi编写hello
qq_41180336的博客
10-06 1314
在消息传递程序中，运行在一个核-内存上的程序 通常称为一个进程。两个进程可以通过调用函数来进行通信。我们将使用消息传递的实现称为消息传递接口(Message-Passing Interface,MPI)。它是一个函数库 #include<mpi.h> 用$开头的，我是用的linux的终端写的，在root模式下 其余的都是普通的.c文件内容 编译（成为可执行文件exe格式） $mpicc...
MPI编译C++程序出现#error "SEEK_SET is #defined but must not be for the C++ binding of MPI" 的解决方法
swanabin的专栏
07-10 1071
转自VC错误:http://www.vcerror.com/?p=1627 问题描述： MPI的mpicc和mpicxx命令编译C++程序时，可能会遇到如下三种错误信息： #error "SEEK_SET is #defined but must not be for the C++ binding of MPI" #error "SEEK_CUR is #defined but must
MPI调用C++程序
CV_2013的专栏
07-16 2163
在运行MPI程序时，有时候需要让并行运算过程中的多个实例共享同一个文件。前提是你安装了MPICH2或者Open MPI，并使用mpicc (for C)和mpicxx（for C++）进行编译。方法如下： //C Version #include "mpi.h" #include *** MPI_Init(&argc,&argv); //MPI初始化 int rank,size;
MPI入门
sttony的专栏
11-20 3063
敢时髦也写点BlogMPI 入门（MPICH2 windows版本）1．什么是MPIMPI是一个在平行计算中传递消息的库的标准，由实现人员和使用人员来遵守。目前的实现版本有MPICH2， Argonne National Laboratory实现，他还有好几个派生子项目。2．安装MPICH2       使用VC7和windows版本的MPICH2。MPICH2 Home
并行作业1：MPI安装，及示例程序运行
八月对半
12-03 7244
运行MPI程序 系统采用vm下ubuntu16.04 一、MPI系统安装 1.1 安装环境 (1)操作系统：Ubuntu 16.04.4 （64位）。 (2)g++ 版本：gcc version 5.4.0 1.2 安装包下载地址及安装包 (1)下载地址：http://www.mpich.org/downloads/ (2)安装包：mpich-3.2.1.tar.gz 1.3安装与测试 (1)将安...
c语言mpi并行程序,高性能计算之并行编程技术MPI并行程序设计(完整版).pdf
weixin_32214771的博客
05-20 1248
高性能计算之并行编程技术MPI并行程序设计(完整版)高性能计算之并行编程技术—— MPI并行程序设计都志辉 编著李三立 审阅陈渝 刘鹏 校对I内容提要本书介绍目前最常见的并行程序— MPI并行程序的设计方法 它适合高校三四年级本科生 非计算机专业研究生作为教材和教学自学参...
水文日常~~Linux命令行下如何编译运行MPI程序(C/C++)
Jing Sir的博客
06-07 8821
Linux命令行下如何编译运行MPI程序 如果是C程序，使用如下命令进行编译： mpicc helloworld.c -o helloworld 如果是C++程序，使用如下命令进行编译： mpicxx hellowprld.cpp -o helloworld 编译完成后，使用如下命令运行 mpirun -np 4 ./helloworld 其中"-np"用来指定进程数 ...
评论 2
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值