1*d99fa3c5SJeremy L Thompson!----------------------------------------------------------------------- 2*d99fa3c5SJeremy L Thompson program test 3*d99fa3c5SJeremy L Thompson 4*d99fa3c5SJeremy L Thompson include 'ceedf.h' 5*d99fa3c5SJeremy L Thompson 6*d99fa3c5SJeremy L Thompson integer ceed,err,i,j 7*d99fa3c5SJeremy L Thompson integer stridesu(3) 8*d99fa3c5SJeremy L Thompson integer erestrictx,erestrictui 9*d99fa3c5SJeremy L Thompson integer erestrictucoarse,erestrictufine 10*d99fa3c5SJeremy L Thompson integer bx,bufine,bucoarse,bctof 11*d99fa3c5SJeremy L Thompson integer qf_setup,qf_mass 12*d99fa3c5SJeremy L Thompson integer op_setup,op_masscoarse,op_massfine 13*d99fa3c5SJeremy L Thompson integer op_prolong,op_restrict 14*d99fa3c5SJeremy L Thompson integer qdata,x,ucoarse,ufine,vcoarse,vfine,pmultfine 15*d99fa3c5SJeremy L Thompson integer nelem,pcoarse,pfine,q 16*d99fa3c5SJeremy L Thompson parameter(nelem=15) 17*d99fa3c5SJeremy L Thompson parameter(pcoarse=3) 18*d99fa3c5SJeremy L Thompson parameter(pfine=5) 19*d99fa3c5SJeremy L Thompson parameter(q=8) 20*d99fa3c5SJeremy L Thompson integer nx,nucoarse,nufine 21*d99fa3c5SJeremy L Thompson parameter(nx=nelem+1) 22*d99fa3c5SJeremy L Thompson parameter(nucoarse=nelem*(pcoarse-1)+1) 23*d99fa3c5SJeremy L Thompson parameter(nufine=nelem*(pfine-1)+1) 24*d99fa3c5SJeremy L Thompson integer indx(nelem*2) 25*d99fa3c5SJeremy L Thompson integer inducoarse(nelem*pcoarse) 26*d99fa3c5SJeremy L Thompson integer indufine(nelem*pfine) 27*d99fa3c5SJeremy L Thompson real*8 arrx(nx) 28*d99fa3c5SJeremy L Thompson integer*8 voffset,xoffset,ioffset 29*d99fa3c5SJeremy L Thompson real*8 val 30*d99fa3c5SJeremy L Thompson integer interpsize 31*d99fa3c5SJeremy L Thompson parameter(interpsize=pcoarse*pfine); 32*d99fa3c5SJeremy L Thompson real*8 interp1d(interpsize),interpctof(interpsize) 33*d99fa3c5SJeremy L Thompson 34*d99fa3c5SJeremy L Thompson real*8 hv(nufine) 35*d99fa3c5SJeremy L Thompson real*8 total 36*d99fa3c5SJeremy L Thompson 37*d99fa3c5SJeremy L Thompson character arg*32 38*d99fa3c5SJeremy L Thompson 39*d99fa3c5SJeremy L Thompson call getarg(1,arg) 40*d99fa3c5SJeremy L Thompson call ceedinit(trim(arg)//char(0),ceed,err) 41*d99fa3c5SJeremy L Thompson 42*d99fa3c5SJeremy L Thompson do i=0,nx-1 43*d99fa3c5SJeremy L Thompson arrx(i+1)=i/(nx-1.d0) 44*d99fa3c5SJeremy L Thompson enddo 45*d99fa3c5SJeremy L Thompson do i=0,nelem-1 46*d99fa3c5SJeremy L Thompson indx(2*i+1)=i 47*d99fa3c5SJeremy L Thompson indx(2*i+2)=i+1 48*d99fa3c5SJeremy L Thompson enddo 49*d99fa3c5SJeremy L Thompson call ceedelemrestrictioncreate(ceed,nelem,2,1,1,nx,ceed_mem_host,& 50*d99fa3c5SJeremy L Thompson & ceed_use_pointer,indx,erestrictx,err) 51*d99fa3c5SJeremy L Thompson 52*d99fa3c5SJeremy L Thompson do i=0,nelem-1 53*d99fa3c5SJeremy L Thompson do j=0,pcoarse-1 54*d99fa3c5SJeremy L Thompson inducoarse(pcoarse*i+j+1)=i*(pcoarse-1)+j 55*d99fa3c5SJeremy L Thompson enddo 56*d99fa3c5SJeremy L Thompson enddo 57*d99fa3c5SJeremy L Thompson call ceedelemrestrictioncreate(ceed,nelem,pcoarse,1,1,nucoarse,& 58*d99fa3c5SJeremy L Thompson & ceed_mem_host,ceed_use_pointer,inducoarse,erestrictucoarse,err) 59*d99fa3c5SJeremy L Thompson 60*d99fa3c5SJeremy L Thompson do i=0,nelem-1 61*d99fa3c5SJeremy L Thompson do j=0,pfine-1 62*d99fa3c5SJeremy L Thompson indufine(pfine*i+j+1)=i*(pfine-1)+j 63*d99fa3c5SJeremy L Thompson enddo 64*d99fa3c5SJeremy L Thompson enddo 65*d99fa3c5SJeremy L Thompson call ceedelemrestrictioncreate(ceed,nelem,pfine,1,1,nufine,& 66*d99fa3c5SJeremy L Thompson & ceed_mem_host,ceed_use_pointer,indufine,erestrictufine,err) 67*d99fa3c5SJeremy L Thompson 68*d99fa3c5SJeremy L Thompson stridesu=[1,q,q] 69*d99fa3c5SJeremy L Thompson call ceedelemrestrictioncreatestrided(ceed,nelem,q,1,q*nelem,stridesu,& 70*d99fa3c5SJeremy L Thompson & erestrictui,err) 71*d99fa3c5SJeremy L Thompson 72*d99fa3c5SJeremy L Thompson call ceedbasiscreatetensorh1lagrange(ceed,1,1,2,q,ceed_gauss,bx,err) 73*d99fa3c5SJeremy L Thompson call ceedbasiscreatetensorh1lagrange(ceed,1,1,pfine,q,ceed_gauss,& 74*d99fa3c5SJeremy L Thompson & bufine,err) 75*d99fa3c5SJeremy L Thompson 76*d99fa3c5SJeremy L Thompson call ceedqfunctioncreateinteriorbyname(ceed,"Mass1DBuild",qf_setup,err) 77*d99fa3c5SJeremy L Thompson call ceedqfunctioncreateinteriorbyname(ceed,"MassApply",qf_mass,err) 78*d99fa3c5SJeremy L Thompson 79*d99fa3c5SJeremy L Thompson call ceedoperatorcreate(ceed,qf_setup,ceed_qfunction_none,& 80*d99fa3c5SJeremy L Thompson & ceed_qfunction_none,op_setup,err) 81*d99fa3c5SJeremy L Thompson call ceedoperatorcreate(ceed,qf_mass,ceed_qfunction_none,& 82*d99fa3c5SJeremy L Thompson & ceed_qfunction_none,op_massfine,err) 83*d99fa3c5SJeremy L Thompson 84*d99fa3c5SJeremy L Thompson call ceedvectorcreate(ceed,nx,x,err) 85*d99fa3c5SJeremy L Thompson xoffset=0 86*d99fa3c5SJeremy L Thompson call ceedvectorsetarray(x,ceed_mem_host,ceed_use_pointer,arrx,xoffset,err) 87*d99fa3c5SJeremy L Thompson call ceedvectorcreate(ceed,nelem*q,qdata,err) 88*d99fa3c5SJeremy L Thompson 89*d99fa3c5SJeremy L Thompson call ceedoperatorsetfield(op_setup,'weights',ceed_elemrestriction_none,& 90*d99fa3c5SJeremy L Thompson & bx,ceed_vector_none,err) 91*d99fa3c5SJeremy L Thompson call ceedoperatorsetfield(op_setup,'dx',erestrictx,bx,& 92*d99fa3c5SJeremy L Thompson & ceed_vector_active,err) 93*d99fa3c5SJeremy L Thompson call ceedoperatorsetfield(op_setup,'qdata',erestrictui,& 94*d99fa3c5SJeremy L Thompson & ceed_basis_collocated,ceed_vector_active,err) 95*d99fa3c5SJeremy L Thompson call ceedoperatorsetfield(op_massfine,'qdata',erestrictui,& 96*d99fa3c5SJeremy L Thompson & ceed_basis_collocated,qdata,err) 97*d99fa3c5SJeremy L Thompson call ceedoperatorsetfield(op_massfine,'u',erestrictufine,bufine,& 98*d99fa3c5SJeremy L Thompson & ceed_vector_active,err) 99*d99fa3c5SJeremy L Thompson call ceedoperatorsetfield(op_massfine,'v',erestrictufine,bufine,& 100*d99fa3c5SJeremy L Thompson & ceed_vector_active,err) 101*d99fa3c5SJeremy L Thompson 102*d99fa3c5SJeremy L Thompson call ceedoperatorapply(op_setup,x,qdata,ceed_request_immediate,err) 103*d99fa3c5SJeremy L Thompson 104*d99fa3c5SJeremy L Thompson! Create multigrid level 105*d99fa3c5SJeremy L Thompson call ceedvectorcreate(ceed,nufine,pmultfine,err) 106*d99fa3c5SJeremy L Thompson val=1.0 107*d99fa3c5SJeremy L Thompson call ceedvectorsetvalue(pmultfine,val,err) 108*d99fa3c5SJeremy L Thompson call ceedbasiscreatetensorh1lagrange(ceed,1,1,pcoarse,q,ceed_gauss,& 109*d99fa3c5SJeremy L Thompson & bucoarse,err) 110*d99fa3c5SJeremy L Thompson call ceedbasiscreatetensorh1lagrange(ceed,1,1,pcoarse,pfine,& 111*d99fa3c5SJeremy L Thompson & ceed_gauss_lobatto,bctof,err) 112*d99fa3c5SJeremy L Thompson call ceedbasisgetinterp1d(bctof,interp1d,ioffset,err) 113*d99fa3c5SJeremy L Thompson do i=1,interpsize 114*d99fa3c5SJeremy L Thompson interpctof(i)=interp1d(ioffset+i) 115*d99fa3c5SJeremy L Thompson enddo 116*d99fa3c5SJeremy L Thompson call ceedoperatormultigridlevelcreateh1(op_massfine,pmultfine,& 117*d99fa3c5SJeremy L Thompson & erestrictucoarse,bucoarse,interpctof,op_masscoarse,& 118*d99fa3c5SJeremy L Thompson & op_prolong,op_restrict,err) 119*d99fa3c5SJeremy L Thompson 120*d99fa3c5SJeremy L Thompson! Coarse problem 121*d99fa3c5SJeremy L Thompson call ceedvectorcreate(ceed,nucoarse,ucoarse,err) 122*d99fa3c5SJeremy L Thompson val=1.0 123*d99fa3c5SJeremy L Thompson call ceedvectorsetvalue(ucoarse,val,err) 124*d99fa3c5SJeremy L Thompson call ceedvectorcreate(ceed,nucoarse,vcoarse,err) 125*d99fa3c5SJeremy L Thompson call ceedoperatorapply(op_masscoarse,ucoarse,vcoarse,& 126*d99fa3c5SJeremy L Thompson & ceed_request_immediate,err) 127*d99fa3c5SJeremy L Thompson 128*d99fa3c5SJeremy L Thompson! Check output 129*d99fa3c5SJeremy L Thompson call ceedvectorgetarrayread(vcoarse,ceed_mem_host,hv,voffset,err) 130*d99fa3c5SJeremy L Thompson total=0. 131*d99fa3c5SJeremy L Thompson do i=1,nucoarse 132*d99fa3c5SJeremy L Thompson total=total+hv(voffset+i) 133*d99fa3c5SJeremy L Thompson enddo 134*d99fa3c5SJeremy L Thompson if (abs(total-1.)>1.0d-10) then 135*d99fa3c5SJeremy L Thompson! LCOV_EXCL_START 136*d99fa3c5SJeremy L Thompson write(*,*) 'Computed Area: ',total,' != True Area: 1.0' 137*d99fa3c5SJeremy L Thompson! LCOV_EXCL_STOP 138*d99fa3c5SJeremy L Thompson endif 139*d99fa3c5SJeremy L Thompson call ceedvectorrestorearrayread(vcoarse,hv,voffset,err) 140*d99fa3c5SJeremy L Thompson 141*d99fa3c5SJeremy L Thompson! Prolong coarse u 142*d99fa3c5SJeremy L Thompson call ceedvectorcreate(ceed,nufine,ufine,err) 143*d99fa3c5SJeremy L Thompson call ceedoperatorapply(op_prolong,ucoarse,ufine,& 144*d99fa3c5SJeremy L Thompson & ceed_request_immediate,err) 145*d99fa3c5SJeremy L Thompson 146*d99fa3c5SJeremy L Thompson! Fine problem 147*d99fa3c5SJeremy L Thompson call ceedvectorcreate(ceed,nufine,vfine,err) 148*d99fa3c5SJeremy L Thompson call ceedoperatorapply(op_massfine,ufine,vfine,& 149*d99fa3c5SJeremy L Thompson & ceed_request_immediate,err) 150*d99fa3c5SJeremy L Thompson 151*d99fa3c5SJeremy L Thompson! Check output 152*d99fa3c5SJeremy L Thompson call ceedvectorgetarrayread(vfine,ceed_mem_host,hv,voffset,err) 153*d99fa3c5SJeremy L Thompson total=0. 154*d99fa3c5SJeremy L Thompson do i=1,nufine 155*d99fa3c5SJeremy L Thompson total=total+hv(voffset+i) 156*d99fa3c5SJeremy L Thompson enddo 157*d99fa3c5SJeremy L Thompson if (abs(total-1.)>1.0d-10) then 158*d99fa3c5SJeremy L Thompson! LCOV_EXCL_START 159*d99fa3c5SJeremy L Thompson write(*,*) 'Computed Area Fine Grid: ',total,' != True Area: 1.0' 160*d99fa3c5SJeremy L Thompson! LCOV_EXCL_STOP 161*d99fa3c5SJeremy L Thompson endif 162*d99fa3c5SJeremy L Thompson call ceedvectorrestorearrayread(vfine,hv,voffset,err) 163*d99fa3c5SJeremy L Thompson 164*d99fa3c5SJeremy L Thompson! Restrict state to coarse grid 165*d99fa3c5SJeremy L Thompson call ceedoperatorapply(op_restrict,vfine,vcoarse,& 166*d99fa3c5SJeremy L Thompson & ceed_request_immediate,err) 167*d99fa3c5SJeremy L Thompson 168*d99fa3c5SJeremy L Thompson! Check output 169*d99fa3c5SJeremy L Thompson call ceedvectorgetarrayread(vcoarse,ceed_mem_host,hv,voffset,err) 170*d99fa3c5SJeremy L Thompson total=0. 171*d99fa3c5SJeremy L Thompson do i=1,nucoarse 172*d99fa3c5SJeremy L Thompson total=total+hv(voffset+i) 173*d99fa3c5SJeremy L Thompson enddo 174*d99fa3c5SJeremy L Thompson if (abs(total-1.)>1.0d-10) then 175*d99fa3c5SJeremy L Thompson! LCOV_EXCL_START 176*d99fa3c5SJeremy L Thompson write(*,*) 'Computed Area Coarse Grid: ',total,' != True Area: 1.0' 177*d99fa3c5SJeremy L Thompson! LCOV_EXCL_STOP 178*d99fa3c5SJeremy L Thompson endif 179*d99fa3c5SJeremy L Thompson call ceedvectorrestorearrayread(vcoarse,hv,voffset,err) 180*d99fa3c5SJeremy L Thompson 181*d99fa3c5SJeremy L Thompson call ceedvectordestroy(qdata,err) 182*d99fa3c5SJeremy L Thompson call ceedvectordestroy(x,err) 183*d99fa3c5SJeremy L Thompson call ceedvectordestroy(ucoarse,err) 184*d99fa3c5SJeremy L Thompson call ceedvectordestroy(ufine,err) 185*d99fa3c5SJeremy L Thompson call ceedvectordestroy(vcoarse,err) 186*d99fa3c5SJeremy L Thompson call ceedvectordestroy(vfine,err) 187*d99fa3c5SJeremy L Thompson call ceedvectordestroy(pmultfine,err) 188*d99fa3c5SJeremy L Thompson call ceedoperatordestroy(op_masscoarse,err) 189*d99fa3c5SJeremy L Thompson call ceedoperatordestroy(op_massfine,err) 190*d99fa3c5SJeremy L Thompson call ceedoperatordestroy(op_prolong,err) 191*d99fa3c5SJeremy L Thompson call ceedoperatordestroy(op_restrict,err) 192*d99fa3c5SJeremy L Thompson call ceedoperatordestroy(op_setup,err) 193*d99fa3c5SJeremy L Thompson call ceedqfunctiondestroy(qf_mass,err) 194*d99fa3c5SJeremy L Thompson call ceedqfunctiondestroy(qf_setup,err) 195*d99fa3c5SJeremy L Thompson call ceedbasisdestroy(bucoarse,err) 196*d99fa3c5SJeremy L Thompson call ceedbasisdestroy(bufine,err) 197*d99fa3c5SJeremy L Thompson call ceedbasisdestroy(bctof,err) 198*d99fa3c5SJeremy L Thompson call ceedbasisdestroy(bx,err) 199*d99fa3c5SJeremy L Thompson call ceedelemrestrictiondestroy(erestrictucoarse,err) 200*d99fa3c5SJeremy L Thompson call ceedelemrestrictiondestroy(erestrictufine,err) 201*d99fa3c5SJeremy L Thompson call ceedelemrestrictiondestroy(erestrictx,err) 202*d99fa3c5SJeremy L Thompson call ceedelemrestrictiondestroy(erestrictui,err) 203*d99fa3c5SJeremy L Thompson call ceeddestroy(ceed,err) 204*d99fa3c5SJeremy L Thompson end 205*d99fa3c5SJeremy L Thompson!----------------------------------------------------------------------- 206