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