exit
(1);
}
intit=0,hk=20,k=hk/2,ext=1;//相关变量的声明和初始化
realx0,y0,x1,y1,Unorm,errU,nu=1,h=1.0/hk,extension=ext*h,itu=1.0;
real[int]relerrU
(2);
mesh[int]Th(mpisize);//子区域网格
if(mpirank==0){
x0=0;x1=1;y0=0;y1=0.5+extension;//定义子区域1
Th[0]=square(hk,k+ext,[x0+(x1-x0)*x,y0+(y1-y0)*y]);//生成子区域1的网格
}
else{
x0=0;x1=1;y0=0.5-extension;y1=1.0;//定义子区域2
Th[1]=square(hk,k+ext,[x0+(x1-x0)*x,y0+(y1-y0)*y]);//生成子区域2的网格
}
broadcast(processor(0),Th[0]);//广播网格
broadcast(processor
(1),Th[1]);
fespaceUh01(Th[0],P2);Uh01uh11,uh12,vh11,vh12,olduh11,olduh12;//有限元空间
fespacePh01(Th[0],P1);Ph01ph1,qh1,oldph1;
fespaceUh02(Th[1],P2);Uh02uh21,uh22,vh21,vh22,olduh21,olduh22;
fespacePh02(Th[1],P1);Ph02ph2,qh2,oldph2;
problemstokes1([uh11,uh12,ph1],[vh11,vh12,qh1],init=it)=//定义子区域1的问题
int2d(Th[0])(nu*(dx(uh11)*dx(vh11)+dy(uh11)*dy(vh11)
+dx(uh12)*dx(vh12)+dy(uh12)*dy(vh12))
-ph1*(dx(vh11)+dy(vh12))+qh1*(dx(uh11)+dy(uh12)))
+on(1,2,4,uh11=0,uh12=0)+on(3,uh11=olduh21,uh12=olduh22,ph1=oldph2);
problemstokes2([uh21,uh22,ph2],[vh21,vh22,qh2],init=it)=//定义子区域2的问题
int2d(Th[1])(nu*(dx(uh21)*dx(vh21)+dy(uh21)*dy(vh21)
+dx(uh22)*dx(vh22)+dy(uh22)*dy(vh22))
-ph2*(dx(vh21)+dy(vh22))+qh2*(dx(uh21)+dy(uh22)))
+on(2,4,uh21=0,uh22=0)+on(3,uh21=1,uh22=0)
+on(1,uh21=olduh11,uh22=olduh12,ph2=oldph1);
uh11=0;uh12=0;ph1=0;uh21=0;uh22=0;ph2=0;//迭代初值
while(it<100&&itu>=1.0e-3){//Schwarz迭代
it++;
if(mpirank==0){olduh11=uh11;olduh12=uh12;oldph1=ph1;}//进程0存储计算结果
if(mpirank==1){olduh21=uh21;olduh22=uh22;oldph2=ph2;}//进程1存储计算结果
broadcast(processor(0),olduh11[]);//广播计算结果
broadcast(processor(0),olduh12[]);
broadcast(processor(0),oldph1[]);
broadcast(processor
(1),olduh21[]);
broadcast(processor
(1),olduh22[]);
broadcast(processor
(1),oldph2[]);
if(mpirank==0){
stokes1;//子区域1上的有限元计算
Unorm=sqrt(int2d(Th[0])(uh11^2+uh12^2));//计算相邻两次迭代解的误差
errU=sqrt(int2d(Th[0])((uh11-olduh11)^2+(uh12-olduh12)^2));
if(Unorm<1.0e-30)Unorm=1.0e-30;relerrU[0]=errU/Unorm;
}
if(mpirank==1){
stokes2;//子区域2上的有限元计算
Unorm=sqrt(int2d(Th[1])(uh21^2+uh22^2));//计算相邻两次迭代解的误差
errU=sqrt(int2d(Th[1])((uh21-olduh21)^2+(uh22-olduh22)^2));
if(Unorm<1.0e-30)Unorm=1.0e-30;relerrU[1]=errU/Unorm;
}
broadcast(processor(0),relerrU[0]);//广播计算所得误差
broadcast(processor
(1),relerrU[1]);
itu=relerrU.max;
}
图2是上述并行程序的计算结果插值到互不重叠子区域的等值线图,其中Schwarz迭代15次。
由以上并行程序可以看出,基于MPI+FreeFem++的有限元并行编程简单而高效,程序员无需将太多精力花在程序代码的编写上,从而集中于算法的设计与分析,实现相关问题的快速、高效有限元并行计算。
5结束语
自20世纪50年代提出有限元法以来,有限元法已经成为当今工程分析和计算中的不可缺少的最重要的工具之一。
本文讨论了基于MPI+FreeFem++的有限元并行计算环境的构建,以及在该环境下有限元并行程序的编写、编译与运行等过程,并通过具体的编程实例表明该并行编程环境的简单和高效。
程序员在该并行编程环境下,无需编写冗长的程序代码,从而将主要精力用于并行算法的设计与分析,快速而高效地实现相关问题的有限元并行计算。
参考文献:
[1]李开泰,黄艾香,黄庆怀.有限元方法及其应用[M].北京:
科学出版社,2006.
[2]HECHTF,PIRONNEAUO,LEHYARICA,OHTSUKAK.FreeFem++[EB/OL].http:
//www.freefem.org
[3]SHANGYUEQIANG,HEYINNIAN.FourieranalysisofSchwarzdomaindecompositionmethodsforthebiharmonicequation[J].AppliedMathematicsandMechanics,2009(9).
[4]TOSELLIA,WIDLUNDO.DomainDecompositionMethods:
AlgorithmsandTheory[M].Springer,Berlin,2005.
ParallelFiniteElementComputationsBasedonMPIandFreeFem++
Abstract:
Finiteelementmethodsareatypeofflexibleandhighlyefficientmethodsfornumericalsolutionofpartialdifferentialequations.Theyareoneoftheimportantandindispensabletoolsinengineeringanalysisandcomputations.Asthepriceofparallelcomputersquicklydropswiththedevelopmentofcomputertechnology,parallelfiniteelementmethodsnowadaysattractpopularattentionsinbothacademicandengineeringfields.ThispaperdiscussestheconstructionofparallelenvironmentbasedonMPIandFreeFem++forfiniteelementcomputations,statesthewriting,compilingandrunningofaparallelfiniteelementprogramintheenvironment.Finally,anexampleisgiventodemonstratethesimplicityandhighefficiencyofparallelfiniteelementcomputationsbasedonMPIandFreeFem++.
KeyWords:
FiniteElementMethod;ParallelComputing;MPI;FreeFem++