xref: /libCEED/tests/t313-basis-f.f90 (revision 52bfb9bbf17f17edbcd45876cdc8689a879bc683) 
18980d4a7Sjeremylt!-----------------------------------------------------------------------
2*52bfb9bbSJeremy L Thompson      subroutine eval(dimn,x,rslt)
3*52bfb9bbSJeremy L Thompson      integer dimn
4*52bfb9bbSJeremy L Thompson      real*8 x(1)
5*52bfb9bbSJeremy L Thompson      real*8 rslt
6*52bfb9bbSJeremy L Thompson      real*8 center
78980d4a7Sjeremylt
8*52bfb9bbSJeremy L Thompson      integer d
98980d4a7Sjeremylt
10*52bfb9bbSJeremy L Thompson      rslt=1
11*52bfb9bbSJeremy L Thompson      center=0.1
128980d4a7Sjeremylt
13*52bfb9bbSJeremy L Thompson      do d=1,dimn
14*52bfb9bbSJeremy L Thompson        rslt=rslt*tanh(x(d)-center)
15*52bfb9bbSJeremy L Thompson        center=center+0.1
16*52bfb9bbSJeremy L Thompson      enddo
178980d4a7Sjeremylt
188980d4a7Sjeremylt      end
198980d4a7Sjeremylt!-----------------------------------------------------------------------
208980d4a7Sjeremylt      program test
218980d4a7Sjeremylt
228980d4a7Sjeremylt      include 'ceedf.h'
238980d4a7Sjeremylt
248980d4a7Sjeremylt      integer ceed,err
25*52bfb9bbSJeremy L Thompson      integer x,xq,u,uq
26*52bfb9bbSJeremy L Thompson      integer bxl,bul,bxg,bug
27*52bfb9bbSJeremy L Thompson      integer dimn,d
28*52bfb9bbSJeremy L Thompson      integer i
29*52bfb9bbSJeremy L Thompson      integer q
30*52bfb9bbSJeremy L Thompson      parameter(q=10)
31*52bfb9bbSJeremy L Thompson      integer maxdim
32*52bfb9bbSJeremy L Thompson      parameter(maxdim=3)
33*52bfb9bbSJeremy L Thompson      integer qdimmax
34*52bfb9bbSJeremy L Thompson      parameter(qdimmax=q**maxdim)
35*52bfb9bbSJeremy L Thompson      integer xdimmax
36*52bfb9bbSJeremy L Thompson      parameter(xdimmax=2**maxdim)
37*52bfb9bbSJeremy L Thompson      integer qdim,xdim
388980d4a7Sjeremylt
39*52bfb9bbSJeremy L Thompson      real*8 xx(xdimmax*maxdim)
40*52bfb9bbSJeremy L Thompson      real*8 xxx(maxdim)
41*52bfb9bbSJeremy L Thompson      real*8 xxq(qdimmax*maxdim)
42*52bfb9bbSJeremy L Thompson      real*8 uuq(qdimmax)
43*52bfb9bbSJeremy L Thompson      real*8 fx
44*52bfb9bbSJeremy L Thompson      integer*8 uqoffset,xoffset,offset1,offset2
458980d4a7Sjeremylt
468980d4a7Sjeremylt      character arg*32
478980d4a7Sjeremylt
488980d4a7Sjeremylt      call getarg(1,arg)
498980d4a7Sjeremylt      call ceedinit(trim(arg)//char(0),ceed,err)
508980d4a7Sjeremylt
51*52bfb9bbSJeremy L Thompson      do dimn=1,maxdim
52*52bfb9bbSJeremy L Thompson        qdim=q**dimn
53*52bfb9bbSJeremy L Thompson        xdim=2**dimn
548980d4a7Sjeremylt
55*52bfb9bbSJeremy L Thompson        do d=0,dimn-1
56*52bfb9bbSJeremy L Thompson          do i=1,xdim
57*52bfb9bbSJeremy L Thompson            if ((mod(i-1,2**(dimn-d))/(2**(dimn-d-1))).ne.0) then
58*52bfb9bbSJeremy L Thompson              xx(d*xdim+i)=1
59*52bfb9bbSJeremy L Thompson            else
60*52bfb9bbSJeremy L Thompson              xx(d*xdim+i)=-1
61*52bfb9bbSJeremy L Thompson            endif
62*52bfb9bbSJeremy L Thompson          enddo
638980d4a7Sjeremylt        enddo
648980d4a7Sjeremylt
65*52bfb9bbSJeremy L Thompson        call ceedvectorcreate(ceed,xdim*dimn,x,err)
66*52bfb9bbSJeremy L Thompson        xoffset=0
67*52bfb9bbSJeremy L Thompson        call ceedvectorsetarray(x,ceed_mem_host,ceed_use_pointer,xx,xoffset,err)
68*52bfb9bbSJeremy L Thompson        call ceedvectorcreate(ceed,qdim*dimn,xq,err)
69*52bfb9bbSJeremy L Thompson        call ceedvectorsetvalue(xq,0.d0,err)
70*52bfb9bbSJeremy L Thompson        call ceedvectorcreate(ceed,qdim,u,err)
71*52bfb9bbSJeremy L Thompson        call ceedvectorsetvalue(u,0.d0,err)
72*52bfb9bbSJeremy L Thompson        call ceedvectorcreate(ceed,qdim,uq,err)
738980d4a7Sjeremylt
74*52bfb9bbSJeremy L Thompson        call ceedbasiscreatetensorh1lagrange(ceed,dimn,dimn,2,q,&
75*52bfb9bbSJeremy L Thompson     &   ceed_gauss_lobatto,bxl,err)
76*52bfb9bbSJeremy L Thompson        call ceedbasiscreatetensorh1lagrange(ceed,dimn,1,q,q,&
77*52bfb9bbSJeremy L Thompson     &   ceed_gauss_lobatto,bul,err)
788980d4a7Sjeremylt
79*52bfb9bbSJeremy L Thompson        call ceedbasisapply(bxl,1,ceed_notranspose,ceed_eval_interp,x,xq,err)
80*52bfb9bbSJeremy L Thompson
81*52bfb9bbSJeremy L Thompson        call ceedvectorgetarrayread(xq,ceed_mem_host,xxq,offset1,err)
82*52bfb9bbSJeremy L Thompson        do i=1,qdim
83*52bfb9bbSJeremy L Thompson          do d=0,dimn-1
84*52bfb9bbSJeremy L Thompson            xxx(d+1)=xxq(d*qdim+i+offset1)
85*52bfb9bbSJeremy L Thompson          enddo
86*52bfb9bbSJeremy L Thompson          call eval(dimn,xxx,uuq(i))
87*52bfb9bbSJeremy L Thompson        enddo
88*52bfb9bbSJeremy L Thompson        call ceedvectorrestorearrayread(xq,xxq,offset1,err)
89*52bfb9bbSJeremy L Thompson        uqoffset=0
90*52bfb9bbSJeremy L Thompson        call ceedvectorsetarray(uq,ceed_mem_host,ceed_use_pointer,uuq,uqoffset,&
91*52bfb9bbSJeremy L Thompson     &   err)
92*52bfb9bbSJeremy L Thompson
93*52bfb9bbSJeremy L Thompson        call ceedbasisapply(bul,1,ceed_transpose,ceed_eval_interp,uq,u,err)
94*52bfb9bbSJeremy L Thompson
95*52bfb9bbSJeremy L Thompson        call ceedbasiscreatetensorh1lagrange(ceed,dimn,dimn,2,q,ceed_gauss,bxg,&
96*52bfb9bbSJeremy L Thompson     &   err)
97*52bfb9bbSJeremy L Thompson        call ceedbasiscreatetensorh1lagrange(ceed,dimn,1,q,q,ceed_gauss,bug,err)
98*52bfb9bbSJeremy L Thompson
99*52bfb9bbSJeremy L Thompson        call ceedbasisapply(bxg,1,ceed_notranspose,ceed_eval_interp,x,xq,err)
100*52bfb9bbSJeremy L Thompson        call ceedbasisapply(bug,1,ceed_notranspose,ceed_eval_interp,u,uq,err)
101*52bfb9bbSJeremy L Thompson
102*52bfb9bbSJeremy L Thompson        call ceedvectorgetarrayread(xq,ceed_mem_host,xxq,offset1,err)
103*52bfb9bbSJeremy L Thompson        call ceedvectorgetarrayread(uq,ceed_mem_host,uuq,offset2,err)
104*52bfb9bbSJeremy L Thompson        do i=1,qdim
105*52bfb9bbSJeremy L Thompson          do d=0,dimn-1
106*52bfb9bbSJeremy L Thompson            xxx(d+1)=xxq(d*qdim+i+offset1)
107*52bfb9bbSJeremy L Thompson          enddo
108*52bfb9bbSJeremy L Thompson          call eval(dimn,xxx,fx)
109*52bfb9bbSJeremy L Thompson
110*52bfb9bbSJeremy L Thompson          if(dabs(uuq(i+offset2)-fx) > 1.0D-4) then
111a2546046Sjeremylt! LCOV_EXCL_START
112*52bfb9bbSJeremy L Thompson          write(*,*) uuq(i+offset2),' not equal to ',fx,dimn,i
113de996c55Sjeremylt! LCOV_EXCL_STOP
1148980d4a7Sjeremylt          endif
1158980d4a7Sjeremylt        enddo
116*52bfb9bbSJeremy L Thompson        call ceedvectorrestorearrayread(xq,xxq,offset1,err)
117*52bfb9bbSJeremy L Thompson        call ceedvectorrestorearrayread(uq,uuq,offest2,err)
1188980d4a7Sjeremylt
119*52bfb9bbSJeremy L Thompson        call ceedvectordestroy(x,err)
120*52bfb9bbSJeremy L Thompson        call ceedvectordestroy(xq,err)
121*52bfb9bbSJeremy L Thompson        call ceedvectordestroy(u,err)
122*52bfb9bbSJeremy L Thompson        call ceedvectordestroy(uq,err)
123*52bfb9bbSJeremy L Thompson        call ceedbasisdestroy(bxl,err)
124*52bfb9bbSJeremy L Thompson        call ceedbasisdestroy(bul,err)
125*52bfb9bbSJeremy L Thompson        call ceedbasisdestroy(bxg,err)
126*52bfb9bbSJeremy L Thompson        call ceedbasisdestroy(bug,err)
127*52bfb9bbSJeremy L Thompson      enddo
128*52bfb9bbSJeremy L Thompson
1298980d4a7Sjeremylt      call ceeddestroy(ceed,err)
1308980d4a7Sjeremylt      end
1318980d4a7Sjeremylt!-----------------------------------------------------------------------
132