1 static const char help[] = "Tests for Plex transforms, including regular refinement"; 2 3 #include <petscdmplex.h> 4 #include <petscsf.h> 5 6 #include <petsc/private/dmpleximpl.h> 7 8 static PetscErrorCode LabelPoints(DM dm) { 9 DMLabel label; 10 PetscInt pStart, pEnd, p; 11 PetscBool flg = PETSC_FALSE; 12 13 PetscFunctionBegin; 14 PetscCall(PetscOptionsGetBool(NULL, NULL, "-label_mesh", &flg, NULL)); 15 if (!flg) PetscFunctionReturn(0); 16 PetscCall(DMCreateLabel(dm, "test")); 17 PetscCall(DMGetLabel(dm, "test", &label)); 18 PetscCall(DMPlexGetChart(dm, &pStart, &pEnd)); 19 for (p = pStart; p < pEnd; ++p) PetscCall(DMLabelSetValue(label, p, p)); 20 PetscFunctionReturn(0); 21 } 22 23 static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm) { 24 PetscFunctionBegin; 25 PetscCall(DMCreate(comm, dm)); 26 PetscCall(DMSetType(*dm, DMPLEX)); 27 PetscCall(DMSetFromOptions(*dm)); 28 PetscCall(LabelPoints(*dm)); 29 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, "post_label_")); 30 PetscCall(DMSetFromOptions(*dm)); 31 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, NULL)); 32 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 33 PetscFunctionReturn(0); 34 } 35 36 int main(int argc, char **argv) { 37 DM dm; 38 39 PetscFunctionBeginUser; 40 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 41 PetscCall(CreateMesh(PETSC_COMM_WORLD, &dm)); 42 PetscCall(DMDestroy(&dm)); 43 PetscCall(PetscFinalize()); 44 return 0; 45 } 46 47 /*TEST 48 test: 49 suffix: ref_seg 50 args: -dm_plex_reference_cell_domain -dm_plex_cell segment -dm_refine 1 -dm_plex_check_all 51 52 test: 53 suffix: ref_tri 54 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_refine 2 -dm_plex_check_all 55 56 test: 57 suffix: box_tri 58 requires: triangle parmetis 59 nsize: {{1 3 5}} 60 args: -dm_plex_box_faces 3,3 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis 61 62 test: 63 suffix: ref_quad 64 args: -dm_plex_reference_cell_domain -dm_plex_cell quadrilateral -dm_refine 2 -dm_plex_check_all 65 66 test: 67 suffix: box_quad 68 nsize: {{1 3 5}} 69 requires: parmetis 70 args: -dm_plex_box_faces 3,3 -dm_plex_simplex 0 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis 71 72 test: 73 suffix: ref_tet 74 args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron -dm_refine 2 -dm_plex_check_all 75 76 test: 77 suffix: box_tet 78 requires: ctetgen parmetis 79 nsize: {{1 3 5}} 80 args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis 81 82 test: 83 suffix: ref_hex 84 args: -dm_plex_reference_cell_domain -dm_plex_cell hexahedron -dm_refine 2 -dm_plex_check_all 85 86 test: 87 suffix: box_hex 88 requires: parmetis 89 nsize: {{1 3 5}} 90 args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_plex_simplex 0 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis 91 92 test: 93 suffix: ref_trip 94 args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism -dm_refine 2 -dm_plex_check_all 95 96 test: 97 suffix: ref_tquad 98 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad -dm_refine 2 -dm_plex_check_all 99 100 test: 101 suffix: ref_ttrip 102 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism -dm_refine 2 -dm_plex_check_all 103 104 test: 105 suffix: ref_tquadp 106 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism -dm_refine 2 -dm_plex_check_all 107 108 test: 109 suffix: ref_pyramid 110 args: -dm_plex_reference_cell_domain -dm_plex_cell pyramid -dm_refine 2 -dm_plex_check_all 111 112 testset: 113 args: -dm_coord_space 0 -dm_plex_transform_type refine_tobox -dm_plex_check_all 114 115 test: 116 suffix: ref_tri_tobox 117 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_refine 2 118 119 test: 120 suffix: box_tri_tobox 121 requires: triangle parmetis 122 nsize: {{1 3 5}} 123 args: -dm_plex_box_faces 3,3 -dm_refine 2 -petscpartitioner_type parmetis 124 125 test: 126 suffix: ref_tet_tobox 127 args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron -dm_refine 2 128 129 test: 130 suffix: box_tet_tobox 131 requires: ctetgen parmetis 132 nsize: {{1 3 5}} 133 args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_refine 2 -petscpartitioner_type parmetis 134 135 test: 136 suffix: ref_trip_tobox 137 args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism -dm_refine 2 138 139 test: 140 suffix: ref_ttrip_tobox 141 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism -dm_refine 2 142 143 test: 144 suffix: ref_tquadp_tobox 145 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism -dm_refine 2 146 147 testset: 148 args: -dm_coord_space 0 -label_mesh -post_label_dm_extrude 2 -post_label_dm_plex_check_all -dm_view ::ascii_info_detail 149 150 test: 151 suffix: extrude_quad 152 args: -dm_plex_simplex 0 153 154 TEST*/ 155