xref: /libCEED/tests/t536-operator-f.f90 (revision 2bba3ffaf13cc90ae8de830effb3b11a0650a2fc)
1b7ec98d8SJeremy L Thompson!-----------------------------------------------------------------------
2b7ec98d8SJeremy L Thompson!
3b7ec98d8SJeremy L Thompson! Header with common subroutine
4b7ec98d8SJeremy L Thompson!
5b7ec98d8SJeremy L Thompson      include 't320-basis-f.h'
6b7ec98d8SJeremy L Thompson!-----------------------------------------------------------------------
7b7ec98d8SJeremy L Thompson      subroutine setup_mass(ctx,q,u1,u2,u3,u4,u5,u6,u7,u8,u9,u10,u11,u12,u13,&
8b7ec98d8SJeremy L Thompson&           u14,u15,u16,v1,v2,v3,v4,v5,v6,v7,v8,v9,v10,v11,v12,v13,v14,v15,&
9b7ec98d8SJeremy L Thompson&           v16,ierr)
10b7ec98d8SJeremy L Thompson      real*8 ctx
11b7ec98d8SJeremy L Thompson      real*8 u1(1)
12b7ec98d8SJeremy L Thompson      real*8 u2(1)
13b7ec98d8SJeremy L Thompson      real*8 v1(1)
14b7ec98d8SJeremy L Thompson      integer q,ierr
15b7ec98d8SJeremy L Thompson
16b7ec98d8SJeremy L Thompson      do i=1,q
17b7ec98d8SJeremy L Thompson        v1(i)=u2(i)*(u1(i+q*0)*u1(i+q*3)-u1(i+q*1)*u1(i+q*2))
18b7ec98d8SJeremy L Thompson      enddo
19b7ec98d8SJeremy L Thompson
20b7ec98d8SJeremy L Thompson      ierr=0
21b7ec98d8SJeremy L Thompson      end
22b7ec98d8SJeremy L Thompson!-----------------------------------------------------------------------
23b7ec98d8SJeremy L Thompson      subroutine setup_diff(ctx,q,u1,u2,u3,u4,u5,u6,u7,u8,u9,u10,u11,u12,u13,&
24b7ec98d8SJeremy L Thompson&           u14,u15,u16,v1,v2,v3,v4,v5,v6,v7,v8,v9,v10,v11,v12,v13,v14,v15,&
25b7ec98d8SJeremy L Thompson&           v16,ierr)
26b7ec98d8SJeremy L Thompson      real*8 ctx
27b7ec98d8SJeremy L Thompson      real*8 u1(1)
28b7ec98d8SJeremy L Thompson      real*8 u2(1)
29b7ec98d8SJeremy L Thompson      real*8 v1(1)
30b7ec98d8SJeremy L Thompson      real*8 w
31b7ec98d8SJeremy L Thompson      integer q,ierr
32b7ec98d8SJeremy L Thompson
33b7ec98d8SJeremy L Thompson      do i=1,q
34b7ec98d8SJeremy L Thompson        w=u2(i)/(u1(i+q*0)*u1(i+q*3)-u1(i+q*1)*u1(i+q*2))
35b7ec98d8SJeremy L Thompson        v1(i+q*0)=w*(u1(i+q*2)*u1(i+q*2)+u1(i+q*3)*u1(i+q*3))
36b7ec98d8SJeremy L Thompson        v1(i+q*1)=w*(u1(i+q*0)*u1(i+q*0)+u1(i+q*1)*u1(i+q*1))
37b7ec98d8SJeremy L Thompson        v1(i+q*2)=-w*(u1(i+q*0)*u1(i+q*2)+u1(i+q*2)*u1(i+q*3))
38b7ec98d8SJeremy L Thompson      enddo
39b7ec98d8SJeremy L Thompson
40b7ec98d8SJeremy L Thompson      ierr=0
41b7ec98d8SJeremy L Thompson      end
42b7ec98d8SJeremy L Thompson!-----------------------------------------------------------------------
43b7ec98d8SJeremy L Thompson      subroutine apply(ctx,q,u1,u2,u3,u4,u5,u6,u7,u8,u9,u10,u11,u12,u13,u14,&
44b7ec98d8SJeremy L Thompson&           u15,u16,v1,v2,v3,v4,v5,v6,v7,v8,v9,v10,v11,v12,v13,v14,v15,v16,ierr)
45b7ec98d8SJeremy L Thompson      real*8 ctx
46b7ec98d8SJeremy L Thompson      real*8 u1(1)
47b7ec98d8SJeremy L Thompson      real*8 u2(1)
48b7ec98d8SJeremy L Thompson      real*8 u3(1)
49b7ec98d8SJeremy L Thompson      real*8 u4(1)
50b7ec98d8SJeremy L Thompson      real*8 v1(1)
51b7ec98d8SJeremy L Thompson      real*8 v2(1)
52b7ec98d8SJeremy L Thompson      real*8 du0,du1
53b7ec98d8SJeremy L Thompson      integer q,ierr
54b7ec98d8SJeremy L Thompson
55b7ec98d8SJeremy L Thompson      do i=1,q
56b7ec98d8SJeremy L Thompson!       mass
57b7ec98d8SJeremy L Thompson        v1(i) = u2(i)*u4(i)
58b7ec98d8SJeremy L Thompson!       diff
59b7ec98d8SJeremy L Thompson        du0=u1(i+q*0)
60b7ec98d8SJeremy L Thompson        du1=u1(i+q*1)
61b7ec98d8SJeremy L Thompson        v2(i+q*0)=u3(i+q*0)*du0+u3(i+q*2)*du1
62b7ec98d8SJeremy L Thompson        v2(i+q*1)=u3(i+q*2)*du0+u3(i+q*1)*du1
63b7ec98d8SJeremy L Thompson      enddo
64b7ec98d8SJeremy L Thompson
65b7ec98d8SJeremy L Thompson      ierr=0
66b7ec98d8SJeremy L Thompson      end
67b7ec98d8SJeremy L Thompson!-----------------------------------------------------------------------
68b7ec98d8SJeremy L Thompson      program test
691f9a83abSJed Brown      implicit none
70b7ec98d8SJeremy L Thompson      include 'ceedf.h'
71b7ec98d8SJeremy L Thompson
72b7ec98d8SJeremy L Thompson      integer ceed,err,i
7315910d16Sjeremylt      integer stridesu(3),stridesqd(3)
7415910d16Sjeremylt      integer erestrictx,erestrictu,erestrictui,erestrictqi
75b7ec98d8SJeremy L Thompson      integer bx,bu
76b7ec98d8SJeremy L Thompson      integer qf_setup_mass,qf_setup_diff,qf_apply
77b7ec98d8SJeremy L Thompson      integer op_setup_mass,op_setup_diff,op_apply
78b7ec98d8SJeremy L Thompson      integer qdata_mass,qdata_diff,x,a,u,v
79b7ec98d8SJeremy L Thompson      integer nelem,p,q,d
80b7ec98d8SJeremy L Thompson      integer row,col,offset
81b7ec98d8SJeremy L Thompson      parameter(nelem=12)
82b7ec98d8SJeremy L Thompson      parameter(p=6)
83b7ec98d8SJeremy L Thompson      parameter(q=4)
84b7ec98d8SJeremy L Thompson      parameter(d=2)
85b7ec98d8SJeremy L Thompson      integer ndofs,nqpts,nx,ny
86b7ec98d8SJeremy L Thompson      parameter(nx=3)
87b7ec98d8SJeremy L Thompson      parameter(ny=2)
88b7ec98d8SJeremy L Thompson      parameter(ndofs=(nx*2+1)*(ny*2+1))
89a8c028e3SNatalie Beams      parameter(nqpts=nelem*q)
90b7ec98d8SJeremy L Thompson      integer indx(nelem*p*p)
91b7ec98d8SJeremy L Thompson      real*8 arrx(d*ndofs),aa(nqpts),uu(ndofs),vv(ndofs),atrue(ndofs)
92b7ec98d8SJeremy L Thompson      integer*8 xoffset,aoffset,uoffset,voffset
93b7ec98d8SJeremy L Thompson
94b7ec98d8SJeremy L Thompson      real*8 qref(d*q)
95b7ec98d8SJeremy L Thompson      real*8 qweight(q)
96b7ec98d8SJeremy L Thompson      real*8 interp(p*q)
97b7ec98d8SJeremy L Thompson      real*8 grad(d*p*q)
981f9a83abSJed Brown      real*8 val
99b7ec98d8SJeremy L Thompson      character arg*32
100b7ec98d8SJeremy L Thompson
1013bd813ffSjeremylt      external setup_mass,setup_diff,apply
102b7ec98d8SJeremy L Thompson
103b7ec98d8SJeremy L Thompson      call getarg(1,arg)
104b7ec98d8SJeremy L Thompson
105b7ec98d8SJeremy L Thompson      call ceedinit(trim(arg)//char(0),ceed,err)
106b7ec98d8SJeremy L Thompson
107b7ec98d8SJeremy L Thompson! DoF Coordinates
108b7ec98d8SJeremy L Thompson      do i=0,ndofs-1
109b7ec98d8SJeremy L Thompson        arrx(i+1)=mod(i,(nx*2+1))
110b7ec98d8SJeremy L Thompson        arrx(i+1)=arrx(i+1)*(1.d0/(nx*2.d0))
111b7ec98d8SJeremy L Thompson        val=(i/(nx*2+1))
112b7ec98d8SJeremy L Thompson        arrx(i+1+ndofs)=val*(1.d0/(ny*2.d0))
113b7ec98d8SJeremy L Thompson      enddo
114b7ec98d8SJeremy L Thompson      call ceedvectorcreate(ceed,d*ndofs,x,err)
115b7ec98d8SJeremy L Thompson      xoffset=0
116b7ec98d8SJeremy L Thompson      call ceedvectorsetarray(x,ceed_mem_host,ceed_use_pointer,arrx,xoffset,err)
117b7ec98d8SJeremy L Thompson
118b7ec98d8SJeremy L Thompson! Qdata Vector
119b7ec98d8SJeremy L Thompson      call ceedvectorcreate(ceed,nqpts,qdata_mass,err)
120b7ec98d8SJeremy L Thompson      call ceedvectorcreate(ceed,nqpts*d*(d+1)/2,qdata_diff,err)
121b7ec98d8SJeremy L Thompson
122b7ec98d8SJeremy L Thompson! Element Setup
123b7ec98d8SJeremy L Thompson      do i=0,5
124b7ec98d8SJeremy L Thompson        col=mod(i,nx)
125b7ec98d8SJeremy L Thompson        row=i/nx
126b7ec98d8SJeremy L Thompson        offset=col*2+row*(nx*2+1)*2
127b7ec98d8SJeremy L Thompson
128b7ec98d8SJeremy L Thompson        indx(i*2*p+1)=2+offset
129b7ec98d8SJeremy L Thompson        indx(i*2*p+2)=9+offset
130b7ec98d8SJeremy L Thompson        indx(i*2*p+3)=16+offset
131b7ec98d8SJeremy L Thompson        indx(i*2*p+4)=1+offset
132b7ec98d8SJeremy L Thompson        indx(i*2*p+5)=8+offset
133b7ec98d8SJeremy L Thompson        indx(i*2*p+6)=0+offset
134b7ec98d8SJeremy L Thompson
135b7ec98d8SJeremy L Thompson        indx(i*2*p+7)=14+offset
136b7ec98d8SJeremy L Thompson        indx(i*2*p+8)=7+offset
137b7ec98d8SJeremy L Thompson        indx(i*2*p+9)=0+offset
138b7ec98d8SJeremy L Thompson        indx(i*2*p+10)=15+offset
139b7ec98d8SJeremy L Thompson        indx(i*2*p+11)=8+offset
140b7ec98d8SJeremy L Thompson        indx(i*2*p+12)=16+offset
141b7ec98d8SJeremy L Thompson      enddo
142b7ec98d8SJeremy L Thompson
143b7ec98d8SJeremy L Thompson! Restrictions
144d979a051Sjeremylt      call ceedelemrestrictioncreate(ceed,nelem,p,d,ndofs,d*ndofs,&
1457509a596Sjeremylt     & ceed_mem_host,ceed_use_pointer,indx,erestrictx,err)
146b7ec98d8SJeremy L Thompson
147d979a051Sjeremylt      call ceedelemrestrictioncreate(ceed,nelem,p,1,1,ndofs,&
1487509a596Sjeremylt     & ceed_mem_host,ceed_use_pointer,indx,erestrictu,err)
1497509a596Sjeremylt      stridesu=[1,q,q]
150d979a051Sjeremylt      call ceedelemrestrictioncreatestrided(ceed,nelem,q,1,nqpts,&
151d979a051Sjeremylt     & stridesu,erestrictui,err)
152b7ec98d8SJeremy L Thompson
1537509a596Sjeremylt      stridesqd=[1,q,q*d*(d+1)/2]
154d979a051Sjeremylt      call ceedelemrestrictioncreatestrided(ceed,nelem,q,d*(d+1)/2,&
155d979a051Sjeremylt     & d*(d+1)/2*nqpts,stridesqd,erestrictqi,err)
156b7ec98d8SJeremy L Thompson
157b7ec98d8SJeremy L Thompson! Bases
158b7ec98d8SJeremy L Thompson      call buildmats(qref,qweight,interp,grad)
159b7ec98d8SJeremy L Thompson      call ceedbasiscreateh1(ceed,ceed_triangle,d,p,q,interp,grad,qref,qweight,&
160b7ec98d8SJeremy L Thompson     & bx,err)
161b7ec98d8SJeremy L Thompson      call buildmats(qref,qweight,interp,grad)
162b7ec98d8SJeremy L Thompson      call ceedbasiscreateh1(ceed,ceed_triangle,1,p,q,interp,grad,qref,qweight,&
163b7ec98d8SJeremy L Thompson     & bu,err)
164b7ec98d8SJeremy L Thompson
165b7ec98d8SJeremy L Thompson! QFunction - setup mass
166b7ec98d8SJeremy L Thompson      call ceedqfunctioncreateinterior(ceed,1,setup_mass,&
167b7ec98d8SJeremy L Thompson     &SOURCE_DIR&
168b7ec98d8SJeremy L Thompson     &//'t532-operator.h:setup_mass'//char(0),qf_setup_mass,err)
169b7ec98d8SJeremy L Thompson      call ceedqfunctionaddinput(qf_setup_mass,'dx',d*d,ceed_eval_grad,err)
170b7ec98d8SJeremy L Thompson      call ceedqfunctionaddinput(qf_setup_mass,'_weight',1,ceed_eval_weight,err)
171b7ec98d8SJeremy L Thompson      call ceedqfunctionaddoutput(qf_setup_mass,'qdata',1,ceed_eval_none,err)
172b7ec98d8SJeremy L Thompson
173b7ec98d8SJeremy L Thompson! Operator - setup mass
174442e7f0bSjeremylt      call ceedoperatorcreate(ceed,qf_setup_mass,ceed_qfunction_none,&
175442e7f0bSjeremylt     & ceed_qfunction_none,op_setup_mass,err)
176b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_setup_mass,'dx',erestrictx,&
177a8d32208Sjeremylt     & bx,ceed_vector_active,err)
17815910d16Sjeremylt      call ceedoperatorsetfield(op_setup_mass,'_weight',&
17915910d16Sjeremylt     & ceed_elemrestriction_none,bx,ceed_vector_none,err)
180b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_setup_mass,'qdata',erestrictui,&
181a8d32208Sjeremylt     & ceed_basis_collocated,ceed_vector_active,err)
182b7ec98d8SJeremy L Thompson
183b7ec98d8SJeremy L Thompson! QFunction - setup diff
184b7ec98d8SJeremy L Thompson      call ceedqfunctioncreateinterior(ceed,1,setup_diff,&
185b7ec98d8SJeremy L Thompson     &SOURCE_DIR&
186b7ec98d8SJeremy L Thompson     &//'t532-operator.h:setup_diff'//char(0),qf_setup_diff,err)
187b7ec98d8SJeremy L Thompson      call ceedqfunctionaddinput(qf_setup_diff,'dx',d*d,ceed_eval_grad,err)
188b7ec98d8SJeremy L Thompson      call ceedqfunctionaddinput(qf_setup_diff,'_weight',1,ceed_eval_weight,err)
189b7ec98d8SJeremy L Thompson      call ceedqfunctionaddoutput(qf_setup_diff,'qdata',&
190b7ec98d8SJeremy L Thompson     & d*(d+1)/2,ceed_eval_none,err)
191b7ec98d8SJeremy L Thompson
192b7ec98d8SJeremy L Thompson! Operator - setup diff
193442e7f0bSjeremylt      call ceedoperatorcreate(ceed,qf_setup_diff,ceed_qfunction_none,&
194442e7f0bSjeremylt     & ceed_qfunction_none,op_setup_diff,err)
195b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_setup_diff,'dx',erestrictx,&
196a8d32208Sjeremylt     & bx,ceed_vector_active,err)
19715910d16Sjeremylt      call ceedoperatorsetfield(op_setup_diff,'_weight',&
19815910d16Sjeremylt     & ceed_elemrestriction_none,bx,ceed_vector_none,err)
199b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_setup_diff,'qdata',erestrictqi,&
200a8d32208Sjeremylt     & ceed_basis_collocated,ceed_vector_active,err)
201b7ec98d8SJeremy L Thompson
202b7ec98d8SJeremy L Thompson! Apply Setup Operators
203b7ec98d8SJeremy L Thompson      call ceedoperatorapply(op_setup_mass,x,qdata_mass,&
204b7ec98d8SJeremy L Thompson     & ceed_request_immediate,err)
205b7ec98d8SJeremy L Thompson      call ceedoperatorapply(op_setup_diff,x,qdata_diff,&
206b7ec98d8SJeremy L Thompson     & ceed_request_immediate,err)
207b7ec98d8SJeremy L Thompson
208b7ec98d8SJeremy L Thompson! QFunction - apply
209b7ec98d8SJeremy L Thompson      call ceedqfunctioncreateinterior(ceed,1,apply,&
210b7ec98d8SJeremy L Thompson     &SOURCE_DIR&
211b7ec98d8SJeremy L Thompson     &//'t532-operator.h:apply'//char(0),qf_apply,err)
212b7ec98d8SJeremy L Thompson      call ceedqfunctionaddinput(qf_apply,'du',d,ceed_eval_grad,err)
213b7ec98d8SJeremy L Thompson      call ceedqfunctionaddinput(qf_apply,'qdata_mass',1,ceed_eval_none,err)
214b7ec98d8SJeremy L Thompson      call ceedqfunctionaddinput(qf_apply,'qdata_diff',&
215b7ec98d8SJeremy L Thompson     & d*(d+1)/2,ceed_eval_none,err)
216b7ec98d8SJeremy L Thompson      call ceedqfunctionaddinput(qf_apply,'u',1,ceed_eval_interp,err)
217b7ec98d8SJeremy L Thompson      call ceedqfunctionaddoutput(qf_apply,'v',1,ceed_eval_interp,err)
218b7ec98d8SJeremy L Thompson      call ceedqfunctionaddoutput(qf_apply,'dv',d,ceed_eval_grad,err)
219b7ec98d8SJeremy L Thompson
220b7ec98d8SJeremy L Thompson! Operator - apply
221442e7f0bSjeremylt      call ceedoperatorcreate(ceed,qf_apply,ceed_qfunction_none,&
222442e7f0bSjeremylt     & ceed_qfunction_none,op_apply,err)
223b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_apply,'du',erestrictu,&
224a8d32208Sjeremylt     & bu,ceed_vector_active,err)
225b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_apply,'qdata_mass',erestrictui,&
226a8d32208Sjeremylt     & ceed_basis_collocated,qdata_mass,err)
227b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_apply,'qdata_diff',erestrictqi,&
228a8d32208Sjeremylt     & ceed_basis_collocated,qdata_diff,err)
229b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_apply,'u',erestrictu,&
230a8d32208Sjeremylt     & bu,ceed_vector_active,err)
231b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_apply,'v',erestrictu,&
232a8d32208Sjeremylt     & bu,ceed_vector_active,err)
233b7ec98d8SJeremy L Thompson      call ceedoperatorsetfield(op_apply,'dv',erestrictu,&
234a8d32208Sjeremylt     & bu,ceed_vector_active,err)
235b7ec98d8SJeremy L Thompson
236b7ec98d8SJeremy L Thompson! Assemble Diagonal
237*2bba3ffaSJeremy L Thompson      call ceedvectorcreate(ceed,ndofs,a,err)
23880ac2e43SJeremy L Thompson      call ceedoperatorlinearassemblediagonal(op_apply,a,&
239b7ec98d8SJeremy L Thompson     & ceed_request_immediate,err)
240b7ec98d8SJeremy L Thompson
241b7ec98d8SJeremy L Thompson! Manually assemble diagonal
242b7ec98d8SJeremy L Thompson      call ceedvectorcreate(ceed,ndofs,u,err)
243b7ec98d8SJeremy L Thompson      call ceedvectorsetvalue(u,0.d0,err)
244b7ec98d8SJeremy L Thompson      call ceedvectorcreate(ceed,ndofs,v,err)
245b7ec98d8SJeremy L Thompson      do i=1,ndofs
246b7ec98d8SJeremy L Thompson        call ceedvectorgetarray(u,ceed_mem_host,uu,uoffset,err)
247b7ec98d8SJeremy L Thompson        uu(i+uoffset)=1.d0
248b7ec98d8SJeremy L Thompson        if (i>1) then
249b7ec98d8SJeremy L Thompson          uu(i-1+uoffset)=0.d0
250b7ec98d8SJeremy L Thompson        endif
251b7ec98d8SJeremy L Thompson        call ceedvectorrestorearray(u,uu,uoffset,err)
252b7ec98d8SJeremy L Thompson
253b7ec98d8SJeremy L Thompson        call ceedoperatorapply(op_apply,u,v,ceed_request_immediate,err)
254b7ec98d8SJeremy L Thompson
255b7ec98d8SJeremy L Thompson        call ceedvectorgetarrayread(v,ceed_mem_host,vv,voffset,err)
256b7ec98d8SJeremy L Thompson        atrue(i)=vv(voffset+i)
257b7ec98d8SJeremy L Thompson        call ceedvectorrestorearrayread(v,vv,voffset,err)
258b7ec98d8SJeremy L Thompson      enddo
259b7ec98d8SJeremy L Thompson
260b7ec98d8SJeremy L Thompson! Check Output
261b7ec98d8SJeremy L Thompson      call ceedvectorgetarrayread(a,ceed_mem_host,aa,aoffset,err)
262b7ec98d8SJeremy L Thompson      do i=1,ndofs
263b7ec98d8SJeremy L Thompson        if (abs(aa(aoffset+i)-atrue(i))>1.0d-14) then
264b7ec98d8SJeremy L Thompson! LCOV_EXCL_START
265b7ec98d8SJeremy L Thompson          write(*,*) '[',i,'] Error in assembly: ',aa(aoffset+i),' != ',&
266b7ec98d8SJeremy L Thompson     &      atrue(i)
267b7ec98d8SJeremy L Thompson! LCOV_EXCL_STOP
268b7ec98d8SJeremy L Thompson        endif
269b7ec98d8SJeremy L Thompson      enddo
270b7ec98d8SJeremy L Thompson      call ceedvectorrestorearrayread(a,aa,aoffset,err)
271b7ec98d8SJeremy L Thompson
272b7ec98d8SJeremy L Thompson! Cleanup
273b7ec98d8SJeremy L Thompson      call ceedqfunctiondestroy(qf_setup_mass,err)
274b7ec98d8SJeremy L Thompson      call ceedqfunctiondestroy(qf_setup_diff,err)
275b7ec98d8SJeremy L Thompson      call ceedqfunctiondestroy(qf_apply,err)
276b7ec98d8SJeremy L Thompson      call ceedoperatordestroy(op_setup_mass,err)
277b7ec98d8SJeremy L Thompson      call ceedoperatordestroy(op_setup_diff,err)
278b7ec98d8SJeremy L Thompson      call ceedoperatordestroy(op_apply,err)
279b7ec98d8SJeremy L Thompson      call ceedelemrestrictiondestroy(erestrictu,err)
280b7ec98d8SJeremy L Thompson      call ceedelemrestrictiondestroy(erestrictx,err)
281b7ec98d8SJeremy L Thompson      call ceedelemrestrictiondestroy(erestrictui,err)
282b7ec98d8SJeremy L Thompson      call ceedelemrestrictiondestroy(erestrictqi,err)
283b7ec98d8SJeremy L Thompson      call ceedbasisdestroy(bu,err)
284b7ec98d8SJeremy L Thompson      call ceedbasisdestroy(bx,err)
285b7ec98d8SJeremy L Thompson      call ceedvectordestroy(x,err)
286b7ec98d8SJeremy L Thompson      call ceedvectordestroy(a,err)
287b7ec98d8SJeremy L Thompson      call ceedvectordestroy(u,err)
288b7ec98d8SJeremy L Thompson      call ceedvectordestroy(v,err)
289b7ec98d8SJeremy L Thompson      call ceedvectordestroy(qdata_mass,err)
290b7ec98d8SJeremy L Thompson      call ceedvectordestroy(qdata_diff,err)
291b7ec98d8SJeremy L Thompson      call ceeddestroy(ceed,err)
292b7ec98d8SJeremy L Thompson      end
293b7ec98d8SJeremy L Thompson!-----------------------------------------------------------------------
294