1 const char help[] = "Test DMPlexInsertBoundaryValues with DMPlexSetClosurePermutationTensor.\n"; 2 3 #include <petscdmplex.h> 4 5 static PetscErrorCode bc_func(PetscInt dim, PetscReal time, const PetscReal coords[], PetscInt num_comp_u, PetscScalar *u, void *ctx) 6 { 7 PetscFunctionBeginUser; 8 for (PetscInt i = 0; i < num_comp_u; i++) u[i] = coords[i]; 9 PetscFunctionReturn(PETSC_SUCCESS); 10 } 11 12 int main(int argc, char **argv) 13 { 14 DM dm; 15 PetscFE fe; 16 Vec U_loc; 17 PetscInt dim, order = 1; 18 PetscBool tensorCoords = PETSC_TRUE; 19 20 /* Initialize PETSc */ 21 PetscFunctionBeginUser; 22 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 23 PetscCall(PetscOptionsGetBool(NULL, NULL, "-tensor_coords", &tensorCoords, NULL)); 24 PetscCall(DMCreate(PETSC_COMM_WORLD, &dm)); 25 PetscCall(DMSetType(dm, DMPLEX)); 26 PetscCall(DMSetFromOptions(dm)); 27 28 PetscCall(DMGetDimension(dm, &dim)); 29 PetscCall(PetscFECreateLagrange(PETSC_COMM_WORLD, dim, dim, PETSC_FALSE, order, order, &fe)); 30 PetscCall(DMAddField(dm, NULL, (PetscObject)fe)); 31 PetscCall(PetscFEDestroy(&fe)); 32 PetscCall(DMCreateDS(dm)); 33 34 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view")); 35 36 DMLabel label; 37 PetscInt marker_ids[] = {1}; 38 PetscCall(DMGetLabel(dm, "marker", &label)); 39 PetscCall(DMAddBoundary(dm, DM_BC_ESSENTIAL, "mms", label, 1, marker_ids, 0, 0, NULL, (PetscVoidFn *)bc_func, NULL, NULL, NULL)); 40 PetscCall(DMPlexSetClosurePermutationTensor(dm, PETSC_DETERMINE, NULL)); 41 { 42 DM cdm; 43 PetscCall(DMGetCoordinateDM(dm, &cdm)); 44 if (tensorCoords) PetscCall(DMPlexSetClosurePermutationTensor(cdm, PETSC_DETERMINE, NULL)); 45 } 46 47 PetscCall(DMCreateLocalVector(dm, &U_loc)); 48 PetscCall(DMPlexInsertBoundaryValues(dm, PETSC_TRUE, U_loc, 1., NULL, NULL, NULL)); 49 PetscCall(VecViewFromOptions(U_loc, NULL, "-u_loc_vec_view")); 50 PetscCall(VecDestroy(&U_loc)); 51 PetscCall(DMDestroy(&dm)); 52 PetscCall(PetscFinalize()); 53 return 0; 54 } 55 56 /*TEST 57 test: 58 suffix: 2d 59 args: -dm_plex_simplex 0 -dm_plex_dim 2 -dm_plex_box_faces 3,3 -u_loc_vec_view 60 test: 61 suffix: 3d 62 args: -dm_plex_simplex 0 -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -u_loc_vec_view 63 TEST*/ 64