1 static char help[] = "Tests for high order geometry\n\n"; 2 3 #include <petscdmplex.h> 4 #include <petscds.h> 5 6 typedef struct { 7 PetscReal volume; /* Analytical volume of the mesh */ 8 PetscReal tol; /* Tolerance for volume check */ 9 } AppCtx; 10 11 PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 12 { 13 PetscFunctionBegin; 14 options->volume = -1.0; 15 options->tol = PETSC_SMALL; 16 17 PetscOptionsBegin(comm, "", "Meshing Interpolation Test Options", "DMPLEX"); 18 PetscCall(PetscOptionsReal("-volume", "The analytical volume of the mesh", "ex33.c", options->volume, &options->volume, NULL)); 19 PetscCall(PetscOptionsReal("-tol", "The tolerance for the volume check", "ex33.c", options->tol, &options->tol, NULL)); 20 PetscOptionsEnd(); 21 PetscFunctionReturn(PETSC_SUCCESS); 22 } 23 24 PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *ctx, DM *dm) 25 { 26 PetscFunctionBegin; 27 PetscCall(DMCreate(comm, dm)); 28 PetscCall(DMSetType(*dm, DMPLEX)); 29 PetscCall(DMSetFromOptions(*dm)); 30 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 31 PetscFunctionReturn(PETSC_SUCCESS); 32 } 33 34 static void volume(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar vol[]) 35 { 36 vol[0] = 1.; 37 } 38 39 static PetscErrorCode CreateDiscretization(DM dm, AppCtx *ctx) 40 { 41 PetscDS ds; 42 PetscFE fe; 43 DMPolytopeType ct; 44 PetscInt dim, cStart; 45 46 PetscFunctionBeginUser; 47 PetscCall(DMGetDimension(dm, &dim)); 48 PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, NULL)); 49 PetscCall(DMPlexGetCellType(dm, cStart, &ct)); 50 PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim, 1, ct, NULL, PETSC_DETERMINE, &fe)); 51 PetscCall(PetscFESetName(fe, "scalar")); 52 PetscCall(DMAddField(dm, NULL, (PetscObject)fe)); 53 PetscCall(PetscFEDestroy(&fe)); 54 PetscCall(DMCreateDS(dm)); 55 PetscCall(DMGetDS(dm, &ds)); 56 PetscCall(PetscDSSetObjective(ds, 0, volume)); 57 PetscFunctionReturn(PETSC_SUCCESS); 58 } 59 60 static PetscErrorCode CheckVolume(DM dm, AppCtx *ctx) 61 { 62 Vec u; 63 PetscScalar result; 64 PetscReal vol, tol = ctx->tol; 65 66 PetscFunctionBeginUser; 67 PetscCall(DMGetGlobalVector(dm, &u)); 68 PetscCall(DMPlexComputeIntegralFEM(dm, u, &result, ctx)); 69 vol = PetscRealPart(result); 70 PetscCall(DMRestoreGlobalVector(dm, &u)); 71 PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "Volume: %g\n", (double)vol)); 72 PetscCheck(ctx->volume <= 0.0 || PetscAbsReal(ctx->volume - vol) <= tol, PetscObjectComm((PetscObject)dm), PETSC_ERR_PLIB, "Calculated volume %g != %g actual volume (error %g > %g tol)", (double)vol, (double)ctx->volume, (double)PetscAbsReal(ctx->volume - vol), (double)tol); 73 PetscFunctionReturn(PETSC_SUCCESS); 74 } 75 76 int main(int argc, char **argv) 77 { 78 DM dm; 79 AppCtx user; 80 81 PetscFunctionBeginUser; 82 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 83 PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user)); 84 PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm)); 85 PetscCall(CreateDiscretization(dm, &user)); 86 PetscCall(CheckVolume(dm, &user)); 87 PetscCall(DMDestroy(&dm)); 88 PetscCall(PetscFinalize()); 89 return 0; 90 } 91 92 /*TEST 93 94 testset: 95 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower -1.,-1. -dm_plex_box_upper 1.,1. -volume 4. 96 97 test: 98 suffix: square_0 99 args: -dm_coord_petscspace_degree 1 100 101 test: 102 suffix: square_1 103 args: -dm_coord_petscspace_degree 2 104 105 test: 106 suffix: square_2 107 args: -dm_refine 1 -dm_coord_petscspace_degree 1 108 109 test: 110 suffix: square_3 111 args: -dm_refine 1 -dm_coord_petscspace_degree 2 112 113 testset: 114 args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 1,1,1 -dm_plex_box_lower -1.,-1.,-1. -dm_plex_box_upper 1.,1.,1. -volume 8. 115 116 test: 117 suffix: cube_0 118 args: -dm_coord_petscspace_degree 1 119 120 test: 121 suffix: cube_1 122 args: -dm_coord_petscspace_degree 2 123 124 test: 125 suffix: cube_2 126 args: -dm_refine 1 -dm_coord_petscspace_degree 1 127 128 test: 129 suffix: cube_3 130 args: -dm_refine 1 -dm_coord_petscspace_degree 2 131 132 testset: 133 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower -1.,-1. -dm_plex_box_upper 1.,1. \ 134 -volume 4. -dm_coord_remap -dm_coord_map shear -dm_coord_map_params 0.0,0.0,3.0 135 136 test: 137 suffix: shear_0 138 args: -dm_coord_petscspace_degree 1 139 140 test: 141 suffix: shear_1 142 args: -dm_coord_petscspace_degree 2 143 144 test: 145 suffix: shear_2 146 args: -dm_refine 1 -dm_coord_petscspace_degree 1 147 148 test: 149 suffix: shear_3 150 args: -dm_refine 1 -dm_coord_petscspace_degree 2 151 152 testset: 153 args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 1,1,1 -dm_plex_box_lower -1.,-1.,-1. -dm_plex_box_upper 1.,1.,1. \ 154 -volume 8. -dm_coord_remap -dm_coord_map shear -dm_coord_map_params 0.0,0.0,3.0,4.0 155 156 test: 157 suffix: shear_4 158 args: -dm_coord_petscspace_degree 1 159 160 test: 161 suffix: shear_5 162 args: -dm_coord_petscspace_degree 2 163 164 test: 165 suffix: shear_6 166 args: -dm_refine 1 -dm_coord_petscspace_degree 1 167 168 test: 169 suffix: shear_7 170 args: -dm_refine 1 -dm_coord_petscspace_degree 2 171 172 testset: 173 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower 1.,-1. -dm_plex_box_upper 3.,1. \ 174 -dm_coord_petscspace_degree 1 -volume 8. -dm_coord_remap -dm_coord_map flare 175 176 test: 177 suffix: flare_0 178 args: 179 180 test: 181 suffix: flare_1 182 args: -dm_refine 2 183 184 testset: 185 # Area: 3/4 \pi = 2.3562 186 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -petscfe_default_quadrature_order 2 \ 187 -volume 2.35619449019235 -dm_coord_remap -dm_coord_map annulus 188 189 test: 190 # Area: (a+b)/2 h = 3/\sqrt{2} (sqrt{2} - 1/\sqrt{2}) = 3/2 191 suffix: annulus_0 192 requires: double 193 args: -dm_coord_petscspace_degree 1 -volume 1.5 194 195 test: 196 suffix: annulus_1 197 requires: double 198 args: -dm_refine 3 -dm_coord_petscspace_degree 1 -tol .016 199 200 test: 201 suffix: annulus_2 202 requires: double 203 args: -dm_refine 3 -dm_coord_petscspace_degree 2 -tol .0038 204 205 test: 206 suffix: annulus_3 207 requires: double 208 args: -dm_refine 3 -dm_coord_petscspace_degree 3 -tol 2.2e-6 209 210 test: 211 suffix: annulus_4 212 requires: double 213 args: -dm_refine 2 -dm_coord_petscspace_degree 2 -tol .00012 214 215 test: 216 suffix: annulus_5 217 requires: double 218 args: -dm_refine 2 -dm_coord_petscspace_degree 3 -petscfe_default_quadrature_order 3 -tol 1.2e-7 219 220 testset: 221 # Volume: 4/3 \pi (8 - 1)/2 = 14/3 \pi = 14.66076571675238 222 args: -dm_plex_dim 3 -dm_plex_simplex 0 -dm_plex_box_faces 1,1,1 -dm_plex_box_lower -1.,-1.,-1. -dm_plex_box_upper 1.,1.,1. \ 223 -volume 14.66076571675238 -dm_coord_remap -dm_coord_map shell 224 225 test: 226 suffix: shell_0 227 requires: double 228 args: -dm_refine 1 -dm_coord_petscspace_degree 1 -petscfe_default_quadrature_order 1 -volume 5.633164922 -tol 1.0e-7 229 230 test: 231 suffix: shell_1 232 requires: double 233 args: -dm_refine 2 -dm_coord_petscspace_degree 1 -petscfe_default_quadrature_order 1 -tol 3.1 234 235 test: 236 suffix: shell_2 237 requires: double 238 args: -dm_refine 2 -dm_coord_petscspace_degree 2 -petscfe_default_quadrature_order 2 -tol .1 239 240 test: 241 suffix: shell_3 242 requires: double 243 args: -dm_refine 2 -dm_coord_petscspace_degree 3 -petscfe_default_quadrature_order 3 -tol .02 244 245 test: 246 suffix: shell_4 247 requires: double 248 args: -dm_refine 2 -dm_coord_petscspace_degree 4 -petscfe_default_quadrature_order 4 -tol .006 249 250 test: 251 # Volume: 1.0 252 suffix: gmsh_q2 253 requires: double 254 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q2.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6 255 256 test: 257 # Volume: 1.0 258 suffix: gmsh_q3 259 requires: double 260 nsize: {{1 2}} 261 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q3.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6 262 263 test: 264 # Volume: 1.0 265 suffix: gmsh_3d_q2 266 requires: double 267 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_q2.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6 268 269 test: 270 # Volume: 1.0 271 suffix: gmsh_3d_q3 272 requires: double 273 nsize: {{1 2}} 274 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_q3.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6 275 276 TEST*/ 277