1 const char help[] = "Construct and set a Lagrange dual space from options, then view it to\n" 2 "understand the effects of different parameters."; 3 4 #include <petscfe.h> 5 #include <petscdmplex.h> 6 7 int main(int argc, char **argv) 8 { 9 PetscInt dim; 10 PetscBool tensorCell; 11 DM K; 12 PetscDualSpace dsp; 13 14 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 15 dim = 2; 16 tensorCell = PETSC_FALSE; 17 PetscOptionsBegin(PETSC_COMM_WORLD,"","Options for PETSCDUALSPACELAGRANGE test","none"); 18 PetscCall(PetscOptionsRangeInt("-dim", "The spatial dimension","ex1.c",dim,&dim,NULL,0,3)); 19 PetscCall(PetscOptionsBool("-tensor", "Whether the cell is a tensor product cell or a simplex","ex1.c",tensorCell,&tensorCell,NULL)); 20 PetscOptionsEnd(); 21 22 PetscCall(PetscDualSpaceCreate(PETSC_COMM_WORLD, &dsp)); 23 PetscCall(PetscObjectSetName((PetscObject)dsp, "Lagrange dual space")); 24 PetscCall(PetscDualSpaceSetType(dsp, PETSCDUALSPACELAGRANGE)); 25 /* While Lagrange nodes don't require the existence of a reference cell to 26 * be refined, when we construct finite element dual spaces we have to be 27 * careful about what kind of continuity is maintained when cells are glued 28 * together to make a mesh. The PetscDualSpace object is responsible for 29 * conveying continuity requirements to a finite element assembly routines, 30 * so a PetscDualSpace needs a reference element: a single element mesh, 31 * whose boundary points are the interstitial points in a mesh */ 32 PetscCall(DMPlexCreateReferenceCell(PETSC_COMM_WORLD, DMPolytopeTypeSimpleShape(dim, (PetscBool) !tensorCell), &K)); 33 PetscCall(PetscDualSpaceSetDM(dsp, K)); 34 /* This gives us the opportunity to change the parameters of the dual space 35 * from the command line, as we do in the tests below. When 36 * PetscDualSpaceSetFromOptions() is called, it also enables other optional 37 * behavior (see the next step) */ 38 PetscCall(PetscDualSpaceSetFromOptions(dsp)); 39 /* This step parses the parameters of the dual space into 40 * sets of functionals that are assigned to each of the mesh points in K. 41 * 42 * The functionals can be accessed individually by 43 * PetscDualSpaceGetFunctional(), or more efficiently all at once by 44 * PetscDualSpaceGetAllData(), which returns a set of quadrature points 45 * at which to evaluate a function, and a matrix that takes those 46 * evaluations and turns them into the evaluation of the dual space's 47 * functionals on the function. 48 * 49 * (TODO: tutorial for PetscDualSpaceGetAllData() and 50 * PetscDualSpaceGetInteriorData().) 51 * 52 * Because we called PetscDualSpaceSetFromOptions(), we have the opportunity 53 * to inspect the results of PetscDualSpaceSetUp() from the command line 54 * with "-petscdualspace_view", followed by an optional description of how 55 * we would like to see the dual space (see examples in the tests below). 56 * */ 57 PetscCall(PetscDualSpaceSetUp(dsp)); 58 PetscCall(DMDestroy(&K)); 59 PetscCall(PetscDualSpaceDestroy(&dsp)); 60 PetscCall(PetscFinalize()); 61 return 0; 62 } 63 64 /*TEST 65 66 # quadratic nodes on the triangle 67 test: 68 suffix: 0 69 filter: sed -E "s/\(\+0, \+0\)/(+0., +0.)/g" 70 args: -dim 2 -tensor 0 -petscdualspace_order 2 -petscdualspace_view ascii::ascii_info_detail 71 72 # linear nodes on the quadrilateral 73 test: 74 suffix: 1 75 args: -dim 2 -tensor 1 -petscdualspace_order 1 -petscdualspace_lagrange_tensor 1 -petscdualspace_view ascii::ascii_info_detail 76 77 # lowest order Raviart-Thomas / Nedelec edge nodes on the hexahedron 78 test: 79 suffix: 2 80 args: -dim 3 -tensor 1 -petscdualspace_order 1 -petscdualspace_components 3 -petscdualspace_form_degree 1 -petscdualspace_lagrange_trimmed 1 -petscdualspace_lagrange_tensor 1 -petscdualspace_view ascii::ascii_info_detail 81 82 # first order Nedelec second type face nodes on the tetrahedron 83 test: 84 suffix: 3 85 args: -dim 3 -tensor 0 -petscdualspace_order 1 -petscdualspace_components 3 -petscdualspace_form_degree -2 -petscdualspace_view ascii::ascii_info_detail 86 87 ## Comparing different node types 88 89 test: 90 suffix: 4 91 args: -dim 2 -tensor 0 -petscdualspace_order 3 -petscdualspace_lagrange_continuity 0 -petscdualspace_lagrange_node_type equispaced -petscdualspace_lagrange_node_endpoints 0 -petscdualspace_view ascii::ascii_info_detail 92 93 test: 94 suffix: 5 95 args: -dim 2 -tensor 0 -petscdualspace_order 3 -petscdualspace_lagrange_continuity 0 -petscdualspace_lagrange_node_type equispaced -petscdualspace_lagrange_node_endpoints 1 -petscdualspace_view ascii::ascii_info_detail 96 97 test: 98 suffix: 6 99 args: -dim 2 -tensor 0 -petscdualspace_order 3 -petscdualspace_lagrange_continuity 0 -petscdualspace_lagrange_node_type gaussjacobi -petscdualspace_lagrange_node_endpoints 0 -petscdualspace_view ascii::ascii_info_detail 100 101 test: 102 suffix: 7 103 args: -dim 2 -tensor 0 -petscdualspace_order 3 -petscdualspace_lagrange_continuity 0 -petscdualspace_lagrange_node_type gaussjacobi -petscdualspace_lagrange_node_endpoints 1 -petscdualspace_view ascii::ascii_info_detail 104 105 TEST*/ 106