1!----------------------------------------------------------------------- 2! 3! Header with QFunctions 4! 5 include 't502-operator-f.h' 6!----------------------------------------------------------------------- 7 program test 8 implicit none 9 include 'ceed/fortran.h' 10 11 integer ceed,err,i,j 12 integer stridesu_small(3),stridesu_large(3) 13 integer erestrictx,erestrictu 14 integer erestrictui_small,erestrictui_large 15 integer bx_small,bu_small,bx_large,bu_large 16 integer qf_setup,qf_mass 17 integer op_setup_small,op_mass_small,op_setup_large,op_mass_large 18 integer qdata_small,qdata_large,x,u,v 19 integer nelem,p,q,scale 20 parameter(nelem=15) 21 parameter(p=5) 22 parameter(q=8) 23 parameter(scale=3) 24 integer nx,nu 25 parameter(nx=nelem+1) 26 parameter(nu=nelem*(p-1)+1) 27 integer indx(nelem*2) 28 integer indu(nelem*p) 29 real*8 arrx(nx) 30 integer*8 voffset,xoffset 31 32 real*8 hu(nu*2),hv(nu*2) 33 real*8 total1,total2 34 35 character arg*32 36 37 external setup,mass 38 39 call getarg(1,arg) 40 call ceedinit(trim(arg)//char(0),ceed,err) 41 42 do i=0,nx-1 43 arrx(i+1)=i/(nx-1.d0) 44 enddo 45 do i=0,nelem-1 46 indx(2*i+1)=i 47 indx(2*i+2)=i+1 48 enddo 49 50 call ceedelemrestrictioncreate(ceed,nelem,2,1,1,nx,ceed_mem_host,& 51 & ceed_use_pointer,indx,erestrictx,err) 52 53 do i=0,nelem-1 54 do j=0,p-1 55 indu(p*i+j+1)=2*(i*(p-1)+j) 56 enddo 57 enddo 58 59 call ceedelemrestrictioncreate(ceed,nelem,p,2,1,2*nu,ceed_mem_host,& 60 & ceed_use_pointer,indu,erestrictu,err) 61 stridesu_small=[1,q,q] 62 call ceedelemrestrictioncreatestrided(ceed,nelem,q,1,q*nelem,& 63 & stridesu_small,erestrictui_small,err) 64 stridesu_large=[1,q*scale,q*scale] 65 call ceedelemrestrictioncreatestrided(ceed,nelem,q*scale,1,& 66 & q*nelem*scale,stridesu_large,erestrictui_large,err) 67 68 call ceedbasiscreatetensorh1lagrange(ceed,1,1,2,q,ceed_gauss,bx_small,err) 69 call ceedbasiscreatetensorh1lagrange(ceed,1,2,p,q,ceed_gauss,bu_small,err) 70 call ceedbasiscreatetensorh1lagrange(ceed,1,1,2,q*scale,& 71 & ceed_gauss,bx_large,err) 72 call ceedbasiscreatetensorh1lagrange(ceed,1,2,p,q*scale,& 73 & ceed_gauss,bu_large,err) 74 75! Common QFunctions 76 77 call ceedqfunctioncreateinterior(ceed,1,setup,& 78 &SOURCE_DIR& 79 &//'t502-operator.h:setup'//char(0),qf_setup,err) 80 call ceedqfunctionaddinput(qf_setup,'weight',1,ceed_eval_weight,err) 81 call ceedqfunctionaddinput(qf_setup,'x',1,ceed_eval_grad,err) 82 call ceedqfunctionaddoutput(qf_setup,'rho',1,ceed_eval_none,err) 83 84 call ceedqfunctioncreateinterior(ceed,1,mass,& 85 &SOURCE_DIR& 86 &//'t502-operator.h:mass'//char(0),qf_mass,err) 87 call ceedqfunctionaddinput(qf_mass,'rho',1,ceed_eval_none,err) 88 call ceedqfunctionaddinput(qf_mass,'u',2,ceed_eval_interp,err) 89 call ceedqfunctionaddoutput(qf_mass,'v',2,ceed_eval_interp,err) 90 91 call ceedvectorcreate(ceed,nx,x,err) 92 xoffset=0 93 call ceedvectorsetarray(x,ceed_mem_host,ceed_use_pointer,arrx,xoffset,err) 94 95! Small operator 96 97 call ceedoperatorcreate(ceed,qf_setup,ceed_qfunction_none,& 98 & ceed_qfunction_none,op_setup_small,err) 99 call ceedoperatorcreate(ceed,qf_mass,ceed_qfunction_none,& 100 & ceed_qfunction_none,op_mass_small,err) 101 102 call ceedvectorcreate(ceed,nelem*q,qdata_small,err) 103 104 call ceedoperatorsetfield(op_setup_small,'weight',& 105 & ceed_elemrestriction_none,bx_small,ceed_vector_none,err) 106 call ceedoperatorsetfield(op_setup_small,'x',erestrictx,& 107 & bx_small,ceed_vector_active,err) 108 call ceedoperatorsetfield(op_setup_small,'rho',erestrictui_small,& 109 ceed_basis_none,ceed_vector_active,err) 110 call ceedoperatorsetfield(op_mass_small,'rho',erestrictui_small,& 111 ceed_basis_none,qdata_small,err) 112 call ceedoperatorsetfield(op_mass_small,'u',erestrictu,& 113 & bu_small,ceed_vector_active,err) 114 call ceedoperatorsetfield(op_mass_small,'v',erestrictu,& 115 & bu_small,ceed_vector_active,err) 116 117! Large operator 118 119 call ceedoperatorcreate(ceed,qf_setup,ceed_qfunction_none,& 120 & ceed_qfunction_none,op_setup_large,err) 121 call ceedoperatorcreate(ceed,qf_mass,ceed_qfunction_none,& 122 & ceed_qfunction_none,op_mass_large,err) 123 124 call ceedvectorcreate(ceed,nelem*q*scale,qdata_large,err) 125 126 call ceedoperatorsetfield(op_setup_large,'weight',& 127 & ceed_elemrestriction_none,bx_large,ceed_vector_none,err) 128 call ceedoperatorsetfield(op_setup_large,'x',erestrictx,& 129 & bx_large,ceed_vector_active,err) 130 call ceedoperatorsetfield(op_setup_large,'rho',erestrictui_large,& 131 ceed_basis_none,ceed_vector_active,err) 132 call ceedoperatorsetfield(op_mass_large,'rho',erestrictui_large,& 133 ceed_basis_none,qdata_large,err) 134 call ceedoperatorsetfield(op_mass_large,'u',erestrictu,& 135 & bu_large,ceed_vector_active,err) 136 call ceedoperatorsetfield(op_mass_large,'v',erestrictu,& 137 & bu_large,ceed_vector_active,err) 138 139! Setup U, V 140 141 call ceedvectorcreate(ceed,2*nu,u,err) 142 call ceedvectorgetarraywrite(u,ceed_mem_host,hu,voffset,err) 143 do i=1,nu 144 hu(voffset+2*i-1)=1. 145 hu(voffset+2*i)=2. 146 enddo 147 call ceedvectorrestorearray(u,hu,voffset,err) 148 call ceedvectorcreate(ceed,2*nu,v,err) 149 150! Small apply 151 152 call ceedoperatorapply(op_setup_small,x,qdata_small,& 153 & ceed_request_immediate,err) 154 call ceedoperatorapply(op_mass_small,u,v,ceed_request_immediate,err) 155 156 call ceedvectorgetarrayread(v,ceed_mem_host,hv,voffset,err) 157 total1=0. 158 total2=0. 159 do i=1,nu 160 total1=total1+hv(voffset+2*i-1) 161 total2=total2+hv(voffset+2*i) 162 enddo 163 if (abs(total1-1.)>1.0d-10) then 164! LCOV_EXCL_START 165 write(*,*) 'Computed Area: ',total1,' != True Area: 1.0' 166! LCOV_EXCL_STOP 167 endif 168 if (abs(total2-2.)>1.0d-10) then 169! LCOV_EXCL_START 170 write(*,*) 'Computed Area: ',total2,' != True Area: 2.0' 171! LCOV_EXCL_STOP 172 endif 173 call ceedvectorrestorearrayread(v,hv,voffset,err) 174 175! Large apply 176 177 call ceedoperatorapply(op_setup_large,x,qdata_large,& 178 & ceed_request_immediate,err) 179 call ceedoperatorapply(op_mass_large,u,v,ceed_request_immediate,err) 180 181 call ceedvectorgetarrayread(v,ceed_mem_host,hv,voffset,err) 182 total1=0. 183 total2=0. 184 do i=1,nu 185 total1=total1+hv(voffset+2*i-1) 186 total2=total2+hv(voffset+2*i) 187 enddo 188 if (abs(total1-1.)>1.0d-10) then 189! LCOV_EXCL_START 190 write(*,*) 'Computed Area: ',total1,' != True Area: 1.0' 191! LCOV_EXCL_STOP 192 endif 193 if (abs(total2-2.)>1.0d-10) then 194! LCOV_EXCL_START 195 write(*,*) 'Computed Area: ',total2,' != True Area: 2.0' 196! LCOV_EXCL_STOP 197 endif 198 call ceedvectorrestorearrayread(v,hv,voffset,err) 199 200 call ceedvectordestroy(qdata_small,err) 201 call ceedvectordestroy(qdata_large,err) 202 call ceedvectordestroy(x,err) 203 call ceedvectordestroy(u,err) 204 call ceedvectordestroy(v,err) 205 call ceedoperatordestroy(op_mass_small,err) 206 call ceedoperatordestroy(op_setup_small,err) 207 call ceedoperatordestroy(op_mass_large,err) 208 call ceedoperatordestroy(op_setup_large,err) 209 call ceedqfunctiondestroy(qf_mass,err) 210 call ceedqfunctiondestroy(qf_setup,err) 211 call ceedbasisdestroy(bu_small,err) 212 call ceedbasisdestroy(bx_small,err) 213 call ceedbasisdestroy(bu_large,err) 214 call ceedbasisdestroy(bx_large,err) 215 call ceedelemrestrictiondestroy(erestrictu,err) 216 call ceedelemrestrictiondestroy(erestrictx,err) 217 call ceedelemrestrictiondestroy(erestrictui_small,err) 218 call ceedelemrestrictiondestroy(erestrictui_large,err) 219 call ceeddestroy(ceed,err) 220 end 221!----------------------------------------------------------------------- 222 223