13c3287b0SJed Brown const char help[] = "Test DMPlexInsertBoundaryValues with DMPlexSetClosurePermutationTensor.\n";
23c3287b0SJed Brown
33c3287b0SJed Brown #include <petscdmplex.h>
43c3287b0SJed Brown
bc_func(PetscInt dim,PetscReal time,const PetscReal coords[],PetscInt num_comp_u,PetscScalar * u,PetscCtx ctx)5*2a8381b2SBarry Smith static PetscErrorCode bc_func(PetscInt dim, PetscReal time, const PetscReal coords[], PetscInt num_comp_u, PetscScalar *u, PetscCtx ctx)
6d71ae5a4SJacob Faibussowitsch {
73c3287b0SJed Brown PetscFunctionBeginUser;
83c3287b0SJed Brown for (PetscInt i = 0; i < num_comp_u; i++) u[i] = coords[i];
93ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
103c3287b0SJed Brown }
113c3287b0SJed Brown
main(int argc,char ** argv)12d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv)
13d71ae5a4SJacob Faibussowitsch {
143c3287b0SJed Brown DM dm;
153c3287b0SJed Brown PetscFE fe;
163c3287b0SJed Brown Vec U_loc;
173c3287b0SJed Brown PetscInt dim, order = 1;
183c3287b0SJed Brown PetscBool tensorCoords = PETSC_TRUE;
193c3287b0SJed Brown
203c3287b0SJed Brown /* Initialize PETSc */
21327415f7SBarry Smith PetscFunctionBeginUser;
223c3287b0SJed Brown PetscCall(PetscInitialize(&argc, &argv, NULL, help));
233c3287b0SJed Brown PetscCall(PetscOptionsGetBool(NULL, NULL, "-tensor_coords", &tensorCoords, NULL));
243c3287b0SJed Brown PetscCall(DMCreate(PETSC_COMM_WORLD, &dm));
253c3287b0SJed Brown PetscCall(DMSetType(dm, DMPLEX));
263c3287b0SJed Brown PetscCall(DMSetFromOptions(dm));
273c3287b0SJed Brown
283c3287b0SJed Brown PetscCall(DMGetDimension(dm, &dim));
293c3287b0SJed Brown PetscCall(PetscFECreateLagrange(PETSC_COMM_WORLD, dim, dim, PETSC_FALSE, order, order, &fe));
303c3287b0SJed Brown PetscCall(DMAddField(dm, NULL, (PetscObject)fe));
313c3287b0SJed Brown PetscCall(PetscFEDestroy(&fe));
323c3287b0SJed Brown PetscCall(DMCreateDS(dm));
333c3287b0SJed Brown
343c3287b0SJed Brown PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
353c3287b0SJed Brown
363c3287b0SJed Brown DMLabel label;
373c3287b0SJed Brown PetscInt marker_ids[] = {1};
383c3287b0SJed Brown PetscCall(DMGetLabel(dm, "marker", &label));
3957d50842SBarry Smith PetscCall(DMAddBoundary(dm, DM_BC_ESSENTIAL, "mms", label, 1, marker_ids, 0, 0, NULL, (PetscVoidFn *)bc_func, NULL, NULL, NULL));
403c3287b0SJed Brown PetscCall(DMPlexSetClosurePermutationTensor(dm, PETSC_DETERMINE, NULL));
413c3287b0SJed Brown {
423c3287b0SJed Brown DM cdm;
433c3287b0SJed Brown PetscCall(DMGetCoordinateDM(dm, &cdm));
443c3287b0SJed Brown if (tensorCoords) PetscCall(DMPlexSetClosurePermutationTensor(cdm, PETSC_DETERMINE, NULL));
453c3287b0SJed Brown }
463c3287b0SJed Brown
473c3287b0SJed Brown PetscCall(DMCreateLocalVector(dm, &U_loc));
483c3287b0SJed Brown PetscCall(DMPlexInsertBoundaryValues(dm, PETSC_TRUE, U_loc, 1., NULL, NULL, NULL));
493c3287b0SJed Brown PetscCall(VecViewFromOptions(U_loc, NULL, "-u_loc_vec_view"));
503c3287b0SJed Brown PetscCall(VecDestroy(&U_loc));
513c3287b0SJed Brown PetscCall(DMDestroy(&dm));
523c3287b0SJed Brown PetscCall(PetscFinalize());
533c3287b0SJed Brown return 0;
543c3287b0SJed Brown }
553c3287b0SJed Brown
563c3287b0SJed Brown /*TEST
573c3287b0SJed Brown test:
583c3287b0SJed Brown suffix: 2d
593c3287b0SJed Brown args: -dm_plex_simplex 0 -dm_plex_dim 2 -dm_plex_box_faces 3,3 -u_loc_vec_view
603c3287b0SJed Brown test:
613c3287b0SJed Brown suffix: 3d
623c3287b0SJed Brown args: -dm_plex_simplex 0 -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -u_loc_vec_view
633c3287b0SJed Brown TEST*/
64