152d6035fSJeremy L Thompson!----------------------------------------------------------------------- 252d6035fSJeremy L Thompson! 352d6035fSJeremy L Thompson! Header with common subroutine 452d6035fSJeremy L Thompson! 552bfb9bbSJeremy L Thompson include 't320-basis-f.h' 6752c3701SJeremy L Thompson! 7752c3701SJeremy L Thompson! Header with QFunctions 8752c3701SJeremy L Thompson! 9752c3701SJeremy L Thompson include 't510-operator-f.h' 1052d6035fSJeremy L Thompson!----------------------------------------------------------------------- 1152d6035fSJeremy L Thompson program test 121f9a83abSJed Brown implicit none 13ec3da8bcSJed Brown include 'ceed/fortran.h' 1452d6035fSJeremy L Thompson 1552d6035fSJeremy L Thompson integer ceed,err,i,j,k 1615910d16Sjeremylt integer stridesutet(3),stridesuhex(3) 1715910d16Sjeremylt integer erestrictxtet,erestrictutet,erestrictuitet,& 1815910d16Sjeremylt& erestrictxhex,erestrictuhex,erestrictuihex 1952d6035fSJeremy L Thompson integer bxtet,butet,bxhex,buhex 2052d6035fSJeremy L Thompson integer qf_setuptet,qf_masstet,qf_setuphex,qf_masshex 2152d6035fSJeremy L Thompson integer op_setuptet,op_masstet,op_setuphex,op_masshex,op_setup,op_mass 2252d6035fSJeremy L Thompson integer qdatatet,qdatahex,x,u,v 2352d6035fSJeremy L Thompson integer nelemtet,nelemhex,ptet,phex,qtet,qhex,d 2452d6035fSJeremy L Thompson integer row,col,offset 2552d6035fSJeremy L Thompson parameter(nelemtet=6) 2652d6035fSJeremy L Thompson parameter(ptet=6) 2752d6035fSJeremy L Thompson parameter(qtet=4) 2852d6035fSJeremy L Thompson parameter(nelemhex=6) 2952d6035fSJeremy L Thompson parameter(phex=3) 3052d6035fSJeremy L Thompson parameter(qhex=4) 3152d6035fSJeremy L Thompson parameter(d=2) 3252d6035fSJeremy L Thompson integer ndofs,nqptstet,nqptshex,nqpts,nx,ny,nxtet,nytet,nxhex 3352d6035fSJeremy L Thompson parameter(nx=3) 3452d6035fSJeremy L Thompson parameter(ny=3) 3552d6035fSJeremy L Thompson parameter(nxtet=3) 3652d6035fSJeremy L Thompson parameter(nytet=1) 3752d6035fSJeremy L Thompson parameter(nxhex=3) 3852d6035fSJeremy L Thompson parameter(ndofs=(nx*2+1)*(ny*2+1)) 3952d6035fSJeremy L Thompson parameter(nqptstet=nelemtet*qtet) 4052d6035fSJeremy L Thompson parameter(nqptshex=nelemhex*qhex*qhex) 4152d6035fSJeremy L Thompson parameter(nqpts=nqptstet+nqptshex) 4252d6035fSJeremy L Thompson integer indxtet(nelemtet*ptet),indxhex(nelemhex*phex*phex) 4352d6035fSJeremy L Thompson real*8 arrx(d*ndofs) 44c8b9fe72Sjeremylt integer*8 voffset,xoffset 4552d6035fSJeremy L Thompson 4652d6035fSJeremy L Thompson real*8 qref(d*qtet) 4752d6035fSJeremy L Thompson real*8 qweight(qtet) 4852d6035fSJeremy L Thompson real*8 interp(ptet*qtet) 4952d6035fSJeremy L Thompson real*8 grad(d*ptet*qtet) 5052d6035fSJeremy L Thompson 5152d6035fSJeremy L Thompson real*8 hv(ndofs) 5252d6035fSJeremy L Thompson real*8 total 5352d6035fSJeremy L Thompson 5452d6035fSJeremy L Thompson character arg*32 5552d6035fSJeremy L Thompson 5652d6035fSJeremy L Thompson external setup,mass 5752d6035fSJeremy L Thompson 5852d6035fSJeremy L Thompson call getarg(1,arg) 5952d6035fSJeremy L Thompson 6052d6035fSJeremy L Thompson call ceedinit(trim(arg)//char(0),ceed,err) 6152d6035fSJeremy L Thompson 6252d6035fSJeremy L Thompson! DoF Coordinates 6352d6035fSJeremy L Thompson do i=0,ny*2 6452d6035fSJeremy L Thompson do j=0,nx*2 6552d6035fSJeremy L Thompson arrx(i+j*(ny*2+1)+0*ndofs+1)=1.d0*i/(2*ny) 6652d6035fSJeremy L Thompson arrx(i+j*(ny*2+1)+1*ndofs+1)=1.d0*j/(2*nx) 6752d6035fSJeremy L Thompson enddo 6852d6035fSJeremy L Thompson enddo 6952d6035fSJeremy L Thompson 7052d6035fSJeremy L Thompson call ceedvectorcreate(ceed,d*ndofs,x,err) 71c8b9fe72Sjeremylt xoffset=0 72c8b9fe72Sjeremylt call ceedvectorsetarray(x,ceed_mem_host,ceed_use_pointer,arrx,xoffset,err) 7352d6035fSJeremy L Thompson 7452d6035fSJeremy L Thompson! Qdata Vectors 7552d6035fSJeremy L Thompson call ceedvectorcreate(ceed,nqptstet,qdatatet,err) 7652d6035fSJeremy L Thompson call ceedvectorcreate(ceed,nqptshex,qdatahex,err) 7752d6035fSJeremy L Thompson 7852d6035fSJeremy L Thompson! Tet Elements 7952d6035fSJeremy L Thompson do i=0,2 8052d6035fSJeremy L Thompson col=mod(i,nx) 8152d6035fSJeremy L Thompson row=i/nx 8252d6035fSJeremy L Thompson offset=col*2+row*(nx*2+1)*2 8352d6035fSJeremy L Thompson 8452d6035fSJeremy L Thompson indxtet(i*2*ptet+1)=2+offset 8552d6035fSJeremy L Thompson indxtet(i*2*ptet+2)=9+offset 8652d6035fSJeremy L Thompson indxtet(i*2*ptet+3)=16+offset 8752d6035fSJeremy L Thompson indxtet(i*2*ptet+4)=1+offset 8852d6035fSJeremy L Thompson indxtet(i*2*ptet+5)=8+offset 8952d6035fSJeremy L Thompson indxtet(i*2*ptet+6)=0+offset 9052d6035fSJeremy L Thompson 9152d6035fSJeremy L Thompson indxtet(i*2*ptet+7)=14+offset 9252d6035fSJeremy L Thompson indxtet(i*2*ptet+8)=7+offset 9352d6035fSJeremy L Thompson indxtet(i*2*ptet+9)=0+offset 9452d6035fSJeremy L Thompson indxtet(i*2*ptet+10)=15+offset 9552d6035fSJeremy L Thompson indxtet(i*2*ptet+11)=8+offset 9652d6035fSJeremy L Thompson indxtet(i*2*ptet+12)=16+offset 9752d6035fSJeremy L Thompson enddo 9852d6035fSJeremy L Thompson 9952d6035fSJeremy L Thompson! -- Restrictions 100d979a051Sjeremylt call ceedelemrestrictioncreate(ceed,nelemtet,ptet,d,ndofs,d*ndofs,& 101a8d32208Sjeremylt & ceed_mem_host,ceed_use_pointer,indxtet,erestrictxtet,err) 10252d6035fSJeremy L Thompson 103d979a051Sjeremylt call ceedelemrestrictioncreate(ceed,nelemtet,ptet,1,1,ndofs,& 104a8d32208Sjeremylt & ceed_mem_host,ceed_use_pointer,indxtet,erestrictutet,err) 1057509a596Sjeremylt stridesutet=[1,qtet,qtet] 106d979a051Sjeremylt call ceedelemrestrictioncreatestrided(ceed,nelemtet,qtet,1,& 107d979a051Sjeremylt & nqptstet,stridesutet,erestrictuitet,err) 10852d6035fSJeremy L Thompson 10952d6035fSJeremy L Thompson! -- Bases 11052d6035fSJeremy L Thompson call buildmats(qref,qweight,interp,grad) 11152d6035fSJeremy L Thompson call ceedbasiscreateh1(ceed,ceed_triangle,d,ptet,qtet,interp,grad,qref,& 11252d6035fSJeremy L Thompson & qweight,bxtet,err) 11352d6035fSJeremy L Thompson call buildmats(qref,qweight,interp,grad) 11452d6035fSJeremy L Thompson call ceedbasiscreateh1(ceed,ceed_triangle,1,ptet,qtet,interp,grad,qref,& 11552d6035fSJeremy L Thompson & qweight,butet,err) 11652d6035fSJeremy L Thompson 11752d6035fSJeremy L Thompson! -- QFunctions 11852d6035fSJeremy L Thompson call ceedqfunctioncreateinterior(ceed,1,setup,& 1192d50dd3dSjeremylt &SOURCE_DIR& 120a05f9790Sjeremylt &//'t510-operator.h:setup'//char(0),qf_setuptet,err) 121a61c78d6SJeremy L Thompson call ceedqfunctionaddinput(qf_setuptet,'weight',1,ceed_eval_weight,err) 1224d537eeaSYohann call ceedqfunctionaddinput(qf_setuptet,'dx',d*d,ceed_eval_grad,err) 12352d6035fSJeremy L Thompson call ceedqfunctionaddoutput(qf_setuptet,'rho',1,ceed_eval_none,err) 12452d6035fSJeremy L Thompson 12552d6035fSJeremy L Thompson call ceedqfunctioncreateinterior(ceed,1,mass,& 1262d50dd3dSjeremylt &SOURCE_DIR& 127a05f9790Sjeremylt &//'t510-operator.h:mass'//char(0),qf_masstet,err) 12852d6035fSJeremy L Thompson call ceedqfunctionaddinput(qf_masstet,'rho',1,ceed_eval_none,err) 12952d6035fSJeremy L Thompson call ceedqfunctionaddinput(qf_masstet,'u',1,ceed_eval_interp,err) 13052d6035fSJeremy L Thompson call ceedqfunctionaddoutput(qf_masstet,'v',1,ceed_eval_interp,err) 13152d6035fSJeremy L Thompson 13252d6035fSJeremy L Thompson! -- Operators 13352d6035fSJeremy L Thompson! ---- Setup Tet 134442e7f0bSjeremylt call ceedoperatorcreate(ceed,qf_setuptet,ceed_qfunction_none,& 135442e7f0bSjeremylt & ceed_qfunction_none,op_setuptet,err) 136a61c78d6SJeremy L Thompson call ceedoperatorsetfield(op_setuptet,'weight',& 13715910d16Sjeremylt & ceed_elemrestriction_none,bxtet,ceed_vector_none,err) 13852d6035fSJeremy L Thompson call ceedoperatorsetfield(op_setuptet,'dx',erestrictxtet,& 139a8d32208Sjeremylt & bxtet,ceed_vector_active,err) 14052d6035fSJeremy L Thompson call ceedoperatorsetfield(op_setuptet,'rho',erestrictuitet,& 141356036faSJeremy L Thompson ceed_basis_none,qdatatet,err) 14252d6035fSJeremy L Thompson! ---- Mass Tet 143442e7f0bSjeremylt call ceedoperatorcreate(ceed,qf_masstet,ceed_qfunction_none,& 144442e7f0bSjeremylt & ceed_qfunction_none,op_masstet,err) 14552d6035fSJeremy L Thompson call ceedoperatorsetfield(op_masstet,'rho',erestrictuitet,& 146356036faSJeremy L Thompson ceed_basis_none,qdatatet,err) 14752d6035fSJeremy L Thompson call ceedoperatorsetfield(op_masstet,'u',erestrictutet,& 148a8d32208Sjeremylt & butet,ceed_vector_active,err) 14952d6035fSJeremy L Thompson call ceedoperatorsetfield(op_masstet,'v',erestrictutet,& 150a8d32208Sjeremylt & butet,ceed_vector_active,err) 15152d6035fSJeremy L Thompson 15252d6035fSJeremy L Thompson! Hex Elements 15352d6035fSJeremy L Thompson do i=0,nelemhex-1 15452d6035fSJeremy L Thompson col=mod(i,nx) 15552d6035fSJeremy L Thompson row=i/nx 15652d6035fSJeremy L Thompson offset=(nxtet*2+1)*(nytet*2)*(1+row)+col*2 15752d6035fSJeremy L Thompson do j=0,phex-1 15852d6035fSJeremy L Thompson do k=0,phex-1 15952d6035fSJeremy L Thompson indxhex(phex*(phex*i+k)+j+1)=offset+k*(nxhex*2+1)+j 16052d6035fSJeremy L Thompson enddo 16152d6035fSJeremy L Thompson enddo 16252d6035fSJeremy L Thompson enddo 16352d6035fSJeremy L Thompson 16452d6035fSJeremy L Thompson! -- Restrictions 165d979a051Sjeremylt call ceedelemrestrictioncreate(ceed,nelemhex,phex*phex,d,ndofs,& 166d979a051Sjeremylt & d*ndofs,ceed_mem_host,ceed_use_pointer,indxhex,erestrictxhex,err) 16752d6035fSJeremy L Thompson 168d979a051Sjeremylt call ceedelemrestrictioncreate(ceed,nelemhex,phex*phex,1,1,ndofs,& 16952d6035fSJeremy L Thompson & ceed_mem_host,ceed_use_pointer,indxhex,erestrictuhex,err) 1707509a596Sjeremylt stridesuhex=[1,qhex*qhex,qhex*qhex] 1717509a596Sjeremylt call ceedelemrestrictioncreatestrided(ceed,nelemhex,qhex*qhex,& 172d979a051Sjeremylt & 1,nqptshex,stridesuhex,erestrictuihex,err) 17352d6035fSJeremy L Thompson 17452d6035fSJeremy L Thompson! -- Bases 17552d6035fSJeremy L Thompson call ceedbasiscreatetensorh1lagrange(ceed,d,d,phex,qhex,ceed_gauss,& 17652d6035fSJeremy L Thompson & bxhex,err) 17752d6035fSJeremy L Thompson call ceedbasiscreatetensorh1lagrange(ceed,d,1,phex,qhex,ceed_gauss,& 17852d6035fSJeremy L Thompson & buhex,err) 17952d6035fSJeremy L Thompson 18052d6035fSJeremy L Thompson! -- QFunctions 18152d6035fSJeremy L Thompson call ceedqfunctioncreateinterior(ceed,1,setup,& 1822d50dd3dSjeremylt &SOURCE_DIR& 183872c4ebbSjeremylt &//'t510-operator.h:setup'//char(0),qf_setuphex,err) 18452d6035fSJeremy L Thompson call ceedqfunctionaddinput(qf_setuphex,'_weight',1,ceed_eval_weight,err) 1854d537eeaSYohann call ceedqfunctionaddinput(qf_setuphex,'dx',d*d,ceed_eval_grad,err) 18652d6035fSJeremy L Thompson call ceedqfunctionaddoutput(qf_setuphex,'rho',1,ceed_eval_none,err) 18752d6035fSJeremy L Thompson 18852d6035fSJeremy L Thompson call ceedqfunctioncreateinterior(ceed,1,mass,& 1892d50dd3dSjeremylt &SOURCE_DIR& 190872c4ebbSjeremylt &//'t510-operator.h:mass'//char(0),qf_masshex,err) 19152d6035fSJeremy L Thompson call ceedqfunctionaddinput(qf_masshex,'rho',1,ceed_eval_none,err) 19252d6035fSJeremy L Thompson call ceedqfunctionaddinput(qf_masshex,'u',1,ceed_eval_interp,err) 19352d6035fSJeremy L Thompson call ceedqfunctionaddoutput(qf_masshex,'v',1,ceed_eval_interp,err) 19452d6035fSJeremy L Thompson 19552d6035fSJeremy L Thompson! -- Operators 19652d6035fSJeremy L Thompson! ---- Setup Hex 197442e7f0bSjeremylt call ceedoperatorcreate(ceed,qf_setuphex,ceed_qfunction_none,& 198442e7f0bSjeremylt & ceed_qfunction_none,op_setuphex,err) 19915910d16Sjeremylt call ceedoperatorsetfield(op_setuphex,'_weight',& 20015910d16Sjeremylt & ceed_elemrestriction_none,bxhex,ceed_vector_none,err) 20152d6035fSJeremy L Thompson call ceedoperatorsetfield(op_setuphex,'dx',erestrictxhex,& 202a8d32208Sjeremylt & bxhex,ceed_vector_active,err) 20352d6035fSJeremy L Thompson call ceedoperatorsetfield(op_setuphex,'rho',erestrictuihex,& 204356036faSJeremy L Thompson ceed_basis_none,qdatahex,err) 20552d6035fSJeremy L Thompson! ---- Mass Hex 206442e7f0bSjeremylt call ceedoperatorcreate(ceed,qf_masshex,ceed_qfunction_none,& 207442e7f0bSjeremylt & ceed_qfunction_none,op_masshex,err) 20852d6035fSJeremy L Thompson call ceedoperatorsetfield(op_masshex,'rho',erestrictuihex,& 209356036faSJeremy L Thompson ceed_basis_none,qdatahex,err) 21052d6035fSJeremy L Thompson call ceedoperatorsetfield(op_masshex,'u',erestrictuhex,& 211a8d32208Sjeremylt & buhex,ceed_vector_active,err) 21252d6035fSJeremy L Thompson call ceedoperatorsetfield(op_masshex,'v',erestrictuhex,& 213a8d32208Sjeremylt & buhex,ceed_vector_active,err) 21452d6035fSJeremy L Thompson 21552d6035fSJeremy L Thompson! Composite Operators 216*ed094490SJeremy L Thompson call ceedoperatorcreatecomposite(ceed,op_setup,err) 217*ed094490SJeremy L Thompson call ceedoperatorcompositeaddsub(op_setup,op_setuptet,err) 218*ed094490SJeremy L Thompson call ceedoperatorcompositeaddsub(op_setup,op_setuphex,err) 21952d6035fSJeremy L Thompson 220*ed094490SJeremy L Thompson call ceedoperatorcreatecomposite(ceed,op_mass,err) 221*ed094490SJeremy L Thompson call ceedoperatorcompositeaddsub(op_mass,op_masstet,err) 222*ed094490SJeremy L Thompson call ceedoperatorcompositeaddsub(op_mass,op_masshex,err) 22352d6035fSJeremy L Thompson 22452d6035fSJeremy L Thompson! Apply Setup Operator 225e97ff134Sjeremylt call ceedoperatorapply(op_setup,x,ceed_vector_none,& 226e97ff134Sjeremylt & ceed_request_immediate,err) 22752d6035fSJeremy L Thompson 22852d6035fSJeremy L Thompson! Apply Mass Operator 22952d6035fSJeremy L Thompson call ceedvectorcreate(ceed,ndofs,u,err) 23052d6035fSJeremy L Thompson call ceedvectorsetvalue(u,1.d0,err) 23152d6035fSJeremy L Thompson call ceedvectorcreate(ceed,ndofs,v,err) 23252d6035fSJeremy L Thompson 23352d6035fSJeremy L Thompson call ceedoperatorapply(op_mass,u,v,ceed_request_immediate,err) 23452d6035fSJeremy L Thompson 23552d6035fSJeremy L Thompson! Check Output 23652d6035fSJeremy L Thompson call ceedvectorgetarrayread(v,ceed_mem_host,hv,voffset,err) 23752d6035fSJeremy L Thompson total=0. 23852d6035fSJeremy L Thompson do i=1,ndofs 23952d6035fSJeremy L Thompson total=total+hv(voffset+i) 24052d6035fSJeremy L Thompson enddo 24152d6035fSJeremy L Thompson if (abs(total-1.)>1.0d-10) then 242a2546046Sjeremylt! LCOV_EXCL_START 24352d6035fSJeremy L Thompson write(*,*) 'Computed Area: ',total,' != True Area: 1.0' 244de996c55Sjeremylt! LCOV_EXCL_STOP 24552d6035fSJeremy L Thompson endif 24652d6035fSJeremy L Thompson call ceedvectorrestorearrayread(v,hv,voffset,err) 24752d6035fSJeremy L Thompson 24852d6035fSJeremy L Thompson! Cleanup 24952d6035fSJeremy L Thompson call ceedqfunctiondestroy(qf_setuptet,err) 25052d6035fSJeremy L Thompson call ceedqfunctiondestroy(qf_masstet,err) 25152d6035fSJeremy L Thompson call ceedoperatordestroy(op_setuptet,err) 25252d6035fSJeremy L Thompson call ceedoperatordestroy(op_masstet,err) 25352d6035fSJeremy L Thompson call ceedqfunctiondestroy(qf_setuphex,err) 25452d6035fSJeremy L Thompson call ceedqfunctiondestroy(qf_masshex,err) 25552d6035fSJeremy L Thompson call ceedoperatordestroy(op_setuphex,err) 25652d6035fSJeremy L Thompson call ceedoperatordestroy(op_masshex,err) 25752d6035fSJeremy L Thompson call ceedoperatordestroy(op_setup,err) 25852d6035fSJeremy L Thompson call ceedoperatordestroy(op_mass,err) 25952d6035fSJeremy L Thompson call ceedelemrestrictiondestroy(erestrictutet,err) 26052d6035fSJeremy L Thompson call ceedelemrestrictiondestroy(erestrictxtet,err) 26152d6035fSJeremy L Thompson call ceedelemrestrictiondestroy(erestrictuitet,err) 26252d6035fSJeremy L Thompson call ceedelemrestrictiondestroy(erestrictuhex,err) 26352d6035fSJeremy L Thompson call ceedelemrestrictiondestroy(erestrictxhex,err) 26452d6035fSJeremy L Thompson call ceedelemrestrictiondestroy(erestrictuihex,err) 26552d6035fSJeremy L Thompson call ceedbasisdestroy(butet,err) 26652d6035fSJeremy L Thompson call ceedbasisdestroy(bxtet,err) 26752d6035fSJeremy L Thompson call ceedbasisdestroy(buhex,err) 26852d6035fSJeremy L Thompson call ceedbasisdestroy(bxhex,err) 26952d6035fSJeremy L Thompson call ceedvectordestroy(x,err) 27052d6035fSJeremy L Thompson call ceedvectordestroy(u,err) 27152d6035fSJeremy L Thompson call ceedvectordestroy(v,err) 27252d6035fSJeremy L Thompson call ceedvectordestroy(qdatatet,err) 27352d6035fSJeremy L Thompson call ceedvectordestroy(qdatahex,err) 27452d6035fSJeremy L Thompson call ceeddestroy(ceed,err) 27552d6035fSJeremy L Thompson end 27652d6035fSJeremy L Thompson!----------------------------------------------------------------------- 277