Actual source code: ex17.c

  1: const char help[] = "Test construction of a mixed form degree dual space";

  3: #include <petscdualspace.h>
  4: #include <petscdmplex.h>

  6: int main(int argc, char **argv)
  7: {
  8:   DM              K;
  9:   PetscDualSpace  h_1, h_div, list[2], joint;
 10:   PetscQuadrature quad;
 11:   Mat             mat;

 13:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 14:   PetscCall(DMPlexCreateReferenceCell(PETSC_COMM_SELF, DM_POLYTOPE_QUADRILATERAL, &K));
 15:   PetscCall(PetscDualSpaceCreate(PETSC_COMM_SELF, &h_1));
 16:   PetscCall(PetscObjectSetName((PetscObject)h_1, "H1"));
 17:   PetscCall(PetscDualSpaceSetType(h_1, PETSCDUALSPACELAGRANGE));
 18:   PetscCall(PetscDualSpaceSetDM(h_1, K));
 19:   PetscCall(PetscDualSpaceSetNumComponents(h_1, 2));
 20:   PetscCall(PetscDualSpaceLagrangeSetTensor(h_1, PETSC_TRUE));
 21:   PetscCall(PetscDualSpaceSetFormDegree(h_1, 0));
 22:   PetscCall(PetscDualSpaceSetOrder(h_1, 1));
 23:   PetscCall(PetscDualSpaceSetUp(h_1));
 24:   PetscCall(PetscDualSpaceCreate(PETSC_COMM_SELF, &h_div));
 25:   PetscCall(PetscObjectSetName((PetscObject)h_div, "Hdiv"));
 26:   PetscCall(PetscDualSpaceSetType(h_div, PETSCDUALSPACELAGRANGE));
 27:   PetscCall(PetscDualSpaceSetDM(h_div, K));
 28:   PetscCall(PetscDualSpaceSetNumComponents(h_div, 2));
 29:   PetscCall(PetscDualSpaceLagrangeSetTensor(h_div, PETSC_TRUE));
 30:   PetscCall(PetscDualSpaceSetFormDegree(h_div, -1));
 31:   PetscCall(PetscDualSpaceLagrangeSetTrimmed(h_div, PETSC_TRUE));
 32:   PetscCall(PetscDualSpaceSetOrder(h_div, 1));
 33:   PetscCall(PetscDualSpaceSetUp(h_div));
 34:   list[0] = h_1;
 35:   list[1] = h_div;
 36:   PetscCall(PetscDualSpaceCreateSum(2, list, PETSC_FALSE, &joint));
 37:   PetscCall(PetscDualSpaceSetUp(joint));
 38:   PetscCall(PetscDualSpaceView(joint, PETSC_VIEWER_STDOUT_SELF));
 39:   PetscCall(PetscDualSpaceGetAllData(joint, &quad, &mat));
 40:   PetscCall(PetscQuadratureView(quad, PETSC_VIEWER_STDOUT_SELF));
 41:   PetscCall(PetscObjectSetName((PetscObject)mat, "Joint dual space projection matrix"));
 42:   PetscCall(MatView(mat, PETSC_VIEWER_STDOUT_SELF));
 43:   PetscCall(PetscDualSpaceDestroy(&joint));
 44:   PetscCall(PetscDualSpaceDestroy(&h_div));
 45:   PetscCall(PetscDualSpaceDestroy(&h_1));
 46:   PetscCall(DMDestroy(&K));
 47:   PetscCall(PetscFinalize());
 48:   return 0;
 49: }

 51: /*TEST

 53:   test:
 54:     suffix: 0

 56: TEST*/