1*cae8b89aSjeremylt!----------------------------------------------------------------------- 2*cae8b89aSjeremylt subroutine setup(ctx,q,u1,u2,u3,u4,u5,u6,u7,u8,u9,u10,u11,u12,u13,u14,& 3*cae8b89aSjeremylt& u15,u16,v1,v2,v3,v4,v5,v6,v7,v8,v9,v10,v11,v12,v13,v14,v15,v16,ierr) 4*cae8b89aSjeremylt real*8 ctx 5*cae8b89aSjeremylt real*8 u1(1) 6*cae8b89aSjeremylt real*8 u2(1) 7*cae8b89aSjeremylt real*8 v1(1) 8*cae8b89aSjeremylt integer q,ierr 9*cae8b89aSjeremylt 10*cae8b89aSjeremylt do i=1,q 11*cae8b89aSjeremylt v1(i)=u1(i)*u2(i) 12*cae8b89aSjeremylt enddo 13*cae8b89aSjeremylt 14*cae8b89aSjeremylt ierr=0 15*cae8b89aSjeremylt end 16*cae8b89aSjeremylt!----------------------------------------------------------------------- 17*cae8b89aSjeremylt subroutine mass(ctx,q,u1,u2,u3,u4,u5,u6,u7,u8,u9,u10,u11,u12,u13,u14,& 18*cae8b89aSjeremylt& u15,u16,v1,v2,v3,v4,v5,v6,v7,v8,v9,v10,v11,v12,v13,v14,v15,v16,ierr) 19*cae8b89aSjeremylt real*8 ctx 20*cae8b89aSjeremylt real*8 u1(1) 21*cae8b89aSjeremylt real*8 u2(1) 22*cae8b89aSjeremylt real*8 v1(1) 23*cae8b89aSjeremylt integer q,ierr 24*cae8b89aSjeremylt 25*cae8b89aSjeremylt do i=1,q 26*cae8b89aSjeremylt v1(i)=u2(i)*u1(i) 27*cae8b89aSjeremylt enddo 28*cae8b89aSjeremylt 29*cae8b89aSjeremylt ierr=0 30*cae8b89aSjeremylt end 31*cae8b89aSjeremylt!----------------------------------------------------------------------- 32*cae8b89aSjeremylt program test 33*cae8b89aSjeremylt 34*cae8b89aSjeremylt include 'ceedf.h' 35*cae8b89aSjeremylt 36*cae8b89aSjeremylt integer ceed,err,i,j 37*cae8b89aSjeremylt integer erestrictx,erestrictu,erestrictxi,erestrictui 38*cae8b89aSjeremylt integer bx,bu 39*cae8b89aSjeremylt integer qf_setup,qf_mass 40*cae8b89aSjeremylt integer op_setup,op_mass 41*cae8b89aSjeremylt integer qdata,x,u,v 42*cae8b89aSjeremylt integer nelem,p,q 43*cae8b89aSjeremylt parameter(nelem=15) 44*cae8b89aSjeremylt parameter(p=5) 45*cae8b89aSjeremylt parameter(q=8) 46*cae8b89aSjeremylt integer nx,nu 47*cae8b89aSjeremylt parameter(nx=nelem+1) 48*cae8b89aSjeremylt parameter(nu=nelem*(p-1)+1) 49*cae8b89aSjeremylt integer indx(nelem*2) 50*cae8b89aSjeremylt integer indu(nelem*p) 51*cae8b89aSjeremylt real*8 arrx(nx) 52*cae8b89aSjeremylt integer*8 voffset,xoffset 53*cae8b89aSjeremylt 54*cae8b89aSjeremylt real*8 hv(nu) 55*cae8b89aSjeremylt real*8 total 56*cae8b89aSjeremylt 57*cae8b89aSjeremylt character arg*32 58*cae8b89aSjeremylt 59*cae8b89aSjeremylt external setup,mass 60*cae8b89aSjeremylt 61*cae8b89aSjeremylt call getarg(1,arg) 62*cae8b89aSjeremylt call ceedinit(trim(arg)//char(0),ceed,err) 63*cae8b89aSjeremylt 64*cae8b89aSjeremylt do i=0,nx-1 65*cae8b89aSjeremylt arrx(i+1)=i/(nx-1.d0) 66*cae8b89aSjeremylt enddo 67*cae8b89aSjeremylt do i=0,nelem-1 68*cae8b89aSjeremylt indx(2*i+1)=i 69*cae8b89aSjeremylt indx(2*i+2)=i+1 70*cae8b89aSjeremylt enddo 71*cae8b89aSjeremylt 72*cae8b89aSjeremylt call ceedelemrestrictioncreate(ceed,nelem,2,nx,1,ceed_mem_host,& 73*cae8b89aSjeremylt & ceed_use_pointer,indx,erestrictx,err) 74*cae8b89aSjeremylt call ceedelemrestrictioncreateidentity(ceed,nelem,2,2*nelem,1,& 75*cae8b89aSjeremylt & erestrictxi,err) 76*cae8b89aSjeremylt 77*cae8b89aSjeremylt do i=0,nelem-1 78*cae8b89aSjeremylt do j=0,p-1 79*cae8b89aSjeremylt indu(p*i+j+1)=i*(p-1)+j 80*cae8b89aSjeremylt enddo 81*cae8b89aSjeremylt enddo 82*cae8b89aSjeremylt 83*cae8b89aSjeremylt call ceedelemrestrictioncreate(ceed,nelem,p,nu,1,ceed_mem_host,& 84*cae8b89aSjeremylt & ceed_use_pointer,indu,erestrictu,err) 85*cae8b89aSjeremylt call ceedelemrestrictioncreateidentity(ceed,nelem,q,q*nelem,1,& 86*cae8b89aSjeremylt & erestrictui,err) 87*cae8b89aSjeremylt 88*cae8b89aSjeremylt call ceedbasiscreatetensorh1lagrange(ceed,1,1,2,q,ceed_gauss,bx,err) 89*cae8b89aSjeremylt call ceedbasiscreatetensorh1lagrange(ceed,1,1,p,q,ceed_gauss,bu,err) 90*cae8b89aSjeremylt 91*cae8b89aSjeremylt call ceedqfunctioncreateinterior(ceed,1,setup,& 92*cae8b89aSjeremylt &SOURCE_DIR& 93*cae8b89aSjeremylt &//'t500-operator.h:setup'//char(0),qf_setup,err) 94*cae8b89aSjeremylt call ceedqfunctionaddinput(qf_setup,'_weight',1,ceed_eval_weight,err) 95*cae8b89aSjeremylt call ceedqfunctionaddinput(qf_setup,'dx',1,ceed_eval_grad,err) 96*cae8b89aSjeremylt call ceedqfunctionaddoutput(qf_setup,'rho',1,ceed_eval_none,err) 97*cae8b89aSjeremylt 98*cae8b89aSjeremylt call ceedqfunctioncreateinterior(ceed,1,mass,& 99*cae8b89aSjeremylt &SOURCE_DIR& 100*cae8b89aSjeremylt &//'t500-operator.h:mass'//char(0),qf_mass,err) 101*cae8b89aSjeremylt call ceedqfunctionaddinput(qf_mass,'rho',1,ceed_eval_none,err) 102*cae8b89aSjeremylt call ceedqfunctionaddinput(qf_mass,'u',1,ceed_eval_interp,err) 103*cae8b89aSjeremylt call ceedqfunctionaddoutput(qf_mass,'v',1,ceed_eval_interp,err) 104*cae8b89aSjeremylt 105*cae8b89aSjeremylt call ceedoperatorcreate(ceed,qf_setup,ceed_qfunction_none,& 106*cae8b89aSjeremylt & ceed_qfunction_none,op_setup,err) 107*cae8b89aSjeremylt call ceedoperatorcreate(ceed,qf_mass,ceed_qfunction_none,& 108*cae8b89aSjeremylt & ceed_qfunction_none,op_mass,err) 109*cae8b89aSjeremylt 110*cae8b89aSjeremylt call ceedvectorcreate(ceed,nx,x,err) 111*cae8b89aSjeremylt xoffset=0 112*cae8b89aSjeremylt call ceedvectorsetarray(x,ceed_mem_host,ceed_use_pointer,arrx,xoffset,err) 113*cae8b89aSjeremylt call ceedvectorcreate(ceed,nelem*q,qdata,err) 114*cae8b89aSjeremylt 115*cae8b89aSjeremylt call ceedoperatorsetfield(op_setup,'_weight',erestrictxi,& 116*cae8b89aSjeremylt & ceed_notranspose,bx,ceed_vector_none,err) 117*cae8b89aSjeremylt call ceedoperatorsetfield(op_setup,'dx',erestrictx,& 118*cae8b89aSjeremylt & ceed_notranspose,bx,ceed_vector_active,err) 119*cae8b89aSjeremylt call ceedoperatorsetfield(op_setup,'rho',erestrictui,& 120*cae8b89aSjeremylt & ceed_notranspose,ceed_basis_collocated,ceed_vector_active,err) 121*cae8b89aSjeremylt call ceedoperatorsetfield(op_mass,'rho',erestrictui,& 122*cae8b89aSjeremylt & ceed_notranspose,ceed_basis_collocated,qdata,err) 123*cae8b89aSjeremylt call ceedoperatorsetfield(op_mass,'u',erestrictu,& 124*cae8b89aSjeremylt & ceed_notranspose,bu,ceed_vector_active,err) 125*cae8b89aSjeremylt call ceedoperatorsetfield(op_mass,'v',erestrictu,& 126*cae8b89aSjeremylt & ceed_notranspose,bu,ceed_vector_active,err) 127*cae8b89aSjeremylt 128*cae8b89aSjeremylt call ceedoperatorapply(op_setup,x,qdata,ceed_request_immediate,err) 129*cae8b89aSjeremylt 130*cae8b89aSjeremylt call ceedvectorcreate(ceed,nu,u,err) 131*cae8b89aSjeremylt call ceedvectorsetvalue(u,1.d0,err) 132*cae8b89aSjeremylt call ceedvectorcreate(ceed,nu,v,err) 133*cae8b89aSjeremylt call ceedvectorsetvalue(v,0.d0,err) 134*cae8b89aSjeremylt 135*cae8b89aSjeremylt call ceedoperatorapplyadd(op_mass,u,v,ceed_request_immediate,err) 136*cae8b89aSjeremylt 137*cae8b89aSjeremylt call ceedvectorgetarrayread(v,ceed_mem_host,hv,voffset,err) 138*cae8b89aSjeremylt total=0. 139*cae8b89aSjeremylt do i=1,nu 140*cae8b89aSjeremylt total=total+hv(voffset+i) 141*cae8b89aSjeremylt enddo 142*cae8b89aSjeremylt if (abs(total-1.)>1.0d-10) then 143*cae8b89aSjeremylt! LCOV_EXCL_START 144*cae8b89aSjeremylt write(*,*) 'Computed Area: ',total,' != True Area: 1.0' 145*cae8b89aSjeremylt! LCOV_EXCL_STOP 146*cae8b89aSjeremylt endif 147*cae8b89aSjeremylt call ceedvectorrestorearrayread(v,hv,voffset,err) 148*cae8b89aSjeremylt 149*cae8b89aSjeremylt call ceedvectorsetvalue(v,1.d0,err) 150*cae8b89aSjeremylt call ceedoperatorapplyadd(op_mass,u,v,ceed_request_immediate,err) 151*cae8b89aSjeremylt 152*cae8b89aSjeremylt call ceedvectorgetarrayread(v,ceed_mem_host,hv,voffset,err) 153*cae8b89aSjeremylt total=-nu 154*cae8b89aSjeremylt do i=1,nu 155*cae8b89aSjeremylt total=total+hv(voffset+i) 156*cae8b89aSjeremylt enddo 157*cae8b89aSjeremylt if (abs(total-1.)>1.0d-10) then 158*cae8b89aSjeremylt! LCOV_EXCL_START 159*cae8b89aSjeremylt write(*,*) 'Computed Area: ',total,' != True Area: 1.0' 160*cae8b89aSjeremylt! LCOV_EXCL_STOP 161*cae8b89aSjeremylt endif 162*cae8b89aSjeremylt call ceedvectorrestorearrayread(v,hv,voffset,err) 163*cae8b89aSjeremylt 164*cae8b89aSjeremylt call ceedvectordestroy(qdata,err) 165*cae8b89aSjeremylt call ceedvectordestroy(x,err) 166*cae8b89aSjeremylt call ceedvectordestroy(u,err) 167*cae8b89aSjeremylt call ceedvectordestroy(v,err) 168*cae8b89aSjeremylt call ceedoperatordestroy(op_mass,err) 169*cae8b89aSjeremylt call ceedoperatordestroy(op_setup,err) 170*cae8b89aSjeremylt call ceedqfunctiondestroy(qf_mass,err) 171*cae8b89aSjeremylt call ceedqfunctiondestroy(qf_setup,err) 172*cae8b89aSjeremylt call ceedbasisdestroy(bu,err) 173*cae8b89aSjeremylt call ceedbasisdestroy(bx,err) 174*cae8b89aSjeremylt call ceedelemrestrictiondestroy(erestrictu,err) 175*cae8b89aSjeremylt call ceedelemrestrictiondestroy(erestrictx,err) 176*cae8b89aSjeremylt call ceedelemrestrictiondestroy(erestrictui,err) 177*cae8b89aSjeremylt call ceedelemrestrictiondestroy(erestrictxi,err) 178*cae8b89aSjeremylt call ceeddestroy(ceed,err) 179*cae8b89aSjeremylt end 180*cae8b89aSjeremylt!----------------------------------------------------------------------- 181