1 static char help[] = "Test for function and field projection\n\n"; 2 3 #include <petscdmplex.h> 4 #include <petscds.h> 5 6 typedef struct { 7 PetscBool multifield; /* Different numbers of input and output fields */ 8 PetscBool subdomain; /* Try with a volumetric submesh */ 9 PetscBool submesh; /* Try with a boundary submesh */ 10 PetscBool auxfield; /* Try with auxiliary fields */ 11 } AppCtx; 12 13 /* (x + y)*dim + d */ 14 static PetscErrorCode linear(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx) { 15 PetscInt c; 16 for (c = 0; c < Nc; ++c) u[c] = (x[0] + x[1]) * Nc + c; 17 return 0; 18 } 19 20 /* {x, y, z} */ 21 static PetscErrorCode linear2(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx) { 22 PetscInt c; 23 for (c = 0; c < Nc; ++c) u[c] = x[c]; 24 return 0; 25 } 26 27 /* {u_x, u_y, u_z} */ 28 static void linear_vector(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 f[]) { 29 PetscInt d; 30 for (d = 0; d < uOff[1] - uOff[0]; ++d) f[d] = u[d + uOff[0]]; 31 } 32 33 /* p */ 34 static void linear_scalar(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 f[]) { 35 f[0] = u[uOff[1]]; 36 } 37 38 /* {div u, p^2} */ 39 static void divergence_sq(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 f[]) { 40 PetscInt d; 41 f[0] = 0.0; 42 for (d = 0; d < dim; ++d) f[0] += u_x[uOff_x[0] + d * dim + d]; 43 f[1] = PetscSqr(u[uOff[1]]); 44 } 45 46 static PetscErrorCode ProcessOptions(AppCtx *options) { 47 PetscFunctionBegin; 48 options->multifield = PETSC_FALSE; 49 options->subdomain = PETSC_FALSE; 50 options->submesh = PETSC_FALSE; 51 options->auxfield = PETSC_FALSE; 52 53 PetscOptionsBegin(PETSC_COMM_SELF, "", "Meshing Problem Options", "DMPLEX"); 54 PetscCall(PetscOptionsBool("-multifield", "Flag for trying different numbers of input/output fields", "ex23.c", options->multifield, &options->multifield, NULL)); 55 PetscCall(PetscOptionsBool("-subdomain", "Flag for trying volumetric submesh", "ex23.c", options->subdomain, &options->subdomain, NULL)); 56 PetscCall(PetscOptionsBool("-submesh", "Flag for trying boundary submesh", "ex23.c", options->submesh, &options->submesh, NULL)); 57 PetscCall(PetscOptionsBool("-auxfield", "Flag for trying auxiliary fields", "ex23.c", options->auxfield, &options->auxfield, NULL)); 58 PetscOptionsEnd(); 59 PetscFunctionReturn(0); 60 } 61 62 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) { 63 PetscFunctionBegin; 64 PetscCall(DMCreate(comm, dm)); 65 PetscCall(DMSetType(*dm, DMPLEX)); 66 PetscCall(DMSetFromOptions(*dm)); 67 PetscCall(DMViewFromOptions(*dm, NULL, "-orig_dm_view")); 68 PetscFunctionReturn(0); 69 } 70 71 static PetscErrorCode SetupDiscretization(DM dm, PetscInt dim, PetscBool simplex, AppCtx *user) { 72 PetscFE fe; 73 MPI_Comm comm; 74 75 PetscFunctionBeginUser; 76 PetscCall(PetscObjectGetComm((PetscObject)dm, &comm)); 77 PetscCall(PetscFECreateDefault(comm, dim, dim, simplex, "velocity_", -1, &fe)); 78 PetscCall(PetscObjectSetName((PetscObject)fe, "velocity")); 79 PetscCall(DMSetField(dm, 0, NULL, (PetscObject)fe)); 80 PetscCall(PetscFEDestroy(&fe)); 81 PetscCall(PetscFECreateDefault(comm, dim, 1, simplex, "pressure_", -1, &fe)); 82 PetscCall(PetscObjectSetName((PetscObject)fe, "pressure")); 83 PetscCall(DMSetField(dm, 1, NULL, (PetscObject)fe)); 84 PetscCall(PetscFEDestroy(&fe)); 85 PetscCall(DMCreateDS(dm)); 86 PetscFunctionReturn(0); 87 } 88 89 static PetscErrorCode SetupOutputDiscretization(DM dm, PetscInt dim, PetscBool simplex, AppCtx *user) { 90 PetscFE fe; 91 MPI_Comm comm; 92 93 PetscFunctionBeginUser; 94 PetscCall(PetscObjectGetComm((PetscObject)dm, &comm)); 95 PetscCall(PetscFECreateDefault(comm, dim, dim, simplex, "output_", -1, &fe)); 96 PetscCall(PetscObjectSetName((PetscObject)fe, "output")); 97 PetscCall(DMSetField(dm, 0, NULL, (PetscObject)fe)); 98 PetscCall(PetscFEDestroy(&fe)); 99 PetscCall(DMCreateDS(dm)); 100 PetscFunctionReturn(0); 101 } 102 103 static PetscErrorCode CreateSubdomainMesh(DM dm, DMLabel *domLabel, DM *subdm, AppCtx *user) { 104 DMLabel label; 105 PetscBool simplex; 106 PetscInt dim, cStart, cEnd, c; 107 108 PetscFunctionBeginUser; 109 PetscCall(DMPlexIsSimplex(dm, &simplex)); 110 PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd)); 111 PetscCall(DMLabelCreate(PETSC_COMM_SELF, "subdomain", &label)); 112 for (c = cStart + (cEnd - cStart) / 2; c < cEnd; ++c) PetscCall(DMLabelSetValue(label, c, 1)); 113 PetscCall(DMPlexFilter(dm, label, 1, subdm)); 114 PetscCall(DMGetDimension(*subdm, &dim)); 115 PetscCall(SetupDiscretization(*subdm, dim, simplex, user)); 116 PetscCall(PetscObjectSetName((PetscObject)*subdm, "subdomain")); 117 PetscCall(DMViewFromOptions(*subdm, NULL, "-sub_dm_view")); 118 if (domLabel) *domLabel = label; 119 else PetscCall(DMLabelDestroy(&label)); 120 PetscFunctionReturn(0); 121 } 122 123 static PetscErrorCode CreateBoundaryMesh(DM dm, DMLabel *bdLabel, DM *subdm, AppCtx *user) { 124 DMLabel label; 125 PetscBool simplex; 126 PetscInt dim; 127 128 PetscFunctionBeginUser; 129 PetscCall(DMPlexIsSimplex(dm, &simplex)); 130 PetscCall(DMLabelCreate(PETSC_COMM_SELF, "sub", &label)); 131 PetscCall(DMPlexMarkBoundaryFaces(dm, 1, label)); 132 PetscCall(DMPlexLabelComplete(dm, label)); 133 PetscCall(DMPlexCreateSubmesh(dm, label, 1, PETSC_TRUE, subdm)); 134 PetscCall(DMGetDimension(*subdm, &dim)); 135 PetscCall(SetupDiscretization(*subdm, dim, simplex, user)); 136 PetscCall(PetscObjectSetName((PetscObject)*subdm, "boundary")); 137 PetscCall(DMViewFromOptions(*subdm, NULL, "-sub_dm_view")); 138 if (bdLabel) *bdLabel = label; 139 else PetscCall(DMLabelDestroy(&label)); 140 PetscFunctionReturn(0); 141 } 142 143 static PetscErrorCode CreateAuxiliaryVec(DM dm, DM *auxdm, Vec *la, AppCtx *user) { 144 PetscErrorCode (**afuncs)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar *, void *); 145 PetscBool simplex; 146 PetscInt dim, Nf, f; 147 148 PetscFunctionBeginUser; 149 PetscCall(DMGetDimension(dm, &dim)); 150 PetscCall(DMPlexIsSimplex(dm, &simplex)); 151 PetscCall(DMGetNumFields(dm, &Nf)); 152 PetscCall(PetscMalloc1(Nf, &afuncs)); 153 for (f = 0; f < Nf; ++f) afuncs[f] = linear; 154 PetscCall(DMClone(dm, auxdm)); 155 PetscCall(SetupDiscretization(*auxdm, dim, simplex, user)); 156 PetscCall(DMCreateLocalVector(*auxdm, la)); 157 PetscCall(DMProjectFunctionLocal(dm, 0.0, afuncs, NULL, INSERT_VALUES, *la)); 158 PetscCall(VecViewFromOptions(*la, NULL, "-local_aux_view")); 159 PetscCall(PetscFree(afuncs)); 160 PetscFunctionReturn(0); 161 } 162 163 static PetscErrorCode TestFunctionProjection(DM dm, DM dmAux, DMLabel label, Vec la, const char name[], AppCtx *user) { 164 PetscErrorCode (**funcs)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar *, void *); 165 Vec x, lx; 166 PetscInt Nf, f; 167 PetscInt val[1] = {1}; 168 char lname[PETSC_MAX_PATH_LEN]; 169 170 PetscFunctionBeginUser; 171 if (dmAux) PetscCall(DMSetAuxiliaryVec(dm, NULL, 0, 0, la)); 172 PetscCall(DMGetNumFields(dm, &Nf)); 173 PetscCall(PetscMalloc1(Nf, &funcs)); 174 for (f = 0; f < Nf; ++f) funcs[f] = linear; 175 PetscCall(DMGetGlobalVector(dm, &x)); 176 PetscCall(PetscStrcpy(lname, "Function ")); 177 PetscCall(PetscStrcat(lname, name)); 178 PetscCall(PetscObjectSetName((PetscObject)x, lname)); 179 if (!label) PetscCall(DMProjectFunction(dm, 0.0, funcs, NULL, INSERT_VALUES, x)); 180 else PetscCall(DMProjectFunctionLabel(dm, 0.0, label, 1, val, 0, NULL, funcs, NULL, INSERT_VALUES, x)); 181 PetscCall(VecViewFromOptions(x, NULL, "-func_view")); 182 PetscCall(DMRestoreGlobalVector(dm, &x)); 183 PetscCall(DMGetLocalVector(dm, &lx)); 184 PetscCall(PetscStrcpy(lname, "Local Function ")); 185 PetscCall(PetscStrcat(lname, name)); 186 PetscCall(PetscObjectSetName((PetscObject)lx, lname)); 187 if (!label) PetscCall(DMProjectFunctionLocal(dm, 0.0, funcs, NULL, INSERT_VALUES, lx)); 188 else PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, label, 1, val, 0, NULL, funcs, NULL, INSERT_VALUES, lx)); 189 PetscCall(VecViewFromOptions(lx, NULL, "-local_func_view")); 190 PetscCall(DMRestoreLocalVector(dm, &lx)); 191 PetscCall(PetscFree(funcs)); 192 if (dmAux) PetscCall(DMSetAuxiliaryVec(dm, NULL, 0, 0, NULL)); 193 PetscFunctionReturn(0); 194 } 195 196 static PetscErrorCode TestFieldProjection(DM dm, DM dmAux, DMLabel label, Vec la, const char name[], AppCtx *user) { 197 PetscErrorCode (**afuncs)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar *, void *); 198 void (**funcs)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]); 199 Vec lx, lu; 200 PetscInt Nf, f; 201 PetscInt val[1] = {1}; 202 char lname[PETSC_MAX_PATH_LEN]; 203 204 PetscFunctionBeginUser; 205 if (dmAux) PetscCall(DMSetAuxiliaryVec(dm, NULL, 0, 0, la)); 206 PetscCall(DMGetNumFields(dm, &Nf)); 207 PetscCall(PetscMalloc2(Nf, &funcs, Nf, &afuncs)); 208 for (f = 0; f < Nf; ++f) afuncs[f] = linear; 209 funcs[0] = linear_vector; 210 funcs[1] = linear_scalar; 211 PetscCall(DMGetLocalVector(dm, &lu)); 212 PetscCall(PetscStrcpy(lname, "Local Field Input ")); 213 PetscCall(PetscStrcat(lname, name)); 214 PetscCall(PetscObjectSetName((PetscObject)lu, lname)); 215 if (!label) PetscCall(DMProjectFunctionLocal(dm, 0.0, afuncs, NULL, INSERT_VALUES, lu)); 216 else PetscCall(DMProjectFunctionLabelLocal(dm, 0.0, label, 1, val, 0, NULL, afuncs, NULL, INSERT_VALUES, lu)); 217 PetscCall(VecViewFromOptions(lu, NULL, "-local_input_view")); 218 PetscCall(DMGetLocalVector(dm, &lx)); 219 PetscCall(PetscStrcpy(lname, "Local Field ")); 220 PetscCall(PetscStrcat(lname, name)); 221 PetscCall(PetscObjectSetName((PetscObject)lx, lname)); 222 if (!label) PetscCall(DMProjectFieldLocal(dm, 0.0, lu, funcs, INSERT_VALUES, lx)); 223 else PetscCall(DMProjectFieldLabelLocal(dm, 0.0, label, 1, val, 0, NULL, lu, funcs, INSERT_VALUES, lx)); 224 PetscCall(VecViewFromOptions(lx, NULL, "-local_field_view")); 225 PetscCall(DMRestoreLocalVector(dm, &lx)); 226 PetscCall(DMRestoreLocalVector(dm, &lu)); 227 PetscCall(PetscFree2(funcs, afuncs)); 228 if (dmAux) PetscCall(DMSetAuxiliaryVec(dm, NULL, 0, 0, NULL)); 229 PetscFunctionReturn(0); 230 } 231 232 static PetscErrorCode TestFieldProjectionMultiple(DM dm, DM dmIn, DM dmAux, DMLabel label, Vec la, const char name[], AppCtx *user) { 233 PetscErrorCode (**afuncs)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar *, void *); 234 void (**funcs)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]); 235 Vec lx, lu; 236 PetscInt Nf, NfIn; 237 PetscInt val[1] = {1}; 238 char lname[PETSC_MAX_PATH_LEN]; 239 240 PetscFunctionBeginUser; 241 if (dmAux) PetscCall(DMSetAuxiliaryVec(dm, NULL, 0, 0, la)); 242 PetscCall(DMGetNumFields(dm, &Nf)); 243 PetscCall(DMGetNumFields(dmIn, &NfIn)); 244 PetscCall(PetscMalloc2(Nf, &funcs, NfIn, &afuncs)); 245 funcs[0] = divergence_sq; 246 afuncs[0] = linear2; 247 afuncs[1] = linear; 248 PetscCall(DMGetLocalVector(dmIn, &lu)); 249 PetscCall(PetscStrcpy(lname, "Local MultiField Input ")); 250 PetscCall(PetscStrcat(lname, name)); 251 PetscCall(PetscObjectSetName((PetscObject)lu, lname)); 252 if (!label) PetscCall(DMProjectFunctionLocal(dmIn, 0.0, afuncs, NULL, INSERT_VALUES, lu)); 253 else PetscCall(DMProjectFunctionLabelLocal(dmIn, 0.0, label, 1, val, 0, NULL, afuncs, NULL, INSERT_VALUES, lu)); 254 PetscCall(VecViewFromOptions(lu, NULL, "-local_input_view")); 255 PetscCall(DMGetLocalVector(dm, &lx)); 256 PetscCall(PetscStrcpy(lname, "Local MultiField ")); 257 PetscCall(PetscStrcat(lname, name)); 258 PetscCall(PetscObjectSetName((PetscObject)lx, lname)); 259 if (!label) PetscCall(DMProjectFieldLocal(dm, 0.0, lu, funcs, INSERT_VALUES, lx)); 260 else PetscCall(DMProjectFieldLabelLocal(dm, 0.0, label, 1, val, 0, NULL, lu, funcs, INSERT_VALUES, lx)); 261 PetscCall(VecViewFromOptions(lx, NULL, "-local_field_view")); 262 PetscCall(DMRestoreLocalVector(dm, &lx)); 263 PetscCall(DMRestoreLocalVector(dmIn, &lu)); 264 PetscCall(PetscFree2(funcs, afuncs)); 265 if (dmAux) PetscCall(DMSetAuxiliaryVec(dm, NULL, 0, 0, NULL)); 266 PetscFunctionReturn(0); 267 } 268 269 int main(int argc, char **argv) { 270 DM dm, subdm, auxdm; 271 Vec la; 272 PetscInt dim; 273 PetscBool simplex; 274 AppCtx user; 275 276 PetscFunctionBeginUser; 277 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 278 PetscCall(ProcessOptions(&user)); 279 PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm)); 280 PetscCall(DMGetDimension(dm, &dim)); 281 PetscCall(DMPlexIsSimplex(dm, &simplex)); 282 PetscCall(SetupDiscretization(dm, dim, simplex, &user)); 283 /* Volumetric Mesh Projection */ 284 if (!user.multifield) { 285 PetscCall(TestFunctionProjection(dm, NULL, NULL, NULL, "Volumetric Primary", &user)); 286 PetscCall(TestFieldProjection(dm, NULL, NULL, NULL, "Volumetric Primary", &user)); 287 } else { 288 DM dmOut; 289 290 PetscCall(DMClone(dm, &dmOut)); 291 PetscCall(SetupOutputDiscretization(dmOut, dim, simplex, &user)); 292 PetscCall(TestFieldProjectionMultiple(dmOut, dm, NULL, NULL, NULL, "Volumetric Primary", &user)); 293 PetscCall(DMDestroy(&dmOut)); 294 } 295 if (user.auxfield) { 296 /* Volumetric Mesh Projection with Volumetric Data */ 297 PetscCall(CreateAuxiliaryVec(dm, &auxdm, &la, &user)); 298 PetscCall(TestFunctionProjection(dm, auxdm, NULL, la, "Volumetric Primary and Volumetric Auxiliary", &user)); 299 PetscCall(TestFieldProjection(dm, auxdm, NULL, la, "Volumetric Primary and Volumetric Auxiliary", &user)); 300 PetscCall(VecDestroy(&la)); 301 /* Update of Volumetric Auxiliary Data with primary Volumetric Data */ 302 PetscCall(DMGetLocalVector(dm, &la)); 303 PetscCall(VecSet(la, 1.0)); 304 PetscCall(TestFieldProjection(auxdm, dm, NULL, la, "Volumetric Auxiliary Update with Volumetric Primary", &user)); 305 PetscCall(DMRestoreLocalVector(dm, &la)); 306 PetscCall(DMDestroy(&auxdm)); 307 } 308 if (user.subdomain) { 309 DMLabel domLabel; 310 311 /* Subdomain Mesh Projection */ 312 PetscCall(CreateSubdomainMesh(dm, &domLabel, &subdm, &user)); 313 PetscCall(TestFunctionProjection(subdm, NULL, NULL, NULL, "Subdomain Primary", &user)); 314 PetscCall(TestFieldProjection(subdm, NULL, NULL, NULL, "Subdomain Primary", &user)); 315 if (user.auxfield) { 316 /* Subdomain Mesh Projection with Subdomain Data */ 317 PetscCall(CreateAuxiliaryVec(subdm, &auxdm, &la, &user)); 318 PetscCall(TestFunctionProjection(subdm, auxdm, NULL, la, "Subdomain Primary and Subdomain Auxiliary", &user)); 319 PetscCall(TestFieldProjection(subdm, auxdm, NULL, la, "Subdomain Primary and Subdomain Auxiliary", &user)); 320 PetscCall(VecDestroy(&la)); 321 PetscCall(DMDestroy(&auxdm)); 322 /* Subdomain Mesh Projection with Volumetric Data */ 323 PetscCall(CreateAuxiliaryVec(dm, &auxdm, &la, &user)); 324 PetscCall(TestFunctionProjection(subdm, auxdm, NULL, la, "Subdomain Primary and Volumetric Auxiliary", &user)); 325 PetscCall(TestFieldProjection(subdm, auxdm, NULL, la, "Subdomain Primary and Volumetric Auxiliary", &user)); 326 PetscCall(VecDestroy(&la)); 327 PetscCall(DMDestroy(&auxdm)); 328 /* Volumetric Mesh Projection with Subdomain Data */ 329 PetscCall(CreateAuxiliaryVec(subdm, &auxdm, &la, &user)); 330 PetscCall(TestFunctionProjection(subdm, auxdm, domLabel, la, "Volumetric Primary and Subdomain Auxiliary", &user)); 331 PetscCall(TestFieldProjection(subdm, auxdm, domLabel, la, "Volumetric Primary and Subdomain Auxiliary", &user)); 332 PetscCall(VecDestroy(&la)); 333 PetscCall(DMDestroy(&auxdm)); 334 } 335 PetscCall(DMDestroy(&subdm)); 336 PetscCall(DMLabelDestroy(&domLabel)); 337 } 338 if (user.submesh) { 339 DMLabel bdLabel; 340 341 /* Boundary Mesh Projection */ 342 PetscCall(CreateBoundaryMesh(dm, &bdLabel, &subdm, &user)); 343 PetscCall(TestFunctionProjection(subdm, NULL, NULL, NULL, "Boundary Primary", &user)); 344 PetscCall(TestFieldProjection(subdm, NULL, NULL, NULL, "Boundary Primary", &user)); 345 if (user.auxfield) { 346 /* Boundary Mesh Projection with Boundary Data */ 347 PetscCall(CreateAuxiliaryVec(subdm, &auxdm, &la, &user)); 348 PetscCall(TestFunctionProjection(subdm, auxdm, NULL, la, "Boundary Primary and Boundary Auxiliary", &user)); 349 PetscCall(TestFieldProjection(subdm, auxdm, NULL, la, "Boundary Primary and Boundary Auxiliary", &user)); 350 PetscCall(VecDestroy(&la)); 351 PetscCall(DMDestroy(&auxdm)); 352 /* Volumetric Mesh Projection with Boundary Data */ 353 PetscCall(CreateAuxiliaryVec(subdm, &auxdm, &la, &user)); 354 PetscCall(TestFunctionProjection(dm, auxdm, bdLabel, la, "Volumetric Primary and Boundary Auxiliary", &user)); 355 PetscCall(TestFieldProjection(dm, auxdm, bdLabel, la, "Volumetric Primary and Boundary Auxiliary", &user)); 356 PetscCall(VecDestroy(&la)); 357 PetscCall(DMDestroy(&auxdm)); 358 } 359 PetscCall(DMLabelDestroy(&bdLabel)); 360 PetscCall(DMDestroy(&subdm)); 361 } 362 PetscCall(DMDestroy(&dm)); 363 PetscCall(PetscFinalize()); 364 return 0; 365 } 366 367 /*TEST 368 369 test: 370 suffix: 0 371 requires: triangle 372 args: -dm_plex_box_faces 1,1 -func_view -local_func_view -local_input_view -local_field_view 373 test: 374 suffix: mf_0 375 requires: triangle 376 args: -dm_plex_box_faces 1,1 -velocity_petscspace_degree 1 -velocity_petscfe_default_quadrature_order 2 \ 377 -pressure_petscspace_degree 2 -pressure_petscfe_default_quadrature_order 2 \ 378 -multifield -output_petscspace_degree 1 -output_petscfe_default_quadrature_order 2 \ 379 -local_input_view -local_field_view 380 test: 381 suffix: 1 382 requires: triangle 383 args: -dm_plex_box_faces 1,1 -velocity_petscspace_degree 1 -velocity_petscfe_default_quadrature_order 2 -pressure_petscspace_degree 2 -pressure_petscfe_default_quadrature_order 2 -func_view -local_func_view -local_input_view -local_field_view -submesh -auxfield 384 test: 385 suffix: 2 386 requires: triangle 387 args: -dm_plex_box_faces 1,1 -velocity_petscspace_degree 1 -velocity_petscfe_default_quadrature_order 2 -pressure_petscspace_degree 2 -pressure_petscfe_default_quadrature_order 2 -func_view -local_func_view -local_input_view -local_field_view -subdomain -auxfield 388 389 TEST*/ 390 391 /* 392 Post-processing wants to project a function of the fields into some FE space 393 - This is DMProjectField() 394 - What about changing the number of components of the output, like displacement to stress? Aux vars 395 396 Update of state variables 397 - This is DMProjectField(), but solution must be the aux var 398 */ 399