c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine mpiutl(isw) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c initialize c if(isw==1) then call MPI_INIT(ierr) call MPI_COMM_RANK(MPI_COMM_WORLD,id,ierr) call MPI_COMM_SIZE(MPI_COMM_WORLD,np,ierr) nproc=np myid=id if(nproc<=1) then mpisw=-1 else mpisw=1 endif if(rknd==rsngl) then mpiflt=MPI_REAL else if(rknd==rdble) then mpiflt=MPI_DOUBLE_PRECISION else mpiflt=MPI_LONG_DOUBLE endif if(iknd==isngl) then mpiint=MPI_INTEGER else mpiint=MPI_LONG_LONG_INT endif else c c quit c call MPI_FINALIZE(ierr) endif return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine star0(cmd) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), save :: mpijnl=0,lowera=97,lowerz=122 character(len=6) :: cmdtyp character(len=80) :: list,cmd,cmd0 common /atest3/mode,jnlsw,jnlr,jnlw,ibatch common /atest4/jcmd,cmdtyp,list common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c master c if(myid==0) then cmd0=cmd if(mode==jnlsw) mpijnl=0 if(mode==1) mpijnl=1 if(mpijnl==1.and.cmdtyp/='mpicmd') return if(mpisw<0) then if(cmdtyp=='popup '.or.cmdtyp=='file ' + .or.cmdtyp=='journl') then ii=ichar(cmd0(1:1)) if(ii>=lowera.and.ii<=lowerz) then cmd0(1:1)=char(ii-32) endif endif else if(cmdtyp=='journl') mpijnl=1 endif call MPI_BCAST(cmd0,80,MPI_CHARACTER,0,MPI_COMM_WORLD,ierr) else c c slave c if(jnlsw>=1.and.cmdtyp/='mpicmd') return call MPI_BCAST(cmd,80,MPI_CHARACTER,0,MPI_COMM_WORLD,ierr) endif c return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine expth(ip,ipath) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), dimension(6,*) :: ipath integer(kind=iknd), dimension(100) :: ip integer(kind=iknd), allocatable, dimension(:,:) :: ipath0 integer(kind=iknd), allocatable, dimension(:) :: bsize integer(kind=isngl), allocatable, dimension(:) :: + counts,offset integer(kind=iknd), dimension(4) :: id,jd common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c condense data structure c maxpth=ip(82) allocate(ipath0(6,maxpth),counts(nproc),offset(nproc), + bsize(nproc)) ishift=nproc+1 ipath(3,1)=ipath(3,nproc+2) ipath(4,1)=ipath(4,nproc+2) ipath(1,1)=ipath(1,nproc+2)-ishift ipath(2,1)=ipath(2,nproc+2)-ishift do j=ipath(1,1),ipath(2,1) do k=1,6 ipath(k,j)=ipath(k,j+ishift) enddo if(ipath(2,j)>0) ipath(2,j)=ipath(2,j)-ishift enddo len=ipath(2,1)*6 c call MPI_ALLGATHER(len,1,mpiint,bsize,1, + mpiint,MPI_COMM_WORLD,ierr) c offset(1)=0 do irgn=1,nproc-1 offset(irgn+1)=offset(irgn)+bsize(irgn) counts(irgn)=bsize(irgn) enddo counts(nproc)=bsize(nproc) isum=(offset(nproc)+counts(nproc))/6+nproc+2 if(isum>maxpth) then ip(25)=72 go to 20 endif c call MPI_ALLGATHERV(ipath,len,mpiint, + ipath0,counts,offset,mpiint,MPI_COMM_WORLD,ierr) c ipath(3,nproc+2)=0 ipath(4,nproc+2)=0 ipath(1,nproc+2)=nproc+3 ipath(2,nproc+2)=nproc+2 do irgn=1,nproc ii=offset(irgn)/6 c ipath(3,nproc+2)=max(ipath(3,nproc+2),ipath0(3,1+ii)) c n0=ipath(4,nproc+2) ipath(3,irgn)=n0+1 ipath(4,irgn)=n0+ipath0(4,1+ii) ipath(4,nproc+2)=ipath(4,irgn) c ipath(1,irgn)=ipath(2,nproc+2)+1 jb0=ipath(1,irgn)-ipath0(1,1+ii) ipath(2,irgn)=ipath0(2,1+ii)+jb0 ipath(2,nproc+2)=ipath(2,irgn) c do j=ipath(1,irgn),ipath(2,irgn) do k=1,6 ipath(k,j)=ipath0(k,j-jb0+ii) enddo if(ipath(2,j)>0) then ipath(2,j)=ipath(2,j)+jb0 endif ipath(3,j)=ipath(3,j)+n0 ipath(4,j)=ipath(4,j)+n0 if(ipath(5,j)>0) then ipath(5,j)=ipath(5,j)+n0 else if(ipath(5,j)<0) then ipath(5,j)=ipath(5,j)-n0 endif enddo enddo c ip(71)=ipath(4,nproc+2) ip(72)=ipath(2,nproc+2) c c global dimensions c id(1)=ip(27) id(2)=ip(28) id(3)=ip(29) id(4)=ip(30) call MPI_ALLREDUCE(id,jd,4,mpiint,MPI_SUM,MPI_COMM_WORLD,ierr) ip(37)=jd(1) ip(38)=jd(2) ip(39)=jd(3) ip(40)=jd(4) c 20 deallocate(ipath0,counts,offset,bsize) return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine exbdy(ipath,ir0,map,gf,nn,num) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), dimension(6,*) :: ipath integer(kind=iknd), dimension(*) :: ir0,map integer(kind=isngl), allocatable, dimension(:) :: + counts,offset real(kind=rknd), allocatable, dimension(:) :: vin,vout real(kind=rknd), dimension(nn,*) :: gf common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c communicate interface data c allocate(vin(nn*num),vout(nn*num),counts(nproc),offset(nproc)) do irgn=1,nproc ii=ipath(3,irgn)-1 nvdd=ipath(4,irgn)-ii counts(irgn)=num*nvdd offset(irgn)=num*ii enddo c irgn=myid+1 ii=ipath(3,irgn)-1 nvdd=ipath(4,irgn)-ii len=counts(irgn) do j=1,num jj=(j-1)*nvdd do i=1,nvdd vin(i+jj)=gf(i+ii,j) enddo enddo c call MPI_ALLGATHERV(vin,len,mpiflt, + vout,counts,offset,mpiflt,MPI_COMM_WORLD,ierr) c do irgn=1,nproc ii=ipath(3,irgn)-1 nvdd=ipath(4,irgn)-ii do j=1,num jj=(j-1)*nvdd+offset(irgn) do i=1,nvdd gf(i+ii,j)=vout(i+jj) enddo enddo enddo c c interpolate to mesh on this processor c lipath=ipath(2,nproc+1) call intrpi(lipath,ipath,ir0,map,nn,num,gf,vin) c deallocate(vin,vout,counts,offset) c return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine pl2ip(t,num) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" real(kind=rknd), dimension(*) :: t real(kind=rknd), dimension(100) :: temp common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c compute all inner products/norms c do i=1,num temp(i)=t(i) enddo c call MPI_ALLREDUCE(temp,t,num,mpiflt, + MPI_SUM,MPI_COMM_WORLD,ierr) return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine bcast(vx,vy,sf,ibndry,itnode,itdof, + ip,rp,sp,iu,ru,su,gf,e) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), dimension(7,*) :: ibndry integer(kind=iknd), dimension(5,*) :: itnode integer(kind=iknd), dimension(100) :: ip,iu integer(kind=iknd), dimension(8,*) :: itdof real(kind=rknd), dimension(*) :: vx,vy,gf,sf,e real(kind=rknd), dimension(100) :: rp,ru character(len=1), allocatable, dimension(:) :: buff character(len=80), dimension(100) :: sp,su common /pltmg7/time(3,50),hist(22,30) common /pltmg6/path(101,6) common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c c inital fan-out of coarse grid data c c master node c if(myid==0) then ntf=ip(1) nvf=ip(2) nbf=ip(3) ndf=ip(4) maxt=ip(83) maxd=ip(85) ngf=ip(77) c leni=(200+13*ntf+7*nbf)*iknd lenr=(200+ngf*ndf+2*nvf+2*nbf+2*ntf+606+660+150)*rknd lens=200*80 len=leni+lenr+lens+16 c c integers c allocate(buff(len)) lenz=len len=0 call MPI_PACK(ip,100,mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(iu,100,mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(itnode,5*ntf,mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(ibndry,7*nbf,mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(itdof,8*ntf,mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) c c floats c call MPI_PACK(rp,100,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(ru,100,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(sf,2*nbf,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(vx,nvf,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(vy,nvf,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) do k=1,ngf iptr=1+(k-1)*maxd call MPI_PACK(gf(iptr),ndf,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) enddo do k=1,2 iptr=1+(k-1)*maxt call MPI_PACK(e(iptr),ntf,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) enddo call MPI_PACK(path,606,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(hist,660,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(time,150,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) c c characters c call MPI_PACK(sp,8000,MPI_CHARACTER, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(su,8000,MPI_CHARACTER, + buff,lenz,len,MPI_COMM_WORLD,ierr) c c the big broadcast c call MPI_BCAST(len,1,mpiint,0,MPI_COMM_WORLD,ierr) call MPI_BCAST(buff,len,MPI_PACKED,0,MPI_COMM_WORLD,ierr) deallocate(buff) c c slave node c else call MPI_BCAST(len,1,mpiint,0,MPI_COMM_WORLD,ierr) allocate(buff(len)) call MPI_BCAST(buff,len,MPI_PACKED,0,MPI_COMM_WORLD,ierr) c c integers c lenz=len len=0 c call MPI_UNPACK(buff,lenz,len,ip,100, + mpiint,MPI_COMM_WORLD,ierr) ip(48)=mpisw ip(49)=nproc ip(50)=myid+1 ntf=ip(1) nvf=ip(2) nbf=ip(3) ndf=ip(4) maxt=ip(83) maxd=ip(85) ngf=ip(77) c call MPI_UNPACK(buff,lenz,len,iu,100, + mpiint,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,itnode,5*ntf, + mpiint,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,ibndry,7*nbf, + mpiint,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,itdof,8*ntf, + mpiint,MPI_COMM_WORLD,ierr) c c floats c call MPI_UNPACK(buff,lenz,len,rp,100, + mpiflt,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,ru,100, + mpiflt,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,sf,2*nbf, + mpiflt,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,vx,nvf, + mpiflt,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,vy,nvf, + mpiflt,MPI_COMM_WORLD,ierr) do k=1,ngf iptr=1+(k-1)*maxd call MPI_UNPACK(buff,lenz,len,gf(iptr),ndf, + mpiflt,MPI_COMM_WORLD,ierr) enddo do k=1,2 iptr=1+(k-1)*maxt call MPI_UNPACK(buff,lenz,len,e(iptr),ntf, + mpiflt,MPI_COMM_WORLD,ierr) enddo call MPI_UNPACK(buff,lenz,len,path,606, + mpiflt,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,hist,660, + mpiflt,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,time,150, + mpiflt,MPI_COMM_WORLD,ierr) c c characters c call MPI_UNPACK(buff,lenz,len,sp,8000, + MPI_CHARACTER,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,su,8000, + MPI_CHARACTER,MPI_COMM_WORLD,ierr) deallocate(buff) endif c c this uses many bcasts but no buffer storage c c call MPI_BCAST(ip,100,mpiint,0,MPI_COMM_WORLD,ierr) c ip(48)=mpisw c ip(49)=nproc c ip(50)=myid+1 c ntf=ip(1) c nvf=ip(2) c nbf=ip(3) c ndf=ip(4) c maxt=ip(83) c maxd=ip(85) c ngf=ip(77) c c call MPI_BCAST(iu,100,mpiint,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(itnode,5*ntf,mpiint,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(ibndry,7*nbf,mpiint,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(itdof,8*ntf,mpiint,0,MPI_COMM_WORLD,ierr) c c floats c c call MPI_BCAST(rp,100,mpiflt,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(ru,100,mpiflt,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(sf,2*nbf,mpiflt,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(vx,nvf,mpiflt,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(vy,nvf,mpiflt,0,MPI_COMM_WORLD,ierr) c do k=1,ngf c iptr=1+(k-1)*maxd c call MPI_BCAST(gf(iptr),ndf,mpiflt,0,MPI_COMM_WORLD,ierr) c enddo c do k=1,2 c iptr=1+(k-1)*maxt c call MPI_BCAST(e(iptr),ntf,mpiflt,0,MPI_COMM_WORLD,ierr) c enddo c call MPI_BCAST(path,606,mpiflt,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(hist,660,mpiflt,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(time,150,mpiflt,0,MPI_COMM_WORLD,ierr) c c characters c c call MPI_BCAST(sp,8000,MPI_CHARACTER,0,MPI_COMM_WORLD,ierr) c call MPI_BCAST(su,8000,MPI_CHARACTER,0,MPI_COMM_WORLD,ierr) return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine exstat(rp,pstat) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" real(kind=rknd), dimension(100) :: rp real(kind=rknd), dimension(4,*) :: pstat real(kind=rknd), dimension(4) :: sendb common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c exchange load balance statistics c sendb(1)=rp(95) sendb(2)=rp(96) sendb(3)=rp(97) sendb(4)=rp(98) c call MPI_ALLGATHER(sendb,4,mpiflt, + pstat,4,mpiflt,MPI_COMM_WORLD,ierr) c return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine exbox(bmin,bmax,num) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" real(kind=rknd), dimension(5) :: tin,tout real(kind=rknd), dimension(*) :: bmin,bmax common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c compute the global bounding box c c this routine doesnt work for quad (june 2012). a fix is c make tin, tout double and do the min/max in double precision. c sub exstep and exqual will not work either) c do i=1,num tin(i)=bmin(i) enddo call MPI_ALLREDUCE(tin,tout,num,mpiflt, + MPI_MIN,MPI_COMM_WORLD,ierr) c do i=1,num bmin(i)=tout(i) tin(i)=bmax(i) enddo c call MPI_ALLREDUCE(tin,tout,num,mpiflt, + MPI_MAX,MPI_COMM_WORLD,ierr) do i=1,num bmax(i)=tout(i) enddo c return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine exdist(kdist,num) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), dimension(*) :: kdist integer(kind=iknd), dimension(22) :: itemp common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c compute all inner products/norms c do i=1,num itemp(i)=kdist(i) enddo c call MPI_ALLREDUCE(itemp,kdist,num,mpiint, + MPI_SUM,MPI_COMM_WORLD,ierr) return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine exqual(num,val) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), dimension(4) :: num,itemp real(kind=rknd), dimension(2) :: temp,val common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c compute global mesh quality statistics (binit) c do i=1,4 itemp(i)=num(i) enddo call MPI_ALLREDUCE(itemp,num,4,mpiint, + MPI_SUM,MPI_COMM_WORLD,ierr) c do i=1,2 temp(i)=val(i) enddo call MPI_ALLREDUCE(temp,val,1,mpiflt, + MPI_MIN,MPI_COMM_WORLD,ierr) call MPI_ALLREDUCE(temp(2),val(2),1,mpiflt, + MPI_SUM,MPI_COMM_WORLD,ierr) c return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine exibox(ibmin,ibmax,num) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), dimension(5) :: itemp integer(kind=iknd), dimension(*) :: ibmin,ibmax common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c compute min and max color (binit) c do i=1,num itemp(i)=ibmin(i) enddo call MPI_ALLREDUCE(itemp,ibmin,num,mpiint, + MPI_MIN,MPI_COMM_WORLD,ierr) c do i=1,num itemp(i)=ibmax(i) enddo call MPI_ALLREDUCE(itemp,ibmax,num,mpiint, + MPI_MAX,MPI_COMM_WORLD,ierr) c return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine exflag(iflag) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c global dimensions c itemp=iflag call MPI_ALLREDUCE(itemp,iflag,1,mpiint, + MPI_MAX,MPI_COMM_WORLD,ierr) c return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine exstep(step) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c compute the maximum step (tpickd) c temp=step call MPI_ALLREDUCE(temp,step,1,mpiflt, + MPI_MIN,MPI_COMM_WORLD,ierr) c return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine extim(atime,ptime) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), save :: len=50 real(kind=rknd), dimension(*) :: ptime real(kind=rknd), dimension(3,50) :: temp real(kind=rknd), dimension(3,*) :: atime real(kind=rknd), dimension(2) :: sendb common /pltmg7/time(3,50),hist(22,30) common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c exchange times c sum=0.0e0_rknd do i=1,len sum=sum+time(2,i) temp(1,i)=time(1,i) temp(2,i)=time(2,i) temp(3,i)=time(3,i) enddo c sendb(1)=sum call MPI_ALLGATHER(sendb,1,mpiflt, + ptime,1,mpiflt,MPI_COMM_WORLD,ierr) c num=3*len call MPI_ALLREDUCE(temp,atime,num,mpiflt, + MPI_SUM,MPI_COMM_WORLD,ierr) c do i=1,len atime(1,i)=atime(1,i)/real(nproc,rknd) atime(2,i)=atime(2,i)/real(nproc,rknd) atime(3,i)=atime(3,i)/real(nproc,rknd) enddo return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine exsze(jd,isw) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), dimension(3) :: id,jd common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c global dimensions c id(1)=jd(1) id(2)=jd(2) id(3)=jd(3) if(isw==0) then call MPI_ALLREDUCE(id,jd,3,mpiint, + MPI_MAX,MPI_COMM_WORLD,ierr) else call MPI_ALLREDUCE(id,jd,3,mpiint, + MPI_SUM,MPI_COMM_WORLD,ierr) endif c return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine glbpix(vx,vy,sf,ibndry,itnode,icolor, + iuvptr,ut,vt,jp,isw) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) include "mpif.h" integer(kind=iknd), dimension(7,*) :: ibndry integer(kind=iknd), dimension(5,*) :: itnode integer(kind=iknd), dimension(25) :: jp integer(kind=iknd), dimension(5) :: ib integer(kind=iknd), dimension(*) :: icolor integer(kind=isngl), dimension(MPI_STATUS_SIZE) :: status integer(kind=iknd), dimension(2,*) :: iuvptr real(kind=rknd), dimension(2,*) :: sf real(kind=rknd), dimension(*) :: vx,vy,ut,vt character(len=1), allocatable, dimension(:) :: buff common /atest6/nproc,myid,mpisw,mpiint,mpiflt cy c collect global graph on master process c itag=4 c if(myid==0) then c c integers c do i=1,nproc-1 ntf=jp(1) nvf=jp(2) nbf=jp(3) ndf=jp(24) if(isw/=0) isv=iuvptr(1,ntf+1)-1 c call MPI_RECV(len,1,mpiint,MPI_ANY_SOURCE, + itag,MPI_COMM_WORLD,status,ierr) allocate(buff(len)) idsrc=status(MPI_SOURCE) call MPI_RECV(buff,len,MPI_PACKED, + idsrc,itag,MPI_COMM_WORLD,status,ierr) c lenz=len len=0 call MPI_UNPACK(buff,lenz,len, + ib,5,mpiint,MPI_COMM_WORLD,ierr) c ntf0=ib(1) nvf0=ib(2) nbf0=ib(3) ndf0=ib(4) nuv0=ib(5) c call MPI_UNPACK(buff,lenz,len,itnode(1,ntf+1),5*ntf0, + mpiint,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len,ibndry(1,nbf+1),7*nbf0, + mpiint,MPI_COMM_WORLD,ierr) if(isw==0) then call MPI_UNPACK(buff,lenz,len,icolor(ntf+1),ntf0, + mpiint,MPI_COMM_WORLD,ierr) else call MPI_UNPACK(buff,lenz,len,iuvptr(1,ntf+1), + 2*(ntf0+1),mpiint,MPI_COMM_WORLD,ierr) endif c c update itnode and ibndry c do k=ntf+1,ntf+ntf0 do j=1,3 itnode(j,k)=itnode(j,k)+nvf enddo enddo do k=nbf+1,nbf+nbf0 do j=1,2 ibndry(j,k)=ibndry(j,k)+nvf enddo if(ibndry(4,k)<0) then ibndry(4,k)=ibndry(4,k)-nbf endif enddo c c fixup iuvptr c if(isw/=0) then do k=ntf+1,ntf+ntf0+1 iuvptr(1,k)=iuvptr(1,k)+isv enddo endif c c floats c call MPI_UNPACK(buff,lenz,len, + sf(1,nbf+1),2*nbf0,mpiflt,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len, + vx(nvf+1),nvf0,mpiflt,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len, + vy(nvf+1),nvf0,mpiflt,MPI_COMM_WORLD,ierr) if(isw/=0) then call MPI_UNPACK(buff,lenz,len, + ut(isv+1),nuv0,mpiflt,MPI_COMM_WORLD,ierr) call MPI_UNPACK(buff,lenz,len, + vt(isv+1),nuv0,mpiflt,MPI_COMM_WORLD,ierr) endif c jp(1)=ntf+ntf0 jp(2)=nvf+nvf0 jp(3)=nbf+nbf0 jp(24)=ndf+ndf0 c deallocate(buff) c enddo else c c slave node c c integers c ntf=jp(1) nvf=jp(2) nbf=jp(3) ndf=jp(24) nuv=0 if(isw/=0) nuv=iuvptr(1,ntf+1)-iuvptr(1,1) c lenz=(5*ntf+7*nbf+5)*iknd+2*(nbf+nvf)*rknd if(isw==0) lenz=lenz+ntf*iknd if(isw/=0) lenz=lenz+(ntf+1)*2*iknd+rknd*2*nuv c allocate(buff(lenz)) c ib(1)=ntf ib(2)=nvf ib(3)=nbf ib(4)=ndf ib(5)=nuv len=0 call MPI_PACK(ib,5,mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(itnode,5*ntf,mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(ibndry,7*nbf,mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) if(isw==0) then call MPI_PACK(icolor,ntf,mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) else call MPI_PACK(iuvptr,2*(ntf+1),mpiint, + buff,lenz,len,MPI_COMM_WORLD,ierr) endif c c floats c call MPI_PACK(sf,2*nbf,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(vx,nvf,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(vy,nvf,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) if(isw/=0) then call MPI_PACK(ut,nuv,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) call MPI_PACK(vt,nuv,mpiflt, + buff,lenz,len,MPI_COMM_WORLD,ierr) endif c call MPI_SEND(len,1,mpiint, + 0,itag,MPI_COMM_WORLD,ierr) call MPI_SEND(buff,len,MPI_PACKED, + 0,itag,MPI_COMM_WORLD,ierr) c deallocate(buff) c endif c return end