北航并行计算矩阵转置作业
北航 数值分析第二次大作业(带双步位移的QR方法)

一、算法设计方案:按题目要求,本程序运用带双步位移的QR方法求解给定矩阵的特征值,并对每一实特征值,求解其相应的特征向量。
总体思路:1)初始化矩阵首先需要将需要求解的矩阵输入程序。
为了防止矩阵在后面的计算中被破坏保存A[][]。
2)对给定的矩阵进行拟上三角化为了尽量减少计算量,提高程序的运行效率,在对矩阵进行QR分解之前,先进行拟上三角化。
由于矩阵的QR 分解不改变矩阵的结构,所以具有拟上三角形状的矩阵的QR分解可以减少大量的计算量。
这里用函数void QuasiTriangularization()来实现,函数形参为double型N维方阵double a[][N]。
3)对拟上三角化后的矩阵进行QR分解对拟上三角化的矩阵进行QR分解会大大减小计算量。
用子程序void QR_decomposition()来实现,将Q、R设为形参,然后将计算出来的结果传入Q和R,然后求出RQ乘积。
4)对拟上三角化后的矩阵进行带双步位移的QR分解为了加速收敛,对QR分解引入双步位移,适当选取位移量,可以避免进行复数运算。
为了进一步减少计算量,在每次进行QR分解之前,先判断是否可以直接得到矩阵的一个特征值或者通过简单的运算得到矩阵的一对特征值。
若可以,则得到特征值,同时对矩阵进行降阶处理;若不可以,则进行QR分解。
这里用函数intTwoStepDisplacement_QR()来实现。
这是用来存储计算得到的特征值的二维数组。
考虑到特征值可能为复数,因此将所有特征值均当成复数处理。
此函数中,QR分解部分用子函数void QR_decompositionMk()实现。
这里形参有三个,分别用来传递引入双步位移后的Mk阵,A矩阵,以及当前目标矩阵的维数m。
5)计算特征向量得到特征值后,计算实特征值相应的特征向量。
这里判断所得特征值的虚数部分是否为零。
求实特征值的特征向量采用求解相应的方程组((A-λI)x=0)的方法。
因此先初始化矩阵Array,计算(A-λI),再求解方程组。
北航数值分析计算实习题目二 矩阵QR分解

数值分析实习二院(系)名称航空科学与工程学院专业名称动力工程及工程热物理学号SY0905303学生姓名解立垚1. 题目试用带双步位移QR 的分解法求矩阵A=[a ij ]10*10的全部特征值,并对其中的每一个实特征值求相应的特征向量。
已知()sin 0.50.2,1.5cos 1.2,ij i j i j a i j i j ⎧⎫+≠⎪⎪=⎨⎬+=⎪⎪⎩⎭(),1,2,...,10i j =。
说明:1、求矩阵特征值时,要求迭代的精度水平为1210ε-=。
2、打印以下内容:算法的设计方案;全部源程序(要求注明主程序和每个子程序的功能); 矩阵A 经过拟上三角话之后所得的矩阵()1n A -;对矩阵()1n A-进行QR 分解方法结束后所得的矩阵;矩阵A 的全部特征值()(),1,2,......10i i iR I i λ=,和A 的相应于实特征值的特征向量;其中()(),.i e i m i R R I I λλ==如果i λ是实数,则令0.i I =3、采用e 型输出数据,并且至少显示12位有效数字。
2. 算法设计方案本题采用带双步位移的QR 分解方法。
为了使程序简洁,自定义类Xmatrix ,其中封装了所需要的函数方法。
在Xmatrix 类中封装了运算符重载的函数,即定义了矩阵的加、减、乘、除、数乘运算及转置运算(T())。
同时为了避免传递数组带来的额外内存开销,使用引用(&)代替值传递,以节省内存空间,避免溢出.(1)此程序的主要部分为Xmatrix 中的doubleQR()方法,具体如下:Step1:使用矩阵拟上三角化的算法将A 化为拟上三角阵A (n-1)(此处调用Xmatrix 中的preQR()方法)Step2:令121,,10k m n ε-===, 其中k 为迭代次数。
Step3:如果,1m m a ε-≤,则得到A 的一个特征值,m m a ,令1m m =-,goto Step4;否则goto Step5.Step4: 如果1m =,则得到A 的一个特征值11a ,goto Step11;如果0m =,则goto Step11;如果1m >,则goto Step3;Step5(Step6):如果2m =,则得到A 的两个特征值12s s 和(12s s 和为右下角两阶子阵对应的特征方程21,1,()det 0m m m m a a D λλ---++=的两个根。
矩阵快速转置算法详解

矩阵快速转置算法详解好啦,今天我们来聊聊矩阵快速转置算法,这可是个让人觉得又爱又恨的话题。
矩阵嘛,想必大家都不陌生,毕竟它就像我们生活中的一张大表格,记录着各种数据。
可这转置,听起来就像是魔法一样,把行变成列,列又变回行,这其中可有不少讲究。
想象一下,一个矩阵就像是你的冰箱,里面装着各种各样的食材。
平常你打开冰箱,看到一层层整齐的食物,心里就觉得舒服。
不过,有时候你也想换个角度,看看这些食材在不同位置的组合。
这个时候,转置就登场啦!转置后的矩阵,就像是你重新整理过的冰箱,虽然东西没变,但看上去就是更有条理,更让人一目了然。
转置到底是怎么个操作呢?其实就是把矩阵里的元素交换位置。
比如,原本在第一行第二列的元素,现在要跑到第二行第一列。
听起来简单对吧?但当你的矩阵大得像个足球场,里面的元素多得跟星星一样,那可就得小心翼翼,别让它们跑丢了。
这里就需要一个“快速”的方法,避免慢吞吞的操作,简直就是搬家时的“火速搬迁”,让你几乎感觉不到麻烦。
说到快速转置算法,最经典的就是“分块转置”。
这就像你在打扫房间的时候,通常不会一次性把所有东西都搬出去,而是先把几个小块整理好。
这个算法也是如此,把大矩阵分成小块,分别转置,然后再合并。
这种方法就像是在料理过程中,把大菜分成几道小菜,最后汇聚成一桌美食,既高效又省力。
分块的方法还有个好处,就是它能更好地利用缓存,提高处理速度,简直是个聪明的“小调皮”。
转置过程中,你可能会碰到一些小坑,比如说对称矩阵。
它们就像是一对双胞胎,无论怎么交换,最后看上去都是一样的。
对此,快速转置算法能巧妙地利用这一点,不用额外的空间,就能轻松搞定,简直让人拍手叫好。
就像生活中遇到的那些小聪明,总是能在关键时刻帮你一把,真是妙不可言。
转置算法的实现可不仅仅是写代码那么简单。
想象一下,代码里的每一行都像是一位舞者,它们需要配合默契,才能让整个矩阵的转置变得流畅自如。
这时候,算法的优化就显得尤为重要。
北航数值分析大作业第二题(fortran)

!计算A(r+1) DO I=1,N DO J=1,N A(I,J)=A(I,J)-W(I)*U(J)-U(I)*P(J) ENDDO ENDDO ENDIF ENDDO RETURN END
!***************符号函数子程序*****************! FUNCTION SGN(X) REAL(8) X IF(X>0) THEN SGN=1 ELSE IF(X<0) THEN SGN=-1 ELSE IF(X==0) THEN SGN=0 ENDIF END
DIMENSION A(N,N),A1(N,N),A2(N,N),C(2,N),Q(N,N),R(N,N),CR(N),CM(N)!C为存储特征值的数 组,1为实部,为虚部 REAL(8) A,A1,A2,C,Q,R,CM E=1E-12 L=1000 !精度水平 !迭代最大次数
OPEN(1,FILE='数值分析大作业第二题计算结果.TXT') DO I=1,N DO J=1,N IF(I==J) THEN A(I,J)=1.52*COS(I+1.2*J) ELSE A(I,J)=SIN(0.5*I+0.2*J) ENDIF ENDDO ENDDO A1=A WRITE(*,"('矩阵A为:')") WRITE(1,"('矩阵A为:')") DO I=1,N DO J=1,N WRITE(*,"(2X,E20.13,2X,\)") A(I,J) WRITE(1,"(2X,E20.13,2X,\)") A(I,J) ENDDO WRITE(*,"(' ')") WRITE(1,"(' ')") ENDDO !使用矩阵的拟上三角化的算法将矩阵A化为拟上三角矩阵A(n-1) CALL HESSENBERG(A,N) WRITE(*,"('拟上三角化后矩阵A(n-1)为:')") WRITE(1,"('拟上三角化后矩阵A(n-1)为:')") DO I=1,N DO J=1,N WRITE(*,"(2X,E20.13,2X,\)") A(I,J) WRITE(1,"(2X,E20.13,2X,\)") A(I,J) ENDDO WRITE(*,"('')") WRITE(1,"('')") ENDDO !计算对矩阵A(n-1)实行QR方法迭代结束后所得矩阵 A2=A CALL QRD(A2,N,Q,R)
北航矩阵考题A

an 1 2an 2 (a , a , 1 2 nan n
, an ) (分解不唯一) ( 4 分)
a1 2a2
nan trA
A2 (trA) A 0
A100 ( )99 A (a1 nan )99 A.
(2)可知 f ( A) 的谱公式为 f ( A) f (1 )G1 f (2 )G2 ,(且 G1 , G2 同上), 令 f ( x) sin x f (7) sin 7, f (2) sin 2 得 sin A 的谱分解为
4 4 5 4 1 sin A (sin 7)G1 (sin 2)G2 sin 7 1 sin 2 9 9 5 5 5 4
(3) 由 (e A ) {e7 , e 2 } 谱半径 (e A ) e7 由 (esin A ) {esin 7 , e sin 2 } 行列式 det(esin A ) esin 7 sin 2
2
1 1 1 1 四、 (15 分)设 A 1 2 3 , b 2 , 0 0 0 0
1 n 1 n 2
n (A 2 ) 为收敛
( 4 分)
以下的六、七题中只需任选一题:
4
六、(15 分)(1) 设矩阵 A 最小式 m( x) ( x 2) 2 且 f ( A) 收敛,推导 f ( A) 的广谱 计算公式。解(1)由 ( A 2)2 0 与台乐公式 f ( A) 0 可得公式
(3)求谱半径 (e A ) 与行列式 det(esin A ) .(4):求 ln( I A ) ?
北航并行计算矩阵转置作业

矩阵转置并行实现1、算法描述:若处理器个数为P,且它们的编号依次是0,1, 则将n阶矩阵A分成p 个大小为mxm的子块,m=[n/p]o p个子块组成一个子块阵列,记其中第i行第j列的子块为Aij,它含有第(i-l)m+l至第im行中的第(j-l)m+l至第jm列的所有元素。
对每一处理器按行主方式赋以二维下标,记编号为i的处理器的二维下标为(v,u),其中u=imod A/p ,将A的子块存入下标为(v,u)表示的对应处理器中。
转置分为两步进行:第一步,子块转置;第二步,处理器内部转置。
为了避免对应子块交换数据是处理器发生死锁,可令下三角块先向与之对应的上三角子块发送数据,然后从上三角子块接收数据;上三角子块数据先存放在缓冲区buffe「中,然后从与之对应的下三角子块接收数据,最后再将缓冲区中的数据发送给下三角子块,流程图如下所示:2、程序代码:#include "stdio.h"#include "stdlib.h"#include ”mpi.h”#include "math.h"#define E 0.0001#define a(x,y) a[x*m+y]#define b(x,y) b[x*m+y]#define A(x,y) A[x*size+y]#define B(x,y) B[x*size+y]#define intsize sizeof(int)#define floatsize sizeof(float)#define charsize sizeof(char)int size,N;int m;int t;float *A, *B;double starttime;double timel;double time2;int my_rank;int p;MPI_Status status;FILE *fdA;void Environment_Finalize(float *ajloat *b){free(a);free(b);}int main(int argc, char **argv){int ij/kmy^rank,group_size;float *a,*b;int u,v;float temp;MPIJ nit(&argc,&argv);MPI_Comm_size(MPI_COMM_WORLD/&group_size);M P l_Co m m_ra nk(MPI_COM M_W0 RLD,&my_ra nk); p=group_size;if(my_ra nk==0) starttime=MPI_Wtime(); fdA=fopen("dataln.txt,,;,r H);fscanf(fdA z"%d %d", &size, &N);if(size != N){puts(H The in put is error!"); exit(O);}A=(float*)malloc(floatsize*size*size);B=(float*)malloc(floatsize*size*size);for(i = 0; i < size; i ++){for(j = 0; j < size; j ++) fscanf(fdA, ”%化A+i*size+j);}fclose(fdA);MPI_Bcast(&size,l,MPI JN TAMPI_COMM_WORLD); t=(int)sqrt(p);if (t>size)t=size;i#si ze%t!=0)for(;;){t-;if(size%t==0)break;p=t*t;m=size/t;a=(float *)malloc(floatsize*m*m); b=(float *)malloc(floatsize*m*m); if(a==NULL||b==NULL)printf("allocate space fail!11);if (my_ra nk==0) for(i=0;i<m;i++)for(j=0;j<m;j++)a(i,j)=A(i,j);}if (my_rank==O){for(i=l;i<p;i++){v=i/t;u=i%t;for(j=v*m;j<(v+l)*m;j++)for(k=u*m;k<(u+l)*m;k++)b((j%m),(k%m))=A(j,k);MPI_Send(bm*m,l\/IPI_FLOAl;ijMPI_COIVIM_WORLD);}}else if (my_rank<p)MPI_Recv(am*m,MPI_FLOATQmy_rank,MPI_COI\/IM_WORLD&status);timel=MPI_Wtime();if ((my_rank/t)>(my_rank%t)&&my_rank<p){v= my_rank/t;u= my_ra nk%t;MPLSend(a/m*m/MPLFLOAT, (u*t+v),(u*t+v)/MPLCOMM_WORLD);MPLRecvfa^^^PLFLOATfu^+vhmy^ran^MPLCOMM^WORL^&status);}if ((my_rank/t)<(my_rank%t)&&my_rank<p){v= my_rank/t;u= my_ra nk%t;for(i=0;i<m;i++)for(j=0;j<m;j++)b(ij)=a(ij);MPI_Recv(am*rrLMPI_FLOA7;(屮t+v), my_rank,MPI_COMM_WORLD&status);MPLSend(b/m*m/MPLFLOAl;(u*t+v)/(u*t+v)/MPLCOMM_WORLD);for(i=l;i<m;i++)for(j=0;j<i;j++){temp=a(ij);a(i,j)=a(j,i);a(jj)=temp;}if (my_rank==O){for(i=0;i<m;i++)for(j=0;j<m;j++)B(i,j)=a(i,j);}if (my_rank==O){for(i=l;i<p;i++){MPLRecv(a z m*m z MPLFLOATJJ z MPLCOMM_WORLD,&status); V= i/t; u= i%t;for(j=v*m;j<(v+l)*m;j++) for(k=u*m;k<(u+l)*m;k++)B(j,k)=a((j%m),(k%m));}}else if(my_rank<p)MPLSend(a/m*m/MPLFLOAT,0,my-rank/MPLCOMM_WORLD);if (my_rank==O){printf(H lnput of file V'dataln.txtV'Xn"); printf(,,%d\t%d\n,,/ size, size);for(i=0;i<size;i++)for(j=0;j<size;j++) printf(,,%f\t,,,A(i/j)); printfCV);printf(M\nOutput of Matrix AT\n M);for(i=0;i<size;i++){for(j=0;j<size;j++) printf(l,%f\t,,,B(i/j)); printfCV);}}time2=MPI_Wtime();if (my_rank==0){printf(”\n”);printf(n processor number: %d\n H z p);printf(H Whole running time = %f seconds\n,,/time2-starttime);printf(H Distribute data time = %f sec onds\n'\timel-starttime); prin tf(H Parallel computetime = %f sec on ds\n':time2・time:l);}MPI_Barrier(MPLCOMM_WORLD); MPI_Finalize();Environme nt_Fin alize(a,b); return(O);。
北航数值分析大作业第二题

北航数值分析大作业第二题本页仅作为文档封面,使用时可以删除This document is for reference only-rar21year.March数值分析第二次大作业史立峰SY1505327一、 方案(1)利用循环结构将sin(0.50.2)()1.5cos( 1.2)(){i j i j ij i j i j a +≠+==(i,j=1,2,……,10)进行赋值,得到需要变换的矩阵A ;(2)然后,对矩阵A 利用Householder 矩阵进行相似变换,把A 化为上三角矩阵A (n-1)。
对A 拟上三角化,得到拟上三角矩阵A (n-1),具体算法如下: 记A(1)=A ,并记A(r)的第r 列至第n 列的元素为()n r r j n i a r ij,,1,;,,2,1)( +==。
对于2,,2,1-=n r 执行 1. 若()n r r i a r ir,,3,2)( ++=全为零,则令A(r+1) =A(r),转5;否则转2。
2. 计算()∑+==nr i r irr a d 12)(()()r r r r r r r r r r d c a d a c ==-=++则取,0sgn )(,1)(,1若 )(,12r rr r r r a c c h +-=3. 令()nTr nrr r r r r r r r R a a c a u ∈-=++)()(,2)(,1,,,,0,,0 。
4. 计算r r T r r h u A p /)(= r r r r h u A q /)(=r r Tr r h u p t /=r r r r u t q -=ωT rr T r r r r p u u A A --=+ω)()1(5. 继续。
(3)使用带双步位移的QR 方法计算矩阵A (n-1)的全部特征值,也是A 的全部特征值,具体算法如下:1. 给定精度水平0>ε和迭代最大次数L 。
北航数值分析报告大作业二

数值分析大作业(二)学院名称宇航学院专业名称航空宇航推进理论与工程学生姓名段毓学号SY16153062016年11月5日1 算法设计方案首先将矩阵A 进行拟上三角化,把矩阵A 进行QR 分解,计算出RQ 。
要得出矩阵A 的全部特征值,首先对A 进行QR 的双步位移得出特征值。
最后,采用列主元的高斯消元法求解特征向量。
1.1 A 的拟上三角化因为对矩阵进行QR 分解并不改变矩阵的结构,因此在进行QR 分解前对矩阵A 进行拟上三角化可以大大减少计算机的计算量,提高程序的运行效率。
具体算法如下所示,记A A =)1(,并记)(r A 的第r 列至第n 列的元素为()n r r j n i a r ij,,1,;,,2,1)(ΛΛ+==。
对于2,,2,1-=n r Λ执行 若()n r r i a r ir,,3,2)(Λ++=全为零,则令)()1(r r A A =+,转5;否则转2。
计算()∑+==nri r ir r a d 12)(()()r r r r r r r r r r d c a d a c ==-=++则取,0sgn )(,1)(,1若)(,12r rr r r r a c c h +-=令()nTr nrr r r r r r r r R a a c a u ∈-=++)()(,2)(,1,,,,0,,0ΛΛ。
计算r r T r r h u A p /)(=r r rr r Tr r h u p t /=r r r r u t q -=ωT rr T r r r r p u u A A --=+ω)()1(继续。
1.2 A 的QR 分解具体算法如下所示,记)1(1-=n A A ,并记[]nn r ij r a A ⨯=)(,令I Q =1 对于1,,2,1-=n r Λ执行 1.若()n r r i a r ir ,,3,1)(Λ++=全为零,则令r r Q Q =+1r r A A =+1,转5;否则转2。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
fclose(fdA);
}
MPI_Bcast(&size,1,MPI_INT,0,MPI_COMM_WORLD);
t=(int)sqrt(p);
if (t>size)
t=size;
if(size%t!=0)
for(;;)
{
t--;
if(size%t==0)
break;
}
p=t*t;
m=size/t;
p=group_size;
if(my_rank==0)
{
starttime=MPI_Wtime();
fdA=fopen("dataIn.txt","r");
fscanf(fdA,"%d %d", &size, &N);
if(size != N)
{
puts("The input is error!");
{
v= my_rank/t;
u= my_rank%t;
for(i=0;i<m;i++)
for(j=0;j<m;j++)
b(i,j)=a(i,j);
MPI_Recv(a,m*m,MPI_FLOAT, (u*t+v), my_rank,MPI_COMM_WORLD,&status);
MPI_Send(b,m*m,MPI_FLOAT, (u*t+v),(u*t+v),MPI_COMM_WORLD);
double time1;
double time2;
int my_rank;
int p;
MPI_Status status;
FILE *fdA;
void Environment_Finalize(float *a,float *b)
{
free(a);
free(b);
}
int main(int argc, char **argv)
MPI_Finalize();
Environment_Finalize(a,b);
return(0);
}
exit(0);
}
A=(float*)malloc(floatsize*size*size);
B=(float*)malloc(floatsize*size*size);
for(i = 0; i < size; i ++)
{
for(j = 0; j < size; j ++) fscanf(fdA, "%f", A+i*size+j);
for(j=0;j<m;j++)
a(i,j)=A(i,j);
}
if (my_rank==0)
{
for(i=1;i<p;i++)
{
v=i/t;
u=i%t;
for(j=v*m;j<(v+1)*m;j++)
for(k=u*m;k<(u+1)*m;k++)
b((j%m),(k%m))=A(j,k);
MPI_Send(b,m*m,MPI_FLOAT,i,i,MPI_COMM_WORLD);
矩阵转置并行实现
1、算法描述:
若处理器个数为p,且它们的编号依次是0,1,…,p-1,则将n阶矩阵A分成p个大小为mxm的子块,m=[n/p]。p个子块组成一个子块阵列,记其中第i行第j列的子块为Aij,它含有第(i-1)m+1至第im行中的第(j-1)m+1至第jm列的所有元素。对每一处理器按行主方式赋以二维下标,记编号为i的处理器的二维下标为(v,u),其中v=[i/ ],u=imod ,将A的子块存入下标为(v,u)表示的对应处理器中。转置分为两步进行:第一步,子块转置;第二步,处理器内部转置。为了避免对应子块交换数据是处理器发生死锁,可令下三角块先向与之对应的上三角子块发送数据,然后从上三角子块接收数据;上三角子块数据先存放在缓冲区buffer中,然后从与之对应的下三角子块接收数据,最后再将缓冲区中的数据发送给下三角子块,流程图如下所示:
2、程序代码:
#include "stdio.h"
#include "stdlib.h"
#include "mpi.h"
#include "math.h"
#define E 0.0001
#define a(x,y) a[x*m+y]
#define b(x,y) b[x*m+y]
#define A(x,y) A[x*size+y]
}
}
else if(my_rank<p)
MPI_Send(a,m*m,MPI_FLOAT, 0,my_rank,MPI_COMM_WORLD);
if (my_rank==0)
{
printf("Input of file \"dataIn.txt\"\n");
printf("%d\t%d\n", size, size);
{
int i,j,k,my_rank,group_size;
float *a,*b;
int u,v;
float temp;
MPI_Init(&argc,&argv);
MPI_Comm_size(MPI_COMM_WORLD,&group_size);
MPI_Comm_rank(MPI_COMM_WORLD,&my_rank);
printf("Distribute data time = %f seconds\n",time1-starttime);
printf("Parallel compute time = %f seconds\n",time2-time1);
}
MPI_Barrier(MPI_COMM_WORLD);
for(i=0;i<size;i++)
{
for(j=0;j<size;j++) printf("%f\t",A(i,j));
printf("\n");
}
printf("\nOutput of Matrix AT\n");
for(i=0;i<size;i++)
{
for(j=0;j<size;j++) printf("%f\t",B(i,j));
}
for(i=1;i<m;i++)
for(j=0;j<i;j++)
{
temp=a(i,j);
a(i,j)=a(j,i);
a(j,i)=temp;
}
if (my_rank==0)
{
for(i=0;i<m;i++)
for(j=0;j<m;j++)
B(i,j)=a(i,j);
}
if (my_rank==0)
u= my_rank%t;
MPI_Send(a,m*m,MPI_FLOAT, (u*t+v),(u*t+v),MPI_COMM_WORLD);
MPI_Recv(a,m*m,MPI_FLOAT, (u*t+v),my_rank,MPI_COMM_WORLD,&status);
}
if ((my_rank/t)<(my_rank%t)&&my_rank<p)
{
for(i=1;i<p;i++)
{
MPI_Recv(a,m*m,MPI_FLOAT,i,i,MPI_COMM_WORLD,&status);
v= i/t;
u= i%t;
for(j=v*m;j<(v+1)*m;j++)
for(k=u*m;k<(u+1)*m;k++)
B(j,k)=a((j%m),(k%m));
a=(float *)malloc(floatsize*m*m);
b=(float *)malloc(floatsize*m*m);
if (a==NULL||b==NULL)
printf("allocate space fail!");
if (my_rank==0)
{
for(i=0;i<m;i++)
printf("\n");
}
}
time2=MPI_Wtime();
if (my_rank==0)
{
printf("\n");
printf("processor number:%d\n",p);
printf("Whole running time = %f seconds\n",time2-starttime);
}
}
else if (my_rank<p)
MPI_Recv(a,m*m,MPI_FLOAT,0,my_rank,MPI_COMM_WORLD,&status);
time1=MPI_Wtime();
if ((my_rank/t)>(my_rank%t)&&my_rank<p)