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