xref: /libCEED/tests/t503-operator-f.f90 (revision ccedf6b0eedde5bdb26f1d628a64390ea2c01243)
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      enddo
28
29      ierr=0
30      end
31!-----------------------------------------------------------------------
32      program test
33
34      include 'ceedf.h'
35
36      integer ceed,err,i,j
37      integer imode
38      parameter(imode=ceed_noninterlaced)
39      integer stridesx(3),stridesu(3)
40      integer erestrictx,erestrictu,erestrictxi,erestrictui
41      integer bx,bu
42      integer qf_setup,qf_mass
43      integer op_setup,op_mass
44      integer qdata,x,u,v
45      integer nelem,p,q
46      parameter(nelem=15)
47      parameter(p=5)
48      parameter(q=8)
49      integer nx,nu
50      parameter(nx=nelem+1)
51      parameter(nu=nelem*(p-1)+1)
52      integer indx(nelem*2)
53      integer indu(nelem*p)
54      real*8 arrx(nx)
55      integer*8 voffset,xoffset
56
57      real*8 hv(nu)
58      real*8 total
59
60      character arg*32
61
62      external setup,mass
63
64      call getarg(1,arg)
65      call ceedinit(trim(arg)//char(0),ceed,err)
66
67      call ceedvectorcreate(ceed,nu,u,err)
68      call ceedvectorsetvalue(u,1.d0,err)
69
70      call ceedvectorcreate(ceed,nu,v,err)
71
72      do i=0,nx-1
73        arrx(i+1)=i/(nx-1.d0)
74      enddo
75      call ceedvectorcreate(ceed,nx,x,err)
76      xoffset=0
77      call ceedvectorsetarray(x,ceed_mem_host,ceed_use_pointer,arrx,xoffset,err)
78
79      call ceedvectorcreate(ceed,nelem*q,qdata,err)
80
81      do i=0,nelem-1
82        indx(2*i+1)=i
83        indx(2*i+2)=i+1
84      enddo
85      call ceedelemrestrictioncreate(ceed,imode,nelem,2,nx,1,ceed_mem_host,&
86     & ceed_use_pointer,indx,erestrictx,err)
87      stridesx=[1,2,2]
88      call ceedelemrestrictioncreatestrided(ceed,nelem,2,2*nelem,1,stridesx,&
89     & erestrictxi,err)
90
91      do i=0,nelem-1
92        do j=0,p-1
93          indu(p*i+j+1)=i*(p-1)+j
94        enddo
95      enddo
96
97      call ceedelemrestrictioncreate(ceed,imode,nelem,p,nu,1,ceed_mem_host,&
98     & ceed_use_pointer,indu,erestrictu,err)
99      stridesu=[1,q,q]
100      call ceedelemrestrictioncreatestrided(ceed,nelem,q,q*nelem,1,stridesu,&
101     & erestrictui,err)
102
103      call ceedbasiscreatetensorh1lagrange(ceed,1,1,2,q,ceed_gauss,bx,err)
104      call ceedbasiscreatetensorh1lagrange(ceed,1,1,p,q,ceed_gauss,bu,err)
105
106      call ceedqfunctioncreateinterior(ceed,1,setup,&
107     &SOURCE_DIR&
108     &//'t500-operator.h:setup'//char(0),qf_setup,err)
109      call ceedqfunctionaddinput(qf_setup,'_weight',1,ceed_eval_weight,err)
110      call ceedqfunctionaddinput(qf_setup,'dx',1,ceed_eval_grad,err)
111      call ceedqfunctionaddoutput(qf_setup,'rho',1,ceed_eval_none,err)
112
113      call ceedqfunctioncreateinterior(ceed,1,mass,&
114     &SOURCE_DIR&
115     &//'t500-operator.h:mass'//char(0),qf_mass,err)
116      call ceedqfunctionaddinput(qf_mass,'rho',1,ceed_eval_none,err)
117      call ceedqfunctionaddinput(qf_mass,'u',1,ceed_eval_interp,err)
118      call ceedqfunctionaddoutput(qf_mass,'v',1,ceed_eval_interp,err)
119
120      call ceedoperatorcreate(ceed,qf_setup,ceed_qfunction_none,&
121     & ceed_qfunction_none,op_setup,err)
122      call ceedoperatorcreate(ceed,qf_mass,ceed_qfunction_none,&
123     & ceed_qfunction_none,op_mass,err)
124
125      call ceedoperatorsetfield(op_setup,'_weight',erestrictxi,bx,&
126     & ceed_vector_none,err)
127      call ceedoperatorsetfield(op_setup,'dx',erestrictx,bx,x,err)
128      call ceedoperatorsetfield(op_setup,'rho',erestrictui,&
129     & ceed_basis_collocated,qdata,err)
130      call ceedoperatorsetfield(op_mass,'rho',erestrictui,&
131     & ceed_basis_collocated,qdata,err)
132      call ceedoperatorsetfield(op_mass,'u',erestrictu,bu,u,err)
133      call ceedoperatorsetfield(op_mass,'v',erestrictu,bu,v,err)
134
135      call ceedoperatorapply(op_setup,ceed_vector_none,ceed_vector_none,&
136     & ceed_request_immediate,err)
137      call ceedoperatorapply(op_mass,ceed_vector_none,ceed_vector_none,&
138     & ceed_request_immediate,err)
139
140      call ceedvectorgetarrayread(v,ceed_mem_host,hv,voffset,err)
141      total=0.
142      do i=1,nu
143        total=total+hv(voffset+i)
144      enddo
145      if (abs(total-1.)>1.0d-10) then
146! LCOV_EXCL_START
147        write(*,*) 'Computed Area: ',total,' != True Area: 1.0'
148! LCOV_EXCL_STOP
149      endif
150      call ceedvectorrestorearrayread(v,hv,voffset,err)
151
152      call ceedvectordestroy(qdata,err)
153      call ceedvectordestroy(x,err)
154      call ceedvectordestroy(u,err)
155      call ceedvectordestroy(v,err)
156      call ceedoperatordestroy(op_mass,err)
157      call ceedoperatordestroy(op_setup,err)
158      call ceedqfunctiondestroy(qf_mass,err)
159      call ceedqfunctiondestroy(qf_setup,err)
160      call ceedbasisdestroy(bu,err)
161      call ceedbasisdestroy(bx,err)
162      call ceedelemrestrictiondestroy(erestrictu,err)
163      call ceedelemrestrictiondestroy(erestrictx,err)
164      call ceedelemrestrictiondestroy(erestrictui,err)
165      call ceedelemrestrictiondestroy(erestrictxi,err)
166      call ceeddestroy(ceed,err)
167      end
168!-----------------------------------------------------------------------
169