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