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