1 #include <petsc/private/petscdsimpl.h> /*I "petscds.h" I*/ 2 3 PetscClassId PETSCDS_CLASSID = 0; 4 5 PetscFunctionList PetscDSList = NULL; 6 PetscBool PetscDSRegisterAllCalled = PETSC_FALSE; 7 8 /* A PetscDS (Discrete System) encodes a set of equations posed in a discrete space, which represents a set of 9 nonlinear continuum equations. The equations can have multiple fields, each field having a different 10 discretization. In addition, different pieces of the domain can have different field combinations and equations. 11 12 The DS provides the user a description of the approximation space on any given cell. It also gives pointwise 13 functions representing the equations. 14 15 Each field is associated with a label, marking the cells on which it is supported. Note that a field can be 16 supported on the closure of a cell not in the label due to overlap of the boundary of neighboring cells. The DM 17 then creates a DS for each set of cells with identical approximation spaces. When assembling, the user asks for 18 the space associated with a given cell. DMPlex uses the labels associated with each DS in the default integration loop. 19 */ 20 21 /*@C 22 PetscDSRegister - Adds a new `PetscDS` implementation 23 24 Not Collective; No Fortran Support 25 26 Input Parameters: 27 + sname - The name of a new user-defined creation routine 28 - function - The creation routine itself 29 30 Example Usage: 31 .vb 32 PetscDSRegister("my_ds", MyPetscDSCreate); 33 .ve 34 35 Then, your PetscDS type can be chosen with the procedural interface via 36 .vb 37 PetscDSCreate(MPI_Comm, PetscDS *); 38 PetscDSSetType(PetscDS, "my_ds"); 39 .ve 40 or at runtime via the option 41 .vb 42 -petscds_type my_ds 43 .ve 44 45 Level: advanced 46 47 Note: 48 `PetscDSRegister()` may be called multiple times to add several user-defined `PetscDSs` 49 50 .seealso: `PetscDSType`, `PetscDS`, `PetscDSRegisterAll()`, `PetscDSRegisterDestroy()` 51 @*/ 52 PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS)) 53 { 54 PetscFunctionBegin; 55 PetscCall(PetscFunctionListAdd(&PetscDSList, sname, function)); 56 PetscFunctionReturn(PETSC_SUCCESS); 57 } 58 59 /*@C 60 PetscDSSetType - Builds a particular `PetscDS` 61 62 Collective; No Fortran Support 63 64 Input Parameters: 65 + prob - The `PetscDS` object 66 - name - The `PetscDSType` 67 68 Options Database Key: 69 . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types 70 71 Level: intermediate 72 73 .seealso: `PetscDSType`, `PetscDS`, `PetscDSGetType()`, `PetscDSCreate()` 74 @*/ 75 PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name) 76 { 77 PetscErrorCode (*r)(PetscDS); 78 PetscBool match; 79 80 PetscFunctionBegin; 81 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 82 PetscCall(PetscObjectTypeCompare((PetscObject)prob, name, &match)); 83 if (match) PetscFunctionReturn(PETSC_SUCCESS); 84 85 PetscCall(PetscDSRegisterAll()); 86 PetscCall(PetscFunctionListFind(PetscDSList, name, &r)); 87 PetscCheck(r, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name); 88 89 PetscTryTypeMethod(prob, destroy); 90 prob->ops->destroy = NULL; 91 92 PetscCall((*r)(prob)); 93 PetscCall(PetscObjectChangeTypeName((PetscObject)prob, name)); 94 PetscFunctionReturn(PETSC_SUCCESS); 95 } 96 97 /*@C 98 PetscDSGetType - Gets the `PetscDSType` name (as a string) from the `PetscDS` 99 100 Not Collective; No Fortran Support 101 102 Input Parameter: 103 . prob - The `PetscDS` 104 105 Output Parameter: 106 . name - The `PetscDSType` name 107 108 Level: intermediate 109 110 .seealso: `PetscDSType`, `PetscDS`, `PetscDSSetType()`, `PetscDSCreate()` 111 @*/ 112 PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name) 113 { 114 PetscFunctionBegin; 115 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 116 PetscAssertPointer(name, 2); 117 PetscCall(PetscDSRegisterAll()); 118 *name = ((PetscObject)prob)->type_name; 119 PetscFunctionReturn(PETSC_SUCCESS); 120 } 121 122 static PetscErrorCode PetscDSView_Ascii(PetscDS ds, PetscViewer viewer) 123 { 124 PetscViewerFormat format; 125 const PetscScalar *constants; 126 PetscInt Nf, numConstants, f; 127 128 PetscFunctionBegin; 129 PetscCall(PetscDSGetNumFields(ds, &Nf)); 130 PetscCall(PetscViewerGetFormat(viewer, &format)); 131 PetscCall(PetscViewerASCIIPrintf(viewer, "Discrete System with %" PetscInt_FMT " fields\n", Nf)); 132 PetscCall(PetscViewerASCIIPushTab(viewer)); 133 PetscCall(PetscViewerASCIIPrintf(viewer, " cell total dim %" PetscInt_FMT " total comp %" PetscInt_FMT "\n", ds->totDim, ds->totComp)); 134 if (ds->isCohesive) PetscCall(PetscViewerASCIIPrintf(viewer, " cohesive cell\n")); 135 for (f = 0; f < Nf; ++f) { 136 DSBoundary b; 137 PetscObject obj; 138 PetscClassId id; 139 PetscQuadrature q; 140 const char *name; 141 PetscInt Nc, Nq, Nqc; 142 143 PetscCall(PetscDSGetDiscretization(ds, f, &obj)); 144 PetscCall(PetscObjectGetClassId(obj, &id)); 145 PetscCall(PetscObjectGetName(obj, &name)); 146 PetscCall(PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>")); 147 PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE)); 148 if (id == PETSCFE_CLASSID) { 149 PetscCall(PetscFEGetNumComponents((PetscFE)obj, &Nc)); 150 PetscCall(PetscFEGetQuadrature((PetscFE)obj, &q)); 151 PetscCall(PetscViewerASCIIPrintf(viewer, " FEM")); 152 } else if (id == PETSCFV_CLASSID) { 153 PetscCall(PetscFVGetNumComponents((PetscFV)obj, &Nc)); 154 PetscCall(PetscFVGetQuadrature((PetscFV)obj, &q)); 155 PetscCall(PetscViewerASCIIPrintf(viewer, " FVM")); 156 } else SETERRQ(PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f); 157 if (Nc > 1) PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " components", Nc)); 158 else PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " component ", Nc)); 159 if (ds->implicit[f]) PetscCall(PetscViewerASCIIPrintf(viewer, " (implicit)")); 160 else PetscCall(PetscViewerASCIIPrintf(viewer, " (explicit)")); 161 if (q) { 162 PetscCall(PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL)); 163 PetscCall(PetscViewerASCIIPrintf(viewer, " (Nq %" PetscInt_FMT " Nqc %" PetscInt_FMT ")", Nq, Nqc)); 164 } 165 PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT "-jet", ds->jetDegree[f])); 166 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 167 PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE)); 168 PetscCall(PetscViewerASCIIPushTab(viewer)); 169 if (id == PETSCFE_CLASSID) PetscCall(PetscFEView((PetscFE)obj, viewer)); 170 else if (id == PETSCFV_CLASSID) PetscCall(PetscFVView((PetscFV)obj, viewer)); 171 PetscCall(PetscViewerASCIIPopTab(viewer)); 172 173 for (b = ds->boundary; b; b = b->next) { 174 char *name; 175 PetscInt c, i; 176 177 if (b->field != f) continue; 178 PetscCall(PetscViewerASCIIPushTab(viewer)); 179 PetscCall(PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->lname, DMBoundaryConditionTypes[b->type])); 180 if (!b->Nc) { 181 PetscCall(PetscViewerASCIIPrintf(viewer, " all components\n")); 182 } else { 183 PetscCall(PetscViewerASCIIPrintf(viewer, " components: ")); 184 PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE)); 185 for (c = 0; c < b->Nc; ++c) { 186 if (c > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", ")); 187 PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->comps[c])); 188 } 189 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 190 PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE)); 191 } 192 PetscCall(PetscViewerASCIIPrintf(viewer, " values: ")); 193 PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE)); 194 for (i = 0; i < b->Nv; ++i) { 195 if (i > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", ")); 196 PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->values[i])); 197 } 198 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 199 PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE)); 200 #if defined(__clang__) 201 PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat-pedantic") 202 #elif defined(__GNUC__) || defined(__GNUG__) 203 PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat") 204 #endif 205 if (b->func) { 206 PetscCall(PetscDLAddr(b->func, &name)); 207 if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func: %s\n", name)); 208 else PetscCall(PetscViewerASCIIPrintf(viewer, " func: %p\n", b->func)); 209 PetscCall(PetscFree(name)); 210 } 211 if (b->func_t) { 212 PetscCall(PetscDLAddr(b->func_t, &name)); 213 if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %s\n", name)); 214 else PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %p\n", b->func_t)); 215 PetscCall(PetscFree(name)); 216 } 217 PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END() 218 PetscCall(PetscWeakFormView(b->wf, viewer)); 219 PetscCall(PetscViewerASCIIPopTab(viewer)); 220 } 221 } 222 PetscCall(PetscDSGetConstants(ds, &numConstants, &constants)); 223 if (numConstants) { 224 PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " constants\n", numConstants)); 225 PetscCall(PetscViewerASCIIPushTab(viewer)); 226 for (f = 0; f < numConstants; ++f) PetscCall(PetscViewerASCIIPrintf(viewer, "%g\n", (double)PetscRealPart(constants[f]))); 227 PetscCall(PetscViewerASCIIPopTab(viewer)); 228 } 229 PetscCall(PetscWeakFormView(ds->wf, viewer)); 230 PetscCall(PetscViewerASCIIPopTab(viewer)); 231 PetscFunctionReturn(PETSC_SUCCESS); 232 } 233 234 /*@C 235 PetscDSViewFromOptions - View a `PetscDS` based on values in the options database 236 237 Collective 238 239 Input Parameters: 240 + A - the `PetscDS` object 241 . obj - Optional object that provides the options prefix used in the search 242 - name - command line option 243 244 Level: intermediate 245 246 .seealso: `PetscDSType`, `PetscDS`, `PetscDSView()`, `PetscObjectViewFromOptions()`, `PetscDSCreate()` 247 @*/ 248 PetscErrorCode PetscDSViewFromOptions(PetscDS A, PetscObject obj, const char name[]) 249 { 250 PetscFunctionBegin; 251 PetscValidHeaderSpecific(A, PETSCDS_CLASSID, 1); 252 PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name)); 253 PetscFunctionReturn(PETSC_SUCCESS); 254 } 255 256 /*@C 257 PetscDSView - Views a `PetscDS` 258 259 Collective 260 261 Input Parameters: 262 + prob - the `PetscDS` object to view 263 - v - the viewer 264 265 Level: developer 266 267 .seealso: `PetscDSType`, `PetscDS`, `PetscViewer`, `PetscDSDestroy()`, `PetscDSViewFromOptions()` 268 @*/ 269 PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v) 270 { 271 PetscBool iascii; 272 273 PetscFunctionBegin; 274 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 275 if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)prob), &v)); 276 else PetscValidHeaderSpecific(v, PETSC_VIEWER_CLASSID, 2); 277 PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii)); 278 if (iascii) PetscCall(PetscDSView_Ascii(prob, v)); 279 PetscTryTypeMethod(prob, view, v); 280 PetscFunctionReturn(PETSC_SUCCESS); 281 } 282 283 /*@ 284 PetscDSSetFromOptions - sets parameters in a `PetscDS` from the options database 285 286 Collective 287 288 Input Parameter: 289 . prob - the `PetscDS` object to set options for 290 291 Options Database Keys: 292 + -petscds_type <type> - Set the `PetscDS` type 293 . -petscds_view <view opt> - View the `PetscDS` 294 . -petscds_jac_pre - Turn formation of a separate Jacobian preconditioner on or off 295 . -bc_<name> <ids> - Specify a list of label ids for a boundary condition 296 - -bc_<name>_comp <comps> - Specify a list of field components to constrain for a boundary condition 297 298 Level: intermediate 299 300 .seealso: `PetscDS`, `PetscDSView()` 301 @*/ 302 PetscErrorCode PetscDSSetFromOptions(PetscDS prob) 303 { 304 DSBoundary b; 305 const char *defaultType; 306 char name[256]; 307 PetscBool flg; 308 309 PetscFunctionBegin; 310 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 311 if (!((PetscObject)prob)->type_name) { 312 defaultType = PETSCDSBASIC; 313 } else { 314 defaultType = ((PetscObject)prob)->type_name; 315 } 316 PetscCall(PetscDSRegisterAll()); 317 318 PetscObjectOptionsBegin((PetscObject)prob); 319 for (b = prob->boundary; b; b = b->next) { 320 char optname[1024]; 321 PetscInt ids[1024], len = 1024; 322 PetscBool flg; 323 324 PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name)); 325 PetscCall(PetscMemzero(ids, sizeof(ids))); 326 PetscCall(PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg)); 327 if (flg) { 328 b->Nv = len; 329 PetscCall(PetscFree(b->values)); 330 PetscCall(PetscMalloc1(len, &b->values)); 331 PetscCall(PetscArraycpy(b->values, ids, len)); 332 PetscCall(PetscWeakFormRewriteKeys(b->wf, b->label, len, b->values)); 333 } 334 len = 1024; 335 PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name)); 336 PetscCall(PetscMemzero(ids, sizeof(ids))); 337 PetscCall(PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg)); 338 if (flg) { 339 b->Nc = len; 340 PetscCall(PetscFree(b->comps)); 341 PetscCall(PetscMalloc1(len, &b->comps)); 342 PetscCall(PetscArraycpy(b->comps, ids, len)); 343 } 344 } 345 PetscCall(PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg)); 346 if (flg) { 347 PetscCall(PetscDSSetType(prob, name)); 348 } else if (!((PetscObject)prob)->type_name) { 349 PetscCall(PetscDSSetType(prob, defaultType)); 350 } 351 PetscCall(PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg)); 352 PetscCall(PetscOptionsBool("-petscds_force_quad", "Discrete System", "PetscDSSetForceQuad", prob->forceQuad, &prob->forceQuad, &flg)); 353 PetscTryTypeMethod(prob, setfromoptions); 354 /* process any options handlers added with PetscObjectAddOptionsHandler() */ 355 PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)prob, PetscOptionsObject)); 356 PetscOptionsEnd(); 357 if (prob->Nf) PetscCall(PetscDSViewFromOptions(prob, NULL, "-petscds_view")); 358 PetscFunctionReturn(PETSC_SUCCESS); 359 } 360 361 /*@C 362 PetscDSSetUp - Construct data structures for the `PetscDS` 363 364 Collective 365 366 Input Parameter: 367 . prob - the `PetscDS` object to setup 368 369 Level: developer 370 371 .seealso: `PetscDS`, `PetscDSView()`, `PetscDSDestroy()` 372 @*/ 373 PetscErrorCode PetscDSSetUp(PetscDS prob) 374 { 375 const PetscInt Nf = prob->Nf; 376 PetscBool hasH = PETSC_FALSE; 377 PetscInt maxOrder[4] = {-1, -1, -1, -1}; 378 PetscInt dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f; 379 380 PetscFunctionBegin; 381 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 382 if (prob->setup) PetscFunctionReturn(PETSC_SUCCESS); 383 /* Calculate sizes */ 384 PetscCall(PetscDSGetSpatialDimension(prob, &dim)); 385 PetscCall(PetscDSGetCoordinateDimension(prob, &dimEmbed)); 386 prob->totDim = prob->totComp = 0; 387 PetscCall(PetscMalloc2(Nf, &prob->Nc, Nf, &prob->Nb)); 388 PetscCall(PetscCalloc2(Nf + 1, &prob->off, Nf + 1, &prob->offDer)); 389 PetscCall(PetscCalloc6(Nf + 1, &prob->offCohesive[0], Nf + 1, &prob->offCohesive[1], Nf + 1, &prob->offCohesive[2], Nf + 1, &prob->offDerCohesive[0], Nf + 1, &prob->offDerCohesive[1], Nf + 1, &prob->offDerCohesive[2])); 390 PetscCall(PetscMalloc2(Nf, &prob->T, Nf, &prob->Tf)); 391 if (prob->forceQuad) { 392 // Note: This assumes we have one kind of cell at each dimension. 393 // We can fix this by having quadrature hold the celltype 394 PetscQuadrature maxQuad[4] = {NULL, NULL, NULL, NULL}; 395 396 for (f = 0; f < Nf; ++f) { 397 PetscObject obj; 398 PetscClassId id; 399 PetscQuadrature q = NULL, fq = NULL; 400 PetscInt dim = -1, order = -1, forder = -1; 401 402 PetscCall(PetscDSGetDiscretization(prob, f, &obj)); 403 if (!obj) continue; 404 PetscCall(PetscObjectGetClassId(obj, &id)); 405 if (id == PETSCFE_CLASSID) { 406 PetscFE fe = (PetscFE)obj; 407 408 PetscCall(PetscFEGetQuadrature(fe, &q)); 409 PetscCall(PetscFEGetFaceQuadrature(fe, &fq)); 410 } else if (id == PETSCFV_CLASSID) { 411 PetscFV fv = (PetscFV)obj; 412 413 PetscCall(PetscFVGetQuadrature(fv, &q)); 414 } 415 if (q) { 416 PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL)); 417 PetscCall(PetscQuadratureGetOrder(q, &order)); 418 if (order > maxOrder[dim]) { 419 maxOrder[dim] = order; 420 maxQuad[dim] = q; 421 } 422 } 423 if (fq) { 424 PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL)); 425 PetscCall(PetscQuadratureGetOrder(fq, &forder)); 426 if (forder > maxOrder[dim]) { 427 maxOrder[dim] = forder; 428 maxQuad[dim] = fq; 429 } 430 } 431 } 432 for (f = 0; f < Nf; ++f) { 433 PetscObject obj; 434 PetscClassId id; 435 PetscQuadrature q; 436 PetscInt dim; 437 438 PetscCall(PetscDSGetDiscretization(prob, f, &obj)); 439 if (!obj) continue; 440 PetscCall(PetscObjectGetClassId(obj, &id)); 441 if (id == PETSCFE_CLASSID) { 442 PetscFE fe = (PetscFE)obj; 443 444 PetscCall(PetscFEGetQuadrature(fe, &q)); 445 PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL)); 446 PetscCall(PetscFESetQuadrature(fe, maxQuad[dim])); 447 PetscCall(PetscFESetFaceQuadrature(fe, dim ? maxQuad[dim - 1] : NULL)); 448 } else if (id == PETSCFV_CLASSID) { 449 PetscFV fv = (PetscFV)obj; 450 451 PetscCall(PetscFVGetQuadrature(fv, &q)); 452 PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL)); 453 PetscCall(PetscFVSetQuadrature(fv, maxQuad[dim])); 454 } 455 } 456 } 457 for (f = 0; f < Nf; ++f) { 458 PetscObject obj; 459 PetscClassId id; 460 PetscQuadrature q = NULL; 461 PetscInt Nq = 0, Nb, Nc; 462 463 PetscCall(PetscDSGetDiscretization(prob, f, &obj)); 464 if (prob->jetDegree[f] > 1) hasH = PETSC_TRUE; 465 if (!obj) { 466 /* Empty mesh */ 467 Nb = Nc = 0; 468 prob->T[f] = prob->Tf[f] = NULL; 469 } else { 470 PetscCall(PetscObjectGetClassId(obj, &id)); 471 if (id == PETSCFE_CLASSID) { 472 PetscFE fe = (PetscFE)obj; 473 474 PetscCall(PetscFEGetQuadrature(fe, &q)); 475 { 476 PetscQuadrature fq; 477 PetscInt dim, order; 478 479 PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL)); 480 PetscCall(PetscQuadratureGetOrder(q, &order)); 481 if (maxOrder[dim] < 0) maxOrder[dim] = order; 482 PetscCheck(order == maxOrder[dim], PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Field %" PetscInt_FMT " cell quadrature order %" PetscInt_FMT " != %" PetscInt_FMT " DS cell quadrature order", f, order, maxOrder[dim]); 483 PetscCall(PetscFEGetFaceQuadrature(fe, &fq)); 484 if (fq) { 485 PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL)); 486 PetscCall(PetscQuadratureGetOrder(fq, &order)); 487 if (maxOrder[dim] < 0) maxOrder[dim] = order; 488 PetscCheck(order == maxOrder[dim], PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Field %" PetscInt_FMT " face quadrature order %" PetscInt_FMT " != %" PetscInt_FMT " DS face quadrature order", f, order, maxOrder[dim]); 489 } 490 } 491 PetscCall(PetscFEGetDimension(fe, &Nb)); 492 PetscCall(PetscFEGetNumComponents(fe, &Nc)); 493 PetscCall(PetscFEGetCellTabulation(fe, prob->jetDegree[f], &prob->T[f])); 494 PetscCall(PetscFEGetFaceTabulation(fe, prob->jetDegree[f], &prob->Tf[f])); 495 } else if (id == PETSCFV_CLASSID) { 496 PetscFV fv = (PetscFV)obj; 497 498 PetscCall(PetscFVGetQuadrature(fv, &q)); 499 PetscCall(PetscFVGetNumComponents(fv, &Nc)); 500 Nb = Nc; 501 PetscCall(PetscFVGetCellTabulation(fv, &prob->T[f])); 502 /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */ 503 } else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f); 504 } 505 prob->Nc[f] = Nc; 506 prob->Nb[f] = Nb; 507 prob->off[f + 1] = Nc + prob->off[f]; 508 prob->offDer[f + 1] = Nc * dim + prob->offDer[f]; 509 prob->offCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[0][f]; 510 prob->offDerCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[0][f]; 511 prob->offCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) + prob->offCohesive[0][f]; 512 prob->offDerCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) * dimEmbed + prob->offDerCohesive[0][f]; 513 prob->offCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[2][f]; 514 prob->offDerCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[2][f]; 515 if (q) PetscCall(PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL)); 516 NqMax = PetscMax(NqMax, Nq); 517 NbMax = PetscMax(NbMax, Nb); 518 NcMax = PetscMax(NcMax, Nc); 519 prob->totDim += Nb; 520 prob->totComp += Nc; 521 /* There are two faces for all fields on a cohesive cell, except for cohesive fields */ 522 if (prob->isCohesive && !prob->cohesive[f]) prob->totDim += Nb; 523 } 524 prob->offCohesive[1][Nf] = prob->offCohesive[0][Nf]; 525 prob->offDerCohesive[1][Nf] = prob->offDerCohesive[0][Nf]; 526 /* Allocate works space */ 527 NsMax = 2; /* A non-cohesive discretizations can be used on a cohesive cell, so we need this extra workspace for all DS */ 528 PetscCall(PetscMalloc3(NsMax * prob->totComp, &prob->u, NsMax * prob->totComp, &prob->u_t, NsMax * prob->totComp * dimEmbed + (hasH ? NsMax * prob->totComp * dimEmbed * dimEmbed : 0), &prob->u_x)); 529 PetscCall(PetscMalloc5(dimEmbed, &prob->x, NbMax * NcMax, &prob->basisReal, NbMax * NcMax * dimEmbed, &prob->basisDerReal, NbMax * NcMax, &prob->testReal, NbMax * NcMax * dimEmbed, &prob->testDerReal)); 530 PetscCall(PetscMalloc6(NsMax * NqMax * NcMax, &prob->f0, NsMax * NqMax * NcMax * dimEmbed, &prob->f1, NsMax * NsMax * NqMax * NcMax * NcMax, &prob->g0, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed, &prob->g1, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed, 531 &prob->g2, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed * dimEmbed, &prob->g3)); 532 PetscTryTypeMethod(prob, setup); 533 prob->setup = PETSC_TRUE; 534 PetscFunctionReturn(PETSC_SUCCESS); 535 } 536 537 static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob) 538 { 539 PetscFunctionBegin; 540 PetscCall(PetscFree2(prob->Nc, prob->Nb)); 541 PetscCall(PetscFree2(prob->off, prob->offDer)); 542 PetscCall(PetscFree6(prob->offCohesive[0], prob->offCohesive[1], prob->offCohesive[2], prob->offDerCohesive[0], prob->offDerCohesive[1], prob->offDerCohesive[2])); 543 PetscCall(PetscFree2(prob->T, prob->Tf)); 544 PetscCall(PetscFree3(prob->u, prob->u_t, prob->u_x)); 545 PetscCall(PetscFree5(prob->x, prob->basisReal, prob->basisDerReal, prob->testReal, prob->testDerReal)); 546 PetscCall(PetscFree6(prob->f0, prob->f1, prob->g0, prob->g1, prob->g2, prob->g3)); 547 PetscFunctionReturn(PETSC_SUCCESS); 548 } 549 550 static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew) 551 { 552 PetscObject *tmpd; 553 PetscBool *tmpi; 554 PetscInt *tmpk; 555 PetscBool *tmpc; 556 PetscPointFunc *tmpup; 557 PetscSimplePointFunc *tmpexactSol, *tmpexactSol_t; 558 void **tmpexactCtx, **tmpexactCtx_t; 559 void **tmpctx; 560 PetscInt Nf = prob->Nf, f; 561 562 PetscFunctionBegin; 563 if (Nf >= NfNew) PetscFunctionReturn(PETSC_SUCCESS); 564 prob->setup = PETSC_FALSE; 565 PetscCall(PetscDSDestroyStructs_Static(prob)); 566 PetscCall(PetscMalloc4(NfNew, &tmpd, NfNew, &tmpi, NfNew, &tmpc, NfNew, &tmpk)); 567 for (f = 0; f < Nf; ++f) { 568 tmpd[f] = prob->disc[f]; 569 tmpi[f] = prob->implicit[f]; 570 tmpc[f] = prob->cohesive[f]; 571 tmpk[f] = prob->jetDegree[f]; 572 } 573 for (f = Nf; f < NfNew; ++f) { 574 tmpd[f] = NULL; 575 tmpi[f] = PETSC_TRUE, tmpc[f] = PETSC_FALSE; 576 tmpk[f] = 1; 577 } 578 PetscCall(PetscFree4(prob->disc, prob->implicit, prob->cohesive, prob->jetDegree)); 579 PetscCall(PetscWeakFormSetNumFields(prob->wf, NfNew)); 580 prob->Nf = NfNew; 581 prob->disc = tmpd; 582 prob->implicit = tmpi; 583 prob->cohesive = tmpc; 584 prob->jetDegree = tmpk; 585 PetscCall(PetscCalloc2(NfNew, &tmpup, NfNew, &tmpctx)); 586 for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f]; 587 for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f]; 588 for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL; 589 for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL; 590 PetscCall(PetscFree2(prob->update, prob->ctx)); 591 prob->update = tmpup; 592 prob->ctx = tmpctx; 593 PetscCall(PetscCalloc4(NfNew, &tmpexactSol, NfNew, &tmpexactCtx, NfNew, &tmpexactSol_t, NfNew, &tmpexactCtx_t)); 594 for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f]; 595 for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f]; 596 for (f = 0; f < Nf; ++f) tmpexactSol_t[f] = prob->exactSol_t[f]; 597 for (f = 0; f < Nf; ++f) tmpexactCtx_t[f] = prob->exactCtx_t[f]; 598 for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL; 599 for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL; 600 for (f = Nf; f < NfNew; ++f) tmpexactSol_t[f] = NULL; 601 for (f = Nf; f < NfNew; ++f) tmpexactCtx_t[f] = NULL; 602 PetscCall(PetscFree4(prob->exactSol, prob->exactCtx, prob->exactSol_t, prob->exactCtx_t)); 603 prob->exactSol = tmpexactSol; 604 prob->exactCtx = tmpexactCtx; 605 prob->exactSol_t = tmpexactSol_t; 606 prob->exactCtx_t = tmpexactCtx_t; 607 PetscFunctionReturn(PETSC_SUCCESS); 608 } 609 610 /*@ 611 PetscDSDestroy - Destroys a `PetscDS` object 612 613 Collective 614 615 Input Parameter: 616 . ds - the `PetscDS` object to destroy 617 618 Level: developer 619 620 .seealso: `PetscDSView()` 621 @*/ 622 PetscErrorCode PetscDSDestroy(PetscDS *ds) 623 { 624 PetscInt f; 625 626 PetscFunctionBegin; 627 if (!*ds) PetscFunctionReturn(PETSC_SUCCESS); 628 PetscValidHeaderSpecific((*ds), PETSCDS_CLASSID, 1); 629 630 if (--((PetscObject)(*ds))->refct > 0) { 631 *ds = NULL; 632 PetscFunctionReturn(PETSC_SUCCESS); 633 } 634 ((PetscObject)(*ds))->refct = 0; 635 if ((*ds)->subprobs) { 636 PetscInt dim, d; 637 638 PetscCall(PetscDSGetSpatialDimension(*ds, &dim)); 639 for (d = 0; d < dim; ++d) PetscCall(PetscDSDestroy(&(*ds)->subprobs[d])); 640 } 641 PetscCall(PetscFree((*ds)->subprobs)); 642 PetscCall(PetscDSDestroyStructs_Static(*ds)); 643 for (f = 0; f < (*ds)->Nf; ++f) PetscCall(PetscObjectDereference((*ds)->disc[f])); 644 PetscCall(PetscFree4((*ds)->disc, (*ds)->implicit, (*ds)->cohesive, (*ds)->jetDegree)); 645 PetscCall(PetscWeakFormDestroy(&(*ds)->wf)); 646 PetscCall(PetscFree2((*ds)->update, (*ds)->ctx)); 647 PetscCall(PetscFree4((*ds)->exactSol, (*ds)->exactCtx, (*ds)->exactSol_t, (*ds)->exactCtx_t)); 648 PetscTryTypeMethod((*ds), destroy); 649 PetscCall(PetscDSDestroyBoundary(*ds)); 650 PetscCall(PetscFree((*ds)->constants)); 651 for (PetscInt c = 0; c < DM_NUM_POLYTOPES; ++c) { 652 const PetscInt Na = DMPolytopeTypeGetNumArrangments((DMPolytopeType)c); 653 if ((*ds)->quadPerm[c]) 654 for (PetscInt o = 0; o < Na; ++o) PetscCall(ISDestroy(&(*ds)->quadPerm[c][o])); 655 PetscCall(PetscFree((*ds)->quadPerm[c])); 656 (*ds)->quadPerm[c] = NULL; 657 } 658 PetscCall(PetscHeaderDestroy(ds)); 659 PetscFunctionReturn(PETSC_SUCCESS); 660 } 661 662 /*@ 663 PetscDSCreate - Creates an empty `PetscDS` object. The type can then be set with `PetscDSSetType()`. 664 665 Collective 666 667 Input Parameter: 668 . comm - The communicator for the `PetscDS` object 669 670 Output Parameter: 671 . ds - The `PetscDS` object 672 673 Level: beginner 674 675 .seealso: `PetscDS`, `PetscDSSetType()`, `PETSCDSBASIC`, `PetscDSType` 676 @*/ 677 PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *ds) 678 { 679 PetscDS p; 680 681 PetscFunctionBegin; 682 PetscAssertPointer(ds, 2); 683 *ds = NULL; 684 PetscCall(PetscDSInitializePackage()); 685 686 PetscCall(PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView)); 687 688 p->Nf = 0; 689 p->setup = PETSC_FALSE; 690 p->numConstants = 0; 691 p->constants = NULL; 692 p->dimEmbed = -1; 693 p->useJacPre = PETSC_TRUE; 694 p->forceQuad = PETSC_TRUE; 695 PetscCall(PetscWeakFormCreate(comm, &p->wf)); 696 PetscCall(PetscArrayzero(p->quadPerm, DM_NUM_POLYTOPES)); 697 698 *ds = p; 699 PetscFunctionReturn(PETSC_SUCCESS); 700 } 701 702 /*@ 703 PetscDSGetNumFields - Returns the number of fields in the `PetscDS` 704 705 Not Collective 706 707 Input Parameter: 708 . prob - The `PetscDS` object 709 710 Output Parameter: 711 . Nf - The number of fields 712 713 Level: beginner 714 715 .seealso: `PetscDS`, `PetscDSGetSpatialDimension()`, `PetscDSCreate()` 716 @*/ 717 PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf) 718 { 719 PetscFunctionBegin; 720 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 721 PetscAssertPointer(Nf, 2); 722 *Nf = prob->Nf; 723 PetscFunctionReturn(PETSC_SUCCESS); 724 } 725 726 /*@ 727 PetscDSGetSpatialDimension - Returns the spatial dimension of the `PetscDS`, meaning the topological dimension of the discretizations 728 729 Not Collective 730 731 Input Parameter: 732 . prob - The `PetscDS` object 733 734 Output Parameter: 735 . dim - The spatial dimension 736 737 Level: beginner 738 739 .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 740 @*/ 741 PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim) 742 { 743 PetscFunctionBegin; 744 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 745 PetscAssertPointer(dim, 2); 746 *dim = 0; 747 if (prob->Nf) { 748 PetscObject obj; 749 PetscClassId id; 750 751 PetscCall(PetscDSGetDiscretization(prob, 0, &obj)); 752 if (obj) { 753 PetscCall(PetscObjectGetClassId(obj, &id)); 754 if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetSpatialDimension((PetscFE)obj, dim)); 755 else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetSpatialDimension((PetscFV)obj, dim)); 756 else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0); 757 } 758 } 759 PetscFunctionReturn(PETSC_SUCCESS); 760 } 761 762 /*@ 763 PetscDSGetCoordinateDimension - Returns the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded 764 765 Not Collective 766 767 Input Parameter: 768 . prob - The `PetscDS` object 769 770 Output Parameter: 771 . dimEmbed - The coordinate dimension 772 773 Level: beginner 774 775 .seealso: `PetscDS`, `PetscDSSetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 776 @*/ 777 PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed) 778 { 779 PetscFunctionBegin; 780 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 781 PetscAssertPointer(dimEmbed, 2); 782 PetscCheck(prob->dimEmbed >= 0, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS"); 783 *dimEmbed = prob->dimEmbed; 784 PetscFunctionReturn(PETSC_SUCCESS); 785 } 786 787 /*@ 788 PetscDSSetCoordinateDimension - Set the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded 789 790 Logically Collective 791 792 Input Parameters: 793 + prob - The `PetscDS` object 794 - dimEmbed - The coordinate dimension 795 796 Level: beginner 797 798 .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 799 @*/ 800 PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed) 801 { 802 PetscFunctionBegin; 803 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 804 PetscCheck(dimEmbed >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %" PetscInt_FMT, dimEmbed); 805 prob->dimEmbed = dimEmbed; 806 PetscFunctionReturn(PETSC_SUCCESS); 807 } 808 809 /*@ 810 PetscDSGetForceQuad - Returns the flag to force matching quadratures among the field discretizations 811 812 Not collective 813 814 Input Parameter: 815 . ds - The `PetscDS` object 816 817 Output Parameter: 818 . forceQuad - The flag 819 820 Level: intermediate 821 822 .seealso: `PetscDS`, `PetscDSSetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 823 @*/ 824 PetscErrorCode PetscDSGetForceQuad(PetscDS ds, PetscBool *forceQuad) 825 { 826 PetscFunctionBegin; 827 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 828 PetscAssertPointer(forceQuad, 2); 829 *forceQuad = ds->forceQuad; 830 PetscFunctionReturn(PETSC_SUCCESS); 831 } 832 833 /*@ 834 PetscDSSetForceQuad - Set the flag to force matching quadratures among the field discretizations 835 836 Logically collective on ds 837 838 Input Parameters: 839 + ds - The `PetscDS` object 840 - forceQuad - The flag 841 842 Level: intermediate 843 844 .seealso: `PetscDS`, `PetscDSGetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 845 @*/ 846 PetscErrorCode PetscDSSetForceQuad(PetscDS ds, PetscBool forceQuad) 847 { 848 PetscFunctionBegin; 849 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 850 ds->forceQuad = forceQuad; 851 PetscFunctionReturn(PETSC_SUCCESS); 852 } 853 854 /*@ 855 PetscDSIsCohesive - Returns the flag indicating that this `PetscDS` is for a cohesive cell 856 857 Not Collective 858 859 Input Parameter: 860 . ds - The `PetscDS` object 861 862 Output Parameter: 863 . isCohesive - The flag 864 865 Level: developer 866 867 .seealso: `PetscDS`, `PetscDSGetNumCohesive()`, `PetscDSGetCohesive()`, `PetscDSSetCohesive()`, `PetscDSCreate()` 868 @*/ 869 PetscErrorCode PetscDSIsCohesive(PetscDS ds, PetscBool *isCohesive) 870 { 871 PetscFunctionBegin; 872 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 873 PetscAssertPointer(isCohesive, 2); 874 *isCohesive = ds->isCohesive; 875 PetscFunctionReturn(PETSC_SUCCESS); 876 } 877 878 /*@ 879 PetscDSGetNumCohesive - Returns the number of cohesive fields, meaning those defined on the interior of a cohesive cell 880 881 Not Collective 882 883 Input Parameter: 884 . ds - The `PetscDS` object 885 886 Output Parameter: 887 . numCohesive - The number of cohesive fields 888 889 Level: developer 890 891 .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSCreate()` 892 @*/ 893 PetscErrorCode PetscDSGetNumCohesive(PetscDS ds, PetscInt *numCohesive) 894 { 895 PetscInt f; 896 897 PetscFunctionBegin; 898 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 899 PetscAssertPointer(numCohesive, 2); 900 *numCohesive = 0; 901 for (f = 0; f < ds->Nf; ++f) *numCohesive += ds->cohesive[f] ? 1 : 0; 902 PetscFunctionReturn(PETSC_SUCCESS); 903 } 904 905 /*@ 906 PetscDSGetCohesive - Returns the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell 907 908 Not Collective 909 910 Input Parameters: 911 + ds - The `PetscDS` object 912 - f - The field index 913 914 Output Parameter: 915 . isCohesive - The flag 916 917 Level: developer 918 919 .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()` 920 @*/ 921 PetscErrorCode PetscDSGetCohesive(PetscDS ds, PetscInt f, PetscBool *isCohesive) 922 { 923 PetscFunctionBegin; 924 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 925 PetscAssertPointer(isCohesive, 3); 926 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 927 *isCohesive = ds->cohesive[f]; 928 PetscFunctionReturn(PETSC_SUCCESS); 929 } 930 931 /*@ 932 PetscDSSetCohesive - Set the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell 933 934 Not Collective 935 936 Input Parameters: 937 + ds - The `PetscDS` object 938 . f - The field index 939 - isCohesive - The flag for a cohesive field 940 941 Level: developer 942 943 .seealso: `PetscDS`, `PetscDSGetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()` 944 @*/ 945 PetscErrorCode PetscDSSetCohesive(PetscDS ds, PetscInt f, PetscBool isCohesive) 946 { 947 PetscInt i; 948 949 PetscFunctionBegin; 950 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 951 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 952 ds->cohesive[f] = isCohesive; 953 ds->isCohesive = PETSC_FALSE; 954 for (i = 0; i < ds->Nf; ++i) ds->isCohesive = ds->isCohesive || ds->cohesive[f] ? PETSC_TRUE : PETSC_FALSE; 955 PetscFunctionReturn(PETSC_SUCCESS); 956 } 957 958 /*@ 959 PetscDSGetTotalDimension - Returns the total size of the approximation space for this system 960 961 Not Collective 962 963 Input Parameter: 964 . prob - The `PetscDS` object 965 966 Output Parameter: 967 . dim - The total problem dimension 968 969 Level: beginner 970 971 .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()` 972 @*/ 973 PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim) 974 { 975 PetscFunctionBegin; 976 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 977 PetscCall(PetscDSSetUp(prob)); 978 PetscAssertPointer(dim, 2); 979 *dim = prob->totDim; 980 PetscFunctionReturn(PETSC_SUCCESS); 981 } 982 983 /*@ 984 PetscDSGetTotalComponents - Returns the total number of components in this system 985 986 Not Collective 987 988 Input Parameter: 989 . prob - The `PetscDS` object 990 991 Output Parameter: 992 . Nc - The total number of components 993 994 Level: beginner 995 996 .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()` 997 @*/ 998 PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc) 999 { 1000 PetscFunctionBegin; 1001 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1002 PetscCall(PetscDSSetUp(prob)); 1003 PetscAssertPointer(Nc, 2); 1004 *Nc = prob->totComp; 1005 PetscFunctionReturn(PETSC_SUCCESS); 1006 } 1007 1008 /*@ 1009 PetscDSGetDiscretization - Returns the discretization object for the given field 1010 1011 Not Collective 1012 1013 Input Parameters: 1014 + prob - The `PetscDS` object 1015 - f - The field number 1016 1017 Output Parameter: 1018 . disc - The discretization object 1019 1020 Level: beginner 1021 1022 .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1023 @*/ 1024 PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc) 1025 { 1026 PetscFunctionBeginHot; 1027 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1028 PetscAssertPointer(disc, 3); 1029 PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf); 1030 *disc = prob->disc[f]; 1031 PetscFunctionReturn(PETSC_SUCCESS); 1032 } 1033 1034 /*@ 1035 PetscDSSetDiscretization - Sets the discretization object for the given field 1036 1037 Not Collective 1038 1039 Input Parameters: 1040 + prob - The `PetscDS` object 1041 . f - The field number 1042 - disc - The discretization object 1043 1044 Level: beginner 1045 1046 .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1047 @*/ 1048 PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc) 1049 { 1050 PetscFunctionBegin; 1051 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1052 if (disc) PetscAssertPointer(disc, 3); 1053 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 1054 PetscCall(PetscDSEnlarge_Static(prob, f + 1)); 1055 PetscCall(PetscObjectDereference(prob->disc[f])); 1056 prob->disc[f] = disc; 1057 PetscCall(PetscObjectReference(disc)); 1058 if (disc) { 1059 PetscClassId id; 1060 1061 PetscCall(PetscObjectGetClassId(disc, &id)); 1062 if (id == PETSCFE_CLASSID) { 1063 PetscCall(PetscDSSetImplicit(prob, f, PETSC_TRUE)); 1064 } else if (id == PETSCFV_CLASSID) { 1065 PetscCall(PetscDSSetImplicit(prob, f, PETSC_FALSE)); 1066 } 1067 PetscCall(PetscDSSetJetDegree(prob, f, 1)); 1068 } 1069 PetscFunctionReturn(PETSC_SUCCESS); 1070 } 1071 1072 /*@ 1073 PetscDSGetWeakForm - Returns the weak form object 1074 1075 Not Collective 1076 1077 Input Parameter: 1078 . ds - The `PetscDS` object 1079 1080 Output Parameter: 1081 . wf - The weak form object 1082 1083 Level: beginner 1084 1085 .seealso: `PetscWeakForm`, `PetscDSSetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1086 @*/ 1087 PetscErrorCode PetscDSGetWeakForm(PetscDS ds, PetscWeakForm *wf) 1088 { 1089 PetscFunctionBegin; 1090 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1091 PetscAssertPointer(wf, 2); 1092 *wf = ds->wf; 1093 PetscFunctionReturn(PETSC_SUCCESS); 1094 } 1095 1096 /*@ 1097 PetscDSSetWeakForm - Sets the weak form object 1098 1099 Not Collective 1100 1101 Input Parameters: 1102 + ds - The `PetscDS` object 1103 - wf - The weak form object 1104 1105 Level: beginner 1106 1107 .seealso: `PetscWeakForm`, `PetscDSGetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1108 @*/ 1109 PetscErrorCode PetscDSSetWeakForm(PetscDS ds, PetscWeakForm wf) 1110 { 1111 PetscFunctionBegin; 1112 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1113 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 2); 1114 PetscCall(PetscObjectDereference((PetscObject)ds->wf)); 1115 ds->wf = wf; 1116 PetscCall(PetscObjectReference((PetscObject)wf)); 1117 PetscCall(PetscWeakFormSetNumFields(wf, ds->Nf)); 1118 PetscFunctionReturn(PETSC_SUCCESS); 1119 } 1120 1121 /*@ 1122 PetscDSAddDiscretization - Adds a discretization object 1123 1124 Not Collective 1125 1126 Input Parameters: 1127 + prob - The `PetscDS` object 1128 - disc - The boundary discretization object 1129 1130 Level: beginner 1131 1132 .seealso: `PetscWeakForm`, `PetscDSGetDiscretization()`, `PetscDSSetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1133 @*/ 1134 PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc) 1135 { 1136 PetscFunctionBegin; 1137 PetscCall(PetscDSSetDiscretization(prob, prob->Nf, disc)); 1138 PetscFunctionReturn(PETSC_SUCCESS); 1139 } 1140 1141 /*@ 1142 PetscDSGetQuadrature - Returns the quadrature, which must agree for all fields in the `PetscDS` 1143 1144 Not Collective 1145 1146 Input Parameter: 1147 . prob - The `PetscDS` object 1148 1149 Output Parameter: 1150 . q - The quadrature object 1151 1152 Level: intermediate 1153 1154 .seealso: `PetscDS`, `PetscQuadrature`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1155 @*/ 1156 PetscErrorCode PetscDSGetQuadrature(PetscDS prob, PetscQuadrature *q) 1157 { 1158 PetscObject obj; 1159 PetscClassId id; 1160 1161 PetscFunctionBegin; 1162 *q = NULL; 1163 if (!prob->Nf) PetscFunctionReturn(PETSC_SUCCESS); 1164 PetscCall(PetscDSGetDiscretization(prob, 0, &obj)); 1165 PetscCall(PetscObjectGetClassId(obj, &id)); 1166 if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetQuadrature((PetscFE)obj, q)); 1167 else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetQuadrature((PetscFV)obj, q)); 1168 else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0); 1169 PetscFunctionReturn(PETSC_SUCCESS); 1170 } 1171 1172 /*@ 1173 PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for `TSIMEX` 1174 1175 Not Collective 1176 1177 Input Parameters: 1178 + prob - The `PetscDS` object 1179 - f - The field number 1180 1181 Output Parameter: 1182 . implicit - The flag indicating what kind of solve to use for this field 1183 1184 Level: developer 1185 1186 .seealso: `TSIMEX`, `PetscDS`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1187 @*/ 1188 PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit) 1189 { 1190 PetscFunctionBegin; 1191 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1192 PetscAssertPointer(implicit, 3); 1193 PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf); 1194 *implicit = prob->implicit[f]; 1195 PetscFunctionReturn(PETSC_SUCCESS); 1196 } 1197 1198 /*@ 1199 PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for `TSIMEX` 1200 1201 Not Collective 1202 1203 Input Parameters: 1204 + prob - The `PetscDS` object 1205 . f - The field number 1206 - implicit - The flag indicating what kind of solve to use for this field 1207 1208 Level: developer 1209 1210 .seealso: `TSIMEX`, `PetscDSGetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1211 @*/ 1212 PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit) 1213 { 1214 PetscFunctionBegin; 1215 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1216 PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf); 1217 prob->implicit[f] = implicit; 1218 PetscFunctionReturn(PETSC_SUCCESS); 1219 } 1220 1221 /*@ 1222 PetscDSGetJetDegree - Returns the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate. 1223 1224 Not Collective 1225 1226 Input Parameters: 1227 + ds - The `PetscDS` object 1228 - f - The field number 1229 1230 Output Parameter: 1231 . k - The highest derivative we need to tabulate 1232 1233 Level: developer 1234 1235 .seealso: `PetscDS`, `PetscDSSetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1236 @*/ 1237 PetscErrorCode PetscDSGetJetDegree(PetscDS ds, PetscInt f, PetscInt *k) 1238 { 1239 PetscFunctionBegin; 1240 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1241 PetscAssertPointer(k, 3); 1242 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 1243 *k = ds->jetDegree[f]; 1244 PetscFunctionReturn(PETSC_SUCCESS); 1245 } 1246 1247 /*@ 1248 PetscDSSetJetDegree - Set the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate. 1249 1250 Not Collective 1251 1252 Input Parameters: 1253 + ds - The `PetscDS` object 1254 . f - The field number 1255 - k - The highest derivative we need to tabulate 1256 1257 Level: developer 1258 1259 .seealso: ``PetscDS`, `PetscDSGetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 1260 @*/ 1261 PetscErrorCode PetscDSSetJetDegree(PetscDS ds, PetscInt f, PetscInt k) 1262 { 1263 PetscFunctionBegin; 1264 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1265 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 1266 ds->jetDegree[f] = k; 1267 PetscFunctionReturn(PETSC_SUCCESS); 1268 } 1269 1270 /*@C 1271 PetscDSGetObjective - Get the pointwise objective function for a given test field 1272 1273 Not Collective 1274 1275 Input Parameters: 1276 + ds - The `PetscDS` 1277 - f - The test field number 1278 1279 Output Parameter: 1280 . obj - integrand for the test function term 1281 1282 Calling sequence of `obj`: 1283 + dim - the spatial dimension 1284 . Nf - the number of fields 1285 . NfAux - the number of auxiliary fields 1286 . uOff - the offset into u[] and u_t[] for each field 1287 . uOff_x - the offset into u_x[] for each field 1288 . u - each field evaluated at the current point 1289 . u_t - the time derivative of each field evaluated at the current point 1290 . u_x - the gradient of each field evaluated at the current point 1291 . aOff - the offset into a[] and a_t[] for each auxiliary field 1292 . aOff_x - the offset into a_x[] for each auxiliary field 1293 . a - each auxiliary field evaluated at the current point 1294 . a_t - the time derivative of each auxiliary field evaluated at the current point 1295 . a_x - the gradient of auxiliary each field evaluated at the current point 1296 . t - current time 1297 . x - coordinates of the current point 1298 . numConstants - number of constant parameters 1299 . constants - constant parameters 1300 - obj - output values at the current point 1301 1302 Level: intermediate 1303 1304 Note: 1305 We are using a first order FEM model for the weak form: \int_\Omega \phi obj(u, u_t, \nabla u, x, t) 1306 1307 .seealso: `PetscDS`, `PetscDSSetObjective()`, `PetscDSGetResidual()` 1308 @*/ 1309 PetscErrorCode PetscDSGetObjective(PetscDS ds, PetscInt f, void (**obj)(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 obj[])) 1310 { 1311 PetscPointFunc *tmp; 1312 PetscInt n; 1313 1314 PetscFunctionBegin; 1315 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1316 PetscAssertPointer(obj, 3); 1317 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 1318 PetscCall(PetscWeakFormGetObjective(ds->wf, NULL, 0, f, 0, &n, &tmp)); 1319 *obj = tmp ? tmp[0] : NULL; 1320 PetscFunctionReturn(PETSC_SUCCESS); 1321 } 1322 1323 /*@C 1324 PetscDSSetObjective - Set the pointwise objective function for a given test field 1325 1326 Not Collective 1327 1328 Input Parameters: 1329 + ds - The `PetscDS` 1330 . f - The test field number 1331 - obj - integrand for the test function term 1332 1333 Calling sequence of `obj`: 1334 + dim - the spatial dimension 1335 . Nf - the number of fields 1336 . NfAux - the number of auxiliary fields 1337 . uOff - the offset into u[] and u_t[] for each field 1338 . uOff_x - the offset into u_x[] for each field 1339 . u - each field evaluated at the current point 1340 . u_t - the time derivative of each field evaluated at the current point 1341 . u_x - the gradient of each field evaluated at the current point 1342 . aOff - the offset into a[] and a_t[] for each auxiliary field 1343 . aOff_x - the offset into a_x[] for each auxiliary field 1344 . a - each auxiliary field evaluated at the current point 1345 . a_t - the time derivative of each auxiliary field evaluated at the current point 1346 . a_x - the gradient of auxiliary each field evaluated at the current point 1347 . t - current time 1348 . x - coordinates of the current point 1349 . numConstants - number of constant parameters 1350 . constants - constant parameters 1351 - obj - output values at the current point 1352 1353 Level: intermediate 1354 1355 Note: 1356 We are using a first order FEM model for the weak form: \int_\Omega \phi obj(u, u_t, \nabla u, x, t) 1357 1358 .seealso: `PetscDS`, `PetscDSGetObjective()`, `PetscDSSetResidual()` 1359 @*/ 1360 PetscErrorCode PetscDSSetObjective(PetscDS ds, PetscInt f, void (*obj)(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 obj[])) 1361 { 1362 PetscFunctionBegin; 1363 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1364 if (obj) PetscValidFunction(obj, 3); 1365 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 1366 PetscCall(PetscWeakFormSetIndexObjective(ds->wf, NULL, 0, f, 0, 0, obj)); 1367 PetscFunctionReturn(PETSC_SUCCESS); 1368 } 1369 1370 /*@C 1371 PetscDSGetResidual - Get the pointwise residual function for a given test field 1372 1373 Not Collective 1374 1375 Input Parameters: 1376 + ds - The `PetscDS` 1377 - f - The test field number 1378 1379 Output Parameters: 1380 + f0 - integrand for the test function term 1381 - f1 - integrand for the test function gradient term 1382 1383 Calling sequence of `f0`: 1384 + dim - the spatial dimension 1385 . Nf - the number of fields 1386 . NfAux - the number of auxiliary fields 1387 . uOff - the offset into u[] and u_t[] for each field 1388 . uOff_x - the offset into u_x[] for each field 1389 . u - each field evaluated at the current point 1390 . u_t - the time derivative of each field evaluated at the current point 1391 . u_x - the gradient of each field evaluated at the current point 1392 . aOff - the offset into a[] and a_t[] for each auxiliary field 1393 . aOff_x - the offset into a_x[] for each auxiliary field 1394 . a - each auxiliary field evaluated at the current point 1395 . a_t - the time derivative of each auxiliary field evaluated at the current point 1396 . a_x - the gradient of auxiliary each field evaluated at the current point 1397 . t - current time 1398 . x - coordinates of the current point 1399 . numConstants - number of constant parameters 1400 . constants - constant parameters 1401 - f0 - output values at the current point 1402 1403 Level: intermediate 1404 1405 Note: 1406 `f1` has an identical form and is omitted for brevity. 1407 1408 We are using a first order FEM model for the weak form: \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t) 1409 1410 .seealso: `PetscDS`, `PetscDSSetResidual()` 1411 @*/ 1412 PetscErrorCode PetscDSGetResidual(PetscDS ds, PetscInt f, void (**f0)(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 f0[]), void (**f1)(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[])) 1413 { 1414 PetscPointFunc *tmp0, *tmp1; 1415 PetscInt n0, n1; 1416 1417 PetscFunctionBegin; 1418 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1419 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 1420 PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1)); 1421 *f0 = tmp0 ? tmp0[0] : NULL; 1422 *f1 = tmp1 ? tmp1[0] : NULL; 1423 PetscFunctionReturn(PETSC_SUCCESS); 1424 } 1425 1426 /*@C 1427 PetscDSSetResidual - Set the pointwise residual function for a given test field 1428 1429 Not Collective 1430 1431 Input Parameters: 1432 + ds - The `PetscDS` 1433 . f - The test field number 1434 . f0 - integrand for the test function term 1435 - f1 - integrand for the test function gradient term 1436 1437 Calling sequence of `f0`: 1438 + dim - the spatial dimension 1439 . Nf - the number of fields 1440 . NfAux - the number of auxiliary fields 1441 . uOff - the offset into u[] and u_t[] for each field 1442 . uOff_x - the offset into u_x[] for each field 1443 . u - each field evaluated at the current point 1444 . u_t - the time derivative of each field evaluated at the current point 1445 . u_x - the gradient of each field evaluated at the current point 1446 . aOff - the offset into a[] and a_t[] for each auxiliary field 1447 . aOff_x - the offset into a_x[] for each auxiliary field 1448 . a - each auxiliary field evaluated at the current point 1449 . a_t - the time derivative of each auxiliary field evaluated at the current point 1450 . a_x - the gradient of auxiliary each field evaluated at the current point 1451 . t - current time 1452 . x - coordinates of the current point 1453 . numConstants - number of constant parameters 1454 . constants - constant parameters 1455 - f0 - output values at the current point 1456 1457 Level: intermediate 1458 1459 Note: 1460 `f1` has an identical form and is omitted for brevity. 1461 1462 We are using a first order FEM model for the weak form: \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t) 1463 1464 .seealso: `PetscDS`, `PetscDSGetResidual()` 1465 @*/ 1466 PetscErrorCode PetscDSSetResidual(PetscDS ds, PetscInt f, void (*f0)(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 f0[]), void (*f1)(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[])) 1467 { 1468 PetscFunctionBegin; 1469 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1470 if (f0) PetscValidFunction(f0, 3); 1471 if (f1) PetscValidFunction(f1, 4); 1472 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 1473 PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1)); 1474 PetscFunctionReturn(PETSC_SUCCESS); 1475 } 1476 1477 /*@C 1478 PetscDSGetRHSResidual - Get the pointwise RHS residual function for explicit timestepping for a given test field 1479 1480 Not Collective 1481 1482 Input Parameters: 1483 + ds - The `PetscDS` 1484 - f - The test field number 1485 1486 Output Parameters: 1487 + f0 - integrand for the test function term 1488 - f1 - integrand for the test function gradient term 1489 1490 Calling sequence of `f0`: 1491 + dim - the spatial dimension 1492 . Nf - the number of fields 1493 . NfAux - the number of auxiliary fields 1494 . uOff - the offset into u[] and u_t[] for each field 1495 . uOff_x - the offset into u_x[] for each field 1496 . u - each field evaluated at the current point 1497 . u_t - the time derivative of each field evaluated at the current point 1498 . u_x - the gradient of each field evaluated at the current point 1499 . aOff - the offset into a[] and a_t[] for each auxiliary field 1500 . aOff_x - the offset into a_x[] for each auxiliary field 1501 . a - each auxiliary field evaluated at the current point 1502 . a_t - the time derivative of each auxiliary field evaluated at the current point 1503 . a_x - the gradient of auxiliary each field evaluated at the current point 1504 . t - current time 1505 . x - coordinates of the current point 1506 . numConstants - number of constant parameters 1507 . constants - constant parameters 1508 - f0 - output values at the current point 1509 1510 Level: intermediate 1511 1512 Note: 1513 `f1` has an identical form and is omitted for brevity. 1514 1515 We are using a first order FEM model for the weak form: \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t) 1516 1517 .seealso: `PetscDS`, `PetscDSSetRHSResidual()` 1518 @*/ 1519 PetscErrorCode PetscDSGetRHSResidual(PetscDS ds, PetscInt f, void (**f0)(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 f0[]), void (**f1)(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[])) 1520 { 1521 PetscPointFunc *tmp0, *tmp1; 1522 PetscInt n0, n1; 1523 1524 PetscFunctionBegin; 1525 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1526 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 1527 PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 100, &n0, &tmp0, &n1, &tmp1)); 1528 *f0 = tmp0 ? tmp0[0] : NULL; 1529 *f1 = tmp1 ? tmp1[0] : NULL; 1530 PetscFunctionReturn(PETSC_SUCCESS); 1531 } 1532 1533 /*@C 1534 PetscDSSetRHSResidual - Set the pointwise residual function for explicit timestepping for a given test field 1535 1536 Not Collective 1537 1538 Input Parameters: 1539 + ds - The `PetscDS` 1540 . f - The test field number 1541 . f0 - integrand for the test function term 1542 - f1 - integrand for the test function gradient term 1543 1544 Calling sequence for the callbacks `f0`: 1545 + dim - the spatial dimension 1546 . Nf - the number of fields 1547 . NfAux - the number of auxiliary fields 1548 . uOff - the offset into u[] and u_t[] for each field 1549 . uOff_x - the offset into u_x[] for each field 1550 . u - each field evaluated at the current point 1551 . u_t - the time derivative of each field evaluated at the current point 1552 . u_x - the gradient of each field evaluated at the current point 1553 . aOff - the offset into a[] and a_t[] for each auxiliary field 1554 . aOff_x - the offset into a_x[] for each auxiliary field 1555 . a - each auxiliary field evaluated at the current point 1556 . a_t - the time derivative of each auxiliary field evaluated at the current point 1557 . a_x - the gradient of auxiliary each field evaluated at the current point 1558 . t - current time 1559 . x - coordinates of the current point 1560 . numConstants - number of constant parameters 1561 . constants - constant parameters 1562 - f0 - output values at the current point 1563 1564 Level: intermediate 1565 1566 Note: 1567 `f1` has an identical form and is omitted for brevity. 1568 1569 We are using a first order FEM model for the weak form: \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t) 1570 1571 .seealso: `PetscDS`, `PetscDSGetResidual()` 1572 @*/ 1573 PetscErrorCode PetscDSSetRHSResidual(PetscDS ds, PetscInt f, void (*f0)(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 f0[]), void (*f1)(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[])) 1574 { 1575 PetscFunctionBegin; 1576 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1577 if (f0) PetscValidFunction(f0, 3); 1578 if (f1) PetscValidFunction(f1, 4); 1579 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 1580 PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 100, 0, f0, 0, f1)); 1581 PetscFunctionReturn(PETSC_SUCCESS); 1582 } 1583 1584 /*@C 1585 PetscDSHasJacobian - Checks that the Jacobian functions have been set 1586 1587 Not Collective 1588 1589 Input Parameter: 1590 . ds - The `PetscDS` 1591 1592 Output Parameter: 1593 . hasJac - flag that pointwise function for the Jacobian has been set 1594 1595 Level: intermediate 1596 1597 .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()` 1598 @*/ 1599 PetscErrorCode PetscDSHasJacobian(PetscDS ds, PetscBool *hasJac) 1600 { 1601 PetscFunctionBegin; 1602 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1603 PetscCall(PetscWeakFormHasJacobian(ds->wf, hasJac)); 1604 PetscFunctionReturn(PETSC_SUCCESS); 1605 } 1606 1607 /*@C 1608 PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field 1609 1610 Not Collective 1611 1612 Input Parameters: 1613 + ds - The `PetscDS` 1614 . f - The test field number 1615 - g - The field number 1616 1617 Output Parameters: 1618 + g0 - integrand for the test and basis function term 1619 . g1 - integrand for the test function and basis function gradient term 1620 . g2 - integrand for the test function gradient and basis function term 1621 - g3 - integrand for the test function gradient and basis function gradient term 1622 1623 Calling sequence of `g0`: 1624 + dim - the spatial dimension 1625 . Nf - the number of fields 1626 . NfAux - the number of auxiliary fields 1627 . uOff - the offset into u[] and u_t[] for each field 1628 . uOff_x - the offset into u_x[] for each field 1629 . u - each field evaluated at the current point 1630 . u_t - the time derivative of each field evaluated at the current point 1631 . u_x - the gradient of each field evaluated at the current point 1632 . aOff - the offset into a[] and a_t[] for each auxiliary field 1633 . aOff_x - the offset into a_x[] for each auxiliary field 1634 . a - each auxiliary field evaluated at the current point 1635 . a_t - the time derivative of each auxiliary field evaluated at the current point 1636 . a_x - the gradient of auxiliary each field evaluated at the current point 1637 . t - current time 1638 . u_tShift - the multiplier a for dF/dU_t 1639 . x - coordinates of the current point 1640 . numConstants - number of constant parameters 1641 . constants - constant parameters 1642 - g0 - output values at the current point 1643 1644 Level: intermediate 1645 1646 Note: 1647 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 1648 1649 We are using a first order FEM model for the weak form\: 1650 \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi 1651 1652 .seealso: `PetscDS`, `PetscDSSetJacobian()` 1653 @*/ 1654 PetscErrorCode PetscDSGetJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(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, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1655 { 1656 PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3; 1657 PetscInt n0, n1, n2, n3; 1658 1659 PetscFunctionBegin; 1660 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1661 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 1662 PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf); 1663 PetscCall(PetscWeakFormGetJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3)); 1664 *g0 = tmp0 ? tmp0[0] : NULL; 1665 *g1 = tmp1 ? tmp1[0] : NULL; 1666 *g2 = tmp2 ? tmp2[0] : NULL; 1667 *g3 = tmp3 ? tmp3[0] : NULL; 1668 PetscFunctionReturn(PETSC_SUCCESS); 1669 } 1670 1671 /*@C 1672 PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields 1673 1674 Not Collective 1675 1676 Input Parameters: 1677 + ds - The `PetscDS` 1678 . f - The test field number 1679 . g - The field number 1680 . g0 - integrand for the test and basis function term 1681 . g1 - integrand for the test function and basis function gradient term 1682 . g2 - integrand for the test function gradient and basis function term 1683 - g3 - integrand for the test function gradient and basis function gradient term 1684 1685 Calling sequence of `g0`: 1686 + dim - the spatial dimension 1687 . Nf - the number of fields 1688 . NfAux - the number of auxiliary fields 1689 . uOff - the offset into u[] and u_t[] for each field 1690 . uOff_x - the offset into u_x[] for each field 1691 . u - each field evaluated at the current point 1692 . u_t - the time derivative of each field evaluated at the current point 1693 . u_x - the gradient of each field evaluated at the current point 1694 . aOff - the offset into a[] and a_t[] for each auxiliary field 1695 . aOff_x - the offset into a_x[] for each auxiliary field 1696 . a - each auxiliary field evaluated at the current point 1697 . a_t - the time derivative of each auxiliary field evaluated at the current point 1698 . a_x - the gradient of auxiliary each field evaluated at the current point 1699 . t - current time 1700 . u_tShift - the multiplier a for dF/dU_t 1701 . x - coordinates of the current point 1702 . numConstants - number of constant parameters 1703 . constants - constant parameters 1704 - g0 - output values at the current point 1705 1706 Level: intermediate 1707 1708 Note: 1709 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 1710 1711 We are using a first order FEM model for the weak form\: 1712 \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi 1713 1714 .seealso: `PetscDS`, `PetscDSGetJacobian()` 1715 @*/ 1716 PetscErrorCode PetscDSSetJacobian(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(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, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1717 { 1718 PetscFunctionBegin; 1719 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1720 if (g0) PetscValidFunction(g0, 4); 1721 if (g1) PetscValidFunction(g1, 5); 1722 if (g2) PetscValidFunction(g2, 6); 1723 if (g3) PetscValidFunction(g3, 7); 1724 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 1725 PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g); 1726 PetscCall(PetscWeakFormSetIndexJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3)); 1727 PetscFunctionReturn(PETSC_SUCCESS); 1728 } 1729 1730 /*@C 1731 PetscDSUseJacobianPreconditioner - Set whether to construct a Jacobian preconditioner 1732 1733 Not Collective 1734 1735 Input Parameters: 1736 + prob - The `PetscDS` 1737 - useJacPre - flag that enables construction of a Jacobian preconditioner 1738 1739 Level: intermediate 1740 1741 .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()` 1742 @*/ 1743 PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre) 1744 { 1745 PetscFunctionBegin; 1746 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1747 prob->useJacPre = useJacPre; 1748 PetscFunctionReturn(PETSC_SUCCESS); 1749 } 1750 1751 /*@C 1752 PetscDSHasJacobianPreconditioner - Checks if a Jacobian preconditioner matrix has been set 1753 1754 Not Collective 1755 1756 Input Parameter: 1757 . ds - The `PetscDS` 1758 1759 Output Parameter: 1760 . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set 1761 1762 Level: intermediate 1763 1764 .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()` 1765 @*/ 1766 PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS ds, PetscBool *hasJacPre) 1767 { 1768 PetscFunctionBegin; 1769 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1770 *hasJacPre = PETSC_FALSE; 1771 if (!ds->useJacPre) PetscFunctionReturn(PETSC_SUCCESS); 1772 PetscCall(PetscWeakFormHasJacobianPreconditioner(ds->wf, hasJacPre)); 1773 PetscFunctionReturn(PETSC_SUCCESS); 1774 } 1775 1776 /*@C 1777 PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian preconditioner function for given test and basis field. If this is missing, 1778 the system matrix is used to build the preconditioner. 1779 1780 Not Collective 1781 1782 Input Parameters: 1783 + ds - The `PetscDS` 1784 . f - The test field number 1785 - g - The field number 1786 1787 Output Parameters: 1788 + g0 - integrand for the test and basis function term 1789 . g1 - integrand for the test function and basis function gradient term 1790 . g2 - integrand for the test function gradient and basis function term 1791 - g3 - integrand for the test function gradient and basis function gradient term 1792 1793 Calling sequence of `g0`: 1794 + dim - the spatial dimension 1795 . Nf - the number of fields 1796 . NfAux - the number of auxiliary fields 1797 . uOff - the offset into u[] and u_t[] for each field 1798 . uOff_x - the offset into u_x[] for each field 1799 . u - each field evaluated at the current point 1800 . u_t - the time derivative of each field evaluated at the current point 1801 . u_x - the gradient of each field evaluated at the current point 1802 . aOff - the offset into a[] and a_t[] for each auxiliary field 1803 . aOff_x - the offset into a_x[] for each auxiliary field 1804 . a - each auxiliary field evaluated at the current point 1805 . a_t - the time derivative of each auxiliary field evaluated at the current point 1806 . a_x - the gradient of auxiliary each field evaluated at the current point 1807 . t - current time 1808 . u_tShift - the multiplier a for dF/dU_t 1809 . x - coordinates of the current point 1810 . numConstants - number of constant parameters 1811 . constants - constant parameters 1812 - g0 - output values at the current point 1813 1814 Level: intermediate 1815 1816 Note: 1817 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 1818 We are using a first order FEM model for the weak form\: 1819 \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi 1820 1821 .seealso: `PetscDS`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()` 1822 @*/ 1823 PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(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, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1824 { 1825 PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3; 1826 PetscInt n0, n1, n2, n3; 1827 1828 PetscFunctionBegin; 1829 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1830 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 1831 PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf); 1832 PetscCall(PetscWeakFormGetJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3)); 1833 *g0 = tmp0 ? tmp0[0] : NULL; 1834 *g1 = tmp1 ? tmp1[0] : NULL; 1835 *g2 = tmp2 ? tmp2[0] : NULL; 1836 *g3 = tmp3 ? tmp3[0] : NULL; 1837 PetscFunctionReturn(PETSC_SUCCESS); 1838 } 1839 1840 /*@C 1841 PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian preconditioner function for given test and basis fields. 1842 If this is missing, the system matrix is used to build the preconditioner. 1843 1844 Not Collective 1845 1846 Input Parameters: 1847 + ds - The `PetscDS` 1848 . f - The test field number 1849 . g - The field number 1850 . g0 - integrand for the test and basis function term 1851 . g1 - integrand for the test function and basis function gradient term 1852 . g2 - integrand for the test function gradient and basis function term 1853 - g3 - integrand for the test function gradient and basis function gradient term 1854 1855 Calling sequence of `g0`: 1856 + dim - the spatial dimension 1857 . Nf - the number of fields 1858 . NfAux - the number of auxiliary fields 1859 . uOff - the offset into u[] and u_t[] for each field 1860 . uOff_x - the offset into u_x[] for each field 1861 . u - each field evaluated at the current point 1862 . u_t - the time derivative of each field evaluated at the current point 1863 . u_x - the gradient of each field evaluated at the current point 1864 . aOff - the offset into a[] and a_t[] for each auxiliary field 1865 . aOff_x - the offset into a_x[] for each auxiliary field 1866 . a - each auxiliary field evaluated at the current point 1867 . a_t - the time derivative of each auxiliary field evaluated at the current point 1868 . a_x - the gradient of auxiliary each field evaluated at the current point 1869 . t - current time 1870 . u_tShift - the multiplier a for dF/dU_t 1871 . x - coordinates of the current point 1872 . numConstants - number of constant parameters 1873 . constants - constant parameters 1874 - g0 - output values at the current point 1875 1876 Level: intermediate 1877 1878 Note: 1879 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 1880 1881 We are using a first order FEM model for the weak form\: 1882 \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi 1883 1884 .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobian()` 1885 @*/ 1886 PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(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, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1887 { 1888 PetscFunctionBegin; 1889 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1890 if (g0) PetscValidFunction(g0, 4); 1891 if (g1) PetscValidFunction(g1, 5); 1892 if (g2) PetscValidFunction(g2, 6); 1893 if (g3) PetscValidFunction(g3, 7); 1894 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 1895 PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g); 1896 PetscCall(PetscWeakFormSetIndexJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3)); 1897 PetscFunctionReturn(PETSC_SUCCESS); 1898 } 1899 1900 /*@C 1901 PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set 1902 1903 Not Collective 1904 1905 Input Parameter: 1906 . ds - The `PetscDS` 1907 1908 Output Parameter: 1909 . hasDynJac - flag that pointwise function for dynamic Jacobian has been set 1910 1911 Level: intermediate 1912 1913 .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscDSGetJacobian()` 1914 @*/ 1915 PetscErrorCode PetscDSHasDynamicJacobian(PetscDS ds, PetscBool *hasDynJac) 1916 { 1917 PetscFunctionBegin; 1918 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1919 PetscCall(PetscWeakFormHasDynamicJacobian(ds->wf, hasDynJac)); 1920 PetscFunctionReturn(PETSC_SUCCESS); 1921 } 1922 1923 /*@C 1924 PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field 1925 1926 Not Collective 1927 1928 Input Parameters: 1929 + ds - The `PetscDS` 1930 . f - The test field number 1931 - g - The field number 1932 1933 Output Parameters: 1934 + g0 - integrand for the test and basis function term 1935 . g1 - integrand for the test function and basis function gradient term 1936 . g2 - integrand for the test function gradient and basis function term 1937 - g3 - integrand for the test function gradient and basis function gradient term 1938 1939 Calling sequence of `g0`: 1940 + dim - the spatial dimension 1941 . Nf - the number of fields 1942 . NfAux - the number of auxiliary fields 1943 . uOff - the offset into u[] and u_t[] for each field 1944 . uOff_x - the offset into u_x[] for each field 1945 . u - each field evaluated at the current point 1946 . u_t - the time derivative of each field evaluated at the current point 1947 . u_x - the gradient of each field evaluated at the current point 1948 . aOff - the offset into a[] and a_t[] for each auxiliary field 1949 . aOff_x - the offset into a_x[] for each auxiliary field 1950 . a - each auxiliary field evaluated at the current point 1951 . a_t - the time derivative of each auxiliary field evaluated at the current point 1952 . a_x - the gradient of auxiliary each field evaluated at the current point 1953 . t - current time 1954 . u_tShift - the multiplier a for dF/dU_t 1955 . x - coordinates of the current point 1956 . numConstants - number of constant parameters 1957 . constants - constant parameters 1958 - g0 - output values at the current point 1959 1960 Level: intermediate 1961 1962 Note: 1963 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 1964 1965 We are using a first order FEM model for the weak form\: 1966 \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi 1967 1968 .seealso: `PetscDS`, `PetscDSSetJacobian()` 1969 @*/ 1970 PetscErrorCode PetscDSGetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(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, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1971 { 1972 PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3; 1973 PetscInt n0, n1, n2, n3; 1974 1975 PetscFunctionBegin; 1976 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 1977 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 1978 PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf); 1979 PetscCall(PetscWeakFormGetDynamicJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3)); 1980 *g0 = tmp0 ? tmp0[0] : NULL; 1981 *g1 = tmp1 ? tmp1[0] : NULL; 1982 *g2 = tmp2 ? tmp2[0] : NULL; 1983 *g3 = tmp3 ? tmp3[0] : NULL; 1984 PetscFunctionReturn(PETSC_SUCCESS); 1985 } 1986 1987 /*@C 1988 PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields 1989 1990 Not Collective 1991 1992 Input Parameters: 1993 + ds - The `PetscDS` 1994 . f - The test field number 1995 . g - The field number 1996 . g0 - integrand for the test and basis function term 1997 . g1 - integrand for the test function and basis function gradient term 1998 . g2 - integrand for the test function gradient and basis function term 1999 - g3 - integrand for the test function gradient and basis function gradient term 2000 2001 Calling sequence of `g0`: 2002 + dim - the spatial dimension 2003 . Nf - the number of fields 2004 . NfAux - the number of auxiliary fields 2005 . uOff - the offset into u[] and u_t[] for each field 2006 . uOff_x - the offset into u_x[] for each field 2007 . u - each field evaluated at the current point 2008 . u_t - the time derivative of each field evaluated at the current point 2009 . u_x - the gradient of each field evaluated at the current point 2010 . aOff - the offset into a[] and a_t[] for each auxiliary field 2011 . aOff_x - the offset into a_x[] for each auxiliary field 2012 . a - each auxiliary field evaluated at the current point 2013 . a_t - the time derivative of each auxiliary field evaluated at the current point 2014 . a_x - the gradient of auxiliary each field evaluated at the current point 2015 . t - current time 2016 . u_tShift - the multiplier a for dF/dU_t 2017 . x - coordinates of the current point 2018 . numConstants - number of constant parameters 2019 . constants - constant parameters 2020 - g0 - output values at the current point 2021 2022 Level: intermediate 2023 2024 Note: 2025 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 2026 2027 We are using a first order FEM model for the weak form\: 2028 \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi 2029 2030 .seealso: `PetscDS`, `PetscDSGetJacobian()` 2031 @*/ 2032 PetscErrorCode PetscDSSetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(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, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 2033 { 2034 PetscFunctionBegin; 2035 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2036 if (g0) PetscValidFunction(g0, 4); 2037 if (g1) PetscValidFunction(g1, 5); 2038 if (g2) PetscValidFunction(g2, 6); 2039 if (g3) PetscValidFunction(g3, 7); 2040 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 2041 PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g); 2042 PetscCall(PetscWeakFormSetIndexDynamicJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3)); 2043 PetscFunctionReturn(PETSC_SUCCESS); 2044 } 2045 2046 /*@C 2047 PetscDSGetRiemannSolver - Returns the Riemann solver for the given field 2048 2049 Not Collective 2050 2051 Input Parameters: 2052 + ds - The `PetscDS` object 2053 - f - The field number 2054 2055 Output Parameter: 2056 . r - Riemann solver 2057 2058 Calling sequence of `r`: 2059 + dim - The spatial dimension 2060 . Nf - The number of fields 2061 . x - The coordinates at a point on the interface 2062 . n - The normal vector to the interface 2063 . uL - The state vector to the left of the interface 2064 . uR - The state vector to the right of the interface 2065 . flux - output array of flux through the interface 2066 . numConstants - number of constant parameters 2067 . constants - constant parameters 2068 - ctx - optional user context 2069 2070 Level: intermediate 2071 2072 .seealso: `PetscDS`, `PetscDSSetRiemannSolver()` 2073 @*/ 2074 PetscErrorCode PetscDSGetRiemannSolver(PetscDS ds, PetscInt f, void (**r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx)) 2075 { 2076 PetscRiemannFunc *tmp; 2077 PetscInt n; 2078 2079 PetscFunctionBegin; 2080 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2081 PetscAssertPointer(r, 3); 2082 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 2083 PetscCall(PetscWeakFormGetRiemannSolver(ds->wf, NULL, 0, f, 0, &n, &tmp)); 2084 *r = tmp ? tmp[0] : NULL; 2085 PetscFunctionReturn(PETSC_SUCCESS); 2086 } 2087 2088 /*@C 2089 PetscDSSetRiemannSolver - Sets the Riemann solver for the given field 2090 2091 Not Collective 2092 2093 Input Parameters: 2094 + ds - The `PetscDS` object 2095 . f - The field number 2096 - r - Riemann solver 2097 2098 Calling sequence of `r`: 2099 + dim - The spatial dimension 2100 . Nf - The number of fields 2101 . x - The coordinates at a point on the interface 2102 . n - The normal vector to the interface 2103 . uL - The state vector to the left of the interface 2104 . uR - The state vector to the right of the interface 2105 . flux - output array of flux through the interface 2106 . numConstants - number of constant parameters 2107 . constants - constant parameters 2108 - ctx - optional user context 2109 2110 Level: intermediate 2111 2112 .seealso: `PetscDS`, `PetscDSGetRiemannSolver()` 2113 @*/ 2114 PetscErrorCode PetscDSSetRiemannSolver(PetscDS ds, PetscInt f, void (*r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx)) 2115 { 2116 PetscFunctionBegin; 2117 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2118 if (r) PetscValidFunction(r, 3); 2119 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 2120 PetscCall(PetscWeakFormSetIndexRiemannSolver(ds->wf, NULL, 0, f, 0, 0, r)); 2121 PetscFunctionReturn(PETSC_SUCCESS); 2122 } 2123 2124 /*@C 2125 PetscDSGetUpdate - Get the pointwise update function for a given field 2126 2127 Not Collective 2128 2129 Input Parameters: 2130 + ds - The `PetscDS` 2131 - f - The field number 2132 2133 Output Parameter: 2134 . update - update function 2135 2136 Calling sequence of `update`: 2137 + dim - the spatial dimension 2138 . Nf - the number of fields 2139 . NfAux - the number of auxiliary fields 2140 . uOff - the offset into u[] and u_t[] for each field 2141 . uOff_x - the offset into u_x[] for each field 2142 . u - each field evaluated at the current point 2143 . u_t - the time derivative of each field evaluated at the current point 2144 . u_x - the gradient of each field evaluated at the current point 2145 . aOff - the offset into a[] and a_t[] for each auxiliary field 2146 . aOff_x - the offset into a_x[] for each auxiliary field 2147 . a - each auxiliary field evaluated at the current point 2148 . a_t - the time derivative of each auxiliary field evaluated at the current point 2149 . a_x - the gradient of auxiliary each field evaluated at the current point 2150 . t - current time 2151 . x - coordinates of the current point 2152 . numConstants - number of constant parameters 2153 . constants - constant parameters 2154 - uNew - new value for field at the current point 2155 2156 Level: intermediate 2157 2158 .seealso: `PetscDS`, `PetscDSSetUpdate()`, `PetscDSSetResidual()` 2159 @*/ 2160 PetscErrorCode PetscDSGetUpdate(PetscDS ds, PetscInt f, void (**update)(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 uNew[])) 2161 { 2162 PetscFunctionBegin; 2163 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2164 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 2165 if (update) { 2166 PetscAssertPointer(update, 3); 2167 *update = ds->update[f]; 2168 } 2169 PetscFunctionReturn(PETSC_SUCCESS); 2170 } 2171 2172 /*@C 2173 PetscDSSetUpdate - Set the pointwise update function for a given field 2174 2175 Not Collective 2176 2177 Input Parameters: 2178 + ds - The `PetscDS` 2179 . f - The field number 2180 - update - update function 2181 2182 Calling sequence of `update`: 2183 + dim - the spatial dimension 2184 . Nf - the number of fields 2185 . NfAux - the number of auxiliary fields 2186 . uOff - the offset into u[] and u_t[] for each field 2187 . uOff_x - the offset into u_x[] for each field 2188 . u - each field evaluated at the current point 2189 . u_t - the time derivative of each field evaluated at the current point 2190 . u_x - the gradient of each field evaluated at the current point 2191 . aOff - the offset into a[] and a_t[] for each auxiliary field 2192 . aOff_x - the offset into a_x[] for each auxiliary field 2193 . a - each auxiliary field evaluated at the current point 2194 . a_t - the time derivative of each auxiliary field evaluated at the current point 2195 . a_x - the gradient of auxiliary each field evaluated at the current point 2196 . t - current time 2197 . x - coordinates of the current point 2198 . numConstants - number of constant parameters 2199 . constants - constant parameters 2200 - uNew - new field values at the current point 2201 2202 Level: intermediate 2203 2204 .seealso: `PetscDS`, `PetscDSGetResidual()` 2205 @*/ 2206 PetscErrorCode PetscDSSetUpdate(PetscDS ds, PetscInt f, void (*update)(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 uNew[])) 2207 { 2208 PetscFunctionBegin; 2209 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2210 if (update) PetscValidFunction(update, 3); 2211 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 2212 PetscCall(PetscDSEnlarge_Static(ds, f + 1)); 2213 ds->update[f] = update; 2214 PetscFunctionReturn(PETSC_SUCCESS); 2215 } 2216 2217 PetscErrorCode PetscDSGetContext(PetscDS ds, PetscInt f, void *ctx) 2218 { 2219 PetscFunctionBegin; 2220 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2221 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 2222 PetscAssertPointer(ctx, 3); 2223 *(void **)ctx = ds->ctx[f]; 2224 PetscFunctionReturn(PETSC_SUCCESS); 2225 } 2226 2227 PetscErrorCode PetscDSSetContext(PetscDS ds, PetscInt f, void *ctx) 2228 { 2229 PetscFunctionBegin; 2230 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2231 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 2232 PetscCall(PetscDSEnlarge_Static(ds, f + 1)); 2233 ds->ctx[f] = ctx; 2234 PetscFunctionReturn(PETSC_SUCCESS); 2235 } 2236 2237 /*@C 2238 PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field 2239 2240 Not Collective 2241 2242 Input Parameters: 2243 + ds - The PetscDS 2244 - f - The test field number 2245 2246 Output Parameters: 2247 + f0 - boundary integrand for the test function term 2248 - f1 - boundary integrand for the test function gradient term 2249 2250 Calling sequence of `f0`: 2251 + dim - the spatial dimension 2252 . Nf - the number of fields 2253 . NfAux - the number of auxiliary fields 2254 . uOff - the offset into u[] and u_t[] for each field 2255 . uOff_x - the offset into u_x[] for each field 2256 . u - each field evaluated at the current point 2257 . u_t - the time derivative of each field evaluated at the current point 2258 . u_x - the gradient of each field evaluated at the current point 2259 . aOff - the offset into a[] and a_t[] for each auxiliary field 2260 . aOff_x - the offset into a_x[] for each auxiliary field 2261 . a - each auxiliary field evaluated at the current point 2262 . a_t - the time derivative of each auxiliary field evaluated at the current point 2263 . a_x - the gradient of auxiliary each field evaluated at the current point 2264 . t - current time 2265 . x - coordinates of the current point 2266 . n - unit normal at the current point 2267 . numConstants - number of constant parameters 2268 . constants - constant parameters 2269 - f0 - output values at the current point 2270 2271 Level: intermediate 2272 2273 Note: 2274 The calling sequence of `f1` is identical, and therefore omitted for brevity. 2275 2276 We are using a first order FEM model for the weak form\: 2277 \int_\Gamma \phi {\vec f}_0(u, u_t, \nabla u, x, t) \cdot \hat n + \nabla\phi \cdot {\overleftrightarrow f}_1(u, u_t, \nabla u, x, t) \cdot \hat n 2278 2279 .seealso: `PetscDS`, `PetscDSSetBdResidual()` 2280 @*/ 2281 PetscErrorCode PetscDSGetBdResidual(PetscDS ds, PetscInt f, void (**f0)(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[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (**f1)(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[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 2282 { 2283 PetscBdPointFunc *tmp0, *tmp1; 2284 PetscInt n0, n1; 2285 2286 PetscFunctionBegin; 2287 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2288 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 2289 PetscCall(PetscWeakFormGetBdResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1)); 2290 *f0 = tmp0 ? tmp0[0] : NULL; 2291 *f1 = tmp1 ? tmp1[0] : NULL; 2292 PetscFunctionReturn(PETSC_SUCCESS); 2293 } 2294 2295 /*@C 2296 PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field 2297 2298 Not Collective 2299 2300 Input Parameters: 2301 + ds - The `PetscDS` 2302 . f - The test field number 2303 . f0 - boundary integrand for the test function term 2304 - f1 - boundary integrand for the test function gradient term 2305 2306 Calling sequence of `f0`: 2307 + dim - the spatial dimension 2308 . Nf - the number of fields 2309 . NfAux - the number of auxiliary fields 2310 . uOff - the offset into u[] and u_t[] for each field 2311 . uOff_x - the offset into u_x[] for each field 2312 . u - each field evaluated at the current point 2313 . u_t - the time derivative of each field evaluated at the current point 2314 . u_x - the gradient of each field evaluated at the current point 2315 . aOff - the offset into a[] and a_t[] for each auxiliary field 2316 . aOff_x - the offset into a_x[] for each auxiliary field 2317 . a - each auxiliary field evaluated at the current point 2318 . a_t - the time derivative of each auxiliary field evaluated at the current point 2319 . a_x - the gradient of auxiliary each field evaluated at the current point 2320 . t - current time 2321 . x - coordinates of the current point 2322 . n - unit normal at the current point 2323 . numConstants - number of constant parameters 2324 . constants - constant parameters 2325 - f0 - output values at the current point 2326 2327 Level: intermediate 2328 2329 Note: 2330 The calling sequence of `f1` is identical, and therefore omitted for brevity. 2331 2332 We are using a first order FEM model for the weak form\: 2333 \int_\Gamma \phi {\vec f}_0(u, u_t, \nabla u, x, t) \cdot \hat n + \nabla\phi \cdot {\overleftrightarrow f}_1(u, u_t, \nabla u, x, t) \cdot \hat n 2334 2335 .seealso: `PetscDS`, `PetscDSGetBdResidual()` 2336 @*/ 2337 PetscErrorCode PetscDSSetBdResidual(PetscDS ds, PetscInt f, void (*f0)(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[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (*f1)(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[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 2338 { 2339 PetscFunctionBegin; 2340 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2341 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 2342 PetscCall(PetscWeakFormSetIndexBdResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1)); 2343 PetscFunctionReturn(PETSC_SUCCESS); 2344 } 2345 2346 /*@ 2347 PetscDSHasBdJacobian - Indicates that boundary Jacobian functions have been set 2348 2349 Not Collective 2350 2351 Input Parameter: 2352 . ds - The `PetscDS` 2353 2354 Output Parameter: 2355 . hasBdJac - flag that pointwise function for the boundary Jacobian has been set 2356 2357 Level: intermediate 2358 2359 .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()` 2360 @*/ 2361 PetscErrorCode PetscDSHasBdJacobian(PetscDS ds, PetscBool *hasBdJac) 2362 { 2363 PetscFunctionBegin; 2364 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2365 PetscAssertPointer(hasBdJac, 2); 2366 PetscCall(PetscWeakFormHasBdJacobian(ds->wf, hasBdJac)); 2367 PetscFunctionReturn(PETSC_SUCCESS); 2368 } 2369 2370 /*@C 2371 PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field 2372 2373 Not Collective 2374 2375 Input Parameters: 2376 + ds - The `PetscDS` 2377 . f - The test field number 2378 - g - The field number 2379 2380 Output Parameters: 2381 + g0 - integrand for the test and basis function term 2382 . g1 - integrand for the test function and basis function gradient term 2383 . g2 - integrand for the test function gradient and basis function term 2384 - g3 - integrand for the test function gradient and basis function gradient term 2385 2386 Calling sequence of `g0`: 2387 + dim - the spatial dimension 2388 . Nf - the number of fields 2389 . NfAux - the number of auxiliary fields 2390 . uOff - the offset into u[] and u_t[] for each field 2391 . uOff_x - the offset into u_x[] for each field 2392 . u - each field evaluated at the current point 2393 . u_t - the time derivative of each field evaluated at the current point 2394 . u_x - the gradient of each field evaluated at the current point 2395 . aOff - the offset into a[] and a_t[] for each auxiliary field 2396 . aOff_x - the offset into a_x[] for each auxiliary field 2397 . a - each auxiliary field evaluated at the current point 2398 . a_t - the time derivative of each auxiliary field evaluated at the current point 2399 . a_x - the gradient of auxiliary each field evaluated at the current point 2400 . t - current time 2401 . u_tShift - the multiplier a for dF/dU_t 2402 . x - coordinates of the current point 2403 . n - normal at the current point 2404 . numConstants - number of constant parameters 2405 . constants - constant parameters 2406 - g0 - output values at the current point 2407 2408 Level: intermediate 2409 2410 Note: 2411 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 2412 2413 We are using a first order FEM model for the weak form\: 2414 \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi 2415 2416 .seealso: `PetscDS`, `PetscDSSetBdJacobian()` 2417 @*/ 2418 PetscErrorCode PetscDSGetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(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, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 2419 { 2420 PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3; 2421 PetscInt n0, n1, n2, n3; 2422 2423 PetscFunctionBegin; 2424 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2425 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 2426 PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf); 2427 PetscCall(PetscWeakFormGetBdJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3)); 2428 *g0 = tmp0 ? tmp0[0] : NULL; 2429 *g1 = tmp1 ? tmp1[0] : NULL; 2430 *g2 = tmp2 ? tmp2[0] : NULL; 2431 *g3 = tmp3 ? tmp3[0] : NULL; 2432 PetscFunctionReturn(PETSC_SUCCESS); 2433 } 2434 2435 /*@C 2436 PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field 2437 2438 Not Collective 2439 2440 Input Parameters: 2441 + ds - The PetscDS 2442 . f - The test field number 2443 . g - The field number 2444 . g0 - integrand for the test and basis function term 2445 . g1 - integrand for the test function and basis function gradient term 2446 . g2 - integrand for the test function gradient and basis function term 2447 - g3 - integrand for the test function gradient and basis function gradient term 2448 2449 Calling sequence of `g0`: 2450 + dim - the spatial dimension 2451 . Nf - the number of fields 2452 . NfAux - the number of auxiliary fields 2453 . uOff - the offset into u[] and u_t[] for each field 2454 . uOff_x - the offset into u_x[] for each field 2455 . u - each field evaluated at the current point 2456 . u_t - the time derivative of each field evaluated at the current point 2457 . u_x - the gradient of each field evaluated at the current point 2458 . aOff - the offset into a[] and a_t[] for each auxiliary field 2459 . aOff_x - the offset into a_x[] for each auxiliary field 2460 . a - each auxiliary field evaluated at the current point 2461 . a_t - the time derivative of each auxiliary field evaluated at the current point 2462 . a_x - the gradient of auxiliary each field evaluated at the current point 2463 . t - current time 2464 . u_tShift - the multiplier a for dF/dU_t 2465 . x - coordinates of the current point 2466 . n - normal at the current point 2467 . numConstants - number of constant parameters 2468 . constants - constant parameters 2469 - g0 - output values at the current point 2470 2471 Level: intermediate 2472 2473 Note: 2474 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 2475 2476 We are using a first order FEM model for the weak form\: 2477 \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi 2478 2479 .seealso: `PetscDS`, `PetscDSGetBdJacobian()` 2480 @*/ 2481 PetscErrorCode PetscDSSetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(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, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 2482 { 2483 PetscFunctionBegin; 2484 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2485 if (g0) PetscValidFunction(g0, 4); 2486 if (g1) PetscValidFunction(g1, 5); 2487 if (g2) PetscValidFunction(g2, 6); 2488 if (g3) PetscValidFunction(g3, 7); 2489 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 2490 PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g); 2491 PetscCall(PetscWeakFormSetIndexBdJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3)); 2492 PetscFunctionReturn(PETSC_SUCCESS); 2493 } 2494 2495 /*@ 2496 PetscDSHasBdJacobianPreconditioner - Signals that boundary Jacobian preconditioner functions have been set 2497 2498 Not Collective 2499 2500 Input Parameter: 2501 . ds - The `PetscDS` 2502 2503 Output Parameter: 2504 . hasBdJacPre - flag that pointwise function for the boundary Jacobian preconditioner has been set 2505 2506 Level: intermediate 2507 2508 .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()` 2509 @*/ 2510 PetscErrorCode PetscDSHasBdJacobianPreconditioner(PetscDS ds, PetscBool *hasBdJacPre) 2511 { 2512 PetscFunctionBegin; 2513 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2514 PetscAssertPointer(hasBdJacPre, 2); 2515 PetscCall(PetscWeakFormHasBdJacobianPreconditioner(ds->wf, hasBdJacPre)); 2516 PetscFunctionReturn(PETSC_SUCCESS); 2517 } 2518 2519 /*@C 2520 PetscDSGetBdJacobianPreconditioner - Get the pointwise boundary Jacobian preconditioner function for given test and basis field 2521 2522 Not Collective; No Fortran Support 2523 2524 Input Parameters: 2525 + ds - The `PetscDS` 2526 . f - The test field number 2527 - g - The field number 2528 2529 Output Parameters: 2530 + g0 - integrand for the test and basis function term 2531 . g1 - integrand for the test function and basis function gradient term 2532 . g2 - integrand for the test function gradient and basis function term 2533 - g3 - integrand for the test function gradient and basis function gradient term 2534 2535 Calling sequence of `g0`: 2536 + dim - the spatial dimension 2537 . Nf - the number of fields 2538 . NfAux - the number of auxiliary fields 2539 . uOff - the offset into u[] and u_t[] for each field 2540 . uOff_x - the offset into u_x[] for each field 2541 . u - each field evaluated at the current point 2542 . u_t - the time derivative of each field evaluated at the current point 2543 . u_x - the gradient of each field evaluated at the current point 2544 . aOff - the offset into a[] and a_t[] for each auxiliary field 2545 . aOff_x - the offset into a_x[] for each auxiliary field 2546 . a - each auxiliary field evaluated at the current point 2547 . a_t - the time derivative of each auxiliary field evaluated at the current point 2548 . a_x - the gradient of auxiliary each field evaluated at the current point 2549 . t - current time 2550 . u_tShift - the multiplier a for dF/dU_t 2551 . x - coordinates of the current point 2552 . n - normal at the current point 2553 . numConstants - number of constant parameters 2554 . constants - constant parameters 2555 - g0 - output values at the current point 2556 2557 Level: intermediate 2558 2559 Note: 2560 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 2561 2562 We are using a first order FEM model for the weak form\: 2563 \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi 2564 2565 .seealso: `PetscDS`, `PetscDSSetBdJacobianPreconditioner()` 2566 @*/ 2567 PetscErrorCode PetscDSGetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(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, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 2568 { 2569 PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3; 2570 PetscInt n0, n1, n2, n3; 2571 2572 PetscFunctionBegin; 2573 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2574 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 2575 PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf); 2576 PetscCall(PetscWeakFormGetBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3)); 2577 *g0 = tmp0 ? tmp0[0] : NULL; 2578 *g1 = tmp1 ? tmp1[0] : NULL; 2579 *g2 = tmp2 ? tmp2[0] : NULL; 2580 *g3 = tmp3 ? tmp3[0] : NULL; 2581 PetscFunctionReturn(PETSC_SUCCESS); 2582 } 2583 2584 /*@C 2585 PetscDSSetBdJacobianPreconditioner - Set the pointwise boundary Jacobian preconditioner function for given test and basis field 2586 2587 Not Collective; No Fortran Support 2588 2589 Input Parameters: 2590 + ds - The `PetscDS` 2591 . f - The test field number 2592 . g - The field number 2593 . g0 - integrand for the test and basis function term 2594 . g1 - integrand for the test function and basis function gradient term 2595 . g2 - integrand for the test function gradient and basis function term 2596 - g3 - integrand for the test function gradient and basis function gradient term 2597 2598 Calling sequence of `g0': 2599 + dim - the spatial dimension 2600 . Nf - the number of fields 2601 . NfAux - the number of auxiliary fields 2602 . uOff - the offset into u[] and u_t[] for each field 2603 . uOff_x - the offset into u_x[] for each field 2604 . u - each field evaluated at the current point 2605 . u_t - the time derivative of each field evaluated at the current point 2606 . u_x - the gradient of each field evaluated at the current point 2607 . aOff - the offset into a[] and a_t[] for each auxiliary field 2608 . aOff_x - the offset into a_x[] for each auxiliary field 2609 . a - each auxiliary field evaluated at the current point 2610 . a_t - the time derivative of each auxiliary field evaluated at the current point 2611 . a_x - the gradient of auxiliary each field evaluated at the current point 2612 . t - current time 2613 . u_tShift - the multiplier a for dF/dU_t 2614 . x - coordinates of the current point 2615 . n - normal at the current point 2616 . numConstants - number of constant parameters 2617 . constants - constant parameters 2618 - g0 - output values at the current point 2619 2620 Level: intermediate 2621 2622 Note: 2623 `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity. 2624 2625 We are using a first order FEM model for the weak form\: 2626 \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi 2627 2628 .seealso: `PetscDS`, `PetscDSGetBdJacobianPreconditioner()` 2629 @*/ 2630 PetscErrorCode PetscDSSetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(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, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 2631 { 2632 PetscFunctionBegin; 2633 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2634 if (g0) PetscValidFunction(g0, 4); 2635 if (g1) PetscValidFunction(g1, 5); 2636 if (g2) PetscValidFunction(g2, 6); 2637 if (g3) PetscValidFunction(g3, 7); 2638 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 2639 PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g); 2640 PetscCall(PetscWeakFormSetIndexBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3)); 2641 PetscFunctionReturn(PETSC_SUCCESS); 2642 } 2643 2644 /*@C 2645 PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field 2646 2647 Not Collective 2648 2649 Input Parameters: 2650 + prob - The PetscDS 2651 - f - The test field number 2652 2653 Output Parameters: 2654 + sol - exact solution for the test field 2655 - ctx - exact solution context 2656 2657 Calling sequence of `exactSol`: 2658 + dim - the spatial dimension 2659 . t - current time 2660 . x - coordinates of the current point 2661 . Nc - the number of field components 2662 . u - the solution field evaluated at the current point 2663 - ctx - a user context 2664 2665 Level: intermediate 2666 2667 .seealso: `PetscDS`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolutionTimeDerivative()` 2668 @*/ 2669 PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx) 2670 { 2671 PetscFunctionBegin; 2672 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2673 PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf); 2674 if (sol) { 2675 PetscAssertPointer(sol, 3); 2676 *sol = prob->exactSol[f]; 2677 } 2678 if (ctx) { 2679 PetscAssertPointer(ctx, 4); 2680 *ctx = prob->exactCtx[f]; 2681 } 2682 PetscFunctionReturn(PETSC_SUCCESS); 2683 } 2684 2685 /*@C 2686 PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field 2687 2688 Not Collective 2689 2690 Input Parameters: 2691 + prob - The `PetscDS` 2692 . f - The test field number 2693 . sol - solution function for the test fields 2694 - ctx - solution context or `NULL` 2695 2696 Calling sequence of `sol`: 2697 .vb 2698 PetscErrorCode sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) 2699 .ve 2700 + dim - the spatial dimension 2701 . t - current time 2702 . x - coordinates of the current point 2703 . Nc - the number of field components 2704 . u - the solution field evaluated at the current point 2705 - ctx - a user context 2706 2707 Level: intermediate 2708 2709 .seealso: `PetscDS`, `PetscDSGetExactSolution()` 2710 @*/ 2711 PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx) 2712 { 2713 PetscFunctionBegin; 2714 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2715 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 2716 PetscCall(PetscDSEnlarge_Static(prob, f + 1)); 2717 if (sol) { 2718 PetscValidFunction(sol, 3); 2719 prob->exactSol[f] = sol; 2720 } 2721 if (ctx) { 2722 PetscValidFunction(ctx, 4); 2723 prob->exactCtx[f] = ctx; 2724 } 2725 PetscFunctionReturn(PETSC_SUCCESS); 2726 } 2727 2728 /*@C 2729 PetscDSGetExactSolutionTimeDerivative - Get the pointwise time derivative of the exact solution function for a given test field 2730 2731 Not Collective 2732 2733 Input Parameters: 2734 + prob - The `PetscDS` 2735 - f - The test field number 2736 2737 Output Parameters: 2738 + sol - time derivative of the exact solution for the test field 2739 - ctx - time derivative of the exact solution context 2740 2741 Calling sequence of `exactSol`: 2742 + dim - the spatial dimension 2743 . t - current time 2744 . x - coordinates of the current point 2745 . Nc - the number of field components 2746 . u - the solution field evaluated at the current point 2747 - ctx - a user context 2748 2749 Level: intermediate 2750 2751 .seealso: `PetscDS`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolution()` 2752 @*/ 2753 PetscErrorCode PetscDSGetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx) 2754 { 2755 PetscFunctionBegin; 2756 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2757 PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf); 2758 if (sol) { 2759 PetscAssertPointer(sol, 3); 2760 *sol = prob->exactSol_t[f]; 2761 } 2762 if (ctx) { 2763 PetscAssertPointer(ctx, 4); 2764 *ctx = prob->exactCtx_t[f]; 2765 } 2766 PetscFunctionReturn(PETSC_SUCCESS); 2767 } 2768 2769 /*@C 2770 PetscDSSetExactSolutionTimeDerivative - Set the pointwise time derivative of the exact solution function for a given test field 2771 2772 Not Collective 2773 2774 Input Parameters: 2775 + prob - The `PetscDS` 2776 . f - The test field number 2777 . sol - time derivative of the solution function for the test fields 2778 - ctx - time derivative of the solution context or `NULL` 2779 2780 Calling sequence of `sol`: 2781 .vb 2782 PetscErrorCode sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) 2783 .ve 2784 + dim - the spatial dimension 2785 . t - current time 2786 . x - coordinates of the current point 2787 . Nc - the number of field components 2788 . u - the solution field evaluated at the current point 2789 - ctx - a user context 2790 2791 Level: intermediate 2792 2793 .seealso: `PetscDS`, `PetscDSGetExactSolutionTimeDerivative()`, `PetscDSSetExactSolution()` 2794 @*/ 2795 PetscErrorCode PetscDSSetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx) 2796 { 2797 PetscFunctionBegin; 2798 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2799 PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f); 2800 PetscCall(PetscDSEnlarge_Static(prob, f + 1)); 2801 if (sol) { 2802 PetscValidFunction(sol, 3); 2803 prob->exactSol_t[f] = sol; 2804 } 2805 if (ctx) { 2806 PetscValidFunction(ctx, 4); 2807 prob->exactCtx_t[f] = ctx; 2808 } 2809 PetscFunctionReturn(PETSC_SUCCESS); 2810 } 2811 2812 /*@C 2813 PetscDSGetConstants - Returns the array of constants passed to point functions 2814 2815 Not Collective 2816 2817 Input Parameter: 2818 . prob - The `PetscDS` object 2819 2820 Output Parameters: 2821 + numConstants - The number of constants 2822 - constants - The array of constants, NULL if there are none 2823 2824 Level: intermediate 2825 2826 .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSCreate()` 2827 @*/ 2828 PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[]) 2829 { 2830 PetscFunctionBegin; 2831 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2832 if (numConstants) { 2833 PetscAssertPointer(numConstants, 2); 2834 *numConstants = prob->numConstants; 2835 } 2836 if (constants) { 2837 PetscAssertPointer(constants, 3); 2838 *constants = prob->constants; 2839 } 2840 PetscFunctionReturn(PETSC_SUCCESS); 2841 } 2842 2843 /*@C 2844 PetscDSSetConstants - Set the array of constants passed to point functions 2845 2846 Not Collective 2847 2848 Input Parameters: 2849 + prob - The `PetscDS` object 2850 . numConstants - The number of constants 2851 - constants - The array of constants, NULL if there are none 2852 2853 Level: intermediate 2854 2855 .seealso: `PetscDS`, `PetscDSGetConstants()`, `PetscDSCreate()` 2856 @*/ 2857 PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[]) 2858 { 2859 PetscFunctionBegin; 2860 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2861 if (numConstants != prob->numConstants) { 2862 PetscCall(PetscFree(prob->constants)); 2863 prob->numConstants = numConstants; 2864 if (prob->numConstants) { 2865 PetscCall(PetscMalloc1(prob->numConstants, &prob->constants)); 2866 } else { 2867 prob->constants = NULL; 2868 } 2869 } 2870 if (prob->numConstants) { 2871 PetscAssertPointer(constants, 3); 2872 PetscCall(PetscArraycpy(prob->constants, constants, prob->numConstants)); 2873 } 2874 PetscFunctionReturn(PETSC_SUCCESS); 2875 } 2876 2877 /*@ 2878 PetscDSGetFieldIndex - Returns the index of the given field 2879 2880 Not Collective 2881 2882 Input Parameters: 2883 + prob - The `PetscDS` object 2884 - disc - The discretization object 2885 2886 Output Parameter: 2887 . f - The field number 2888 2889 Level: beginner 2890 2891 .seealso: `PetscDS`, `PetscGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 2892 @*/ 2893 PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f) 2894 { 2895 PetscInt g; 2896 2897 PetscFunctionBegin; 2898 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2899 PetscAssertPointer(f, 3); 2900 *f = -1; 2901 for (g = 0; g < prob->Nf; ++g) { 2902 if (disc == prob->disc[g]) break; 2903 } 2904 PetscCheck(g != prob->Nf, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS."); 2905 *f = g; 2906 PetscFunctionReturn(PETSC_SUCCESS); 2907 } 2908 2909 /*@ 2910 PetscDSGetFieldSize - Returns the size of the given field in the full space basis 2911 2912 Not Collective 2913 2914 Input Parameters: 2915 + prob - The `PetscDS` object 2916 - f - The field number 2917 2918 Output Parameter: 2919 . size - The size 2920 2921 Level: beginner 2922 2923 .seealso: `PetscDS`, `PetscDSGetFieldOffset()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 2924 @*/ 2925 PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size) 2926 { 2927 PetscFunctionBegin; 2928 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2929 PetscAssertPointer(size, 3); 2930 PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf); 2931 PetscCall(PetscDSSetUp(prob)); 2932 *size = prob->Nb[f]; 2933 PetscFunctionReturn(PETSC_SUCCESS); 2934 } 2935 2936 /*@ 2937 PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis 2938 2939 Not Collective 2940 2941 Input Parameters: 2942 + prob - The `PetscDS` object 2943 - f - The field number 2944 2945 Output Parameter: 2946 . off - The offset 2947 2948 Level: beginner 2949 2950 .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 2951 @*/ 2952 PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off) 2953 { 2954 PetscInt size, g; 2955 2956 PetscFunctionBegin; 2957 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2958 PetscAssertPointer(off, 3); 2959 PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf); 2960 *off = 0; 2961 for (g = 0; g < f; ++g) { 2962 PetscCall(PetscDSGetFieldSize(prob, g, &size)); 2963 *off += size; 2964 } 2965 PetscFunctionReturn(PETSC_SUCCESS); 2966 } 2967 2968 /*@ 2969 PetscDSGetFieldOffsetCohesive - Returns the offset of the given field in the full space basis on a cohesive cell 2970 2971 Not Collective 2972 2973 Input Parameters: 2974 + ds - The `PetscDS` object 2975 - f - The field number 2976 2977 Output Parameter: 2978 . off - The offset 2979 2980 Level: beginner 2981 2982 .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 2983 @*/ 2984 PetscErrorCode PetscDSGetFieldOffsetCohesive(PetscDS ds, PetscInt f, PetscInt *off) 2985 { 2986 PetscInt size, g; 2987 2988 PetscFunctionBegin; 2989 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2990 PetscAssertPointer(off, 3); 2991 PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf); 2992 *off = 0; 2993 for (g = 0; g < f; ++g) { 2994 PetscBool cohesive; 2995 2996 PetscCall(PetscDSGetCohesive(ds, g, &cohesive)); 2997 PetscCall(PetscDSGetFieldSize(ds, g, &size)); 2998 *off += cohesive ? size : size * 2; 2999 } 3000 PetscFunctionReturn(PETSC_SUCCESS); 3001 } 3002 3003 /*@ 3004 PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point 3005 3006 Not Collective 3007 3008 Input Parameter: 3009 . prob - The `PetscDS` object 3010 3011 Output Parameter: 3012 . dimensions - The number of dimensions 3013 3014 Level: beginner 3015 3016 .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 3017 @*/ 3018 PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[]) 3019 { 3020 PetscFunctionBegin; 3021 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3022 PetscCall(PetscDSSetUp(prob)); 3023 PetscAssertPointer(dimensions, 2); 3024 *dimensions = prob->Nb; 3025 PetscFunctionReturn(PETSC_SUCCESS); 3026 } 3027 3028 /*@ 3029 PetscDSGetComponents - Returns the number of components for each field on an evaluation point 3030 3031 Not Collective 3032 3033 Input Parameter: 3034 . prob - The `PetscDS` object 3035 3036 Output Parameter: 3037 . components - The number of components 3038 3039 Level: beginner 3040 3041 .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()` 3042 @*/ 3043 PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[]) 3044 { 3045 PetscFunctionBegin; 3046 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3047 PetscCall(PetscDSSetUp(prob)); 3048 PetscAssertPointer(components, 2); 3049 *components = prob->Nc; 3050 PetscFunctionReturn(PETSC_SUCCESS); 3051 } 3052 3053 /*@ 3054 PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point 3055 3056 Not Collective 3057 3058 Input Parameters: 3059 + prob - The `PetscDS` object 3060 - f - The field number 3061 3062 Output Parameter: 3063 . off - The offset 3064 3065 Level: beginner 3066 3067 .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()` 3068 @*/ 3069 PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off) 3070 { 3071 PetscFunctionBegin; 3072 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3073 PetscAssertPointer(off, 3); 3074 PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf); 3075 PetscCall(PetscDSSetUp(prob)); 3076 *off = prob->off[f]; 3077 PetscFunctionReturn(PETSC_SUCCESS); 3078 } 3079 3080 /*@ 3081 PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point 3082 3083 Not Collective 3084 3085 Input Parameter: 3086 . prob - The `PetscDS` object 3087 3088 Output Parameter: 3089 . offsets - The offsets 3090 3091 Level: beginner 3092 3093 .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()` 3094 @*/ 3095 PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[]) 3096 { 3097 PetscFunctionBegin; 3098 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3099 PetscAssertPointer(offsets, 2); 3100 PetscCall(PetscDSSetUp(prob)); 3101 *offsets = prob->off; 3102 PetscFunctionReturn(PETSC_SUCCESS); 3103 } 3104 3105 /*@ 3106 PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point 3107 3108 Not Collective 3109 3110 Input Parameter: 3111 . prob - The `PetscDS` object 3112 3113 Output Parameter: 3114 . offsets - The offsets 3115 3116 Level: beginner 3117 3118 .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()` 3119 @*/ 3120 PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[]) 3121 { 3122 PetscFunctionBegin; 3123 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3124 PetscAssertPointer(offsets, 2); 3125 PetscCall(PetscDSSetUp(prob)); 3126 *offsets = prob->offDer; 3127 PetscFunctionReturn(PETSC_SUCCESS); 3128 } 3129 3130 /*@ 3131 PetscDSGetComponentOffsetsCohesive - Returns the offset of each field on an evaluation point 3132 3133 Not Collective 3134 3135 Input Parameters: 3136 + ds - The `PetscDS` object 3137 - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive 3138 3139 Output Parameter: 3140 . offsets - The offsets 3141 3142 Level: beginner 3143 3144 .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()` 3145 @*/ 3146 PetscErrorCode PetscDSGetComponentOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[]) 3147 { 3148 PetscFunctionBegin; 3149 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 3150 PetscAssertPointer(offsets, 3); 3151 PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS"); 3152 PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s); 3153 PetscCall(PetscDSSetUp(ds)); 3154 *offsets = ds->offCohesive[s]; 3155 PetscFunctionReturn(PETSC_SUCCESS); 3156 } 3157 3158 /*@ 3159 PetscDSGetComponentDerivativeOffsetsCohesive - Returns the offset of each field derivative on an evaluation point 3160 3161 Not Collective 3162 3163 Input Parameters: 3164 + ds - The `PetscDS` object 3165 - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive 3166 3167 Output Parameter: 3168 . offsets - The offsets 3169 3170 Level: beginner 3171 3172 .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()` 3173 @*/ 3174 PetscErrorCode PetscDSGetComponentDerivativeOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[]) 3175 { 3176 PetscFunctionBegin; 3177 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 3178 PetscAssertPointer(offsets, 3); 3179 PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS"); 3180 PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s); 3181 PetscCall(PetscDSSetUp(ds)); 3182 *offsets = ds->offDerCohesive[s]; 3183 PetscFunctionReturn(PETSC_SUCCESS); 3184 } 3185 3186 /*@C 3187 PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization 3188 3189 Not Collective 3190 3191 Input Parameter: 3192 . prob - The `PetscDS` object 3193 3194 Output Parameter: 3195 . T - The basis function and derivatives tabulation at quadrature points for each field 3196 3197 Level: intermediate 3198 3199 .seealso: `PetscDS`, `PetscTabulation`, `PetscDSCreate()` 3200 @*/ 3201 PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[]) 3202 { 3203 PetscFunctionBegin; 3204 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3205 PetscAssertPointer(T, 2); 3206 PetscCall(PetscDSSetUp(prob)); 3207 *T = prob->T; 3208 PetscFunctionReturn(PETSC_SUCCESS); 3209 } 3210 3211 /*@C 3212 PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces 3213 3214 Not Collective 3215 3216 Input Parameter: 3217 . prob - The `PetscDS` object 3218 3219 Output Parameter: 3220 . Tf - The basis function and derivative tabulation on each local face at quadrature points for each and field 3221 3222 Level: intermediate 3223 3224 .seealso: `PetscTabulation`, `PetscDS`, `PetscDSGetTabulation()`, `PetscDSCreate()` 3225 @*/ 3226 PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[]) 3227 { 3228 PetscFunctionBegin; 3229 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3230 PetscAssertPointer(Tf, 2); 3231 PetscCall(PetscDSSetUp(prob)); 3232 *Tf = prob->Tf; 3233 PetscFunctionReturn(PETSC_SUCCESS); 3234 } 3235 3236 PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x) 3237 { 3238 PetscFunctionBegin; 3239 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3240 PetscCall(PetscDSSetUp(prob)); 3241 if (u) { 3242 PetscAssertPointer(u, 2); 3243 *u = prob->u; 3244 } 3245 if (u_t) { 3246 PetscAssertPointer(u_t, 3); 3247 *u_t = prob->u_t; 3248 } 3249 if (u_x) { 3250 PetscAssertPointer(u_x, 4); 3251 *u_x = prob->u_x; 3252 } 3253 PetscFunctionReturn(PETSC_SUCCESS); 3254 } 3255 3256 PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3) 3257 { 3258 PetscFunctionBegin; 3259 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3260 PetscCall(PetscDSSetUp(prob)); 3261 if (f0) { 3262 PetscAssertPointer(f0, 2); 3263 *f0 = prob->f0; 3264 } 3265 if (f1) { 3266 PetscAssertPointer(f1, 3); 3267 *f1 = prob->f1; 3268 } 3269 if (g0) { 3270 PetscAssertPointer(g0, 4); 3271 *g0 = prob->g0; 3272 } 3273 if (g1) { 3274 PetscAssertPointer(g1, 5); 3275 *g1 = prob->g1; 3276 } 3277 if (g2) { 3278 PetscAssertPointer(g2, 6); 3279 *g2 = prob->g2; 3280 } 3281 if (g3) { 3282 PetscAssertPointer(g3, 7); 3283 *g3 = prob->g3; 3284 } 3285 PetscFunctionReturn(PETSC_SUCCESS); 3286 } 3287 3288 PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal) 3289 { 3290 PetscFunctionBegin; 3291 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3292 PetscCall(PetscDSSetUp(prob)); 3293 if (x) { 3294 PetscAssertPointer(x, 2); 3295 *x = prob->x; 3296 } 3297 if (basisReal) { 3298 PetscAssertPointer(basisReal, 3); 3299 *basisReal = prob->basisReal; 3300 } 3301 if (basisDerReal) { 3302 PetscAssertPointer(basisDerReal, 4); 3303 *basisDerReal = prob->basisDerReal; 3304 } 3305 if (testReal) { 3306 PetscAssertPointer(testReal, 5); 3307 *testReal = prob->testReal; 3308 } 3309 if (testDerReal) { 3310 PetscAssertPointer(testDerReal, 6); 3311 *testDerReal = prob->testDerReal; 3312 } 3313 PetscFunctionReturn(PETSC_SUCCESS); 3314 } 3315 3316 /*@C 3317 PetscDSAddBoundary - Add a boundary condition to the model. 3318 3319 Collective 3320 3321 Input Parameters: 3322 + ds - The PetscDS object 3323 . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann) 3324 . name - The BC name 3325 . label - The label defining constrained points 3326 . Nv - The number of `DMLabel` values for constrained points 3327 . values - An array of label values for constrained points 3328 . field - The field to constrain 3329 . Nc - The number of constrained field components (0 will constrain all fields) 3330 . comps - An array of constrained component numbers 3331 . bcFunc - A pointwise function giving boundary values 3332 . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL 3333 - ctx - An optional user context for bcFunc 3334 3335 Output Parameter: 3336 . bd - The boundary number 3337 3338 Options Database Keys: 3339 + -bc_<boundary name> <num> - Overrides the boundary ids 3340 - -bc_<boundary name>_comp <num> - Overrides the boundary components 3341 3342 Level: developer 3343 3344 Note: 3345 Both `bcFunc` and `bcFunc_t` will depend on the boundary condition type. If the type if `DM_BC_ESSENTIAL`, then the calling sequence is\: 3346 3347 $ void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[]) 3348 3349 If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, then the calling sequence is\: 3350 .vb 3351 void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux, 3352 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 3353 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 3354 PetscReal time, const PetscReal x[], PetscScalar bcval[]) 3355 .ve 3356 + dim - the spatial dimension 3357 . Nf - the number of fields 3358 . uOff - the offset into u[] and u_t[] for each field 3359 . uOff_x - the offset into u_x[] for each field 3360 . u - each field evaluated at the current point 3361 . u_t - the time derivative of each field evaluated at the current point 3362 . u_x - the gradient of each field evaluated at the current point 3363 . aOff - the offset into a[] and a_t[] for each auxiliary field 3364 . aOff_x - the offset into a_x[] for each auxiliary field 3365 . a - each auxiliary field evaluated at the current point 3366 . a_t - the time derivative of each auxiliary field evaluated at the current point 3367 . a_x - the gradient of auxiliary each field evaluated at the current point 3368 . t - current time 3369 . x - coordinates of the current point 3370 . numConstants - number of constant parameters 3371 . constants - constant parameters 3372 - bcval - output values at the current point 3373 3374 Notes: 3375 The pointwise functions are used to provide boundary values for essential boundary 3376 conditions. In FEM, they are acting upon by dual basis functionals to generate FEM 3377 coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary 3378 integrals should be performed, using the kernels from `PetscDSSetBdResidual()`. 3379 3380 .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundaryByName()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()` 3381 @*/ 3382 PetscErrorCode PetscDSAddBoundary(PetscDS ds, DMBoundaryConditionType type, const char name[], DMLabel label, PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], void (*bcFunc)(void), void (*bcFunc_t)(void), void *ctx, PetscInt *bd) 3383 { 3384 DSBoundary head = ds->boundary, b; 3385 PetscInt n = 0; 3386 const char *lname; 3387 3388 PetscFunctionBegin; 3389 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 3390 PetscValidLogicalCollectiveEnum(ds, type, 2); 3391 PetscAssertPointer(name, 3); 3392 PetscValidHeaderSpecific(label, DMLABEL_CLASSID, 4); 3393 PetscValidLogicalCollectiveInt(ds, Nv, 5); 3394 PetscValidLogicalCollectiveInt(ds, field, 7); 3395 PetscValidLogicalCollectiveInt(ds, Nc, 8); 3396 PetscCheck(field >= 0 && field < ds->Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_OUTOFRANGE, "Field %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", field, ds->Nf); 3397 if (Nc > 0) { 3398 PetscInt *fcomps; 3399 PetscInt c; 3400 3401 PetscCall(PetscDSGetComponents(ds, &fcomps)); 3402 PetscCheck(Nc <= fcomps[field], PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_OUTOFRANGE, "Number of constrained components %" PetscInt_FMT " > %" PetscInt_FMT " components for field %" PetscInt_FMT, Nc, fcomps[field], field); 3403 for (c = 0; c < Nc; ++c) { 3404 PetscCheck(comps[c] >= 0 && comps[c] < fcomps[field], PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_OUTOFRANGE, "Constrained component[%" PetscInt_FMT "] %" PetscInt_FMT " not in [0, %" PetscInt_FMT ") components for field %" PetscInt_FMT, c, comps[c], fcomps[field], field); 3405 } 3406 } 3407 PetscCall(PetscNew(&b)); 3408 PetscCall(PetscStrallocpy(name, (char **)&b->name)); 3409 PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf)); 3410 PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf)); 3411 PetscCall(PetscMalloc1(Nv, &b->values)); 3412 if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv)); 3413 PetscCall(PetscMalloc1(Nc, &b->comps)); 3414 if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc)); 3415 PetscCall(PetscObjectGetName((PetscObject)label, &lname)); 3416 PetscCall(PetscStrallocpy(lname, (char **)&b->lname)); 3417 b->type = type; 3418 b->label = label; 3419 b->Nv = Nv; 3420 b->field = field; 3421 b->Nc = Nc; 3422 b->func = bcFunc; 3423 b->func_t = bcFunc_t; 3424 b->ctx = ctx; 3425 b->next = NULL; 3426 /* Append to linked list so that we can preserve the order */ 3427 if (!head) ds->boundary = b; 3428 while (head) { 3429 if (!head->next) { 3430 head->next = b; 3431 head = b; 3432 } 3433 head = head->next; 3434 ++n; 3435 } 3436 if (bd) { 3437 PetscAssertPointer(bd, 13); 3438 *bd = n; 3439 } 3440 PetscFunctionReturn(PETSC_SUCCESS); 3441 } 3442 3443 // PetscClangLinter pragma ignore: -fdoc-section-header-unknown 3444 /*@C 3445 PetscDSAddBoundaryByName - Add a boundary condition to the model. 3446 3447 Collective 3448 3449 Input Parameters: 3450 + ds - The `PetscDS` object 3451 . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann) 3452 . name - The BC name 3453 . lname - The naem of the label defining constrained points 3454 . Nv - The number of `DMLabel` values for constrained points 3455 . values - An array of label values for constrained points 3456 . field - The field to constrain 3457 . Nc - The number of constrained field components (0 will constrain all fields) 3458 . comps - An array of constrained component numbers 3459 . bcFunc - A pointwise function giving boundary values 3460 . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL 3461 - ctx - An optional user context for bcFunc 3462 3463 Output Parameter: 3464 . bd - The boundary number 3465 3466 Options Database Keys: 3467 + -bc_<boundary name> <num> - Overrides the boundary ids 3468 - -bc_<boundary name>_comp <num> - Overrides the boundary components 3469 3470 Calling Sequence of `bcFunc` and `bcFunc_t`: 3471 If the type is `DM_BC_ESSENTIAL` 3472 .vb 3473 void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[]) 3474 .ve 3475 If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, 3476 .vb 3477 void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux, 3478 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 3479 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 3480 PetscReal time, const PetscReal x[], PetscScalar bcval[]) 3481 .ve 3482 + dim - the spatial dimension 3483 . Nf - the number of fields 3484 . uOff - the offset into u[] and u_t[] for each field 3485 . uOff_x - the offset into u_x[] for each field 3486 . u - each field evaluated at the current point 3487 . u_t - the time derivative of each field evaluated at the current point 3488 . u_x - the gradient of each field evaluated at the current point 3489 . aOff - the offset into a[] and a_t[] for each auxiliary field 3490 . aOff_x - the offset into a_x[] for each auxiliary field 3491 . a - each auxiliary field evaluated at the current point 3492 . a_t - the time derivative of each auxiliary field evaluated at the current point 3493 . a_x - the gradient of auxiliary each field evaluated at the current point 3494 . t - current time 3495 . x - coordinates of the current point 3496 . numConstants - number of constant parameters 3497 . constants - constant parameters 3498 - bcval - output values at the current point 3499 3500 Level: developer 3501 3502 Notes: 3503 The pointwise functions are used to provide boundary values for essential boundary 3504 conditions. In FEM, they are acting upon by dual basis functionals to generate FEM 3505 coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary 3506 integrals should be performed, using the kernels from `PetscDSSetBdResidual()`. 3507 3508 This function should only be used with `DMFOREST` currently, since labels cannot be defined before the underlying `DMPLEX` is built. 3509 3510 .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()` 3511 @*/ 3512 PetscErrorCode PetscDSAddBoundaryByName(PetscDS ds, DMBoundaryConditionType type, const char name[], const char lname[], PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], void (*bcFunc)(void), void (*bcFunc_t)(void), void *ctx, PetscInt *bd) 3513 { 3514 DSBoundary head = ds->boundary, b; 3515 PetscInt n = 0; 3516 3517 PetscFunctionBegin; 3518 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 3519 PetscValidLogicalCollectiveEnum(ds, type, 2); 3520 PetscAssertPointer(name, 3); 3521 PetscAssertPointer(lname, 4); 3522 PetscValidLogicalCollectiveInt(ds, Nv, 5); 3523 PetscValidLogicalCollectiveInt(ds, field, 7); 3524 PetscValidLogicalCollectiveInt(ds, Nc, 8); 3525 PetscCall(PetscNew(&b)); 3526 PetscCall(PetscStrallocpy(name, (char **)&b->name)); 3527 PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf)); 3528 PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf)); 3529 PetscCall(PetscMalloc1(Nv, &b->values)); 3530 if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv)); 3531 PetscCall(PetscMalloc1(Nc, &b->comps)); 3532 if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc)); 3533 PetscCall(PetscStrallocpy(lname, (char **)&b->lname)); 3534 b->type = type; 3535 b->label = NULL; 3536 b->Nv = Nv; 3537 b->field = field; 3538 b->Nc = Nc; 3539 b->func = bcFunc; 3540 b->func_t = bcFunc_t; 3541 b->ctx = ctx; 3542 b->next = NULL; 3543 /* Append to linked list so that we can preserve the order */ 3544 if (!head) ds->boundary = b; 3545 while (head) { 3546 if (!head->next) { 3547 head->next = b; 3548 head = b; 3549 } 3550 head = head->next; 3551 ++n; 3552 } 3553 if (bd) { 3554 PetscAssertPointer(bd, 13); 3555 *bd = n; 3556 } 3557 PetscFunctionReturn(PETSC_SUCCESS); 3558 } 3559 3560 /*@C 3561 PetscDSUpdateBoundary - Change a boundary condition for the model. 3562 3563 Input Parameters: 3564 + ds - The `PetscDS` object 3565 . bd - The boundary condition number 3566 . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann) 3567 . name - The BC name 3568 . label - The label defining constrained points 3569 . Nv - The number of `DMLabel` ids for constrained points 3570 . values - An array of ids for constrained points 3571 . field - The field to constrain 3572 . Nc - The number of constrained field components 3573 . comps - An array of constrained component numbers 3574 . bcFunc - A pointwise function giving boundary values 3575 . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL 3576 - ctx - An optional user context for bcFunc 3577 3578 Level: developer 3579 3580 Notes: 3581 The pointwise functions are used to provide boundary values for essential boundary 3582 conditions. In FEM, they are acting upon by dual basis functionals to generate FEM 3583 coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary 3584 integrals should be performed, using the kernels from `PetscDSSetBdResidual()`. 3585 3586 The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from `PetscDSGetNumBoundary()`. 3587 See `PetscDSAddBoundary()` for a description of the calling sequences for the callbacks. 3588 3589 .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSGetNumBoundary()`, `DMLabel` 3590 @*/ 3591 PetscErrorCode PetscDSUpdateBoundary(PetscDS ds, PetscInt bd, DMBoundaryConditionType type, const char name[], DMLabel label, PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], void (*bcFunc)(void), void (*bcFunc_t)(void), void *ctx) 3592 { 3593 DSBoundary b = ds->boundary; 3594 PetscInt n = 0; 3595 3596 PetscFunctionBegin; 3597 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 3598 while (b) { 3599 if (n == bd) break; 3600 b = b->next; 3601 ++n; 3602 } 3603 PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n); 3604 if (name) { 3605 PetscCall(PetscFree(b->name)); 3606 PetscCall(PetscStrallocpy(name, (char **)&b->name)); 3607 } 3608 b->type = type; 3609 if (label) { 3610 const char *name; 3611 3612 b->label = label; 3613 PetscCall(PetscFree(b->lname)); 3614 PetscCall(PetscObjectGetName((PetscObject)label, &name)); 3615 PetscCall(PetscStrallocpy(name, (char **)&b->lname)); 3616 } 3617 if (Nv >= 0) { 3618 b->Nv = Nv; 3619 PetscCall(PetscFree(b->values)); 3620 PetscCall(PetscMalloc1(Nv, &b->values)); 3621 if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv)); 3622 } 3623 if (field >= 0) b->field = field; 3624 if (Nc >= 0) { 3625 b->Nc = Nc; 3626 PetscCall(PetscFree(b->comps)); 3627 PetscCall(PetscMalloc1(Nc, &b->comps)); 3628 if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc)); 3629 } 3630 if (bcFunc) b->func = bcFunc; 3631 if (bcFunc_t) b->func_t = bcFunc_t; 3632 if (ctx) b->ctx = ctx; 3633 PetscFunctionReturn(PETSC_SUCCESS); 3634 } 3635 3636 /*@ 3637 PetscDSGetNumBoundary - Get the number of registered BC 3638 3639 Input Parameter: 3640 . ds - The `PetscDS` object 3641 3642 Output Parameter: 3643 . numBd - The number of BC 3644 3645 Level: intermediate 3646 3647 .seealso: `PetscDS`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()` 3648 @*/ 3649 PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd) 3650 { 3651 DSBoundary b = ds->boundary; 3652 3653 PetscFunctionBegin; 3654 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 3655 PetscAssertPointer(numBd, 2); 3656 *numBd = 0; 3657 while (b) { 3658 ++(*numBd); 3659 b = b->next; 3660 } 3661 PetscFunctionReturn(PETSC_SUCCESS); 3662 } 3663 3664 /*@C 3665 PetscDSGetBoundary - Gets a boundary condition to the model 3666 3667 Input Parameters: 3668 + ds - The `PetscDS` object 3669 - bd - The BC number 3670 3671 Output Parameters: 3672 + wf - The `PetscWeakForm` holding the pointwise functions 3673 . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann) 3674 . name - The BC name 3675 . label - The label defining constrained points 3676 . Nv - The number of `DMLabel` ids for constrained points 3677 . values - An array of ids for constrained points 3678 . field - The field to constrain 3679 . Nc - The number of constrained field components 3680 . comps - An array of constrained component numbers 3681 . func - A pointwise function giving boundary values 3682 . func_t - A pointwise function giving the time derivative of the boundary values 3683 - ctx - An optional user context for bcFunc 3684 3685 Options Database Keys: 3686 + -bc_<boundary name> <num> - Overrides the boundary ids 3687 - -bc_<boundary name>_comp <num> - Overrides the boundary components 3688 3689 Level: developer 3690 3691 .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `DMLabel` 3692 @*/ 3693 PetscErrorCode PetscDSGetBoundary(PetscDS ds, PetscInt bd, PetscWeakForm *wf, DMBoundaryConditionType *type, const char *name[], DMLabel *label, PetscInt *Nv, const PetscInt *values[], PetscInt *field, PetscInt *Nc, const PetscInt *comps[], void (**func)(void), void (**func_t)(void), void **ctx) 3694 { 3695 DSBoundary b = ds->boundary; 3696 PetscInt n = 0; 3697 3698 PetscFunctionBegin; 3699 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 3700 while (b) { 3701 if (n == bd) break; 3702 b = b->next; 3703 ++n; 3704 } 3705 PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n); 3706 if (wf) { 3707 PetscAssertPointer(wf, 3); 3708 *wf = b->wf; 3709 } 3710 if (type) { 3711 PetscAssertPointer(type, 4); 3712 *type = b->type; 3713 } 3714 if (name) { 3715 PetscAssertPointer(name, 5); 3716 *name = b->name; 3717 } 3718 if (label) { 3719 PetscAssertPointer(label, 6); 3720 *label = b->label; 3721 } 3722 if (Nv) { 3723 PetscAssertPointer(Nv, 7); 3724 *Nv = b->Nv; 3725 } 3726 if (values) { 3727 PetscAssertPointer(values, 8); 3728 *values = b->values; 3729 } 3730 if (field) { 3731 PetscAssertPointer(field, 9); 3732 *field = b->field; 3733 } 3734 if (Nc) { 3735 PetscAssertPointer(Nc, 10); 3736 *Nc = b->Nc; 3737 } 3738 if (comps) { 3739 PetscAssertPointer(comps, 11); 3740 *comps = b->comps; 3741 } 3742 if (func) { 3743 PetscAssertPointer(func, 12); 3744 *func = b->func; 3745 } 3746 if (func_t) { 3747 PetscAssertPointer(func_t, 13); 3748 *func_t = b->func_t; 3749 } 3750 if (ctx) { 3751 PetscAssertPointer(ctx, 14); 3752 *ctx = b->ctx; 3753 } 3754 PetscFunctionReturn(PETSC_SUCCESS); 3755 } 3756 3757 static PetscErrorCode DSBoundaryDuplicate_Internal(DSBoundary b, DSBoundary *bNew) 3758 { 3759 PetscFunctionBegin; 3760 PetscCall(PetscNew(bNew)); 3761 PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &(*bNew)->wf)); 3762 PetscCall(PetscWeakFormCopy(b->wf, (*bNew)->wf)); 3763 PetscCall(PetscStrallocpy(b->name, (char **)&((*bNew)->name))); 3764 PetscCall(PetscStrallocpy(b->lname, (char **)&((*bNew)->lname))); 3765 (*bNew)->type = b->type; 3766 (*bNew)->label = b->label; 3767 (*bNew)->Nv = b->Nv; 3768 PetscCall(PetscMalloc1(b->Nv, &(*bNew)->values)); 3769 PetscCall(PetscArraycpy((*bNew)->values, b->values, b->Nv)); 3770 (*bNew)->field = b->field; 3771 (*bNew)->Nc = b->Nc; 3772 PetscCall(PetscMalloc1(b->Nc, &(*bNew)->comps)); 3773 PetscCall(PetscArraycpy((*bNew)->comps, b->comps, b->Nc)); 3774 (*bNew)->func = b->func; 3775 (*bNew)->func_t = b->func_t; 3776 (*bNew)->ctx = b->ctx; 3777 PetscFunctionReturn(PETSC_SUCCESS); 3778 } 3779 3780 /*@ 3781 PetscDSCopyBoundary - Copy all boundary condition objects to the new problem 3782 3783 Not Collective 3784 3785 Input Parameters: 3786 + ds - The source `PetscDS` object 3787 . numFields - The number of selected fields, or `PETSC_DEFAULT` for all fields 3788 - fields - The selected fields, or NULL for all fields 3789 3790 Output Parameter: 3791 . newds - The target `PetscDS`, now with a copy of the boundary conditions 3792 3793 Level: intermediate 3794 3795 .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()` 3796 @*/ 3797 PetscErrorCode PetscDSCopyBoundary(PetscDS ds, PetscInt numFields, const PetscInt fields[], PetscDS newds) 3798 { 3799 DSBoundary b, *lastnext; 3800 3801 PetscFunctionBegin; 3802 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 3803 PetscValidHeaderSpecific(newds, PETSCDS_CLASSID, 4); 3804 if (ds == newds) PetscFunctionReturn(PETSC_SUCCESS); 3805 PetscCall(PetscDSDestroyBoundary(newds)); 3806 lastnext = &(newds->boundary); 3807 for (b = ds->boundary; b; b = b->next) { 3808 DSBoundary bNew; 3809 PetscInt fieldNew = -1; 3810 3811 if (numFields > 0 && fields) { 3812 PetscInt f; 3813 3814 for (f = 0; f < numFields; ++f) 3815 if (b->field == fields[f]) break; 3816 if (f == numFields) continue; 3817 fieldNew = f; 3818 } 3819 PetscCall(DSBoundaryDuplicate_Internal(b, &bNew)); 3820 bNew->field = fieldNew < 0 ? b->field : fieldNew; 3821 *lastnext = bNew; 3822 lastnext = &(bNew->next); 3823 } 3824 PetscFunctionReturn(PETSC_SUCCESS); 3825 } 3826 3827 /*@ 3828 PetscDSDestroyBoundary - Remove all `DMBoundary` objects from the `PetscDS` 3829 3830 Not Collective 3831 3832 Input Parameter: 3833 . ds - The `PetscDS` object 3834 3835 Level: intermediate 3836 3837 .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()` 3838 @*/ 3839 PetscErrorCode PetscDSDestroyBoundary(PetscDS ds) 3840 { 3841 DSBoundary next = ds->boundary; 3842 3843 PetscFunctionBegin; 3844 while (next) { 3845 DSBoundary b = next; 3846 3847 next = b->next; 3848 PetscCall(PetscWeakFormDestroy(&b->wf)); 3849 PetscCall(PetscFree(b->name)); 3850 PetscCall(PetscFree(b->lname)); 3851 PetscCall(PetscFree(b->values)); 3852 PetscCall(PetscFree(b->comps)); 3853 PetscCall(PetscFree(b)); 3854 } 3855 PetscFunctionReturn(PETSC_SUCCESS); 3856 } 3857 3858 /*@ 3859 PetscDSSelectDiscretizations - Copy discretizations to the new problem with different field layout 3860 3861 Not Collective 3862 3863 Input Parameters: 3864 + prob - The `PetscDS` object 3865 . numFields - Number of new fields 3866 - fields - Old field number for each new field 3867 3868 Output Parameter: 3869 . newprob - The `PetscDS` copy 3870 3871 Level: intermediate 3872 3873 .seealso: `PetscDS`, `PetscDSSelectEquations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()` 3874 @*/ 3875 PetscErrorCode PetscDSSelectDiscretizations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob) 3876 { 3877 PetscInt Nf, Nfn, fn; 3878 3879 PetscFunctionBegin; 3880 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3881 if (fields) PetscAssertPointer(fields, 3); 3882 PetscValidHeaderSpecific(newprob, PETSCDS_CLASSID, 4); 3883 PetscCall(PetscDSGetNumFields(prob, &Nf)); 3884 PetscCall(PetscDSGetNumFields(newprob, &Nfn)); 3885 numFields = numFields < 0 ? Nf : numFields; 3886 for (fn = 0; fn < numFields; ++fn) { 3887 const PetscInt f = fields ? fields[fn] : fn; 3888 PetscObject disc; 3889 3890 if (f >= Nf) continue; 3891 PetscCall(PetscDSGetDiscretization(prob, f, &disc)); 3892 PetscCall(PetscDSSetDiscretization(newprob, fn, disc)); 3893 } 3894 PetscFunctionReturn(PETSC_SUCCESS); 3895 } 3896 3897 /*@ 3898 PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout 3899 3900 Not Collective 3901 3902 Input Parameters: 3903 + prob - The `PetscDS` object 3904 . numFields - Number of new fields 3905 - fields - Old field number for each new field 3906 3907 Output Parameter: 3908 . newprob - The `PetscDS` copy 3909 3910 Level: intermediate 3911 3912 .seealso: `PetscDS`, `PetscDSSelectDiscretizations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()` 3913 @*/ 3914 PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob) 3915 { 3916 PetscInt Nf, Nfn, fn, gn; 3917 3918 PetscFunctionBegin; 3919 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3920 if (fields) PetscAssertPointer(fields, 3); 3921 PetscValidHeaderSpecific(newprob, PETSCDS_CLASSID, 4); 3922 PetscCall(PetscDSGetNumFields(prob, &Nf)); 3923 PetscCall(PetscDSGetNumFields(newprob, &Nfn)); 3924 PetscCheck(numFields <= Nfn, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields %" PetscInt_FMT " to transfer must not be greater then the total number of fields %" PetscInt_FMT, numFields, Nfn); 3925 for (fn = 0; fn < numFields; ++fn) { 3926 const PetscInt f = fields ? fields[fn] : fn; 3927 PetscPointFunc obj; 3928 PetscPointFunc f0, f1; 3929 PetscBdPointFunc f0Bd, f1Bd; 3930 PetscRiemannFunc r; 3931 3932 if (f >= Nf) continue; 3933 PetscCall(PetscDSGetObjective(prob, f, &obj)); 3934 PetscCall(PetscDSGetResidual(prob, f, &f0, &f1)); 3935 PetscCall(PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd)); 3936 PetscCall(PetscDSGetRiemannSolver(prob, f, &r)); 3937 PetscCall(PetscDSSetObjective(newprob, fn, obj)); 3938 PetscCall(PetscDSSetResidual(newprob, fn, f0, f1)); 3939 PetscCall(PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd)); 3940 PetscCall(PetscDSSetRiemannSolver(newprob, fn, r)); 3941 for (gn = 0; gn < numFields; ++gn) { 3942 const PetscInt g = fields ? fields[gn] : gn; 3943 PetscPointJac g0, g1, g2, g3; 3944 PetscPointJac g0p, g1p, g2p, g3p; 3945 PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd; 3946 3947 if (g >= Nf) continue; 3948 PetscCall(PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3)); 3949 PetscCall(PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p)); 3950 PetscCall(PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd)); 3951 PetscCall(PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3)); 3952 PetscCall(PetscDSSetJacobianPreconditioner(newprob, fn, gn, g0p, g1p, g2p, g3p)); 3953 PetscCall(PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd)); 3954 } 3955 } 3956 PetscFunctionReturn(PETSC_SUCCESS); 3957 } 3958 3959 /*@ 3960 PetscDSCopyEquations - Copy all pointwise function pointers to another `PetscDS` 3961 3962 Not Collective 3963 3964 Input Parameter: 3965 . prob - The `PetscDS` object 3966 3967 Output Parameter: 3968 . newprob - The `PetscDS` copy 3969 3970 Level: intermediate 3971 3972 .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()` 3973 @*/ 3974 PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob) 3975 { 3976 PetscWeakForm wf, newwf; 3977 PetscInt Nf, Ng; 3978 3979 PetscFunctionBegin; 3980 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 3981 PetscValidHeaderSpecific(newprob, PETSCDS_CLASSID, 2); 3982 PetscCall(PetscDSGetNumFields(prob, &Nf)); 3983 PetscCall(PetscDSGetNumFields(newprob, &Ng)); 3984 PetscCheck(Nf == Ng, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %" PetscInt_FMT " != %" PetscInt_FMT, Nf, Ng); 3985 PetscCall(PetscDSGetWeakForm(prob, &wf)); 3986 PetscCall(PetscDSGetWeakForm(newprob, &newwf)); 3987 PetscCall(PetscWeakFormCopy(wf, newwf)); 3988 PetscFunctionReturn(PETSC_SUCCESS); 3989 } 3990 3991 /*@ 3992 PetscDSCopyConstants - Copy all constants to another `PetscDS` 3993 3994 Not Collective 3995 3996 Input Parameter: 3997 . prob - The `PetscDS` object 3998 3999 Output Parameter: 4000 . newprob - The `PetscDS` copy 4001 4002 Level: intermediate 4003 4004 .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()` 4005 @*/ 4006 PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob) 4007 { 4008 PetscInt Nc; 4009 const PetscScalar *constants; 4010 4011 PetscFunctionBegin; 4012 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 4013 PetscValidHeaderSpecific(newprob, PETSCDS_CLASSID, 2); 4014 PetscCall(PetscDSGetConstants(prob, &Nc, &constants)); 4015 PetscCall(PetscDSSetConstants(newprob, Nc, (PetscScalar *)constants)); 4016 PetscFunctionReturn(PETSC_SUCCESS); 4017 } 4018 4019 /*@ 4020 PetscDSCopyExactSolutions - Copy all exact solutions to another `PetscDS` 4021 4022 Not Collective 4023 4024 Input Parameter: 4025 . ds - The `PetscDS` object 4026 4027 Output Parameter: 4028 . newds - The `PetscDS` copy 4029 4030 Level: intermediate 4031 4032 .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()` 4033 @*/ 4034 PetscErrorCode PetscDSCopyExactSolutions(PetscDS ds, PetscDS newds) 4035 { 4036 PetscSimplePointFunc sol; 4037 void *ctx; 4038 PetscInt Nf, f; 4039 4040 PetscFunctionBegin; 4041 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 4042 PetscValidHeaderSpecific(newds, PETSCDS_CLASSID, 2); 4043 PetscCall(PetscDSGetNumFields(ds, &Nf)); 4044 for (f = 0; f < Nf; ++f) { 4045 PetscCall(PetscDSGetExactSolution(ds, f, &sol, &ctx)); 4046 PetscCall(PetscDSSetExactSolution(newds, f, sol, ctx)); 4047 PetscCall(PetscDSGetExactSolutionTimeDerivative(ds, f, &sol, &ctx)); 4048 PetscCall(PetscDSSetExactSolutionTimeDerivative(newds, f, sol, ctx)); 4049 } 4050 PetscFunctionReturn(PETSC_SUCCESS); 4051 } 4052 4053 PetscErrorCode PetscDSCopy(PetscDS ds, DM dmNew, PetscDS dsNew) 4054 { 4055 DSBoundary b; 4056 PetscInt cdim, Nf, f, d; 4057 PetscBool isCohesive; 4058 void *ctx; 4059 4060 PetscFunctionBegin; 4061 PetscCall(PetscDSCopyConstants(ds, dsNew)); 4062 PetscCall(PetscDSCopyExactSolutions(ds, dsNew)); 4063 PetscCall(PetscDSSelectDiscretizations(ds, PETSC_DETERMINE, NULL, dsNew)); 4064 PetscCall(PetscDSCopyEquations(ds, dsNew)); 4065 PetscCall(PetscDSGetNumFields(ds, &Nf)); 4066 for (f = 0; f < Nf; ++f) { 4067 PetscCall(PetscDSGetContext(ds, f, &ctx)); 4068 PetscCall(PetscDSSetContext(dsNew, f, ctx)); 4069 PetscCall(PetscDSGetCohesive(ds, f, &isCohesive)); 4070 PetscCall(PetscDSSetCohesive(dsNew, f, isCohesive)); 4071 PetscCall(PetscDSGetJetDegree(ds, f, &d)); 4072 PetscCall(PetscDSSetJetDegree(dsNew, f, d)); 4073 } 4074 if (Nf) { 4075 PetscCall(PetscDSGetCoordinateDimension(ds, &cdim)); 4076 PetscCall(PetscDSSetCoordinateDimension(dsNew, cdim)); 4077 } 4078 PetscCall(PetscDSCopyBoundary(ds, PETSC_DETERMINE, NULL, dsNew)); 4079 for (b = dsNew->boundary; b; b = b->next) { 4080 PetscCall(DMGetLabel(dmNew, b->lname, &b->label)); 4081 /* Do not check if label exists here, since p4est calls this for the reference tree which does not have the labels */ 4082 //PetscCheck(b->label,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Label %s missing in new DM", name); 4083 } 4084 PetscFunctionReturn(PETSC_SUCCESS); 4085 } 4086 4087 PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob) 4088 { 4089 PetscInt dim, Nf, f; 4090 4091 PetscFunctionBegin; 4092 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 4093 PetscAssertPointer(subprob, 3); 4094 if (height == 0) { 4095 *subprob = prob; 4096 PetscFunctionReturn(PETSC_SUCCESS); 4097 } 4098 PetscCall(PetscDSGetNumFields(prob, &Nf)); 4099 PetscCall(PetscDSGetSpatialDimension(prob, &dim)); 4100 PetscCheck(height <= dim, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %" PetscInt_FMT "], not %" PetscInt_FMT, dim, height); 4101 if (!prob->subprobs) PetscCall(PetscCalloc1(dim, &prob->subprobs)); 4102 if (!prob->subprobs[height - 1]) { 4103 PetscInt cdim; 4104 4105 PetscCall(PetscDSCreate(PetscObjectComm((PetscObject)prob), &prob->subprobs[height - 1])); 4106 PetscCall(PetscDSGetCoordinateDimension(prob, &cdim)); 4107 PetscCall(PetscDSSetCoordinateDimension(prob->subprobs[height - 1], cdim)); 4108 for (f = 0; f < Nf; ++f) { 4109 PetscFE subfe; 4110 PetscObject obj; 4111 PetscClassId id; 4112 4113 PetscCall(PetscDSGetDiscretization(prob, f, &obj)); 4114 PetscCall(PetscObjectGetClassId(obj, &id)); 4115 if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetHeightSubspace((PetscFE)obj, height, &subfe)); 4116 else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %" PetscInt_FMT, f); 4117 PetscCall(PetscDSSetDiscretization(prob->subprobs[height - 1], f, (PetscObject)subfe)); 4118 } 4119 } 4120 *subprob = prob->subprobs[height - 1]; 4121 PetscFunctionReturn(PETSC_SUCCESS); 4122 } 4123 4124 PetscErrorCode PetscDSPermuteQuadPoint(PetscDS ds, PetscInt ornt, PetscInt field, PetscInt q, PetscInt *qperm) 4125 { 4126 IS permIS; 4127 PetscQuadrature quad; 4128 DMPolytopeType ct; 4129 const PetscInt *perm; 4130 PetscInt Na, Nq; 4131 4132 PetscFunctionBeginHot; 4133 PetscCall(PetscFEGetQuadrature((PetscFE)ds->disc[field], &quad)); 4134 PetscCall(PetscQuadratureGetData(quad, NULL, NULL, &Nq, NULL, NULL)); 4135 PetscCall(PetscQuadratureGetCellType(quad, &ct)); 4136 PetscCheck(q >= 0 && q < Nq, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Quadrature point %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", q, Nq); 4137 Na = DMPolytopeTypeGetNumArrangments(ct) / 2; 4138 PetscCheck(ornt >= -Na && ornt < Na, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Orientation %" PetscInt_FMT " of %s is not in [%" PetscInt_FMT ", %" PetscInt_FMT ")", ornt, DMPolytopeTypes[ct], -Na, Na); 4139 if (!ds->quadPerm[(PetscInt)ct]) PetscCall(PetscQuadratureComputePermutations(quad, NULL, &ds->quadPerm[(PetscInt)ct])); 4140 permIS = ds->quadPerm[(PetscInt)ct][ornt + Na]; 4141 PetscCall(ISGetIndices(permIS, &perm)); 4142 *qperm = perm[q]; 4143 PetscCall(ISRestoreIndices(permIS, &perm)); 4144 PetscFunctionReturn(PETSC_SUCCESS); 4145 } 4146 4147 PetscErrorCode PetscDSGetDiscType_Internal(PetscDS ds, PetscInt f, PetscDiscType *disctype) 4148 { 4149 PetscObject obj; 4150 PetscClassId id; 4151 PetscInt Nf; 4152 4153 PetscFunctionBegin; 4154 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 4155 PetscAssertPointer(disctype, 3); 4156 *disctype = PETSC_DISC_NONE; 4157 PetscCall(PetscDSGetNumFields(ds, &Nf)); 4158 PetscCheck(f < Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_SIZ, "Field %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, Nf); 4159 PetscCall(PetscDSGetDiscretization(ds, f, &obj)); 4160 if (obj) { 4161 PetscCall(PetscObjectGetClassId(obj, &id)); 4162 if (id == PETSCFE_CLASSID) *disctype = PETSC_DISC_FE; 4163 else *disctype = PETSC_DISC_FV; 4164 } 4165 PetscFunctionReturn(PETSC_SUCCESS); 4166 } 4167 4168 static PetscErrorCode PetscDSDestroy_Basic(PetscDS ds) 4169 { 4170 PetscFunctionBegin; 4171 PetscCall(PetscFree(ds->data)); 4172 PetscFunctionReturn(PETSC_SUCCESS); 4173 } 4174 4175 static PetscErrorCode PetscDSInitialize_Basic(PetscDS ds) 4176 { 4177 PetscFunctionBegin; 4178 ds->ops->setfromoptions = NULL; 4179 ds->ops->setup = NULL; 4180 ds->ops->view = NULL; 4181 ds->ops->destroy = PetscDSDestroy_Basic; 4182 PetscFunctionReturn(PETSC_SUCCESS); 4183 } 4184 4185 /*MC 4186 PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions 4187 4188 Level: intermediate 4189 4190 .seealso: `PetscDSType`, `PetscDSCreate()`, `PetscDSSetType()` 4191 M*/ 4192 4193 PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS ds) 4194 { 4195 PetscDS_Basic *b; 4196 4197 PetscFunctionBegin; 4198 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 4199 PetscCall(PetscNew(&b)); 4200 ds->data = b; 4201 4202 PetscCall(PetscDSInitialize_Basic(ds)); 4203 PetscFunctionReturn(PETSC_SUCCESS); 4204 } 4205