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], grad_1d[j + P * i]); 28 // LCOV_EXCL_START 29 } 30 } 31 } 32 CeedBasisDestroy(&b); 33 34 // Q = P, not already collocated 35 CeedBasisCreateTensorH1Lagrange(ceed, 1, 1, P, P, CEED_GAUSS, &b); 36 CeedBasisGetCollocatedGrad(b, collo_grad_1d); 37 38 CeedBasisGetQRef(b, &q_ref); 39 for (CeedInt i = 0; i < P; i++) x_2[i] = q_ref[i] * q_ref[i]; 40 41 // Verify collo_grad * x^2 = 2x for quadrature points 42 for (CeedInt i = 0; i < P; i++) { 43 sum = 0.0; 44 for (CeedInt j = 0; j < P; j++) sum += collo_grad_1d[j + P * i] * x_2[j]; 45 if (fabs(sum - 2 * q_ref[i]) > 100 * CEED_EPSILON) printf("Error in collocated gradient %f != %f\n", sum, 2 * q_ref[i]); 46 } 47 CeedBasisDestroy(&b); 48 49 // Q = P + 2, not already collocated 50 CeedBasisCreateTensorH1Lagrange(ceed, 1, 1, P, P + 2, CEED_GAUSS, &b); 51 CeedBasisGetCollocatedGrad(b, collo_grad_1d); 52 53 CeedBasisGetQRef(b, &q_ref); 54 for (CeedInt i = 0; i < P + 2; i++) x_2[i] = q_ref[i] * q_ref[i]; 55 56 // Verify collo_grad * x^2 = 2x for quadrature points 57 for (CeedInt i = 0; i < P + 2; i++) { 58 sum = 0.0; 59 for (CeedInt j = 0; j < P + 2; j++) sum += collo_grad_1d[j + (P + 2) * i] * x_2[j]; 60 if (fabs(sum - 2 * q_ref[i]) > 100 * CEED_EPSILON) printf("Error in collocated gradient %f != %f\n", sum, 2 * q_ref[i]); 61 } 62 CeedBasisDestroy(&b); 63 64 CeedDestroy(&ceed); 65 return 0; 66 } 67