xref: /libCEED/tests/t509-operator.c (revision 9df49d7ef0a77c7a3baec2427d8a7274681409b6)
1 /// @file
2 /// Test creation, action, and destruction for identity operator
3 /// \test Test creation, action, and destruction for identity operator
4 #include <ceed.h>
5 #include <stdlib.h>
6 #include <math.h>
7 
8 int main(int argc, char **argv) {
9   Ceed ceed;
10   CeedElemRestriction elem_restr_u, elem_restr_u_i;
11   CeedBasis basis_u;
12   CeedQFunction qf_identity;
13   CeedOperator op_identity;
14   CeedVector U, V;
15   const CeedScalar *hv;
16   CeedInt num_elem = 15, P = 5, Q = 8;
17   CeedInt elem_size = P, num_nodes = num_elem*(P-1)+1;
18   CeedInt ind_u[num_elem*P];
19 
20   CeedInit(argv[1], &ceed);
21 
22   // Restrictions
23   for (CeedInt i=0; i<num_elem; i++) {
24     for (CeedInt j=0; j<P; j++) {
25       ind_u[P*i+j] = i*(P-1) + j;
26     }
27   }
28   CeedElemRestrictionCreate(ceed, num_elem, elem_size, 1, 1, num_nodes,
29                             CEED_MEM_HOST,
30                             CEED_USE_POINTER, ind_u, &elem_restr_u);
31   CeedInt strides_u_i[3] = {1, P, P};
32   CeedElemRestrictionCreateStrided(ceed, num_elem, elem_size, 1,
33                                    elem_size*num_elem,
34                                    strides_u_i, &elem_restr_u_i);
35 
36   // Bases
37   CeedBasisCreateTensorH1Lagrange(ceed, 1, 1, P, Q, CEED_GAUSS, &basis_u);
38 
39   // QFunction
40   CeedQFunctionCreateIdentity(ceed, 1, CEED_EVAL_NONE, CEED_EVAL_NONE,
41                               &qf_identity);
42 
43   // Operators
44   CeedOperatorCreate(ceed, qf_identity, CEED_QFUNCTION_NONE, CEED_QFUNCTION_NONE,
45                      &op_identity);
46   CeedOperatorSetField(op_identity, "input", elem_restr_u, CEED_BASIS_COLLOCATED,
47                        CEED_VECTOR_ACTIVE);
48   CeedOperatorSetField(op_identity, "output", elem_restr_u_i,
49                        CEED_BASIS_COLLOCATED,
50                        CEED_VECTOR_ACTIVE);
51   CeedOperatorSetNumQuadraturePoints(op_identity, elem_size);
52 
53   CeedVectorCreate(ceed, num_nodes, &U);
54   CeedVectorSetValue(U, 3.0);
55   CeedVectorCreate(ceed, elem_size*num_elem, &V);
56   CeedOperatorApply(op_identity, U, V, CEED_REQUEST_IMMEDIATE);
57 
58   // Check output
59   CeedVectorGetArrayRead(V, CEED_MEM_HOST, &hv);
60   for (CeedInt i=0; i<elem_size*num_elem; i++)
61     if (fabs(hv[i]-3.)>1e-14)
62       // LCOV_EXCL_START
63       printf("%d: Computed Value: %f != True Value: 1.0\n", i, hv[i]);
64   // LCOV_EXCL_STOP
65   CeedVectorRestoreArrayRead(V, &hv);
66 
67   CeedQFunctionDestroy(&qf_identity);
68   CeedOperatorDestroy(&op_identity);
69   CeedElemRestrictionDestroy(&elem_restr_u);
70   CeedElemRestrictionDestroy(&elem_restr_u_i);
71   CeedBasisDestroy(&basis_u);
72   CeedVectorDestroy(&U);
73   CeedVectorDestroy(&V);
74   CeedDestroy(&ceed);
75   return 0;
76 }
77