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