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