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 #undef __FUNCT__ 9 #define __FUNCT__ "PetscDSRegister" 10 /*@C 11 PetscDSRegister - Adds a new PetscDS implementation 12 13 Not Collective 14 15 Input Parameters: 16 + name - The name of a new user-defined creation routine 17 - create_func - The creation routine itself 18 19 Notes: 20 PetscDSRegister() may be called multiple times to add several user-defined PetscDSs 21 22 Sample usage: 23 .vb 24 PetscDSRegister("my_ds", MyPetscDSCreate); 25 .ve 26 27 Then, your PetscDS type can be chosen with the procedural interface via 28 .vb 29 PetscDSCreate(MPI_Comm, PetscDS *); 30 PetscDSSetType(PetscDS, "my_ds"); 31 .ve 32 or at runtime via the option 33 .vb 34 -petscds_type my_ds 35 .ve 36 37 Level: advanced 38 39 .keywords: PetscDS, register 40 .seealso: PetscDSRegisterAll(), PetscDSRegisterDestroy() 41 42 @*/ 43 PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS)) 44 { 45 PetscErrorCode ierr; 46 47 PetscFunctionBegin; 48 ierr = PetscFunctionListAdd(&PetscDSList, sname, function);CHKERRQ(ierr); 49 PetscFunctionReturn(0); 50 } 51 52 #undef __FUNCT__ 53 #define __FUNCT__ "PetscDSSetType" 54 /*@C 55 PetscDSSetType - Builds a particular PetscDS 56 57 Collective on PetscDS 58 59 Input Parameters: 60 + prob - The PetscDS object 61 - name - The kind of system 62 63 Options Database Key: 64 . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types 65 66 Level: intermediate 67 68 .keywords: PetscDS, set, type 69 .seealso: PetscDSGetType(), PetscDSCreate() 70 @*/ 71 PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name) 72 { 73 PetscErrorCode (*r)(PetscDS); 74 PetscBool match; 75 PetscErrorCode ierr; 76 77 PetscFunctionBegin; 78 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 79 ierr = PetscObjectTypeCompare((PetscObject) prob, name, &match);CHKERRQ(ierr); 80 if (match) PetscFunctionReturn(0); 81 82 ierr = PetscDSRegisterAll();CHKERRQ(ierr); 83 ierr = PetscFunctionListFind(PetscDSList, name, &r);CHKERRQ(ierr); 84 if (!r) SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name); 85 86 if (prob->ops->destroy) { 87 ierr = (*prob->ops->destroy)(prob);CHKERRQ(ierr); 88 prob->ops->destroy = NULL; 89 } 90 ierr = (*r)(prob);CHKERRQ(ierr); 91 ierr = PetscObjectChangeTypeName((PetscObject) prob, name);CHKERRQ(ierr); 92 PetscFunctionReturn(0); 93 } 94 95 #undef __FUNCT__ 96 #define __FUNCT__ "PetscDSGetType" 97 /*@C 98 PetscDSGetType - Gets the PetscDS type name (as a string) from the object. 99 100 Not Collective 101 102 Input Parameter: 103 . prob - The PetscDS 104 105 Output Parameter: 106 . name - The PetscDS type name 107 108 Level: intermediate 109 110 .keywords: PetscDS, get, type, name 111 .seealso: PetscDSSetType(), PetscDSCreate() 112 @*/ 113 PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name) 114 { 115 PetscErrorCode ierr; 116 117 PetscFunctionBegin; 118 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 119 PetscValidPointer(name, 2); 120 ierr = PetscDSRegisterAll();CHKERRQ(ierr); 121 *name = ((PetscObject) prob)->type_name; 122 PetscFunctionReturn(0); 123 } 124 125 #undef __FUNCT__ 126 #define __FUNCT__ "PetscDSView_Ascii" 127 static PetscErrorCode PetscDSView_Ascii(PetscDS prob, PetscViewer viewer) 128 { 129 PetscViewerFormat format; 130 PetscInt f; 131 PetscErrorCode ierr; 132 133 PetscFunctionBegin; 134 ierr = PetscViewerGetFormat(viewer, &format);CHKERRQ(ierr); 135 ierr = PetscViewerASCIIPrintf(viewer, "Discrete System with %d fields\n", prob->Nf);CHKERRQ(ierr); 136 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 137 for (f = 0; f < prob->Nf; ++f) { 138 PetscObject obj; 139 PetscClassId id; 140 const char *name; 141 PetscInt Nc; 142 143 ierr = PetscDSGetDiscretization(prob, f, &obj);CHKERRQ(ierr); 144 ierr = PetscObjectGetClassId(obj, &id);CHKERRQ(ierr); 145 ierr = PetscObjectGetName(obj, &name);CHKERRQ(ierr); 146 ierr = PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>");CHKERRQ(ierr); 147 if (id == PETSCFE_CLASSID) { 148 ierr = PetscFEGetNumComponents((PetscFE) obj, &Nc);CHKERRQ(ierr); 149 ierr = PetscViewerASCIIPrintf(viewer, " FEM");CHKERRQ(ierr); 150 } else if (id == PETSCFV_CLASSID) { 151 ierr = PetscFVGetNumComponents((PetscFV) obj, &Nc);CHKERRQ(ierr); 152 ierr = PetscViewerASCIIPrintf(viewer, " FVM");CHKERRQ(ierr); 153 } 154 else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f); 155 if (Nc > 1) {ierr = PetscViewerASCIIPrintf(viewer, "%d components", Nc);CHKERRQ(ierr);} 156 else {ierr = PetscViewerASCIIPrintf(viewer, "%d component ", Nc);CHKERRQ(ierr);} 157 if (prob->implicit[f]) {ierr = PetscViewerASCIIPrintf(viewer, " (implicit)");CHKERRQ(ierr);} 158 else {ierr = PetscViewerASCIIPrintf(viewer, " (explicit)");CHKERRQ(ierr);} 159 if (prob->adjacency[f*2+0]) { 160 if (prob->adjacency[f*2+1]) {ierr = PetscViewerASCIIPrintf(viewer, " (adj FVM++)");CHKERRQ(ierr);} 161 else {ierr = PetscViewerASCIIPrintf(viewer, " (adj FVM)");CHKERRQ(ierr);} 162 } else { 163 if (prob->adjacency[f*2+1]) {ierr = PetscViewerASCIIPrintf(viewer, " (adj FEM)");CHKERRQ(ierr);} 164 else {ierr = PetscViewerASCIIPrintf(viewer, " (adj FUNKY)");CHKERRQ(ierr);} 165 } 166 ierr = PetscViewerASCIIPrintf(viewer, "\n");CHKERRQ(ierr); 167 if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) { 168 if (id == PETSCFE_CLASSID) {ierr = PetscFEView((PetscFE) obj, viewer);CHKERRQ(ierr);} 169 else if (id == PETSCFV_CLASSID) {ierr = PetscFVView((PetscFV) obj, viewer);CHKERRQ(ierr);} 170 } 171 } 172 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 173 PetscFunctionReturn(0); 174 } 175 176 #undef __FUNCT__ 177 #define __FUNCT__ "PetscDSView" 178 /*@C 179 PetscDSView - Views a PetscDS 180 181 Collective on PetscDS 182 183 Input Parameter: 184 + prob - the PetscDS object to view 185 - v - the viewer 186 187 Level: developer 188 189 .seealso PetscDSDestroy() 190 @*/ 191 PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v) 192 { 193 PetscBool iascii; 194 PetscErrorCode ierr; 195 196 PetscFunctionBegin; 197 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 198 if (!v) {ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) prob), &v);CHKERRQ(ierr);} 199 else {PetscValidHeaderSpecific(v, PETSC_VIEWER_CLASSID, 2);} 200 ierr = PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 201 if (iascii) {ierr = PetscDSView_Ascii(prob, v);CHKERRQ(ierr);} 202 if (prob->ops->view) {ierr = (*prob->ops->view)(prob, v);CHKERRQ(ierr);} 203 PetscFunctionReturn(0); 204 } 205 206 #undef __FUNCT__ 207 #define __FUNCT__ "PetscDSSetFromOptions" 208 /*@ 209 PetscDSSetFromOptions - sets parameters in a PetscDS from the options database 210 211 Collective on PetscDS 212 213 Input Parameter: 214 . prob - the PetscDS object to set options for 215 216 Options Database: 217 218 Level: developer 219 220 .seealso PetscDSView() 221 @*/ 222 PetscErrorCode PetscDSSetFromOptions(PetscDS prob) 223 { 224 const char *defaultType; 225 char name[256]; 226 PetscBool flg; 227 PetscErrorCode ierr; 228 229 PetscFunctionBegin; 230 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 231 if (!((PetscObject) prob)->type_name) { 232 defaultType = PETSCDSBASIC; 233 } else { 234 defaultType = ((PetscObject) prob)->type_name; 235 } 236 ierr = PetscDSRegisterAll();CHKERRQ(ierr); 237 238 ierr = PetscObjectOptionsBegin((PetscObject) prob);CHKERRQ(ierr); 239 ierr = PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg);CHKERRQ(ierr); 240 if (flg) { 241 ierr = PetscDSSetType(prob, name);CHKERRQ(ierr); 242 } else if (!((PetscObject) prob)->type_name) { 243 ierr = PetscDSSetType(prob, defaultType);CHKERRQ(ierr); 244 } 245 if (prob->ops->setfromoptions) {ierr = (*prob->ops->setfromoptions)(prob);CHKERRQ(ierr);} 246 /* process any options handlers added with PetscObjectAddOptionsHandler() */ 247 ierr = PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject) prob);CHKERRQ(ierr); 248 ierr = PetscOptionsEnd();CHKERRQ(ierr); 249 ierr = PetscDSViewFromOptions(prob, NULL, "-petscds_view");CHKERRQ(ierr); 250 PetscFunctionReturn(0); 251 } 252 253 #undef __FUNCT__ 254 #define __FUNCT__ "PetscDSSetUp" 255 /*@C 256 PetscDSSetUp - Construct data structures for the PetscDS 257 258 Collective on PetscDS 259 260 Input Parameter: 261 . prob - the PetscDS object to setup 262 263 Level: developer 264 265 .seealso PetscDSView(), PetscDSDestroy() 266 @*/ 267 PetscErrorCode PetscDSSetUp(PetscDS prob) 268 { 269 const PetscInt Nf = prob->Nf; 270 PetscInt dim, work, NcMax = 0, NqMax = 0, f; 271 PetscErrorCode ierr; 272 273 PetscFunctionBegin; 274 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 275 if (prob->setup) PetscFunctionReturn(0); 276 /* Calculate sizes */ 277 ierr = PetscDSGetSpatialDimension(prob, &dim);CHKERRQ(ierr); 278 prob->totDim = prob->totDimBd = prob->totComp = 0; 279 ierr = PetscCalloc4(Nf+1,&prob->off,Nf+1,&prob->offDer,Nf+1,&prob->offBd,Nf+1,&prob->offDerBd);CHKERRQ(ierr); 280 ierr = PetscMalloc4(Nf,&prob->basis,Nf,&prob->basisDer,Nf,&prob->basisBd,Nf,&prob->basisDerBd);CHKERRQ(ierr); 281 for (f = 0; f < Nf; ++f) { 282 PetscFE feBd = (PetscFE) prob->discBd[f]; 283 PetscObject obj; 284 PetscClassId id; 285 PetscQuadrature q; 286 PetscInt Nq = 0, Nb, Nc; 287 288 ierr = PetscDSGetDiscretization(prob, f, &obj);CHKERRQ(ierr); 289 ierr = PetscObjectGetClassId(obj, &id);CHKERRQ(ierr); 290 if (id == PETSCFE_CLASSID) { 291 PetscFE fe = (PetscFE) obj; 292 293 ierr = PetscFEGetQuadrature(fe, &q);CHKERRQ(ierr); 294 ierr = PetscFEGetDimension(fe, &Nb);CHKERRQ(ierr); 295 ierr = PetscFEGetNumComponents(fe, &Nc);CHKERRQ(ierr); 296 ierr = PetscFEGetDefaultTabulation(fe, &prob->basis[f], &prob->basisDer[f], NULL);CHKERRQ(ierr); 297 } else if (id == PETSCFV_CLASSID) { 298 PetscFV fv = (PetscFV) obj; 299 300 ierr = PetscFVGetQuadrature(fv, &q);CHKERRQ(ierr); 301 Nb = 1; 302 ierr = PetscFVGetNumComponents(fv, &Nc);CHKERRQ(ierr); 303 ierr = PetscFVGetDefaultTabulation(fv, &prob->basis[f], &prob->basisDer[f], NULL);CHKERRQ(ierr); 304 } else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f); 305 prob->off[f+1] = Nc + prob->off[f]; 306 prob->offDer[f+1] = Nc*dim + prob->offDer[f]; 307 if (q) {ierr = PetscQuadratureGetData(q, NULL, &Nq, NULL, NULL);CHKERRQ(ierr);} 308 NqMax = PetscMax(NqMax, Nq); 309 NcMax = PetscMax(NcMax, Nc); 310 prob->totDim += Nb*Nc; 311 prob->totComp += Nc; 312 if (feBd) { 313 ierr = PetscFEGetDimension(feBd, &Nb);CHKERRQ(ierr); 314 ierr = PetscFEGetNumComponents(feBd, &Nc);CHKERRQ(ierr); 315 ierr = PetscFEGetDefaultTabulation(feBd, &prob->basisBd[f], &prob->basisDerBd[f], NULL);CHKERRQ(ierr); 316 prob->totDimBd += Nb*Nc; 317 prob->offBd[f+1] = Nc + prob->offBd[f]; 318 prob->offDerBd[f+1] = Nc*dim + prob->offDerBd[f]; 319 } 320 } 321 work = PetscMax(prob->totComp*dim, PetscSqr(NcMax*dim)); 322 /* Allocate works space */ 323 ierr = PetscMalloc5(prob->totComp,&prob->u,prob->totComp,&prob->u_t,prob->totComp*dim,&prob->u_x,dim,&prob->x,work,&prob->refSpaceDer);CHKERRQ(ierr); 324 ierr = PetscMalloc6(NqMax*NcMax,&prob->f0,NqMax*NcMax*dim,&prob->f1,NqMax*NcMax*NcMax,&prob->g0,NqMax*NcMax*NcMax*dim,&prob->g1,NqMax*NcMax*NcMax*dim,&prob->g2,NqMax*NcMax*NcMax*dim*dim,&prob->g3);CHKERRQ(ierr); 325 if (prob->ops->setup) {ierr = (*prob->ops->setup)(prob);CHKERRQ(ierr);} 326 prob->setup = PETSC_TRUE; 327 PetscFunctionReturn(0); 328 } 329 330 #undef __FUNCT__ 331 #define __FUNCT__ "PetscDSDestroyStructs_Static" 332 static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob) 333 { 334 PetscErrorCode ierr; 335 336 PetscFunctionBegin; 337 ierr = PetscFree4(prob->off,prob->offDer,prob->offBd,prob->offDerBd);CHKERRQ(ierr); 338 ierr = PetscFree4(prob->basis,prob->basisDer,prob->basisBd,prob->basisDerBd);CHKERRQ(ierr); 339 ierr = PetscFree5(prob->u,prob->u_t,prob->u_x,prob->x,prob->refSpaceDer);CHKERRQ(ierr); 340 ierr = PetscFree6(prob->f0,prob->f1,prob->g0,prob->g1,prob->g2,prob->g3);CHKERRQ(ierr); 341 PetscFunctionReturn(0); 342 } 343 344 #undef __FUNCT__ 345 #define __FUNCT__ "PetscDSEnlarge_Static" 346 static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew) 347 { 348 PetscObject *tmpd, *tmpdbd; 349 PetscBool *tmpi, *tmpa; 350 PetscPointFunc *tmpobj, *tmpf; 351 PetscPointJac *tmpg, *tmpgp; 352 PetscBdPointFunc *tmpfbd; 353 PetscBdPointJac *tmpgbd; 354 PetscRiemannFunc *tmpr; 355 void **tmpctx; 356 PetscInt Nf = prob->Nf, f, i; 357 PetscErrorCode ierr; 358 359 PetscFunctionBegin; 360 if (Nf >= NfNew) PetscFunctionReturn(0); 361 prob->setup = PETSC_FALSE; 362 ierr = PetscDSDestroyStructs_Static(prob);CHKERRQ(ierr); 363 ierr = PetscMalloc4(NfNew, &tmpd, NfNew, &tmpdbd, NfNew, &tmpi, NfNew*2, &tmpa);CHKERRQ(ierr); 364 for (f = 0; f < Nf; ++f) {tmpd[f] = prob->disc[f]; tmpdbd[f] = prob->discBd[f]; tmpi[f] = prob->implicit[f]; for (i = 0; i < 2; ++i) tmpa[f*2+i] = prob->adjacency[f*2+i];} 365 for (f = Nf; f < NfNew; ++f) {ierr = PetscContainerCreate(PetscObjectComm((PetscObject) prob), (PetscContainer *) &tmpd[f]);CHKERRQ(ierr); 366 tmpdbd[f] = NULL; tmpi[f] = PETSC_TRUE; tmpa[f*2+0] = PETSC_FALSE; tmpa[f*2+1] = PETSC_TRUE;} 367 ierr = PetscFree4(prob->disc, prob->discBd, prob->implicit, prob->adjacency);CHKERRQ(ierr); 368 prob->Nf = NfNew; 369 prob->disc = tmpd; 370 prob->discBd = tmpdbd; 371 prob->implicit = tmpi; 372 prob->adjacency = tmpa; 373 ierr = PetscCalloc6(NfNew, &tmpobj, NfNew*2, &tmpf, NfNew*NfNew*4, &tmpg, NfNew*NfNew*4, &tmpgp, NfNew, &tmpr, NfNew, &tmpctx);CHKERRQ(ierr); 374 for (f = 0; f < Nf; ++f) tmpobj[f] = prob->obj[f]; 375 for (f = 0; f < Nf*2; ++f) tmpf[f] = prob->f[f]; 376 for (f = 0; f < Nf*Nf*4; ++f) tmpg[f] = prob->g[f]; 377 for (f = 0; f < Nf*Nf*4; ++f) tmpgp[f] = prob->gp[f]; 378 for (f = 0; f < Nf; ++f) tmpr[f] = prob->r[f]; 379 for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f]; 380 for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL; 381 for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL; 382 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL; 383 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgp[f] = NULL; 384 for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL; 385 for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL; 386 ierr = PetscFree6(prob->obj, prob->f, prob->g, prob->gp, prob->r, prob->ctx);CHKERRQ(ierr); 387 prob->obj = tmpobj; 388 prob->f = tmpf; 389 prob->g = tmpg; 390 prob->gp = tmpgp; 391 prob->r = tmpr; 392 prob->ctx = tmpctx; 393 ierr = PetscCalloc2(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd);CHKERRQ(ierr); 394 for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f]; 395 for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f]; 396 for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL; 397 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL; 398 ierr = PetscFree2(prob->fBd, prob->gBd);CHKERRQ(ierr); 399 prob->fBd = tmpfbd; 400 prob->gBd = tmpgbd; 401 PetscFunctionReturn(0); 402 } 403 404 #undef __FUNCT__ 405 #define __FUNCT__ "PetscDSDestroy" 406 /*@ 407 PetscDSDestroy - Destroys a PetscDS object 408 409 Collective on PetscDS 410 411 Input Parameter: 412 . prob - the PetscDS object to destroy 413 414 Level: developer 415 416 .seealso PetscDSView() 417 @*/ 418 PetscErrorCode PetscDSDestroy(PetscDS *prob) 419 { 420 PetscInt f; 421 PetscErrorCode ierr; 422 423 PetscFunctionBegin; 424 if (!*prob) PetscFunctionReturn(0); 425 PetscValidHeaderSpecific((*prob), PETSCDS_CLASSID, 1); 426 427 if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; PetscFunctionReturn(0);} 428 ((PetscObject) (*prob))->refct = 0; 429 ierr = PetscDSDestroyStructs_Static(*prob);CHKERRQ(ierr); 430 for (f = 0; f < (*prob)->Nf; ++f) { 431 ierr = PetscObjectDereference((*prob)->disc[f]);CHKERRQ(ierr); 432 ierr = PetscObjectDereference((*prob)->discBd[f]);CHKERRQ(ierr); 433 } 434 ierr = PetscFree4((*prob)->disc, (*prob)->discBd, (*prob)->implicit, (*prob)->adjacency);CHKERRQ(ierr); 435 ierr = PetscFree6((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->gp,(*prob)->r,(*prob)->ctx);CHKERRQ(ierr); 436 ierr = PetscFree2((*prob)->fBd,(*prob)->gBd);CHKERRQ(ierr); 437 if ((*prob)->ops->destroy) {ierr = (*(*prob)->ops->destroy)(*prob);CHKERRQ(ierr);} 438 ierr = PetscHeaderDestroy(prob);CHKERRQ(ierr); 439 PetscFunctionReturn(0); 440 } 441 442 #undef __FUNCT__ 443 #define __FUNCT__ "PetscDSCreate" 444 /*@ 445 PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType(). 446 447 Collective on MPI_Comm 448 449 Input Parameter: 450 . comm - The communicator for the PetscDS object 451 452 Output Parameter: 453 . prob - The PetscDS object 454 455 Level: beginner 456 457 .seealso: PetscDSSetType(), PETSCDSBASIC 458 @*/ 459 PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob) 460 { 461 PetscDS p; 462 PetscErrorCode ierr; 463 464 PetscFunctionBegin; 465 PetscValidPointer(prob, 2); 466 *prob = NULL; 467 ierr = PetscDSInitializePackage();CHKERRQ(ierr); 468 469 ierr = PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);CHKERRQ(ierr); 470 471 p->Nf = 0; 472 p->setup = PETSC_FALSE; 473 474 *prob = p; 475 PetscFunctionReturn(0); 476 } 477 478 #undef __FUNCT__ 479 #define __FUNCT__ "PetscDSGetNumFields" 480 /*@ 481 PetscDSGetNumFields - Returns the number of fields in the DS 482 483 Not collective 484 485 Input Parameter: 486 . prob - The PetscDS object 487 488 Output Parameter: 489 . Nf - The number of fields 490 491 Level: beginner 492 493 .seealso: PetscDSGetSpatialDimension(), PetscDSCreate() 494 @*/ 495 PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf) 496 { 497 PetscFunctionBegin; 498 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 499 PetscValidPointer(Nf, 2); 500 *Nf = prob->Nf; 501 PetscFunctionReturn(0); 502 } 503 504 #undef __FUNCT__ 505 #define __FUNCT__ "PetscDSGetSpatialDimension" 506 /*@ 507 PetscDSGetSpatialDimension - Returns the spatial dimension of the DS 508 509 Not collective 510 511 Input Parameter: 512 . prob - The PetscDS object 513 514 Output Parameter: 515 . dim - The spatial dimension 516 517 Level: beginner 518 519 .seealso: PetscDSGetNumFields(), PetscDSCreate() 520 @*/ 521 PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim) 522 { 523 PetscErrorCode ierr; 524 525 PetscFunctionBegin; 526 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 527 PetscValidPointer(dim, 2); 528 *dim = 0; 529 if (prob->Nf) { 530 PetscObject obj; 531 PetscClassId id; 532 533 ierr = PetscDSGetDiscretization(prob, 0, &obj);CHKERRQ(ierr); 534 ierr = PetscObjectGetClassId(obj, &id);CHKERRQ(ierr); 535 if (id == PETSCFE_CLASSID) {ierr = PetscFEGetSpatialDimension((PetscFE) obj, dim);CHKERRQ(ierr);} 536 else if (id == PETSCFV_CLASSID) {ierr = PetscFVGetSpatialDimension((PetscFV) obj, dim);CHKERRQ(ierr);} 537 else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0); 538 } 539 PetscFunctionReturn(0); 540 } 541 542 #undef __FUNCT__ 543 #define __FUNCT__ "PetscDSGetTotalDimension" 544 /*@ 545 PetscDSGetTotalDimension - Returns the total size of the approximation space for this system 546 547 Not collective 548 549 Input Parameter: 550 . prob - The PetscDS object 551 552 Output Parameter: 553 . dim - The total problem dimension 554 555 Level: beginner 556 557 .seealso: PetscDSGetNumFields(), PetscDSCreate() 558 @*/ 559 PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim) 560 { 561 PetscErrorCode ierr; 562 563 PetscFunctionBegin; 564 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 565 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 566 PetscValidPointer(dim, 2); 567 *dim = prob->totDim; 568 PetscFunctionReturn(0); 569 } 570 571 #undef __FUNCT__ 572 #define __FUNCT__ "PetscDSGetTotalBdDimension" 573 /*@ 574 PetscDSGetTotalBdDimension - Returns the total size of the boundary approximation space for this system 575 576 Not collective 577 578 Input Parameter: 579 . prob - The PetscDS object 580 581 Output Parameter: 582 . dim - The total boundary problem dimension 583 584 Level: beginner 585 586 .seealso: PetscDSGetNumFields(), PetscDSCreate() 587 @*/ 588 PetscErrorCode PetscDSGetTotalBdDimension(PetscDS prob, PetscInt *dim) 589 { 590 PetscErrorCode ierr; 591 592 PetscFunctionBegin; 593 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 594 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 595 PetscValidPointer(dim, 2); 596 *dim = prob->totDimBd; 597 PetscFunctionReturn(0); 598 } 599 600 #undef __FUNCT__ 601 #define __FUNCT__ "PetscDSGetTotalComponents" 602 /*@ 603 PetscDSGetTotalComponents - Returns the total number of components in this system 604 605 Not collective 606 607 Input Parameter: 608 . prob - The PetscDS object 609 610 Output Parameter: 611 . dim - The total number of components 612 613 Level: beginner 614 615 .seealso: PetscDSGetNumFields(), PetscDSCreate() 616 @*/ 617 PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc) 618 { 619 PetscErrorCode ierr; 620 621 PetscFunctionBegin; 622 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 623 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 624 PetscValidPointer(Nc, 2); 625 *Nc = prob->totComp; 626 PetscFunctionReturn(0); 627 } 628 629 #undef __FUNCT__ 630 #define __FUNCT__ "PetscDSGetDiscretization" 631 /*@ 632 PetscDSGetDiscretization - Returns the discretization object for the given field 633 634 Not collective 635 636 Input Parameters: 637 + prob - The PetscDS object 638 - f - The field number 639 640 Output Parameter: 641 . disc - The discretization object 642 643 Level: beginner 644 645 .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 646 @*/ 647 PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc) 648 { 649 PetscFunctionBegin; 650 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 651 PetscValidPointer(disc, 3); 652 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); 653 *disc = prob->disc[f]; 654 PetscFunctionReturn(0); 655 } 656 657 #undef __FUNCT__ 658 #define __FUNCT__ "PetscDSGetBdDiscretization" 659 /*@ 660 PetscDSGetBdDiscretization - Returns the boundary discretization object for the given field 661 662 Not collective 663 664 Input Parameters: 665 + prob - The PetscDS object 666 - f - The field number 667 668 Output Parameter: 669 . disc - The boundary discretization object 670 671 Level: beginner 672 673 .seealso: PetscDSSetBdDiscretization(), PetscDSAddBdDiscretization(), PetscDSGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 674 @*/ 675 PetscErrorCode PetscDSGetBdDiscretization(PetscDS prob, PetscInt f, PetscObject *disc) 676 { 677 PetscFunctionBegin; 678 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 679 PetscValidPointer(disc, 3); 680 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); 681 *disc = prob->discBd[f]; 682 PetscFunctionReturn(0); 683 } 684 685 #undef __FUNCT__ 686 #define __FUNCT__ "PetscDSSetDiscretization" 687 /*@ 688 PetscDSSetDiscretization - Sets the discretization object for the given field 689 690 Not collective 691 692 Input Parameters: 693 + prob - The PetscDS object 694 . f - The field number 695 - disc - The discretization object 696 697 Level: beginner 698 699 .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 700 @*/ 701 PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc) 702 { 703 PetscErrorCode ierr; 704 705 PetscFunctionBegin; 706 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 707 PetscValidPointer(disc, 3); 708 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 709 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 710 if (prob->disc[f]) {ierr = PetscObjectDereference(prob->disc[f]);CHKERRQ(ierr);} 711 prob->disc[f] = disc; 712 ierr = PetscObjectReference(disc);CHKERRQ(ierr); 713 { 714 PetscClassId id; 715 716 ierr = PetscObjectGetClassId(disc, &id);CHKERRQ(ierr); 717 if (id == PETSCFV_CLASSID) { 718 prob->implicit[f] = PETSC_FALSE; 719 prob->adjacency[f*2+0] = PETSC_TRUE; 720 prob->adjacency[f*2+1] = PETSC_FALSE; 721 } 722 } 723 PetscFunctionReturn(0); 724 } 725 726 #undef __FUNCT__ 727 #define __FUNCT__ "PetscDSSetBdDiscretization" 728 /*@ 729 PetscDSSetBdDiscretization - Sets the boundary discretization object for the given field 730 731 Not collective 732 733 Input Parameters: 734 + prob - The PetscDS object 735 . f - The field number 736 - disc - The boundary discretization object 737 738 Level: beginner 739 740 .seealso: PetscDSGetBdDiscretization(), PetscDSAddBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 741 @*/ 742 PetscErrorCode PetscDSSetBdDiscretization(PetscDS prob, PetscInt f, PetscObject disc) 743 { 744 PetscErrorCode ierr; 745 746 PetscFunctionBegin; 747 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 748 if (disc) PetscValidPointer(disc, 3); 749 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 750 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 751 if (prob->discBd[f]) {ierr = PetscObjectDereference(prob->discBd[f]);CHKERRQ(ierr);} 752 prob->discBd[f] = disc; 753 ierr = PetscObjectReference(disc);CHKERRQ(ierr); 754 PetscFunctionReturn(0); 755 } 756 757 #undef __FUNCT__ 758 #define __FUNCT__ "PetscDSAddDiscretization" 759 /*@ 760 PetscDSAddDiscretization - Adds a discretization object 761 762 Not collective 763 764 Input Parameters: 765 + prob - The PetscDS object 766 - disc - The boundary discretization object 767 768 Level: beginner 769 770 .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 771 @*/ 772 PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc) 773 { 774 PetscErrorCode ierr; 775 776 PetscFunctionBegin; 777 ierr = PetscDSSetDiscretization(prob, prob->Nf, disc);CHKERRQ(ierr); 778 PetscFunctionReturn(0); 779 } 780 781 #undef __FUNCT__ 782 #define __FUNCT__ "PetscDSAddBdDiscretization" 783 /*@ 784 PetscDSAddBdDiscretization - Adds a boundary discretization object 785 786 Not collective 787 788 Input Parameters: 789 + prob - The PetscDS object 790 - disc - The boundary discretization object 791 792 Level: beginner 793 794 .seealso: PetscDSGetBdDiscretization(), PetscDSSetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 795 @*/ 796 PetscErrorCode PetscDSAddBdDiscretization(PetscDS prob, PetscObject disc) 797 { 798 PetscErrorCode ierr; 799 800 PetscFunctionBegin; 801 ierr = PetscDSSetBdDiscretization(prob, prob->Nf, disc);CHKERRQ(ierr); 802 PetscFunctionReturn(0); 803 } 804 805 #undef __FUNCT__ 806 #define __FUNCT__ "PetscDSGetImplicit" 807 /*@ 808 PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX 809 810 Not collective 811 812 Input Parameters: 813 + prob - The PetscDS object 814 - f - The field number 815 816 Output Parameter: 817 . implicit - The flag indicating what kind of solve to use for this field 818 819 Level: developer 820 821 .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 822 @*/ 823 PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit) 824 { 825 PetscFunctionBegin; 826 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 827 PetscValidPointer(implicit, 3); 828 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); 829 *implicit = prob->implicit[f]; 830 PetscFunctionReturn(0); 831 } 832 833 #undef __FUNCT__ 834 #define __FUNCT__ "PetscDSSetImplicit" 835 /*@ 836 PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for IMEX 837 838 Not collective 839 840 Input Parameters: 841 + prob - The PetscDS object 842 . f - The field number 843 - implicit - The flag indicating what kind of solve to use for this field 844 845 Level: developer 846 847 .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 848 @*/ 849 PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit) 850 { 851 PetscFunctionBegin; 852 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 853 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); 854 prob->implicit[f] = implicit; 855 PetscFunctionReturn(0); 856 } 857 858 #undef __FUNCT__ 859 #define __FUNCT__ "PetscDSGetAdjacency" 860 /*@ 861 PetscDSGetAdjacency - Returns the flags for determining variable influence 862 863 Not collective 864 865 Input Parameters: 866 + prob - The PetscDS object 867 - f - The field number 868 869 Output Parameter: 870 + useCone - Flag for variable influence starting with the cone operation 871 - useClosure - Flag for variable influence using transitive closure 872 873 Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure() 874 875 Level: developer 876 877 .seealso: PetscDSSetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 878 @*/ 879 PetscErrorCode PetscDSGetAdjacency(PetscDS prob, PetscInt f, PetscBool *useCone, PetscBool *useClosure) 880 { 881 PetscFunctionBegin; 882 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 883 PetscValidPointer(useCone, 3); 884 PetscValidPointer(useClosure, 4); 885 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); 886 *useCone = prob->adjacency[f*2+0]; 887 *useClosure = prob->adjacency[f*2+1]; 888 PetscFunctionReturn(0); 889 } 890 891 #undef __FUNCT__ 892 #define __FUNCT__ "PetscDSSetAdjacency" 893 /*@ 894 PetscDSSetAdjacency - Set the flags for determining variable influence 895 896 Not collective 897 898 Input Parameters: 899 + prob - The PetscDS object 900 . f - The field number 901 . useCone - Flag for variable influence starting with the cone operation 902 - useClosure - Flag for variable influence using transitive closure 903 904 Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure() 905 906 Level: developer 907 908 .seealso: PetscDSGetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 909 @*/ 910 PetscErrorCode PetscDSSetAdjacency(PetscDS prob, PetscInt f, PetscBool useCone, PetscBool useClosure) 911 { 912 PetscFunctionBegin; 913 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 914 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); 915 prob->adjacency[f*2+0] = useCone; 916 prob->adjacency[f*2+1] = useClosure; 917 PetscFunctionReturn(0); 918 } 919 920 #undef __FUNCT__ 921 #define __FUNCT__ "PetscDSGetObjective" 922 PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f, 923 void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 924 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 925 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 926 PetscReal t, const PetscReal x[], PetscScalar obj[])) 927 { 928 PetscFunctionBegin; 929 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 930 PetscValidPointer(obj, 2); 931 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); 932 *obj = prob->obj[f]; 933 PetscFunctionReturn(0); 934 } 935 936 #undef __FUNCT__ 937 #define __FUNCT__ "PetscDSSetObjective" 938 PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f, 939 void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 940 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 941 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 942 PetscReal t, const PetscReal x[], PetscScalar obj[])) 943 { 944 PetscErrorCode ierr; 945 946 PetscFunctionBegin; 947 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 948 PetscValidFunction(obj, 2); 949 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 950 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 951 prob->obj[f] = obj; 952 PetscFunctionReturn(0); 953 } 954 955 #undef __FUNCT__ 956 #define __FUNCT__ "PetscDSGetResidual" 957 /*@C 958 PetscDSGetResidual - Get the pointwise residual function for a given test field 959 960 Not collective 961 962 Input Parameters: 963 + prob - The PetscDS 964 - f - The test field number 965 966 Output Parameters: 967 + f0 - integrand for the test function term 968 - f1 - integrand for the test function gradient term 969 970 Note: We are using a first order FEM model for the weak form: 971 972 \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) 973 974 The calling sequence for the callbacks f0 and f1 is given by: 975 976 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 977 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 978 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 979 $ PetscReal t, const PetscReal x[], PetscScalar f0[]) 980 981 + dim - the spatial dimension 982 . Nf - the number of fields 983 . uOff - the offset into u[] and u_t[] for each field 984 . uOff_x - the offset into u_x[] for each field 985 . u - each field evaluated at the current point 986 . u_t - the time derivative of each field evaluated at the current point 987 . u_x - the gradient of each field evaluated at the current point 988 . aOff - the offset into a[] and a_t[] for each auxiliary field 989 . aOff_x - the offset into a_x[] for each auxiliary field 990 . a - each auxiliary field evaluated at the current point 991 . a_t - the time derivative of each auxiliary field evaluated at the current point 992 . a_x - the gradient of auxiliary each field evaluated at the current point 993 . t - current time 994 . x - coordinates of the current point 995 - f0 - output values at the current point 996 997 Level: intermediate 998 999 .seealso: PetscDSSetResidual() 1000 @*/ 1001 PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f, 1002 void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1003 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1004 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1005 PetscReal t, const PetscReal x[], PetscScalar f0[]), 1006 void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1007 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1008 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1009 PetscReal t, const PetscReal x[], PetscScalar f1[])) 1010 { 1011 PetscFunctionBegin; 1012 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1013 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); 1014 if (f0) {PetscValidPointer(f0, 3); *f0 = prob->f[f*2+0];} 1015 if (f1) {PetscValidPointer(f1, 4); *f1 = prob->f[f*2+1];} 1016 PetscFunctionReturn(0); 1017 } 1018 1019 #undef __FUNCT__ 1020 #define __FUNCT__ "PetscDSSetResidual" 1021 /*@C 1022 PetscDSSetResidual - Set the pointwise residual function for a given test field 1023 1024 Not collective 1025 1026 Input Parameters: 1027 + prob - The PetscDS 1028 . f - The test field number 1029 . f0 - integrand for the test function term 1030 - f1 - integrand for the test function gradient term 1031 1032 Note: We are using a first order FEM model for the weak form: 1033 1034 \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) 1035 1036 The calling sequence for the callbacks f0 and f1 is given by: 1037 1038 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1039 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1040 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1041 $ PetscReal t, const PetscReal x[], PetscScalar f0[]) 1042 1043 + dim - the spatial dimension 1044 . Nf - the number of fields 1045 . uOff - the offset into u[] and u_t[] for each field 1046 . uOff_x - the offset into u_x[] for each field 1047 . u - each field evaluated at the current point 1048 . u_t - the time derivative of each field evaluated at the current point 1049 . u_x - the gradient of each field evaluated at the current point 1050 . aOff - the offset into a[] and a_t[] for each auxiliary field 1051 . aOff_x - the offset into a_x[] for each auxiliary field 1052 . a - each auxiliary field evaluated at the current point 1053 . a_t - the time derivative of each auxiliary field evaluated at the current point 1054 . a_x - the gradient of auxiliary each field evaluated at the current point 1055 . t - current time 1056 . x - coordinates of the current point 1057 - f0 - output values at the current point 1058 1059 Level: intermediate 1060 1061 .seealso: PetscDSGetResidual() 1062 @*/ 1063 PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f, 1064 void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1065 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1066 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1067 PetscReal t, const PetscReal x[], PetscScalar f0[]), 1068 void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1069 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1070 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1071 PetscReal t, const PetscReal x[], PetscScalar f1[])) 1072 { 1073 PetscErrorCode ierr; 1074 1075 PetscFunctionBegin; 1076 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1077 if (f0) PetscValidFunction(f0, 3); 1078 if (f1) PetscValidFunction(f1, 4); 1079 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1080 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1081 prob->f[f*2+0] = f0; 1082 prob->f[f*2+1] = f1; 1083 PetscFunctionReturn(0); 1084 } 1085 1086 #undef __FUNCT__ 1087 #define __FUNCT__ "PetscDSGetJacobian" 1088 /*@C 1089 PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field 1090 1091 Not collective 1092 1093 Input Parameters: 1094 + prob - The PetscDS 1095 . f - The test field number 1096 - g - The field number 1097 1098 Output Parameters: 1099 + g0 - integrand for the test and basis function term 1100 . g1 - integrand for the test function and basis function gradient term 1101 . g2 - integrand for the test function gradient and basis function term 1102 - g3 - integrand for the test function gradient and basis function gradient term 1103 1104 Note: We are using a first order FEM model for the weak form: 1105 1106 \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 1107 1108 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1109 1110 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1111 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1112 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1113 $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]) 1114 1115 + dim - the spatial dimension 1116 . Nf - the number of fields 1117 . uOff - the offset into u[] and u_t[] for each field 1118 . uOff_x - the offset into u_x[] for each field 1119 . u - each field evaluated at the current point 1120 . u_t - the time derivative of each field evaluated at the current point 1121 . u_x - the gradient of each field evaluated at the current point 1122 . aOff - the offset into a[] and a_t[] for each auxiliary field 1123 . aOff_x - the offset into a_x[] for each auxiliary field 1124 . a - each auxiliary field evaluated at the current point 1125 . a_t - the time derivative of each auxiliary field evaluated at the current point 1126 . a_x - the gradient of auxiliary each field evaluated at the current point 1127 . t - current time 1128 . u_tShift - the multiplier a for dF/dU_t 1129 . x - coordinates of the current point 1130 - g0 - output values at the current point 1131 1132 Level: intermediate 1133 1134 .seealso: PetscDSSetJacobian() 1135 @*/ 1136 PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g, 1137 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1138 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1139 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1140 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]), 1141 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1142 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1143 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1144 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g1[]), 1145 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1146 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1147 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1148 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g2[]), 1149 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1150 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1151 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1152 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g3[])) 1153 { 1154 PetscFunctionBegin; 1155 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1156 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); 1157 if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf); 1158 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->g[(f*prob->Nf + g)*4+0];} 1159 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->g[(f*prob->Nf + g)*4+1];} 1160 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->g[(f*prob->Nf + g)*4+2];} 1161 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->g[(f*prob->Nf + g)*4+3];} 1162 PetscFunctionReturn(0); 1163 } 1164 1165 #undef __FUNCT__ 1166 #define __FUNCT__ "PetscDSSetJacobian" 1167 /*@C 1168 PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields 1169 1170 Not collective 1171 1172 Input Parameters: 1173 + prob - The PetscDS 1174 . f - The test field number 1175 . g - The field number 1176 . g0 - integrand for the test and basis function term 1177 . g1 - integrand for the test function and basis function gradient term 1178 . g2 - integrand for the test function gradient and basis function term 1179 - g3 - integrand for the test function gradient and basis function gradient term 1180 1181 Note: We are using a first order FEM model for the weak form: 1182 1183 \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 1184 1185 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1186 1187 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1188 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1189 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1190 $ PetscReal t, const PetscReal x[], PetscScalar g0[]) 1191 1192 + dim - the spatial dimension 1193 . Nf - the number of fields 1194 . uOff - the offset into u[] and u_t[] for each field 1195 . uOff_x - the offset into u_x[] for each field 1196 . u - each field evaluated at the current point 1197 . u_t - the time derivative of each field evaluated at the current point 1198 . u_x - the gradient of each field evaluated at the current point 1199 . aOff - the offset into a[] and a_t[] for each auxiliary field 1200 . aOff_x - the offset into a_x[] for each auxiliary field 1201 . a - each auxiliary field evaluated at the current point 1202 . a_t - the time derivative of each auxiliary field evaluated at the current point 1203 . a_x - the gradient of auxiliary each field evaluated at the current point 1204 . t - current time 1205 . u_tShift - the multiplier a for dF/dU_t 1206 . x - coordinates of the current point 1207 - g0 - output values at the current point 1208 1209 Level: intermediate 1210 1211 .seealso: PetscDSGetJacobian() 1212 @*/ 1213 PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g, 1214 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1215 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1216 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1217 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]), 1218 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1219 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1220 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1221 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g1[]), 1222 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1223 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1224 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1225 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g2[]), 1226 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1227 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1228 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1229 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g3[])) 1230 { 1231 PetscErrorCode ierr; 1232 1233 PetscFunctionBegin; 1234 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1235 if (g0) PetscValidFunction(g0, 4); 1236 if (g1) PetscValidFunction(g1, 5); 1237 if (g2) PetscValidFunction(g2, 6); 1238 if (g3) PetscValidFunction(g3, 7); 1239 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1240 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1241 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1242 prob->g[(f*prob->Nf + g)*4+0] = g0; 1243 prob->g[(f*prob->Nf + g)*4+1] = g1; 1244 prob->g[(f*prob->Nf + g)*4+2] = g2; 1245 prob->g[(f*prob->Nf + g)*4+3] = g3; 1246 PetscFunctionReturn(0); 1247 } 1248 1249 #undef __FUNCT__ 1250 #define __FUNCT__ "PetscDSHasJacobianPreconditioner" 1251 /*@C 1252 PetscDSHasJacobianPreconditioner - Signals that a Jacobian preconditioner matrix has been set 1253 1254 Not collective 1255 1256 Input Parameter: 1257 . prob - The PetscDS 1258 1259 Output Parameter: 1260 . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set 1261 1262 Level: intermediate 1263 1264 .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian() 1265 @*/ 1266 PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS prob, PetscBool *hasJacPre) 1267 { 1268 PetscInt f, g, h; 1269 1270 PetscFunctionBegin; 1271 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1272 *hasJacPre = PETSC_FALSE; 1273 for (f = 0; f < prob->Nf; ++f) { 1274 for (g = 0; g < prob->Nf; ++g) { 1275 for (h = 0; h < 4; ++h) { 1276 if (prob->gp[(f*prob->Nf + g)*4+h]) *hasJacPre = PETSC_TRUE; 1277 } 1278 } 1279 } 1280 PetscFunctionReturn(0); 1281 } 1282 1283 #undef __FUNCT__ 1284 #define __FUNCT__ "PetscDSGetJacobianPreconditioner" 1285 /*@C 1286 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. 1287 1288 Not collective 1289 1290 Input Parameters: 1291 + prob - The PetscDS 1292 . f - The test field number 1293 - g - The field number 1294 1295 Output Parameters: 1296 + g0 - integrand for the test and basis function term 1297 . g1 - integrand for the test function and basis function gradient term 1298 . g2 - integrand for the test function gradient and basis function term 1299 - g3 - integrand for the test function gradient and basis function gradient term 1300 1301 Note: We are using a first order FEM model for the weak form: 1302 1303 \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 1304 1305 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1306 1307 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1308 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1309 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1310 $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]) 1311 1312 + dim - the spatial dimension 1313 . Nf - the number of fields 1314 . uOff - the offset into u[] and u_t[] for each field 1315 . uOff_x - the offset into u_x[] for each field 1316 . u - each field evaluated at the current point 1317 . u_t - the time derivative of each field evaluated at the current point 1318 . u_x - the gradient of each field evaluated at the current point 1319 . aOff - the offset into a[] and a_t[] for each auxiliary field 1320 . aOff_x - the offset into a_x[] for each auxiliary field 1321 . a - each auxiliary field evaluated at the current point 1322 . a_t - the time derivative of each auxiliary field evaluated at the current point 1323 . a_x - the gradient of auxiliary each field evaluated at the current point 1324 . t - current time 1325 . u_tShift - the multiplier a for dF/dU_t 1326 . x - coordinates of the current point 1327 - g0 - output values at the current point 1328 1329 Level: intermediate 1330 1331 .seealso: PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian() 1332 @*/ 1333 PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS prob, 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[], 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[], 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[], 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[], PetscScalar g3[])) 1350 { 1351 PetscFunctionBegin; 1352 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1353 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); 1354 if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf); 1355 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->gp[(f*prob->Nf + g)*4+0];} 1356 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->gp[(f*prob->Nf + g)*4+1];} 1357 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->gp[(f*prob->Nf + g)*4+2];} 1358 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->gp[(f*prob->Nf + g)*4+3];} 1359 PetscFunctionReturn(0); 1360 } 1361 1362 #undef __FUNCT__ 1363 #define __FUNCT__ "PetscDSSetJacobianPreconditioner" 1364 /*@C 1365 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. 1366 1367 Not collective 1368 1369 Input Parameters: 1370 + prob - The PetscDS 1371 . f - The test field number 1372 . g - The field number 1373 . g0 - integrand for the test and basis function term 1374 . g1 - integrand for the test function and basis function gradient term 1375 . g2 - integrand for the test function gradient and basis function term 1376 - g3 - integrand for the test function gradient and basis function gradient term 1377 1378 Note: We are using a first order FEM model for the weak form: 1379 1380 \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 1381 1382 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1383 1384 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1385 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1386 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1387 $ PetscReal t, const PetscReal x[], PetscScalar g0[]) 1388 1389 + dim - the spatial dimension 1390 . Nf - the number of fields 1391 . uOff - the offset into u[] and u_t[] for each field 1392 . uOff_x - the offset into u_x[] for each field 1393 . u - each field evaluated at the current point 1394 . u_t - the time derivative of each field evaluated at the current point 1395 . u_x - the gradient of each field evaluated at the current point 1396 . aOff - the offset into a[] and a_t[] for each auxiliary field 1397 . aOff_x - the offset into a_x[] for each auxiliary field 1398 . a - each auxiliary field evaluated at the current point 1399 . a_t - the time derivative of each auxiliary field evaluated at the current point 1400 . a_x - the gradient of auxiliary each field evaluated at the current point 1401 . t - current time 1402 . u_tShift - the multiplier a for dF/dU_t 1403 . x - coordinates of the current point 1404 - g0 - output values at the current point 1405 1406 Level: intermediate 1407 1408 .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobian() 1409 @*/ 1410 PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g, 1411 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1412 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1413 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1414 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]), 1415 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1416 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1417 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1418 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g1[]), 1419 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1420 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1421 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1422 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g2[]), 1423 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1424 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1425 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1426 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g3[])) 1427 { 1428 PetscErrorCode ierr; 1429 1430 PetscFunctionBegin; 1431 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1432 if (g0) PetscValidFunction(g0, 4); 1433 if (g1) PetscValidFunction(g1, 5); 1434 if (g2) PetscValidFunction(g2, 6); 1435 if (g3) PetscValidFunction(g3, 7); 1436 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1437 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1438 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1439 prob->gp[(f*prob->Nf + g)*4+0] = g0; 1440 prob->gp[(f*prob->Nf + g)*4+1] = g1; 1441 prob->gp[(f*prob->Nf + g)*4+2] = g2; 1442 prob->gp[(f*prob->Nf + g)*4+3] = g3; 1443 PetscFunctionReturn(0); 1444 } 1445 1446 #undef __FUNCT__ 1447 #define __FUNCT__ "PetscDSGetRiemannSolver" 1448 /*@C 1449 PetscDSGetRiemannSolver - Returns the Riemann solver for the given field 1450 1451 Not collective 1452 1453 Input Arguments: 1454 + prob - The PetscDS object 1455 - f - The field number 1456 1457 Output Argument: 1458 . r - Riemann solver 1459 1460 Calling sequence for r: 1461 1462 $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx) 1463 1464 + dim - The spatial dimension 1465 . Nf - The number of fields 1466 . x - The coordinates at a point on the interface 1467 . n - The normal vector to the interface 1468 . uL - The state vector to the left of the interface 1469 . uR - The state vector to the right of the interface 1470 . flux - output array of flux through the interface 1471 - ctx - optional user context 1472 1473 Level: intermediate 1474 1475 .seealso: PetscDSSetRiemannSolver() 1476 @*/ 1477 PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f, 1478 void (**r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)) 1479 { 1480 PetscFunctionBegin; 1481 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1482 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); 1483 PetscValidPointer(r, 3); 1484 *r = prob->r[f]; 1485 PetscFunctionReturn(0); 1486 } 1487 1488 #undef __FUNCT__ 1489 #define __FUNCT__ "PetscDSSetRiemannSolver" 1490 /*@C 1491 PetscDSSetRiemannSolver - Sets the Riemann solver for the given field 1492 1493 Not collective 1494 1495 Input Arguments: 1496 + prob - The PetscDS object 1497 . f - The field number 1498 - r - Riemann solver 1499 1500 Calling sequence for r: 1501 1502 $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx) 1503 1504 + dim - The spatial dimension 1505 . Nf - The number of fields 1506 . x - The coordinates at a point on the interface 1507 . n - The normal vector to the interface 1508 . uL - The state vector to the left of the interface 1509 . uR - The state vector to the right of the interface 1510 . flux - output array of flux through the interface 1511 - ctx - optional user context 1512 1513 Level: intermediate 1514 1515 .seealso: PetscDSGetRiemannSolver() 1516 @*/ 1517 PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f, 1518 void (*r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)) 1519 { 1520 PetscErrorCode ierr; 1521 1522 PetscFunctionBegin; 1523 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1524 PetscValidFunction(r, 3); 1525 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1526 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1527 prob->r[f] = r; 1528 PetscFunctionReturn(0); 1529 } 1530 1531 #undef __FUNCT__ 1532 #define __FUNCT__ "PetscDSGetContext" 1533 PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx) 1534 { 1535 PetscFunctionBegin; 1536 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1537 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); 1538 PetscValidPointer(ctx, 3); 1539 *ctx = prob->ctx[f]; 1540 PetscFunctionReturn(0); 1541 } 1542 1543 #undef __FUNCT__ 1544 #define __FUNCT__ "PetscDSSetContext" 1545 PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx) 1546 { 1547 PetscErrorCode ierr; 1548 1549 PetscFunctionBegin; 1550 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1551 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1552 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1553 prob->ctx[f] = ctx; 1554 PetscFunctionReturn(0); 1555 } 1556 1557 #undef __FUNCT__ 1558 #define __FUNCT__ "PetscDSGetBdResidual" 1559 /*@C 1560 PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field 1561 1562 Not collective 1563 1564 Input Parameters: 1565 + prob - The PetscDS 1566 - f - The test field number 1567 1568 Output Parameters: 1569 + f0 - boundary integrand for the test function term 1570 - f1 - boundary integrand for the test function gradient term 1571 1572 Note: We are using a first order FEM model for the weak form: 1573 1574 \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 1575 1576 The calling sequence for the callbacks f0 and f1 is given by: 1577 1578 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1579 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1580 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1581 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]) 1582 1583 + dim - the spatial dimension 1584 . Nf - the number of fields 1585 . uOff - the offset into u[] and u_t[] for each field 1586 . uOff_x - the offset into u_x[] for each field 1587 . u - each field evaluated at the current point 1588 . u_t - the time derivative of each field evaluated at the current point 1589 . u_x - the gradient of each field evaluated at the current point 1590 . aOff - the offset into a[] and a_t[] for each auxiliary field 1591 . aOff_x - the offset into a_x[] for each auxiliary field 1592 . a - each auxiliary field evaluated at the current point 1593 . a_t - the time derivative of each auxiliary field evaluated at the current point 1594 . a_x - the gradient of auxiliary each field evaluated at the current point 1595 . t - current time 1596 . x - coordinates of the current point 1597 . n - unit normal at the current point 1598 - f0 - output values at the current point 1599 1600 Level: intermediate 1601 1602 .seealso: PetscDSSetBdResidual() 1603 @*/ 1604 PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f, 1605 void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1606 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1607 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1608 PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]), 1609 void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1610 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1611 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1612 PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f1[])) 1613 { 1614 PetscFunctionBegin; 1615 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1616 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); 1617 if (f0) {PetscValidPointer(f0, 3); *f0 = prob->fBd[f*2+0];} 1618 if (f1) {PetscValidPointer(f1, 4); *f1 = prob->fBd[f*2+1];} 1619 PetscFunctionReturn(0); 1620 } 1621 1622 #undef __FUNCT__ 1623 #define __FUNCT__ "PetscDSSetBdResidual" 1624 /*@C 1625 PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field 1626 1627 Not collective 1628 1629 Input Parameters: 1630 + prob - The PetscDS 1631 . f - The test field number 1632 . f0 - boundary integrand for the test function term 1633 - f1 - boundary integrand for the test function gradient term 1634 1635 Note: We are using a first order FEM model for the weak form: 1636 1637 \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 1638 1639 The calling sequence for the callbacks f0 and f1 is given by: 1640 1641 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1642 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1643 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1644 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]) 1645 1646 + dim - the spatial dimension 1647 . Nf - the number of fields 1648 . uOff - the offset into u[] and u_t[] for each field 1649 . uOff_x - the offset into u_x[] for each field 1650 . u - each field evaluated at the current point 1651 . u_t - the time derivative of each field evaluated at the current point 1652 . u_x - the gradient of each field evaluated at the current point 1653 . aOff - the offset into a[] and a_t[] for each auxiliary field 1654 . aOff_x - the offset into a_x[] for each auxiliary field 1655 . a - each auxiliary field evaluated at the current point 1656 . a_t - the time derivative of each auxiliary field evaluated at the current point 1657 . a_x - the gradient of auxiliary each field evaluated at the current point 1658 . t - current time 1659 . x - coordinates of the current point 1660 . n - unit normal at the current point 1661 - f0 - output values at the current point 1662 1663 Level: intermediate 1664 1665 .seealso: PetscDSGetBdResidual() 1666 @*/ 1667 PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f, 1668 void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1669 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1670 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1671 PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]), 1672 void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1673 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1674 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1675 PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f1[])) 1676 { 1677 PetscErrorCode ierr; 1678 1679 PetscFunctionBegin; 1680 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1681 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1682 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1683 if (f0) {PetscValidFunction(f0, 3); prob->fBd[f*2+0] = f0;} 1684 if (f1) {PetscValidFunction(f1, 4); prob->fBd[f*2+1] = f1;} 1685 PetscFunctionReturn(0); 1686 } 1687 1688 #undef __FUNCT__ 1689 #define __FUNCT__ "PetscDSGetBdJacobian" 1690 /*@C 1691 PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field 1692 1693 Not collective 1694 1695 Input Parameters: 1696 + prob - The PetscDS 1697 . f - The test field number 1698 - g - The field number 1699 1700 Output Parameters: 1701 + g0 - integrand for the test and basis function term 1702 . g1 - integrand for the test function and basis function gradient term 1703 . g2 - integrand for the test function gradient and basis function term 1704 - g3 - integrand for the test function gradient and basis function gradient term 1705 1706 Note: We are using a first order FEM model for the weak form: 1707 1708 \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 1709 1710 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1711 1712 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1713 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1714 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1715 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[]) 1716 1717 + dim - the spatial dimension 1718 . Nf - the number of fields 1719 . uOff - the offset into u[] and u_t[] for each field 1720 . uOff_x - the offset into u_x[] for each field 1721 . u - each field evaluated at the current point 1722 . u_t - the time derivative of each field evaluated at the current point 1723 . u_x - the gradient of each field evaluated at the current point 1724 . aOff - the offset into a[] and a_t[] for each auxiliary field 1725 . aOff_x - the offset into a_x[] for each auxiliary field 1726 . a - each auxiliary field evaluated at the current point 1727 . a_t - the time derivative of each auxiliary field evaluated at the current point 1728 . a_x - the gradient of auxiliary each field evaluated at the current point 1729 . t - current time 1730 . u_tShift - the multiplier a for dF/dU_t 1731 . x - coordinates of the current point 1732 . n - normal at the current point 1733 - g0 - output values at the current point 1734 1735 Level: intermediate 1736 1737 .seealso: PetscDSSetBdJacobian() 1738 @*/ 1739 PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g, 1740 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1741 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1742 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1743 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g0[]), 1744 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1745 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1746 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1747 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g1[]), 1748 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1749 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1750 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1751 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g2[]), 1752 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1753 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1754 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1755 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g3[])) 1756 { 1757 PetscFunctionBegin; 1758 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1759 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); 1760 if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf); 1761 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->gBd[(f*prob->Nf + g)*4+0];} 1762 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->gBd[(f*prob->Nf + g)*4+1];} 1763 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->gBd[(f*prob->Nf + g)*4+2];} 1764 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->gBd[(f*prob->Nf + g)*4+3];} 1765 PetscFunctionReturn(0); 1766 } 1767 1768 #undef __FUNCT__ 1769 #define __FUNCT__ "PetscDSSetBdJacobian" 1770 /*@C 1771 PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field 1772 1773 Not collective 1774 1775 Input Parameters: 1776 + prob - The PetscDS 1777 . f - The test field number 1778 . g - The field number 1779 . g0 - integrand for the test and basis function term 1780 . g1 - integrand for the test function and basis function gradient term 1781 . g2 - integrand for the test function gradient and basis function term 1782 - g3 - integrand for the test function gradient and basis function gradient term 1783 1784 Note: We are using a first order FEM model for the weak form: 1785 1786 \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 1787 1788 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1789 1790 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1791 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1792 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1793 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[]) 1794 1795 + dim - the spatial dimension 1796 . Nf - the number of fields 1797 . uOff - the offset into u[] and u_t[] for each field 1798 . uOff_x - the offset into u_x[] for each field 1799 . u - each field evaluated at the current point 1800 . u_t - the time derivative of each field evaluated at the current point 1801 . u_x - the gradient of each field evaluated at the current point 1802 . aOff - the offset into a[] and a_t[] for each auxiliary field 1803 . aOff_x - the offset into a_x[] for each auxiliary field 1804 . a - each auxiliary field evaluated at the current point 1805 . a_t - the time derivative of each auxiliary field evaluated at the current point 1806 . a_x - the gradient of auxiliary each field evaluated at the current point 1807 . t - current time 1808 . u_tShift - the multiplier a for dF/dU_t 1809 . x - coordinates of the current point 1810 . n - normal at the current point 1811 - g0 - output values at the current point 1812 1813 Level: intermediate 1814 1815 .seealso: PetscDSGetBdJacobian() 1816 @*/ 1817 PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g, 1818 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1819 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1820 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1821 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g0[]), 1822 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1823 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1824 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1825 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g1[]), 1826 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1827 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1828 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1829 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g2[]), 1830 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1831 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1832 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1833 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g3[])) 1834 { 1835 PetscErrorCode ierr; 1836 1837 PetscFunctionBegin; 1838 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1839 if (g0) PetscValidFunction(g0, 4); 1840 if (g1) PetscValidFunction(g1, 5); 1841 if (g2) PetscValidFunction(g2, 6); 1842 if (g3) PetscValidFunction(g3, 7); 1843 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1844 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1845 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1846 prob->gBd[(f*prob->Nf + g)*4+0] = g0; 1847 prob->gBd[(f*prob->Nf + g)*4+1] = g1; 1848 prob->gBd[(f*prob->Nf + g)*4+2] = g2; 1849 prob->gBd[(f*prob->Nf + g)*4+3] = g3; 1850 PetscFunctionReturn(0); 1851 } 1852 1853 #undef __FUNCT__ 1854 #define __FUNCT__ "PetscDSGetFieldOffset" 1855 /*@ 1856 PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis 1857 1858 Not collective 1859 1860 Input Parameters: 1861 + prob - The PetscDS object 1862 - f - The field number 1863 1864 Output Parameter: 1865 . off - The offset 1866 1867 Level: beginner 1868 1869 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1870 @*/ 1871 PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off) 1872 { 1873 PetscInt g; 1874 PetscErrorCode ierr; 1875 1876 PetscFunctionBegin; 1877 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1878 PetscValidPointer(off, 3); 1879 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); 1880 *off = 0; 1881 for (g = 0; g < f; ++g) { 1882 PetscFE fe = (PetscFE) prob->disc[g]; 1883 PetscInt Nb, Nc; 1884 1885 ierr = PetscFEGetDimension(fe, &Nb);CHKERRQ(ierr); 1886 ierr = PetscFEGetNumComponents(fe, &Nc);CHKERRQ(ierr); 1887 *off += Nb*Nc; 1888 } 1889 PetscFunctionReturn(0); 1890 } 1891 1892 #undef __FUNCT__ 1893 #define __FUNCT__ "PetscDSGetBdFieldOffset" 1894 /*@ 1895 PetscDSGetBdFieldOffset - Returns the offset of the given field in the full space boundary basis 1896 1897 Not collective 1898 1899 Input Parameters: 1900 + prob - The PetscDS object 1901 - f - The field number 1902 1903 Output Parameter: 1904 . off - The boundary offset 1905 1906 Level: beginner 1907 1908 .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1909 @*/ 1910 PetscErrorCode PetscDSGetBdFieldOffset(PetscDS prob, PetscInt f, PetscInt *off) 1911 { 1912 PetscInt g; 1913 PetscErrorCode ierr; 1914 1915 PetscFunctionBegin; 1916 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1917 PetscValidPointer(off, 3); 1918 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); 1919 *off = 0; 1920 for (g = 0; g < f; ++g) { 1921 PetscFE fe = (PetscFE) prob->discBd[g]; 1922 PetscInt Nb, Nc; 1923 1924 ierr = PetscFEGetDimension(fe, &Nb);CHKERRQ(ierr); 1925 ierr = PetscFEGetNumComponents(fe, &Nc);CHKERRQ(ierr); 1926 *off += Nb*Nc; 1927 } 1928 PetscFunctionReturn(0); 1929 } 1930 1931 #undef __FUNCT__ 1932 #define __FUNCT__ "PetscDSGetComponentOffset" 1933 /*@ 1934 PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point 1935 1936 Not collective 1937 1938 Input Parameters: 1939 + prob - The PetscDS object 1940 - f - The field number 1941 1942 Output Parameter: 1943 . off - The offset 1944 1945 Level: beginner 1946 1947 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1948 @*/ 1949 PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off) 1950 { 1951 PetscInt g; 1952 PetscErrorCode ierr; 1953 1954 PetscFunctionBegin; 1955 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1956 PetscValidPointer(off, 3); 1957 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); 1958 *off = 0; 1959 for (g = 0; g < f; ++g) { 1960 PetscFE fe = (PetscFE) prob->disc[g]; 1961 PetscInt Nc; 1962 1963 ierr = PetscFEGetNumComponents(fe, &Nc);CHKERRQ(ierr); 1964 *off += Nc; 1965 } 1966 PetscFunctionReturn(0); 1967 } 1968 1969 #undef __FUNCT__ 1970 #define __FUNCT__ "PetscDSGetComponentOffsets" 1971 /*@ 1972 PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point 1973 1974 Not collective 1975 1976 Input Parameter: 1977 . prob - The PetscDS object 1978 1979 Output Parameter: 1980 . offsets - The offsets 1981 1982 Level: beginner 1983 1984 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1985 @*/ 1986 PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[]) 1987 { 1988 PetscFunctionBegin; 1989 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1990 PetscValidPointer(offsets, 2); 1991 *offsets = prob->off; 1992 PetscFunctionReturn(0); 1993 } 1994 1995 #undef __FUNCT__ 1996 #define __FUNCT__ "PetscDSGetComponentDerivativeOffsets" 1997 /*@ 1998 PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point 1999 2000 Not collective 2001 2002 Input Parameter: 2003 . prob - The PetscDS object 2004 2005 Output Parameter: 2006 . offsets - The offsets 2007 2008 Level: beginner 2009 2010 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 2011 @*/ 2012 PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[]) 2013 { 2014 PetscFunctionBegin; 2015 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2016 PetscValidPointer(offsets, 2); 2017 *offsets = prob->offDer; 2018 PetscFunctionReturn(0); 2019 } 2020 2021 #undef __FUNCT__ 2022 #define __FUNCT__ "PetscDSGetComponentBdOffsets" 2023 /*@ 2024 PetscDSGetComponentBdOffsets - Returns the offset of each field on a boundary evaluation point 2025 2026 Not collective 2027 2028 Input Parameter: 2029 . prob - The PetscDS object 2030 2031 Output Parameter: 2032 . offsets - The offsets 2033 2034 Level: beginner 2035 2036 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 2037 @*/ 2038 PetscErrorCode PetscDSGetComponentBdOffsets(PetscDS prob, PetscInt *offsets[]) 2039 { 2040 PetscFunctionBegin; 2041 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2042 PetscValidPointer(offsets, 2); 2043 *offsets = prob->offBd; 2044 PetscFunctionReturn(0); 2045 } 2046 2047 #undef __FUNCT__ 2048 #define __FUNCT__ "PetscDSGetComponentBdDerivativeOffsets" 2049 /*@ 2050 PetscDSGetComponentBdDerivativeOffsets - Returns the offset of each field derivative on a boundary evaluation point 2051 2052 Not collective 2053 2054 Input Parameter: 2055 . prob - The PetscDS object 2056 2057 Output Parameter: 2058 . offsets - The offsets 2059 2060 Level: beginner 2061 2062 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 2063 @*/ 2064 PetscErrorCode PetscDSGetComponentBdDerivativeOffsets(PetscDS prob, PetscInt *offsets[]) 2065 { 2066 PetscFunctionBegin; 2067 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2068 PetscValidPointer(offsets, 2); 2069 *offsets = prob->offDerBd; 2070 PetscFunctionReturn(0); 2071 } 2072 2073 #undef __FUNCT__ 2074 #define __FUNCT__ "PetscDSGetTabulation" 2075 /*@C 2076 PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization 2077 2078 Not collective 2079 2080 Input Parameter: 2081 . prob - The PetscDS object 2082 2083 Output Parameters: 2084 + basis - The basis function tabulation at quadrature points 2085 - basisDer - The basis function derivative tabulation at quadrature points 2086 2087 Level: intermediate 2088 2089 .seealso: PetscDSGetBdTabulation(), PetscDSCreate() 2090 @*/ 2091 PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer) 2092 { 2093 PetscErrorCode ierr; 2094 2095 PetscFunctionBegin; 2096 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2097 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2098 if (basis) {PetscValidPointer(basis, 2); *basis = prob->basis;} 2099 if (basisDer) {PetscValidPointer(basisDer, 3); *basisDer = prob->basisDer;} 2100 PetscFunctionReturn(0); 2101 } 2102 2103 #undef __FUNCT__ 2104 #define __FUNCT__ "PetscDSGetBdTabulation" 2105 /*@C 2106 PetscDSGetBdTabulation - Return the basis tabulation at quadrature points for the boundary discretization 2107 2108 Not collective 2109 2110 Input Parameter: 2111 . prob - The PetscDS object 2112 2113 Output Parameters: 2114 + basis - The basis function tabulation at quadrature points 2115 - basisDer - The basis function derivative tabulation at quadrature points 2116 2117 Level: intermediate 2118 2119 .seealso: PetscDSGetTabulation(), PetscDSCreate() 2120 @*/ 2121 PetscErrorCode PetscDSGetBdTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer) 2122 { 2123 PetscErrorCode ierr; 2124 2125 PetscFunctionBegin; 2126 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2127 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2128 if (basis) {PetscValidPointer(basis, 2); *basis = prob->basisBd;} 2129 if (basisDer) {PetscValidPointer(basisDer, 3); *basisDer = prob->basisDerBd;} 2130 PetscFunctionReturn(0); 2131 } 2132 2133 #undef __FUNCT__ 2134 #define __FUNCT__ "PetscDSGetEvaluationArrays" 2135 PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x) 2136 { 2137 PetscErrorCode ierr; 2138 2139 PetscFunctionBegin; 2140 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2141 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2142 if (u) {PetscValidPointer(u, 2); *u = prob->u;} 2143 if (u_t) {PetscValidPointer(u_t, 3); *u_t = prob->u_t;} 2144 if (u_x) {PetscValidPointer(u_x, 4); *u_x = prob->u_x;} 2145 PetscFunctionReturn(0); 2146 } 2147 2148 #undef __FUNCT__ 2149 #define __FUNCT__ "PetscDSGetWeakFormArrays" 2150 PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3) 2151 { 2152 PetscErrorCode ierr; 2153 2154 PetscFunctionBegin; 2155 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2156 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2157 if (f0) {PetscValidPointer(f0, 2); *f0 = prob->f0;} 2158 if (f1) {PetscValidPointer(f1, 3); *f1 = prob->f1;} 2159 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->g0;} 2160 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->g1;} 2161 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->g2;} 2162 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->g3;} 2163 PetscFunctionReturn(0); 2164 } 2165 2166 #undef __FUNCT__ 2167 #define __FUNCT__ "PetscDSGetRefCoordArrays" 2168 PetscErrorCode PetscDSGetRefCoordArrays(PetscDS prob, PetscReal **x, PetscScalar **refSpaceDer) 2169 { 2170 PetscErrorCode ierr; 2171 2172 PetscFunctionBegin; 2173 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2174 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2175 if (x) {PetscValidPointer(x, 2); *x = prob->x;} 2176 if (refSpaceDer) {PetscValidPointer(refSpaceDer, 3); *refSpaceDer = prob->refSpaceDer;} 2177 PetscFunctionReturn(0); 2178 } 2179 2180 #undef __FUNCT__ 2181 #define __FUNCT__ "PetscDSCopyEquations" 2182 /*@ 2183 PetscDSCopyEquations - Copy all pointwise function pointers to the new problem 2184 2185 Not collective 2186 2187 Input Parameter: 2188 . prob - The PetscDS object 2189 2190 Output Parameter: 2191 . newprob - The PetscDS copy 2192 2193 Level: intermediate 2194 2195 .seealso: PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate() 2196 @*/ 2197 PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob) 2198 { 2199 PetscInt Nf, Ng, f, g; 2200 PetscErrorCode ierr; 2201 2202 PetscFunctionBegin; 2203 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2204 PetscValidHeaderSpecific(newprob, PETSCDS_CLASSID, 2); 2205 ierr = PetscDSGetNumFields(prob, &Nf);CHKERRQ(ierr); 2206 ierr = PetscDSGetNumFields(newprob, &Ng);CHKERRQ(ierr); 2207 if (Nf != Ng) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %D != %D", Nf, Ng);CHKERRQ(ierr); 2208 for (f = 0; f < Nf; ++f) { 2209 PetscPointFunc obj; 2210 PetscPointFunc f0, f1; 2211 PetscPointJac g0, g1, g2, g3; 2212 PetscBdPointFunc f0Bd, f1Bd; 2213 PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd; 2214 PetscRiemannFunc r; 2215 2216 ierr = PetscDSGetObjective(prob, f, &obj);CHKERRQ(ierr); 2217 ierr = PetscDSGetResidual(prob, f, &f0, &f1);CHKERRQ(ierr); 2218 ierr = PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);CHKERRQ(ierr); 2219 ierr = PetscDSGetRiemannSolver(prob, f, &r);CHKERRQ(ierr); 2220 ierr = PetscDSSetObjective(newprob, f, obj);CHKERRQ(ierr); 2221 ierr = PetscDSSetResidual(newprob, f, f0, f1);CHKERRQ(ierr); 2222 ierr = PetscDSSetBdResidual(newprob, f, f0Bd, f1Bd);CHKERRQ(ierr); 2223 ierr = PetscDSSetRiemannSolver(newprob, f, r);CHKERRQ(ierr); 2224 for (g = 0; g < Nf; ++g) { 2225 ierr = PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);CHKERRQ(ierr); 2226 ierr = PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);CHKERRQ(ierr); 2227 ierr = PetscDSSetJacobian(newprob, f, g, g0, g1, g2, g3);CHKERRQ(ierr); 2228 ierr = PetscDSSetBdJacobian(newprob, f, g, g0Bd, g1Bd, g2Bd, g3Bd);CHKERRQ(ierr); 2229 } 2230 } 2231 PetscFunctionReturn(0); 2232 } 2233 2234 #undef __FUNCT__ 2235 #define __FUNCT__ "PetscDSDestroy_Basic" 2236 static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob) 2237 { 2238 PetscFunctionBegin; 2239 PetscFunctionReturn(0); 2240 } 2241 2242 #undef __FUNCT__ 2243 #define __FUNCT__ "PetscDSInitialize_Basic" 2244 static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob) 2245 { 2246 PetscFunctionBegin; 2247 prob->ops->setfromoptions = NULL; 2248 prob->ops->setup = NULL; 2249 prob->ops->view = NULL; 2250 prob->ops->destroy = PetscDSDestroy_Basic; 2251 PetscFunctionReturn(0); 2252 } 2253 2254 /*MC 2255 PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions 2256 2257 Level: intermediate 2258 2259 .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType() 2260 M*/ 2261 2262 #undef __FUNCT__ 2263 #define __FUNCT__ "PetscDSCreate_Basic" 2264 PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob) 2265 { 2266 PetscDS_Basic *b; 2267 PetscErrorCode ierr; 2268 2269 PetscFunctionBegin; 2270 PetscValidHeaderSpecific(prob, PETSCSPACE_CLASSID, 1); 2271 ierr = PetscNewLog(prob, &b);CHKERRQ(ierr); 2272 prob->data = b; 2273 2274 ierr = PetscDSInitialize_Basic(prob);CHKERRQ(ierr); 2275 PetscFunctionReturn(0); 2276 } 2277