xref: /libCEED/tests/t324-basis.c (revision a9e65696a8c8214eb82d2dcf9ed1f28a32d2c94e)
1 /// @file
2 /// Test grad transpose with a 2D Simplex non-tensor H1 basis
3 /// \test Test grad transposewith a 2D Simplex non-tensor H1 basis
4 #include <ceed.h>
5 #include <math.h>
6 
7 #include "t320-basis.h"
8 
9 int main(int argc, char **argv) {
10   Ceed              ceed;
11   CeedVector        In, Out;
12   const CeedInt     P = 6, Q = 4, dim = 2;
13   CeedBasis         b;
14   CeedScalar        q_ref[dim * Q], q_weight[Q];
15   CeedScalar        interp[P * Q], grad[dim * P * Q];
16   const CeedScalar *out;
17   CeedScalar        colsum[P];
18 
19   buildmats(q_ref, q_weight, interp, grad);
20 
21   CeedInit(argv[1], &ceed);
22 
23   for (int i = 0; i < P; i++) {
24     colsum[i] = 0;
25     for (int j = 0; j < Q * dim; j++) {
26       colsum[i] += grad[i + j * P];
27     }
28   }
29 
30   CeedBasisCreateH1(ceed, CEED_TOPOLOGY_TRIANGLE, 1, P, Q, interp, grad, q_ref, q_weight, &b);
31 
32   CeedVectorCreate(ceed, Q * dim, &In);
33   CeedVectorSetValue(In, 1);
34   CeedVectorCreate(ceed, P, &Out);
35   CeedVectorSetValue(Out, 0);
36 
37   CeedBasisApply(b, 1, CEED_TRANSPOSE, CEED_EVAL_GRAD, In, Out);
38 
39   // Check values at quadrature points
40   CeedVectorGetArrayRead(Out, CEED_MEM_HOST, &out);
41   for (int i = 0; i < P; i++) {
42     if (fabs(colsum[i] - out[i]) > 100. * CEED_EPSILON) printf("[%" CeedInt_FMT "] %f != %f\n", i, out[i], colsum[i]);
43   }
44   CeedVectorRestoreArrayRead(Out, &out);
45 
46   CeedVectorDestroy(&In);
47   CeedVectorDestroy(&Out);
48   CeedBasisDestroy(&b);
49   CeedDestroy(&ceed);
50   return 0;
51 }
52