c*********************** problem name: square ************************ c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine a1xy(x,y,u,ux,uy,rl,itag,values) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) real(kind=rknd), dimension(*) :: values character(len=80) :: su common /val0/k0,ku,kx,ky,kl,kuu,kux,kuy,kul, + kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll common /atest2/ibc(4),icont,iu(95),a1x,a1y,a1u,a2x,a2y,a2u, + bux,buy,cu0,cu1,cu2,cahn,cexp,cir,csin,du0,du1,eu0,f0, 1 ru(81),su(100) cy call setrl(rl) values(k0)=a1x*ux+a1y*uy+a1u*u values(ku)=a1u values(kx)=a1x values(ky)=a1y if(icont==1) then values(kl)=ux values(klx)=1.0e0_rknd else if(icont==2) then values(kl)=uy values(kly)=1.0e0_rknd else if(icont==3) then values(kl)=u values(klu)=1.0e0_rknd endif return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine a2xy(x,y,u,ux,uy,rl,itag,values) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) real(kind=rknd), dimension(*) :: values character(len=80) :: su common /val0/k0,ku,kx,ky,kl,kuu,kux,kuy,kul, + kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll common /atest2/ibc(4),icont,iu(95),a1x,a1y,a1u,a2x,a2y,a2u, + bux,buy,cu0,cu1,cu2,cahn,cexp,cir,csin,du0,du1,eu0,f0, 1 ru(81),su(100) cy call setrl(rl) values(k0)=a2x*ux+a2y*uy+a2u*u values(ku)=a2u values(kx)=a2x values(ky)=a2y if(icont==4) then values(kl)=ux values(klx)=1.0e0_rknd else if(icont==5) then values(kl)=uy values(kly)=1.0e0_rknd else if(icont==6) then values(kl)=u values(klu)=1.0e0_rknd endif return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine fxy(x,y,u,ux,uy,rl,itag,values) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) real(kind=rknd), dimension(*) :: values character(len=80) :: su common /val0/k0,ku,kx,ky,kl,kuu,kux,kuy,kul, + kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll common /atest2/ibc(4),icont,iu(95),a1x,a1y,a1u,a2x,a2y,a2u, + bux,buy,cu0,cu1,cu2,cahn,cexp,cir,csin,du0,du1,eu0,f0, 1 ru(81),su(100) cy call setrl(rl) values(k0)= - bux*ux - buy*uy + - cir*(ux*(y-0.5e0_rknd)-uy*(x-0.5e0_rknd)) 1 - cu0 - cu1*u - cu2*u**2 2 - cahn*u*(1.0e0_rknd-u**2) -f0*(y-x) values(ku)= - cu1 - cu2*u*2.0e0_rknd + - cahn*(1.0e0_rknd-3.0e0_rknd*u**2) values(kuu)= - cu2*2.0e0_rknd + cu3*6.0e0_rknd*u values(kx) = - bux - cir*(y-0.5e0_rknd) values(ky) = - buy + cir*(x-0.5e0_rknd) c if(cexp/=0.0e0_rknd.or.icont==13) then ee=exp(u) values(k0)=values(k0) - cexp*ee values(ku)=values(ku) - cexp*ee values(kuu)=values(kuu) - cexp*ee endif if(csin/=0.0e0_rknd.or.icont==15) then ss=sin(u) cc=cos(u) values(k0)=values(k0) - csin*ss values(ku)=values(ku) - csin*cc values(kuu)=values(kuu) + csin*ss endif c if(icont==7) then values(kl) = - ux values(klx) = - 1.0e0_rknd else if(icont==8) then values(kl) = - uy values(kly) = - 1.0e0_rknd else if(icont==9) then values(kl) = - 1.0e0_rknd else if(icont==10) then values(kl) = - u values(klu) = - 1.0e0_rknd else if(icont==11) then values(kl) = - u**2 values(klu) = - 2.0e0_rknd*u else if(icont==12) then values(kl) = - u*(1.0e0_rknd-u**2) values(klu) = - (1.0e0_rknd-3.0e0_rknd*u**2) else if(icont==13) then values(kl) = - ee values(klu)= - ee else if(icont==14) then values(kl) = - (ux*(y-0.5e0_rknd)-uy*(x-0.5e0_rknd)) values(klx) = - (y-0.5e0_rknd) values(kly) = - (x-0.5e0_rknd) else if(icont==15) then values(kl) = - ss values(klu) = - cc endif return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine gnxy(x,y,u,rl,itag,values) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) real(kind=rknd), dimension(*) :: values character(len=80) :: su common /val1/k0,ku,kl,kuu,kul,klu,kll common /atest2/ibc(4),icont,iu(95),a1x,a1y,a1u,a2x,a2y,a2u, + bux,buy,cu0,cu1,cu2,cahn,cexp,cir,csin,du0,du1,eu0,f0, 1 ru(81),su(100) cy call setrl(rl) values(k0)= - du0 - du1*u values(ku)= - du1 if(icont==16) then values(kl)=-1.0e0_rknd else if(icont==17) then values(kl)=-u values(klu)=-1.0e0_rknd endif return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine gdxy(x,y,rl,itag,values) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) real(kind=rknd), dimension(*) :: values character(len=80) :: su common /val2/k0,kl,kll,klb,kub,kic,kim,kil common /atest2/ibc(4),icont,iu(95),a1x,a1y,a1u,a2x,a2y,a2u, + bux,buy,cu0,cu1,cu2,cahn,cexp,cir,csin,du0,du1,eu0,f0, 1 ru(81),su(100) cy call setrl(rl) values(k0)=-eu0 if(icont==18) values(kl)=-1.0e0_rknd return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine p1xy(x,y,u,ux,uy,rl,itag,values) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) real(kind=rknd), dimension(*) :: values common /val0/k0,ku,kx,ky,kl,kuu,kux,kuy,kul, + kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll cy values(k0)=u*u values(ku)=2.0e0_rknd*u values(kuu)=2.0e0_rknd return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine p2xy(x,y,dx,dy,u,ux,uy,rl,itag,jtag,values) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) real(kind=rknd), dimension(*) :: values common /val0/k0,ku,kx,ky,kl,kuu,kux,kuy,kul, + kxu,kxx,kxy,kxl,kyu,kyx,kyy,kyl,klu,klx,kly,kll cy return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine qxy(x,y,u,ux,uy,rl,itag,values) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) real(kind=rknd), dimension(*) :: values common /val3/kf,kf1,kf2,kad cy r=(x-0.5e0_rknd)**2+(y-0.5e0_rknd)**2 s=64.0e0_rknd d=1.0e0_rknd/16.0e0_rknd t=tanh(s*(d-r)) values(kf)=t values(kad)=t return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine sxy(rl,s,itag,values) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) real(kind=rknd), dimension(*) :: values common /val4/jx,jy,jxs,jys,jxl,jyl,jxss,jyss,jxll,jyll, + jxsl,jysl,jxls,jyls cy return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine setrl(rl) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) character(len=80) :: su common /atest2/ibc(4),icont,iu(95),ru(100),su(100) cy if(icont>0) ru(icont)=rl return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine usrcmd(vx,vy,sf,itnode,ibndry,ip,rp,sp,iu,ru,su) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) integer(kind=iknd), dimension(5,*) :: itnode integer(kind=iknd), dimension(7,*) :: ibndry integer(kind=iknd), dimension(100) :: ip,iu integer(kind=iknd), dimension(5) :: isv integer(kind=iknd), save :: len real(kind=rknd), dimension(*) :: vx,vy real(kind=rknd), dimension(2,*) :: sf real(kind=rknd), dimension(100) :: rp,ru character(len=80), dimension(100) :: sp,su character(len=80), save, dimension(80) :: file cy data len/56/ data (file(i),i= 1, 10)/ + 'n i= 1,n=left, a= l,t=i', 1 'n i= 2,n=top, a= t,t=i', 2 'n i= 3,n=right, a= r,t=i', 3 'n i= 4,n=bottom,a= b,t=i', 4 'n i= 5,n=icont, a= i,t=i', 5 'n i= 1,n=a1x, a=x1,t=r', 6 'n i= 2,n=a1y, a=y1,t=r', 7 'n i= 3,n=a1u, a=u1,t=r', 8 'n i= 4,n=a2x, a=x2,t=r', 9 'n i= 5,n=a2y, a=y2,t=r'/ data (file(i),i= 11, 20)/ + 'n i= 6,n=a2u, a=u2,t=r', 1 'n i= 7,n=bux, a=bx,t=r', 2 'n i= 8,n=buy, a=by,t=r', 3 'n i= 9,n=cu0, a=c0,t=r', 4 'n i=10,n=cu1, a=c1,t=r', 5 'n i=11,n=cu2, a=c2,t=r', 6 'n i=12,n=cahn, a=ca,t=r', 7 'n i=13,n=cexp, a=cx,t=r', 8 'n i=14,n=cir, a=cr,t=r', 9 'n i=15,n=csin, a=cs,t=r'/ data (file(i),i= 21, 30)/ + 'n i=16,n=du0, a=d0,t=r', 1 'n i=17,n=du1, a=d1,t=r', 2 'n i=18,n=eu0, a=e0,t=r', 3 'n i=19,n=f0, a=f0,t=r', 4 's n=left ,v=1,l="neumann"', 5 's n=left ,v=2,l="dirichlet"', 6 's n=left ,v=0,l="periodic"', 7 's n=top ,v=1,l="neumann"', 8 's n=top ,v=2,l="dirichlet"', 9 's n=top ,v=0,l="periodic"'/ data (file(i),i= 31, 40)/ + 's n=right ,v=1,l="neumann"', 1 's n=right ,v=2,l="dirichlet"', 2 's n=right ,v=0,l="periodic"', 3 's n=bottom,v=1,l="neumann"', 4 's n=bottom,v=2,l="dirichlet"', 5 's n=bottom,v=0,l="periodic"', 6 's n=icont ,v= 0,l="none"', 7 's n=icont ,v= 1,l="a1x"', 8 's n=icont ,v= 2,l="a1y"', 9 's n=icont ,v= 3,l="a1u"'/ data (file(i),i= 41, 50)/ + 's n=icont ,v= 4,l="a2x"', 1 's n=icont ,v= 5,l="a2y"', 2 's n=icont ,v= 6,l="a2u"', 3 's n=icont ,v= 7,l="bux"', 4 's n=icont ,v= 8,l="buy"', 5 's n=icont ,v= 9,l="cu0"', 6 's n=icont ,v=10,l="cu1"', 7 's n=icont ,v=11,l="cu2"', 8 's n=icont ,v=12,l="cahn"', 9 's n=icont ,v=13,l="cexp"'/ data (file(i),i= 51, 56)/ + 's n=icont ,v=14,l="cir"', 1 's n=icont ,v=15,l="csin"', 2 's n=icont ,v=16,l="du0"', 3 's n=icont ,v=17,l="du1"', 4 's n=icont ,v=18,l="eu0"', 5 's n=icont ,v=19,l="f0"'/ c c c save integer parameters c do i=1,5 isv(i)=iu(i) enddo c c enter input mode c call usrset(file,len,iu,ru,su) c c if any of the integer parameters have changed, call gdata c do i=1,4 if(iu(i)<0.or.iu(i)>2) iu(i)=2 if(iu(i)/=isv(i)) ip(41)=-1 enddo if(iu(5)==0) then if(isv(5)/=0) ip(41)=-1 else if(isv(5)==0) ip(41)=-1 rp(1)=ru(iu(5)) endif return end c----------------------------------------------------------------------- c c piecewise lagrange triangle multi grid package c c edition 11.0 - - - june, 2012 c c----------------------------------------------------------------------- subroutine gdata(vx,vy,sf,itnode,ibndry,ip,rp,sp,iu,ru,su,sxy) cx use mthdef implicit real(kind=rknd) (a-h,o-z) implicit integer(kind=iknd) (i-n) integer(kind=iknd), dimension(5,*) :: itnode integer(kind=iknd), dimension(7,*) :: ibndry integer(kind=iknd), dimension(100) :: ip,iu integer(kind=iknd), save, dimension(9) :: ix,iy integer(kind=iknd), save :: ntf,nvf,nbf,ispd,ifirst real(kind=rknd), dimension(*) :: vx,vy real(kind=rknd), dimension(2,*) :: sf real(kind=rknd), dimension(100) :: rp,ru character(len=80), dimension(100) :: sp,su cy external sxy data ix/0,0,1,2,2,2,1,0,1/ data iy/1,2,2,2,1,0,0,0,1/ data ntf,nvf,nbf,ispd,ifirst/8,9,8,0,1/ c c common /atest2/ibc(4),icont,iu(95),a1x,a1y,a1u,a2x,a2y,a2u, c + bux,buy,cu0,cu1,cu2,cahn,cexp,cir,csin,du0,du1,eu0,f0, c 1 ru(81),su(100) c if(ip(41)==1) then sp(2)='square' sp(1)='square' sp(3)='square' sp(6)='square_mpixxx.rw' sp(7)='square.jnl' sp(9)='square_mpixxx.out' c c initialize as laplacian with dirichlet b.c. c do i=1,19 ru(i)=0.0e0_rknd enddo ru(1)=1.0e0_rknd ru(5)=1.0e0_rknd ru(9)=1.0e0_rknd do i=1,4 iu(i)=2 enddo iu(5)=0 endif c ip(1)=ntf ip(2)=nvf ip(3)=nbf ip(5)=max(ip(5),ifirst) ip(6)=1 if(iu(5)/=0) ip(6)=3 ip(7)=0 ip(8)=ispd ip(22)=100000 do i=1,ntf itnode(1,i)=9 itnode(2,i)=i itnode(3,i)=i-1 itnode(4,i)=0 itnode(5,i)=i ccc if(i>4) itnode(5,i)=1 ibndry(1,i)=i ibndry(2,i)=i-1 ibndry(3,i)=0 k=(i+1)/2 ibndry(4,i)=iu(k) ibndry(5,i)=0 ibndry(6,i)=0 ibndry(7,i)=(i+1)/2 sf(1,i)=0.0e0_rknd sf(2,i)=0.0e0_rknd enddo itnode(3,1)=8 ibndry(2,1)=8 c do i=1,nvf vx(i)=real(ix(i),rknd)/2.0e0_rknd vy(i)=real(iy(i),rknd)/2.0e0_rknd enddo c c check for periodic boundary conditions c if(iu(1)==0.and.iu(3)==0) then ibndry(4,1)=-6 ibndry(4,6)=-1 ibndry(4,2)=-5 ibndry(4,5)=-2 endif if(iu(2)==0.and.iu(4)==0) then ibndry(4,3)=-8 ibndry(4,8)=-3 ibndry(4,4)=-7 ibndry(4,7)=-4 endif c return end