1 /// @file
2 /// Test interpolation ApplyAdd in multiple dimensions
3 /// \test Test interpolation ApplyAdd in multiple dimensions
4 #include <ceed.h>
5 #include <math.h>
6 #include <stdio.h>
7
main(int argc,char ** argv)8 int main(int argc, char **argv) {
9 Ceed ceed;
10
11 CeedInit(argv[1], &ceed);
12
13 for (CeedInt dim = 1; dim <= 3; dim++) {
14 CeedVector u, u_q, v, v_q, w_q;
15 CeedBasis basis;
16 CeedInt p = 4, q = 5, p_dim = CeedIntPow(p, dim), q_dim = CeedIntPow(q, dim);
17
18 CeedVectorCreate(ceed, p_dim, &u);
19 CeedVectorCreate(ceed, p_dim, &v);
20 CeedVectorSetValue(u, 1.0);
21 CeedVectorSetValue(v, 0.0);
22 CeedVectorCreate(ceed, q_dim, &u_q);
23 CeedVectorCreate(ceed, q_dim, &v_q);
24 CeedVectorCreate(ceed, q_dim, &w_q);
25
26 CeedBasisCreateTensorH1Lagrange(ceed, dim, 1, p, q, CEED_GAUSS, &basis);
27
28 // Compute area
29 CeedBasisApply(basis, 1, CEED_NOTRANSPOSE, CEED_EVAL_INTERP, u, u_q);
30 CeedBasisApply(basis, 1, CEED_NOTRANSPOSE, CEED_EVAL_WEIGHT, CEED_VECTOR_NONE, w_q);
31 CeedVectorPointwiseMult(v_q, u_q, w_q);
32 CeedBasisApply(basis, 1, CEED_TRANSPOSE, CEED_EVAL_INTERP, v_q, v);
33 // Double area computed
34 CeedBasisApplyAdd(basis, 1, CEED_TRANSPOSE, CEED_EVAL_INTERP, v_q, v);
35
36 // Check area
37 {
38 const CeedScalar *v_array;
39 CeedScalar area = 0.0;
40
41 CeedVectorGetArrayRead(v, CEED_MEM_HOST, &v_array);
42 for (CeedInt i = 0; i < p_dim; i++) area += v_array[i];
43 if (fabs(area - 2.0 * CeedIntPow(2, dim)) > 5E-6) printf("Incorrect area computed %f != %f\n", area, 2.0 * CeedIntPow(2, dim));
44 CeedVectorRestoreArrayRead(v, &v_array);
45 }
46
47 CeedVectorDestroy(&u);
48 CeedVectorDestroy(&v);
49 CeedVectorDestroy(&u_q);
50 CeedVectorDestroy(&v_q);
51 CeedVectorDestroy(&w_q);
52 CeedBasisDestroy(&basis);
53 }
54 CeedDestroy(&ceed);
55 return 0;
56 }
57