xref: /libCEED/tests/t302-basis.c (revision ea6b58218a3c4883c2efd66165b4d6b684f89f5a)
1 /// @file
2 /// Test Collocated Grad calculated matches basis with Lobatto points
3 /// \test Test Collocated Grad calculated matches basis with Lobatto points
4 #include <ceed.h>
5 #include <ceed/backend.h>
6 #include <math.h>
7 
8 int main(int argc, char **argv) {
9   Ceed ceed;
10   CeedInt P = 4;
11   CeedScalar collo_grad_1d[(P+2)*(P+2)], x_2[P+2];
12   const CeedScalar *grad_1d, *q_ref;
13   CeedScalar sum = 0.0;
14   CeedBasis b;
15 
16   CeedInit(argv[1], &ceed);
17 
18   // Already collocated, GetCollocatedGrad will return grad_1d
19   CeedBasisCreateTensorH1Lagrange(ceed, 1, 1, P, P, CEED_GAUSS_LOBATTO, &b);
20   CeedBasisGetCollocatedGrad(b, collo_grad_1d);
21   CeedBasisGetGrad(b, &grad_1d);
22 
23   for (CeedInt i=0; i<P; i++)
24     for (CeedInt j=0; j<P; j++)
25       if (fabs(collo_grad_1d[j+P*i] - grad_1d[j+P*i]) > 100*CEED_EPSILON)
26         // LCOV_EXCL_START
27         printf("Error in collocated gradient %f != %f\n", collo_grad_1d[j+P*i],
28                grad_1d[j+P*i]);
29   // LCOV_EXCL_START
30   CeedBasisDestroy(&b);
31 
32   // Q = P, not already collocated
33   CeedBasisCreateTensorH1Lagrange(ceed, 1,  1, P, P, CEED_GAUSS, &b);
34   CeedBasisGetCollocatedGrad(b, collo_grad_1d);
35 
36   CeedBasisGetQRef(b, &q_ref);
37   for (CeedInt i=0; i<P; i++)
38     x_2[i] = q_ref[i]*q_ref[i];
39 
40   // Verify collo_grad * x^2 = 2x for quadrature points
41   for (CeedInt i=0; i<P; i++) {
42     sum = 0.0;
43     for (CeedInt j=0; j<P; j++)
44       sum += collo_grad_1d[j+P*i]*x_2[j];
45     if (fabs(sum - 2*q_ref[i]) > 100*CEED_EPSILON)
46       // LCOV_EXCL_START
47       printf("Error in collocated gradient %f != %f\n", sum, 2*q_ref[i]);
48     // LCOV_EXCL_STOP
49   }
50   CeedBasisDestroy(&b);
51 
52   // Q = P + 2, not already collocated
53   CeedBasisCreateTensorH1Lagrange(ceed, 1,  1, P, P+2, CEED_GAUSS, &b);
54   CeedBasisGetCollocatedGrad(b, collo_grad_1d);
55 
56   CeedBasisGetQRef(b, &q_ref);
57   for (CeedInt i=0; i<P+2; i++)
58     x_2[i] = q_ref[i]*q_ref[i];
59 
60   // Verify collo_grad * x^2 = 2x for quadrature points
61   for (CeedInt i=0; i<P+2; i++) {
62     sum = 0.0;
63     for (CeedInt j=0; j<P+2; j++)
64       sum += collo_grad_1d[j+(P+2)*i]*x_2[j];
65     if (fabs(sum - 2*q_ref[i]) > 100*CEED_EPSILON)
66       // LCOV_EXCL_START
67       printf("Error in collocated gradient %f != %f\n", sum, 2*q_ref[i]);
68     // LCOV_EXCL_STOP
69   }
70   CeedBasisDestroy(&b);
71 
72   CeedDestroy(&ceed);
73   return 0;
74 }
75