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 typedef struct { 9 PetscInt test; 10 } AppCtx; 11 12 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 13 { 14 PetscFunctionBegin; 15 options->test = 0; 16 PetscOptionsBegin(comm, "", "Zero-sized DMPlexGetCellFields Test Options", "DMPLEX"); 17 PetscCall(PetscOptionsBoundedInt("-test", "Test to run", FILENAME, options->test, &options->test, NULL,0)); 18 PetscOptionsEnd(); 19 PetscFunctionReturn(0); 20 } 21 22 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *options, DM *dm) 23 { 24 PetscFunctionBegin; 25 PetscCall(DMCreate(comm, dm)); 26 PetscCall(DMSetType(*dm, DMPLEX)); 27 PetscCall(DMSetFromOptions(*dm)); 28 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 29 PetscFunctionReturn(0); 30 } 31 32 /* no discretization is given so DMGetNumFields yields 0 */ 33 static PetscErrorCode test0(DM dm, AppCtx *options) 34 { 35 Vec locX; 36 37 PetscFunctionBegin; 38 PetscCall(DMGetLocalVector(dm, &locX)); 39 PetscCall(DMRestoreLocalVector(dm, &locX)); 40 PetscFunctionReturn(0); 41 } 42 43 /* no discretization is given so DMGetNumFields and PetscDSGetTotalDimension yield 0 */ 44 static PetscErrorCode test1(DM dm, AppCtx *options) 45 { 46 IS cells; 47 Vec locX, locX_t, locA; 48 PetscScalar *u, *u_t, *a; 49 50 PetscFunctionBegin; 51 PetscCall(ISCreateStride(PETSC_COMM_SELF, 0, 0, 1, &cells)); 52 PetscCall(DMGetLocalVector(dm, &locX)); 53 PetscCall(DMGetLocalVector(dm, &locX_t)); 54 PetscCall(DMGetLocalVector(dm, &locA)); 55 PetscCall(DMPlexGetCellFields( dm, cells, locX, locX_t, locA, &u, &u_t, &a)); 56 PetscCall(DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a)); 57 PetscCall(DMRestoreLocalVector(dm, &locX)); 58 PetscCall(DMRestoreLocalVector(dm, &locX_t)); 59 PetscCall(DMRestoreLocalVector(dm, &locA)); 60 PetscCall(ISDestroy(&cells)); 61 PetscFunctionReturn(0); 62 } 63 64 /* no discretization is given so DMGetNumFields and PetscDSGetTotalDimension yield 0 */ 65 static PetscErrorCode test2(DM dm, AppCtx *options) 66 { 67 IS cells; 68 Vec locX, locX_t, locA; 69 PetscScalar *u, *u_t, *a; 70 PetscMPIInt rank; 71 72 PetscFunctionBegin; 73 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank)); 74 PetscCall(ISCreateStride(PETSC_COMM_SELF, rank ? 0 : 1, 0, 1, &cells)); 75 PetscCall(DMGetLocalVector(dm, &locX)); 76 PetscCall(DMGetLocalVector(dm, &locX_t)); 77 PetscCall(DMGetLocalVector(dm, &locA)); 78 PetscCall(DMPlexGetCellFields( dm, cells, locX, locX_t, locA, &u, &u_t, &a)); 79 PetscCall(DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a)); 80 PetscCall(DMRestoreLocalVector(dm, &locX)); 81 PetscCall(DMRestoreLocalVector(dm, &locX_t)); 82 PetscCall(DMRestoreLocalVector(dm, &locA)); 83 PetscCall(ISDestroy(&cells)); 84 PetscFunctionReturn(0); 85 } 86 87 static PetscErrorCode test3(DM dm, AppCtx *options) 88 { 89 PetscDS ds; 90 PetscFE fe; 91 PetscInt dim; 92 PetscBool simplex; 93 94 PetscFunctionBegin; 95 PetscCall(DMGetDimension(dm, &dim)); 96 PetscCall(DMPlexIsSimplex(dm, &simplex)); 97 PetscCall(DMGetDS(dm, &ds)); 98 PetscCall(PetscFECreateDefault(PetscObjectComm((PetscObject) dm), dim, 1, simplex, NULL, -1, &fe)); 99 PetscCall(PetscDSSetDiscretization(ds, 0, (PetscObject)fe)); 100 PetscCall(PetscFEDestroy(&fe)); 101 PetscCall(test1(dm, options)); 102 PetscFunctionReturn(0); 103 } 104 105 static PetscErrorCode test4(DM dm, AppCtx *options) 106 { 107 PetscFE fe; 108 PetscInt dim; 109 PetscBool simplex; 110 111 PetscFunctionBegin; 112 PetscCall(DMGetDimension(dm, &dim)); 113 PetscCall(DMPlexIsSimplex(dm, &simplex)); 114 PetscCall(PetscFECreateDefault(PetscObjectComm((PetscObject) dm), dim, 1, simplex, NULL, -1, &fe)); 115 PetscCall(DMSetField(dm, 0, NULL, (PetscObject)fe)); 116 PetscCall(PetscFEDestroy(&fe)); 117 PetscCall(DMCreateDS(dm)); 118 PetscCall(test2(dm, options)); 119 PetscFunctionReturn(0); 120 } 121 122 static PetscErrorCode test5(DM dm, AppCtx *options) 123 { 124 IS cells; 125 Vec locX, locX_t, locA; 126 PetscScalar *u, *u_t, *a; 127 128 PetscFunctionBegin; 129 locX_t = NULL; 130 locA = NULL; 131 PetscCall(ISCreateStride(PETSC_COMM_SELF, 0, 0, 1, &cells)); 132 PetscCall(DMGetLocalVector(dm, &locX)); 133 PetscCall(DMPlexGetCellFields( dm, cells, locX, locX_t, locA, &u, &u_t, &a)); 134 PetscCall(DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a)); 135 PetscCall(DMRestoreLocalVector(dm, &locX)); 136 PetscCall(ISDestroy(&cells)); 137 PetscFunctionReturn(0); 138 } 139 140 static PetscErrorCode test6(DM dm, AppCtx *options) 141 { 142 IS cells; 143 Vec locX, locX_t, locA; 144 PetscScalar *u, *u_t, *a; 145 PetscMPIInt rank; 146 147 PetscFunctionBegin; 148 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank)); 149 locX_t = NULL; 150 locA = NULL; 151 PetscCall(ISCreateStride(PETSC_COMM_SELF, rank ? 0 : 1, 0, 1, &cells)); 152 PetscCall(DMGetLocalVector(dm, &locX)); 153 PetscCall(DMPlexGetCellFields( dm, cells, locX, locX_t, locA, &u, &u_t, &a)); 154 PetscCall(DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a)); 155 PetscCall(DMRestoreLocalVector(dm, &locX)); 156 PetscCall(ISDestroy(&cells)); 157 PetscFunctionReturn(0); 158 } 159 160 static PetscErrorCode test7(DM dm, AppCtx *options) 161 { 162 PetscFE fe; 163 PetscInt dim; 164 PetscBool simplex; 165 166 PetscFunctionBegin; 167 PetscCall(DMGetDimension(dm, &dim)); 168 PetscCall(DMPlexIsSimplex(dm, &simplex)); 169 PetscCall(PetscFECreateDefault(PetscObjectComm((PetscObject) dm), dim, 1, simplex, NULL, -1, &fe)); 170 PetscCall(DMSetField(dm, 0, NULL, (PetscObject)fe)); 171 PetscCall(PetscFEDestroy(&fe)); 172 PetscCall(DMCreateDS(dm)); 173 PetscCall(test5(dm, options)); 174 PetscFunctionReturn(0); 175 } 176 177 static PetscErrorCode test8(DM dm, AppCtx *options) 178 { 179 PetscFE fe; 180 PetscInt dim; 181 PetscBool simplex; 182 183 PetscFunctionBegin; 184 PetscCall(DMGetDimension(dm, &dim)); 185 PetscCall(DMPlexIsSimplex(dm, &simplex)); 186 PetscCall(PetscFECreateDefault(PetscObjectComm((PetscObject) dm), dim, 1, simplex, NULL, -1, &fe)); 187 PetscCall(DMSetField(dm, 0, NULL, (PetscObject)fe)); 188 PetscCall(PetscFEDestroy(&fe)); 189 PetscCall(DMCreateDS(dm)); 190 PetscCall(test6(dm, options)); 191 PetscFunctionReturn(0); 192 } 193 194 int main(int argc, char **argv) 195 { 196 MPI_Comm comm; 197 DM dm; 198 AppCtx options; 199 200 PetscCall(PetscInitialize(&argc, &argv, NULL,help)); 201 comm = PETSC_COMM_WORLD; 202 PetscCall(ProcessOptions(comm, &options)); 203 PetscCall(CreateMesh(comm, &options, &dm)); 204 205 switch (options.test) { 206 case 0: PetscCall(test0(dm, &options)); break; 207 case 1: PetscCall(test1(dm, &options)); break; 208 case 2: PetscCall(test2(dm, &options)); break; 209 case 3: PetscCall(test3(dm, &options)); break; 210 case 4: PetscCall(test4(dm, &options)); break; 211 case 5: PetscCall(test5(dm, &options)); break; 212 case 6: PetscCall(test6(dm, &options)); break; 213 case 7: PetscCall(test7(dm, &options)); break; 214 case 8: PetscCall(test8(dm, &options)); break; 215 default: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No such test: %" PetscInt_FMT, options.test); 216 } 217 218 PetscCall(DMDestroy(&dm)); 219 PetscCall(PetscFinalize()); 220 return 0; 221 } 222 223 /*TEST 224 225 testset: 226 args: -dm_plex_dim 3 -dm_plex_interpolate 0 227 228 test: 229 suffix: 0 230 requires: ctetgen 231 args: -test 0 232 test: 233 suffix: 1 234 requires: ctetgen 235 args: -test 1 236 test: 237 suffix: 2 238 requires: ctetgen 239 args: -test 2 240 test: 241 suffix: 3 242 requires: ctetgen 243 args: -test 3 244 test: 245 suffix: 4 246 requires: ctetgen 247 args: -test 4 248 test: 249 suffix: 5 250 requires: ctetgen 251 args: -test 5 252 test: 253 suffix: 6 254 requires: ctetgen 255 args: -test 6 256 test: 257 suffix: 7 258 requires: ctetgen 259 args: -test 7 260 test: 261 suffix: 8 262 requires: ctetgen 263 args: -test 8 264 test: 265 suffix: 9 266 requires: ctetgen 267 nsize: 2 268 args: -test 1 269 test: 270 suffix: 10 271 requires: ctetgen 272 nsize: 2 273 args: -test 2 274 test: 275 suffix: 11 276 requires: ctetgen 277 nsize: 2 278 args: -test 3 279 test: 280 suffix: 12 281 requires: ctetgen 282 nsize: 2 283 args: -test 4 284 285 TEST*/ 286