1 static char help[] = "Test DMCreateLocalVector_Plex, DMPlexGetCellFields and DMPlexRestoreCellFields work properly for 0 fields/cells/DS dimension\n\n"; 2 static char FILENAME[] = "ex25.c"; 3 4 #include <petscdmplex.h> 5 #include <petscds.h> 6 #include <petscsnes.h> 7 8 9 typedef struct { 10 PetscInt test; 11 PetscInt dim; /* The topological mesh dimension */ 12 PetscInt faces[3]; /* Number of faces per dimension */ 13 PetscBool simplex; /* Use simplices or hexes */ 14 PetscBool interpolate; /* Interpolate mesh */ 15 char filename[PETSC_MAX_PATH_LEN]; /* Import mesh from file */ 16 } AppCtx; 17 18 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 19 { 20 PetscInt dim; 21 PetscErrorCode ierr; 22 23 PetscFunctionBegin; 24 options->test = 0; 25 options->dim = 3; 26 options->simplex = PETSC_TRUE; 27 options->interpolate = PETSC_FALSE; 28 options->filename[0] = '\0'; 29 ierr = PetscOptionsBegin(comm, "", "Zero-sized DMPlexGetCellFields Test Options", "DMPLEX");CHKERRQ(ierr); 30 ierr = PetscOptionsBoundedInt("-test", "Test to run", FILENAME, options->test, &options->test, NULL,0);CHKERRQ(ierr); 31 ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", FILENAME, options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr); 32 if (options->dim > 3) SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "dimension set to %d, must be <= 3", options->dim); 33 ierr = PetscOptionsBool("-simplex", "Use simplices if true, otherwise hexes", FILENAME, options->simplex, &options->simplex, NULL);CHKERRQ(ierr); 34 ierr = PetscOptionsBool("-interpolate", "Interpolate the mesh", FILENAME, options->interpolate, &options->interpolate, NULL);CHKERRQ(ierr); 35 ierr = PetscOptionsString("-filename", "The mesh file", FILENAME, options->filename, options->filename, PETSC_MAX_PATH_LEN, NULL);CHKERRQ(ierr); 36 options->faces[0] = 1; options->faces[1] = 1; options->faces[2] = 1; 37 dim = options->dim; 38 ierr = PetscOptionsIntArray("-faces", "Number of faces per dimension", FILENAME, options->faces, &dim, NULL);CHKERRQ(ierr); 39 if (dim) options->dim = dim; 40 ierr = PetscOptionsEnd(); 41 PetscFunctionReturn(0); 42 } 43 44 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *options, DM *dm) 45 { 46 PetscInt dim = options->dim; 47 PetscInt *faces = options->faces; 48 PetscBool simplex = options->simplex; 49 PetscBool interpolate = options->interpolate; 50 const char *filename = options->filename; 51 size_t len; 52 PetscMPIInt rank; 53 PetscErrorCode ierr; 54 55 PetscFunctionBegin; 56 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 57 ierr = PetscStrlen(filename, &len);CHKERRQ(ierr); 58 if (len) { 59 ierr = DMPlexCreateFromFile(comm, filename, interpolate, dm);CHKERRQ(ierr); 60 ierr = DMGetDimension(*dm, &options->dim);CHKERRQ(ierr); 61 } else { 62 ierr = DMPlexCreateBoxMesh(comm, dim, simplex, faces, NULL, NULL, NULL, interpolate, dm);CHKERRQ(ierr); 63 } 64 PetscFunctionReturn(0); 65 } 66 67 /* no discretization is given so DMGetNumFields yields 0 */ 68 static PetscErrorCode test0(DM dm, AppCtx *options) 69 { 70 Vec locX; 71 PetscErrorCode ierr; 72 73 PetscFunctionBegin; 74 ierr = DMGetLocalVector(dm, &locX);CHKERRQ(ierr); 75 ierr = DMRestoreLocalVector(dm, &locX);CHKERRQ(ierr); 76 PetscFunctionReturn(0); 77 } 78 79 /* no discretization is given so DMGetNumFields and PetscDSGetTotalDimension yield 0 */ 80 static PetscErrorCode test1(DM dm, AppCtx *options) 81 { 82 IS cells; 83 Vec locX, locX_t, locA; 84 PetscScalar *u, *u_t, *a; 85 PetscErrorCode ierr; 86 87 PetscFunctionBegin; 88 ierr = ISCreateStride(PETSC_COMM_SELF, 0, 0, 1, &cells);CHKERRQ(ierr); 89 ierr = DMGetLocalVector(dm, &locX);CHKERRQ(ierr); 90 ierr = DMGetLocalVector(dm, &locX_t);CHKERRQ(ierr); 91 ierr = DMGetLocalVector(dm, &locA);CHKERRQ(ierr); 92 ierr = DMPlexGetCellFields( dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr); 93 ierr = DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr); 94 ierr = DMRestoreLocalVector(dm, &locX);CHKERRQ(ierr); 95 ierr = DMRestoreLocalVector(dm, &locX_t);CHKERRQ(ierr); 96 ierr = DMRestoreLocalVector(dm, &locA);CHKERRQ(ierr); 97 ierr = ISDestroy(&cells);CHKERRQ(ierr); 98 PetscFunctionReturn(0); 99 } 100 101 /* no discretization is given so DMGetNumFields and PetscDSGetTotalDimension yield 0 */ 102 static PetscErrorCode test2(DM dm, AppCtx *options) 103 { 104 IS cells; 105 Vec locX, locX_t, locA; 106 PetscScalar *u, *u_t, *a; 107 PetscMPIInt rank; 108 PetscErrorCode ierr; 109 110 PetscFunctionBegin; 111 ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank);CHKERRQ(ierr); 112 ierr = ISCreateStride(PETSC_COMM_SELF, rank ? 0 : 1, 0, 1, &cells);CHKERRQ(ierr); 113 ierr = DMGetLocalVector(dm, &locX);CHKERRQ(ierr); 114 ierr = DMGetLocalVector(dm, &locX_t);CHKERRQ(ierr); 115 ierr = DMGetLocalVector(dm, &locA);CHKERRQ(ierr); 116 ierr = DMPlexGetCellFields( dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr); 117 ierr = DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr); 118 ierr = DMRestoreLocalVector(dm, &locX);CHKERRQ(ierr); 119 ierr = DMRestoreLocalVector(dm, &locX_t);CHKERRQ(ierr); 120 ierr = DMRestoreLocalVector(dm, &locA);CHKERRQ(ierr); 121 ierr = ISDestroy(&cells);CHKERRQ(ierr); 122 PetscFunctionReturn(0); 123 } 124 125 static PetscErrorCode test3(DM dm, AppCtx *options) 126 { 127 PetscDS ds; 128 PetscFE fe; 129 PetscErrorCode ierr; 130 131 PetscFunctionBegin; 132 ierr = DMGetDS(dm, &ds);CHKERRQ(ierr); 133 ierr = PetscFECreateDefault(PetscObjectComm((PetscObject) dm), options->dim, 1, options->simplex, NULL, -1, &fe);CHKERRQ(ierr); 134 ierr = PetscDSSetDiscretization(ds, 0, (PetscObject)fe);CHKERRQ(ierr); 135 ierr = PetscFEDestroy(&fe);CHKERRQ(ierr); 136 ierr = test1(dm, options);CHKERRQ(ierr); 137 PetscFunctionReturn(0); 138 } 139 140 static PetscErrorCode test4(DM dm, AppCtx *options) 141 { 142 PetscFE fe; 143 PetscErrorCode ierr; 144 145 PetscFunctionBegin; 146 ierr = PetscFECreateDefault(PetscObjectComm((PetscObject) dm), options->dim, 1, options->simplex, NULL, -1, &fe);CHKERRQ(ierr); 147 ierr = DMSetField(dm, 0, NULL, (PetscObject)fe);CHKERRQ(ierr); 148 ierr = PetscFEDestroy(&fe);CHKERRQ(ierr); 149 ierr = DMCreateDS(dm);CHKERRQ(ierr); 150 ierr = test2(dm, options);CHKERRQ(ierr); 151 PetscFunctionReturn(0); 152 } 153 154 static PetscErrorCode test5(DM dm, AppCtx *options) 155 { 156 IS cells; 157 Vec locX, locX_t, locA; 158 PetscScalar *u, *u_t, *a; 159 PetscErrorCode ierr; 160 161 PetscFunctionBegin; 162 locX_t = NULL; 163 locA = NULL; 164 ierr = ISCreateStride(PETSC_COMM_SELF, 0, 0, 1, &cells);CHKERRQ(ierr); 165 ierr = DMGetLocalVector(dm, &locX);CHKERRQ(ierr); 166 ierr = DMPlexGetCellFields( dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr); 167 ierr = DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr); 168 ierr = DMRestoreLocalVector(dm, &locX);CHKERRQ(ierr); 169 ierr = ISDestroy(&cells);CHKERRQ(ierr); 170 PetscFunctionReturn(0); 171 } 172 173 static PetscErrorCode test6(DM dm, AppCtx *options) 174 { 175 IS cells; 176 Vec locX, locX_t, locA; 177 PetscScalar *u, *u_t, *a; 178 PetscMPIInt rank; 179 PetscErrorCode ierr; 180 181 PetscFunctionBegin; 182 ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank);CHKERRQ(ierr); 183 locX_t = NULL; 184 locA = NULL; 185 ierr = ISCreateStride(PETSC_COMM_SELF, rank ? 0 : 1, 0, 1, &cells);CHKERRQ(ierr); 186 ierr = DMGetLocalVector(dm, &locX);CHKERRQ(ierr); 187 ierr = DMPlexGetCellFields( dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr); 188 ierr = DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr); 189 ierr = DMRestoreLocalVector(dm, &locX);CHKERRQ(ierr); 190 ierr = ISDestroy(&cells);CHKERRQ(ierr); 191 PetscFunctionReturn(0); 192 } 193 194 static PetscErrorCode test7(DM dm, AppCtx *options) 195 { 196 PetscFE fe; 197 PetscErrorCode ierr; 198 199 PetscFunctionBegin; 200 ierr = PetscFECreateDefault(PetscObjectComm((PetscObject) dm), options->dim, 1, options->simplex, NULL, -1, &fe);CHKERRQ(ierr); 201 ierr = DMSetField(dm, 0, NULL, (PetscObject)fe);CHKERRQ(ierr); 202 ierr = PetscFEDestroy(&fe);CHKERRQ(ierr); 203 ierr = DMCreateDS(dm);CHKERRQ(ierr); 204 ierr = test5(dm, options);CHKERRQ(ierr); 205 PetscFunctionReturn(0); 206 } 207 208 static PetscErrorCode test8(DM dm, AppCtx *options) 209 { 210 PetscFE fe; 211 PetscErrorCode ierr; 212 213 PetscFunctionBegin; 214 ierr = PetscFECreateDefault(PetscObjectComm((PetscObject) dm), options->dim, 1, options->simplex, NULL, -1, &fe);CHKERRQ(ierr); 215 ierr = DMSetField(dm, 0, NULL, (PetscObject)fe);CHKERRQ(ierr); 216 ierr = PetscFEDestroy(&fe);CHKERRQ(ierr); 217 ierr = DMCreateDS(dm);CHKERRQ(ierr); 218 ierr = test6(dm, options);CHKERRQ(ierr); 219 PetscFunctionReturn(0); 220 } 221 222 int main(int argc, char **argv) 223 { 224 MPI_Comm comm; 225 DM dm; 226 AppCtx options; 227 PetscErrorCode ierr; 228 229 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 230 comm = PETSC_COMM_WORLD; 231 ierr = ProcessOptions(comm, &options);CHKERRQ(ierr); 232 ierr = CreateMesh(comm, &options, &dm);CHKERRQ(ierr); 233 234 switch (options.test) { 235 case 0: ierr = test0(dm, &options);CHKERRQ(ierr); break; 236 case 1: ierr = test1(dm, &options);CHKERRQ(ierr); break; 237 case 2: ierr = test2(dm, &options);CHKERRQ(ierr); break; 238 case 3: ierr = test3(dm, &options);CHKERRQ(ierr); break; 239 case 4: ierr = test4(dm, &options);CHKERRQ(ierr); break; 240 case 5: ierr = test5(dm, &options);CHKERRQ(ierr); break; 241 case 6: ierr = test6(dm, &options);CHKERRQ(ierr); break; 242 case 7: ierr = test7(dm, &options);CHKERRQ(ierr); break; 243 case 8: ierr = test8(dm, &options);CHKERRQ(ierr); break; 244 default: SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No such test: %D", options.test); 245 } 246 247 ierr = DMDestroy(&dm);CHKERRQ(ierr); 248 ierr = PetscFinalize(); 249 return ierr; 250 } 251 252 /*TEST 253 254 test: 255 suffix: 0 256 requires: ctetgen 257 args: -test 0 258 test: 259 suffix: 1 260 requires: ctetgen 261 args: -test 1 262 test: 263 suffix: 2 264 requires: ctetgen 265 args: -test 2 266 test: 267 suffix: 3 268 requires: ctetgen 269 args: -test 3 270 test: 271 suffix: 4 272 requires: ctetgen 273 args: -test 4 274 test: 275 suffix: 5 276 requires: ctetgen 277 args: -test 5 278 test: 279 suffix: 6 280 requires: ctetgen 281 args: -test 6 282 test: 283 suffix: 7 284 requires: ctetgen 285 args: -test 7 286 test: 287 suffix: 8 288 requires: ctetgen 289 args: -test 8 290 test: 291 suffix: 9 292 requires: ctetgen 293 nsize: 2 294 args: -test 1 295 test: 296 suffix: 10 297 requires: ctetgen 298 nsize: 2 299 args: -test 2 300 test: 301 suffix: 11 302 requires: ctetgen 303 nsize: 2 304 args: -test 3 305 test: 306 suffix: 12 307 requires: ctetgen 308 nsize: 2 309 args: -test 4 310 311 TEST*/ 312 313