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 if (ctx->volume > 0.0 && PetscAbsReal(ctx->volume - vol) > tol) { 73 SETERRQ(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); 74 } 75 PetscFunctionReturn(PETSC_SUCCESS); 76 } 77 78 int main(int argc, char **argv) 79 { 80 DM dm; 81 AppCtx user; 82 83 PetscFunctionBeginUser; 84 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 85 PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user)); 86 PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm)); 87 PetscCall(CreateDiscretization(dm, &user)); 88 PetscCall(CheckVolume(dm, &user)); 89 PetscCall(DMDestroy(&dm)); 90 PetscCall(PetscFinalize()); 91 return 0; 92 } 93 94 /*TEST 95 96 testset: 97 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. 98 99 test: 100 suffix: square_0 101 args: -dm_coord_petscspace_degree 1 102 103 test: 104 suffix: square_1 105 args: -dm_coord_petscspace_degree 2 106 107 test: 108 suffix: square_2 109 args: -dm_refine 1 -dm_coord_petscspace_degree 1 110 111 test: 112 suffix: square_3 113 args: -dm_refine 1 -dm_coord_petscspace_degree 2 114 115 testset: 116 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. 117 118 test: 119 suffix: cube_0 120 args: -dm_coord_petscspace_degree 1 121 122 test: 123 suffix: cube_1 124 args: -dm_coord_petscspace_degree 2 125 126 test: 127 suffix: cube_2 128 args: -dm_refine 1 -dm_coord_petscspace_degree 1 129 130 test: 131 suffix: cube_3 132 args: -dm_refine 1 -dm_coord_petscspace_degree 2 133 134 testset: 135 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower -1.,-1. -dm_plex_box_upper 1.,1. \ 136 -volume 4. -dm_coord_remap -dm_coord_map shear -dm_coord_map_params 0.0,0.0,3.0 137 138 test: 139 suffix: shear_0 140 args: -dm_coord_petscspace_degree 1 141 142 test: 143 suffix: shear_1 144 args: -dm_coord_petscspace_degree 2 145 146 test: 147 suffix: shear_2 148 args: -dm_refine 1 -dm_coord_petscspace_degree 1 149 150 test: 151 suffix: shear_3 152 args: -dm_refine 1 -dm_coord_petscspace_degree 2 153 154 testset: 155 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. \ 156 -volume 8. -dm_coord_remap -dm_coord_map shear -dm_coord_map_params 0.0,0.0,3.0,4.0 157 158 test: 159 suffix: shear_4 160 args: -dm_coord_petscspace_degree 1 161 162 test: 163 suffix: shear_5 164 args: -dm_coord_petscspace_degree 2 165 166 test: 167 suffix: shear_6 168 args: -dm_refine 1 -dm_coord_petscspace_degree 1 169 170 test: 171 suffix: shear_7 172 args: -dm_refine 1 -dm_coord_petscspace_degree 2 173 174 testset: 175 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -dm_plex_box_lower 1.,-1. -dm_plex_box_upper 3.,1. \ 176 -dm_coord_petscspace_degree 1 -volume 8. -dm_coord_remap -dm_coord_map flare 177 178 test: 179 suffix: flare_0 180 args: 181 182 test: 183 suffix: flare_1 184 args: -dm_refine 2 185 186 testset: 187 # Area: 3/4 \pi = 2.3562 188 args: -dm_plex_simplex 0 -dm_plex_box_faces 1,1 -petscfe_default_quadrature_order 2 \ 189 -volume 2.35619449019235 -dm_coord_remap -dm_coord_map annulus 190 191 test: 192 # Area: (a+b)/2 h = 3/\sqrt{2} (sqrt{2} - 1/\sqrt{2}) = 3/2 193 suffix: annulus_0 194 requires: double 195 args: -dm_coord_petscspace_degree 1 -volume 1.5 196 197 test: 198 suffix: annulus_1 199 requires: double 200 args: -dm_refine 3 -dm_coord_petscspace_degree 1 -tol .016 201 202 test: 203 suffix: annulus_2 204 requires: double 205 args: -dm_refine 3 -dm_coord_petscspace_degree 2 -tol .0038 206 207 test: 208 suffix: annulus_3 209 requires: double 210 args: -dm_refine 3 -dm_coord_petscspace_degree 3 -tol 2.2e-6 211 212 test: 213 suffix: annulus_4 214 requires: double 215 args: -dm_refine 2 -dm_coord_petscspace_degree 2 -tol .00012 216 217 test: 218 suffix: annulus_5 219 requires: double 220 args: -dm_refine 2 -dm_coord_petscspace_degree 3 -petscfe_default_quadrature_order 3 -tol 1.2e-7 221 222 testset: 223 # Volume: 4/3 \pi (8 - 1)/2 = 14/3 \pi = 14.66076571675238 224 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. \ 225 -volume 14.66076571675238 -dm_coord_remap -dm_coord_map shell 226 227 test: 228 suffix: shell_0 229 requires: double 230 args: -dm_refine 1 -dm_coord_petscspace_degree 1 -petscfe_default_quadrature_order 1 -volume 5.633164922 -tol 1.0e-7 231 232 test: 233 suffix: shell_1 234 requires: double 235 args: -dm_refine 2 -dm_coord_petscspace_degree 1 -petscfe_default_quadrature_order 1 -tol 3.1 236 237 test: 238 suffix: shell_2 239 requires: double 240 args: -dm_refine 2 -dm_coord_petscspace_degree 2 -petscfe_default_quadrature_order 2 -tol .1 241 242 test: 243 suffix: shell_3 244 requires: double 245 args: -dm_refine 2 -dm_coord_petscspace_degree 3 -petscfe_default_quadrature_order 3 -tol .02 246 247 test: 248 suffix: shell_4 249 requires: double 250 args: -dm_refine 2 -dm_coord_petscspace_degree 4 -petscfe_default_quadrature_order 4 -tol .006 251 252 test: 253 # Volume: 1.0 254 suffix: gmsh_q2 255 requires: double 256 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q2.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6 257 258 test: 259 # Volume: 1.0 260 suffix: gmsh_q3 261 requires: double 262 nsize: {{1 2}} 263 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/quads-q3.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6 264 265 test: 266 # Volume: 1.0 267 suffix: gmsh_3d_q2 268 requires: double 269 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_q2.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6 270 271 test: 272 # Volume: 1.0 273 suffix: gmsh_3d_q3 274 requires: double 275 nsize: {{1 2}} 276 args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_q3.msh -dm_plex_gmsh_project -volume 1.0 -tol 1e-6 277 278 TEST*/ 279