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