150c301a5SRezgar Shakeri /// @file 250c301a5SRezgar Shakeri /// Test GetInterp and BasisApply for a 2D Quad non-tensor H(div) basis 350c301a5SRezgar Shakeri /// \test Test GetInterp and BasisApply for a 2D Quad non-tensor H(div) basis 450c301a5SRezgar Shakeri #include <ceed.h> 550c301a5SRezgar Shakeri #include <math.h> 6*2b730f8bSJeremy L Thompson 750c301a5SRezgar Shakeri #include "t330-basis.h" 850c301a5SRezgar Shakeri 950c301a5SRezgar Shakeri int main(int argc, char **argv) { 1050c301a5SRezgar Shakeri Ceed ceed; 1150c301a5SRezgar Shakeri const CeedInt num_nodes = 4, Q = 3, dim = 2, num_qpts = Q * Q; 1250c301a5SRezgar Shakeri CeedInt num_comp = 1; // one vector componenet 1350c301a5SRezgar Shakeri CeedInt P = dim * num_nodes; // dof per element! 1450c301a5SRezgar Shakeri CeedBasis b; 1550c301a5SRezgar Shakeri CeedScalar q_ref[dim * num_qpts], q_weights[num_qpts]; 1650c301a5SRezgar Shakeri CeedScalar div[P * num_qpts], interp[P * dim * num_qpts]; 1750c301a5SRezgar Shakeri CeedVector X, Y; 1850c301a5SRezgar Shakeri const CeedScalar *y, *x; 1950c301a5SRezgar Shakeri 2050c301a5SRezgar Shakeri CeedInit(argv[1], &ceed); 2150c301a5SRezgar Shakeri 2250c301a5SRezgar Shakeri HdivBasisQuad(Q, q_ref, q_weights, interp, div, CEED_GAUSS); 23*2b730f8bSJeremy L Thompson CeedBasisCreateHdiv(ceed, CEED_TOPOLOGY_QUAD, num_comp, P, num_qpts, interp, div, q_ref, q_weights, &b); 2450c301a5SRezgar Shakeri 2550c301a5SRezgar Shakeri // Test GetInterp for H(div) 2650c301a5SRezgar Shakeri const CeedScalar *interp2; 2750c301a5SRezgar Shakeri CeedBasisGetInterp(b, &interp2); 2850c301a5SRezgar Shakeri for (CeedInt i = 0; i < P * dim * num_qpts; i++) { 29*2b730f8bSJeremy L Thompson if (fabs(interp[i] - interp2[i]) > 100. * CEED_EPSILON) { 3050c301a5SRezgar Shakeri // LCOV_EXCL_START 3150c301a5SRezgar Shakeri printf("%f != %f\n", interp[i], interp2[i]); 3250c301a5SRezgar Shakeri // LCOV_EXCL_STOP 3350c301a5SRezgar Shakeri } 34*2b730f8bSJeremy L Thompson } 3550c301a5SRezgar Shakeri 3650c301a5SRezgar Shakeri CeedVectorCreate(ceed, P, &X); 3750c301a5SRezgar Shakeri CeedVectorSetValue(X, 1.0); 3850c301a5SRezgar Shakeri CeedVectorCreate(ceed, num_qpts * dim, &Y); 3950c301a5SRezgar Shakeri CeedVectorSetValue(Y, 0.); 4050c301a5SRezgar Shakeri // BasisApply for H(div): CEED_EVAL_INTERP, NOTRANSPOSE case 4150c301a5SRezgar Shakeri CeedBasisApply(b, 1, CEED_NOTRANSPOSE, CEED_EVAL_INTERP, X, Y); 4250c301a5SRezgar Shakeri 4350c301a5SRezgar Shakeri CeedVectorGetArrayRead(Y, CEED_MEM_HOST, &y); 4450c301a5SRezgar Shakeri for (CeedInt i = 0; i < dim * num_qpts; i++) { 45*2b730f8bSJeremy L Thompson if (fabs(q_ref[i] - y[i]) > 100. * CEED_EPSILON) { 4650c301a5SRezgar Shakeri // LCOV_EXCL_START 4750c301a5SRezgar Shakeri printf("%f != %f\n", q_ref[i], y[i]); 4850c301a5SRezgar Shakeri // LCOV_EXCL_STOP 4950c301a5SRezgar Shakeri } 50*2b730f8bSJeremy L Thompson } 5150c301a5SRezgar Shakeri CeedVectorRestoreArrayRead(Y, &y); 5250c301a5SRezgar Shakeri 5350c301a5SRezgar Shakeri CeedVectorSetValue(Y, 1.0); 5450c301a5SRezgar Shakeri CeedVectorSetValue(X, 0.0); 5550c301a5SRezgar Shakeri // BasisApply for Hdiv: CEED_EVAL_INTERP, TRANSPOSE case 5650c301a5SRezgar Shakeri CeedBasisApply(b, 1, CEED_TRANSPOSE, CEED_EVAL_INTERP, Y, X); 5750c301a5SRezgar Shakeri 5850c301a5SRezgar Shakeri CeedVectorGetArrayRead(X, CEED_MEM_HOST, &x); 5950c301a5SRezgar Shakeri CeedScalar sum = 0.; 6050c301a5SRezgar Shakeri for (CeedInt i = 0; i < P; i++) { 6150c301a5SRezgar Shakeri sum += x[i]; 6250c301a5SRezgar Shakeri } 63*2b730f8bSJeremy L Thompson if (fabs(sum) > 100. * CEED_EPSILON) printf("sum of array %f != %f\n", sum, 0.0); 6450c301a5SRezgar Shakeri CeedVectorRestoreArrayRead(X, &x); 6550c301a5SRezgar Shakeri 6650c301a5SRezgar Shakeri CeedBasisDestroy(&b); 6750c301a5SRezgar Shakeri CeedVectorDestroy(&X); 6850c301a5SRezgar Shakeri CeedVectorDestroy(&Y); 6950c301a5SRezgar Shakeri CeedDestroy(&ceed); 7050c301a5SRezgar Shakeri return 0; 7150c301a5SRezgar Shakeri } 72