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 /*@C 9 PetscDSRegister - Adds a new PetscDS implementation 10 11 Not Collective 12 13 Input Parameters: 14 + name - The name of a new user-defined creation routine 15 - create_func - The creation routine itself 16 17 Notes: 18 PetscDSRegister() may be called multiple times to add several user-defined PetscDSs 19 20 Sample usage: 21 .vb 22 PetscDSRegister("my_ds", MyPetscDSCreate); 23 .ve 24 25 Then, your PetscDS type can be chosen with the procedural interface via 26 .vb 27 PetscDSCreate(MPI_Comm, PetscDS *); 28 PetscDSSetType(PetscDS, "my_ds"); 29 .ve 30 or at runtime via the option 31 .vb 32 -petscds_type my_ds 33 .ve 34 35 Level: advanced 36 37 .keywords: PetscDS, register 38 .seealso: PetscDSRegisterAll(), PetscDSRegisterDestroy() 39 40 @*/ 41 PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS)) 42 { 43 PetscErrorCode ierr; 44 45 PetscFunctionBegin; 46 ierr = PetscFunctionListAdd(&PetscDSList, sname, function);CHKERRQ(ierr); 47 PetscFunctionReturn(0); 48 } 49 50 /*@C 51 PetscDSSetType - Builds a particular PetscDS 52 53 Collective on PetscDS 54 55 Input Parameters: 56 + prob - The PetscDS object 57 - name - The kind of system 58 59 Options Database Key: 60 . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types 61 62 Level: intermediate 63 64 .keywords: PetscDS, set, type 65 .seealso: PetscDSGetType(), PetscDSCreate() 66 @*/ 67 PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name) 68 { 69 PetscErrorCode (*r)(PetscDS); 70 PetscBool match; 71 PetscErrorCode ierr; 72 73 PetscFunctionBegin; 74 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 75 ierr = PetscObjectTypeCompare((PetscObject) prob, name, &match);CHKERRQ(ierr); 76 if (match) PetscFunctionReturn(0); 77 78 ierr = PetscDSRegisterAll();CHKERRQ(ierr); 79 ierr = PetscFunctionListFind(PetscDSList, name, &r);CHKERRQ(ierr); 80 if (!r) SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name); 81 82 if (prob->ops->destroy) { 83 ierr = (*prob->ops->destroy)(prob);CHKERRQ(ierr); 84 prob->ops->destroy = NULL; 85 } 86 ierr = (*r)(prob);CHKERRQ(ierr); 87 ierr = PetscObjectChangeTypeName((PetscObject) prob, name);CHKERRQ(ierr); 88 PetscFunctionReturn(0); 89 } 90 91 /*@C 92 PetscDSGetType - Gets the PetscDS type name (as a string) from the object. 93 94 Not Collective 95 96 Input Parameter: 97 . prob - The PetscDS 98 99 Output Parameter: 100 . name - The PetscDS type name 101 102 Level: intermediate 103 104 .keywords: PetscDS, get, type, name 105 .seealso: PetscDSSetType(), PetscDSCreate() 106 @*/ 107 PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name) 108 { 109 PetscErrorCode ierr; 110 111 PetscFunctionBegin; 112 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 113 PetscValidPointer(name, 2); 114 ierr = PetscDSRegisterAll();CHKERRQ(ierr); 115 *name = ((PetscObject) prob)->type_name; 116 PetscFunctionReturn(0); 117 } 118 119 static PetscErrorCode PetscDSView_Ascii(PetscDS prob, PetscViewer viewer) 120 { 121 PetscViewerFormat format; 122 const PetscScalar *constants; 123 PetscInt numConstants, f; 124 PetscErrorCode ierr; 125 126 PetscFunctionBegin; 127 ierr = PetscViewerGetFormat(viewer, &format);CHKERRQ(ierr); 128 ierr = PetscViewerASCIIPrintf(viewer, "Discrete System with %d fields\n", prob->Nf);CHKERRQ(ierr); 129 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 130 for (f = 0; f < prob->Nf; ++f) { 131 PetscObject obj; 132 PetscClassId id; 133 const char *name; 134 PetscInt Nc; 135 136 ierr = PetscDSGetDiscretization(prob, f, &obj);CHKERRQ(ierr); 137 ierr = PetscObjectGetClassId(obj, &id);CHKERRQ(ierr); 138 ierr = PetscObjectGetName(obj, &name);CHKERRQ(ierr); 139 ierr = PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>");CHKERRQ(ierr); 140 if (id == PETSCFE_CLASSID) { 141 ierr = PetscFEGetNumComponents((PetscFE) obj, &Nc);CHKERRQ(ierr); 142 ierr = PetscViewerASCIIPrintf(viewer, " FEM");CHKERRQ(ierr); 143 } else if (id == PETSCFV_CLASSID) { 144 ierr = PetscFVGetNumComponents((PetscFV) obj, &Nc);CHKERRQ(ierr); 145 ierr = PetscViewerASCIIPrintf(viewer, " FVM");CHKERRQ(ierr); 146 } 147 else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %D", f); 148 if (Nc > 1) {ierr = PetscViewerASCIIPrintf(viewer, "%D components", Nc);CHKERRQ(ierr);} 149 else {ierr = PetscViewerASCIIPrintf(viewer, "%D component ", Nc);CHKERRQ(ierr);} 150 if (prob->implicit[f]) {ierr = PetscViewerASCIIPrintf(viewer, " (implicit)");CHKERRQ(ierr);} 151 else {ierr = PetscViewerASCIIPrintf(viewer, " (explicit)");CHKERRQ(ierr);} 152 if (prob->adjacency[f*2+0]) { 153 if (prob->adjacency[f*2+1]) {ierr = PetscViewerASCIIPrintf(viewer, " (adj FVM++)");CHKERRQ(ierr);} 154 else {ierr = PetscViewerASCIIPrintf(viewer, " (adj FVM)");CHKERRQ(ierr);} 155 } else { 156 if (prob->adjacency[f*2+1]) {ierr = PetscViewerASCIIPrintf(viewer, " (adj FEM)");CHKERRQ(ierr);} 157 else {ierr = PetscViewerASCIIPrintf(viewer, " (adj FUNKY)");CHKERRQ(ierr);} 158 } 159 ierr = PetscViewerASCIIPrintf(viewer, "\n");CHKERRQ(ierr); 160 if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) { 161 if (id == PETSCFE_CLASSID) {ierr = PetscFEView((PetscFE) obj, viewer);CHKERRQ(ierr);} 162 else if (id == PETSCFV_CLASSID) {ierr = PetscFVView((PetscFV) obj, viewer);CHKERRQ(ierr);} 163 } 164 } 165 ierr = PetscDSGetConstants(prob, &numConstants, &constants);CHKERRQ(ierr); 166 if (numConstants) { 167 ierr = PetscViewerASCIIPrintf(viewer, "%D constants\n", numConstants);CHKERRQ(ierr); 168 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 169 for (f = 0; f < numConstants; ++f) {ierr = PetscViewerASCIIPrintf(viewer, "%g\n", (double) PetscRealPart(constants[f]));CHKERRQ(ierr);} 170 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 171 } 172 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 173 PetscFunctionReturn(0); 174 } 175 176 /*@C 177 PetscDSView - Views a PetscDS 178 179 Collective on PetscDS 180 181 Input Parameter: 182 + prob - the PetscDS object to view 183 - v - the viewer 184 185 Level: developer 186 187 .seealso PetscDSDestroy() 188 @*/ 189 PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v) 190 { 191 PetscBool iascii; 192 PetscErrorCode ierr; 193 194 PetscFunctionBegin; 195 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 196 if (!v) {ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) prob), &v);CHKERRQ(ierr);} 197 else {PetscValidHeaderSpecific(v, PETSC_VIEWER_CLASSID, 2);} 198 ierr = PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 199 if (iascii) {ierr = PetscDSView_Ascii(prob, v);CHKERRQ(ierr);} 200 if (prob->ops->view) {ierr = (*prob->ops->view)(prob, v);CHKERRQ(ierr);} 201 PetscFunctionReturn(0); 202 } 203 204 /*@ 205 PetscDSSetFromOptions - sets parameters in a PetscDS from the options database 206 207 Collective on PetscDS 208 209 Input Parameter: 210 . prob - the PetscDS object to set options for 211 212 Options Database: 213 214 Level: developer 215 216 .seealso PetscDSView() 217 @*/ 218 PetscErrorCode PetscDSSetFromOptions(PetscDS prob) 219 { 220 DSBoundary b; 221 const char *defaultType; 222 char name[256]; 223 PetscBool flg; 224 PetscErrorCode ierr; 225 226 PetscFunctionBegin; 227 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 228 if (!((PetscObject) prob)->type_name) { 229 defaultType = PETSCDSBASIC; 230 } else { 231 defaultType = ((PetscObject) prob)->type_name; 232 } 233 ierr = PetscDSRegisterAll();CHKERRQ(ierr); 234 235 ierr = PetscObjectOptionsBegin((PetscObject) prob);CHKERRQ(ierr); 236 for (b = prob->boundary; b; b = b->next) { 237 char optname[1024]; 238 PetscInt ids[1024], len = 1024; 239 PetscBool flg; 240 241 ierr = PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name);CHKERRQ(ierr); 242 ierr = PetscMemzero(ids, sizeof(ids));CHKERRQ(ierr); 243 ierr = PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg);CHKERRQ(ierr); 244 if (flg) { 245 b->numids = len; 246 ierr = PetscFree(b->ids);CHKERRQ(ierr); 247 ierr = PetscMalloc1(len, &b->ids);CHKERRQ(ierr); 248 ierr = PetscMemcpy(b->ids, ids, len*sizeof(PetscInt));CHKERRQ(ierr); 249 } 250 len = 1024; 251 ierr = PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name);CHKERRQ(ierr); 252 ierr = PetscMemzero(ids, sizeof(ids));CHKERRQ(ierr); 253 ierr = PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg);CHKERRQ(ierr); 254 if (flg) { 255 b->numcomps = len; 256 ierr = PetscFree(b->comps);CHKERRQ(ierr); 257 ierr = PetscMalloc1(len, &b->comps);CHKERRQ(ierr); 258 ierr = PetscMemcpy(b->comps, ids, len*sizeof(PetscInt));CHKERRQ(ierr); 259 } 260 } 261 ierr = PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg);CHKERRQ(ierr); 262 if (flg) { 263 ierr = PetscDSSetType(prob, name);CHKERRQ(ierr); 264 } else if (!((PetscObject) prob)->type_name) { 265 ierr = PetscDSSetType(prob, defaultType);CHKERRQ(ierr); 266 } 267 if (prob->ops->setfromoptions) {ierr = (*prob->ops->setfromoptions)(prob);CHKERRQ(ierr);} 268 /* process any options handlers added with PetscObjectAddOptionsHandler() */ 269 ierr = PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject) prob);CHKERRQ(ierr); 270 ierr = PetscOptionsEnd();CHKERRQ(ierr); 271 ierr = PetscDSViewFromOptions(prob, NULL, "-petscds_view");CHKERRQ(ierr); 272 PetscFunctionReturn(0); 273 } 274 275 /*@C 276 PetscDSSetUp - Construct data structures for the PetscDS 277 278 Collective on PetscDS 279 280 Input Parameter: 281 . prob - the PetscDS object to setup 282 283 Level: developer 284 285 .seealso PetscDSView(), PetscDSDestroy() 286 @*/ 287 PetscErrorCode PetscDSSetUp(PetscDS prob) 288 { 289 const PetscInt Nf = prob->Nf; 290 PetscInt dim, work, NcMax = 0, NqMax = 0, f; 291 PetscErrorCode ierr; 292 293 PetscFunctionBegin; 294 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 295 if (prob->setup) PetscFunctionReturn(0); 296 /* Calculate sizes */ 297 ierr = PetscDSGetSpatialDimension(prob, &dim);CHKERRQ(ierr); 298 prob->totDim = prob->totComp = 0; 299 ierr = PetscMalloc2(Nf,&prob->Nc,Nf,&prob->Nb);CHKERRQ(ierr); 300 ierr = PetscCalloc2(Nf+1,&prob->off,Nf+1,&prob->offDer);CHKERRQ(ierr); 301 ierr = PetscMalloc4(Nf,&prob->basis,Nf,&prob->basisDer,Nf,&prob->basisFace,Nf,&prob->basisDerFace);CHKERRQ(ierr); 302 for (f = 0; f < Nf; ++f) { 303 PetscObject obj; 304 PetscClassId id; 305 PetscQuadrature q; 306 PetscInt Nq = 0, Nb, Nc; 307 308 ierr = PetscDSGetDiscretization(prob, f, &obj);CHKERRQ(ierr); 309 ierr = PetscObjectGetClassId(obj, &id);CHKERRQ(ierr); 310 if (id == PETSCFE_CLASSID) { 311 PetscFE fe = (PetscFE) obj; 312 313 ierr = PetscFEGetQuadrature(fe, &q);CHKERRQ(ierr); 314 ierr = PetscFEGetDimension(fe, &Nb);CHKERRQ(ierr); 315 ierr = PetscFEGetNumComponents(fe, &Nc);CHKERRQ(ierr); 316 ierr = PetscFEGetDefaultTabulation(fe, &prob->basis[f], &prob->basisDer[f], NULL);CHKERRQ(ierr); 317 ierr = PetscFEGetFaceTabulation(fe, &prob->basisFace[f], &prob->basisDerFace[f], NULL);CHKERRQ(ierr); 318 } else if (id == PETSCFV_CLASSID) { 319 PetscFV fv = (PetscFV) obj; 320 321 ierr = PetscFVGetQuadrature(fv, &q);CHKERRQ(ierr); 322 ierr = PetscFVGetNumComponents(fv, &Nc);CHKERRQ(ierr); 323 Nb = Nc; 324 ierr = PetscFVGetDefaultTabulation(fv, &prob->basis[f], &prob->basisDer[f], NULL);CHKERRQ(ierr); 325 /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */ 326 } else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f); 327 prob->Nc[f] = Nc; 328 prob->Nb[f] = Nb; 329 prob->off[f+1] = Nc + prob->off[f]; 330 prob->offDer[f+1] = Nc*dim + prob->offDer[f]; 331 if (q) {ierr = PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL);CHKERRQ(ierr);} 332 NqMax = PetscMax(NqMax, Nq); 333 NcMax = PetscMax(NcMax, Nc); 334 prob->totDim += Nb; 335 prob->totComp += Nc; 336 } 337 work = PetscMax(prob->totComp*dim, PetscSqr(NcMax*dim)); 338 /* Allocate works space */ 339 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); 340 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); 341 if (prob->ops->setup) {ierr = (*prob->ops->setup)(prob);CHKERRQ(ierr);} 342 prob->setup = PETSC_TRUE; 343 PetscFunctionReturn(0); 344 } 345 346 static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob) 347 { 348 PetscErrorCode ierr; 349 350 PetscFunctionBegin; 351 ierr = PetscFree2(prob->Nc,prob->Nb);CHKERRQ(ierr); 352 ierr = PetscFree2(prob->off,prob->offDer);CHKERRQ(ierr); 353 ierr = PetscFree4(prob->basis,prob->basisDer,prob->basisFace,prob->basisDerFace);CHKERRQ(ierr); 354 ierr = PetscFree5(prob->u,prob->u_t,prob->u_x,prob->x,prob->refSpaceDer);CHKERRQ(ierr); 355 ierr = PetscFree6(prob->f0,prob->f1,prob->g0,prob->g1,prob->g2,prob->g3);CHKERRQ(ierr); 356 PetscFunctionReturn(0); 357 } 358 359 static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew) 360 { 361 PetscObject *tmpd; 362 PetscBool *tmpi, *tmpa; 363 PetscPointFunc *tmpobj, *tmpf, *tmpup; 364 PetscPointJac *tmpg, *tmpgp, *tmpgt; 365 PetscBdPointFunc *tmpfbd; 366 PetscBdPointJac *tmpgbd; 367 PetscRiemannFunc *tmpr; 368 PetscSimplePointFunc *tmpexactSol; 369 void **tmpctx; 370 PetscInt Nf = prob->Nf, f, i; 371 PetscErrorCode ierr; 372 373 PetscFunctionBegin; 374 if (Nf >= NfNew) PetscFunctionReturn(0); 375 prob->setup = PETSC_FALSE; 376 ierr = PetscDSDestroyStructs_Static(prob);CHKERRQ(ierr); 377 ierr = PetscMalloc3(NfNew, &tmpd, NfNew, &tmpi, NfNew*2, &tmpa);CHKERRQ(ierr); 378 for (f = 0; f < Nf; ++f) {tmpd[f] = prob->disc[f]; tmpi[f] = prob->implicit[f]; for (i = 0; i < 2; ++i) tmpa[f*2+i] = prob->adjacency[f*2+i];} 379 for (f = Nf; f < NfNew; ++f) {tmpd[f] = NULL; tmpi[f] = PETSC_TRUE; tmpa[f*2+0] = PETSC_FALSE; tmpa[f*2+1] = PETSC_TRUE;} 380 ierr = PetscFree3(prob->disc, prob->implicit, prob->adjacency);CHKERRQ(ierr); 381 prob->Nf = NfNew; 382 prob->disc = tmpd; 383 prob->implicit = tmpi; 384 prob->adjacency = tmpa; 385 ierr = PetscCalloc7(NfNew, &tmpobj, NfNew*2, &tmpf, NfNew*NfNew*4, &tmpg, NfNew*NfNew*4, &tmpgp, NfNew*NfNew*4, &tmpgt, NfNew, &tmpr, NfNew, &tmpctx);CHKERRQ(ierr); 386 ierr = PetscCalloc1(NfNew, &tmpup);CHKERRQ(ierr); 387 for (f = 0; f < Nf; ++f) tmpobj[f] = prob->obj[f]; 388 for (f = 0; f < Nf*2; ++f) tmpf[f] = prob->f[f]; 389 for (f = 0; f < Nf*Nf*4; ++f) tmpg[f] = prob->g[f]; 390 for (f = 0; f < Nf*Nf*4; ++f) tmpgp[f] = prob->gp[f]; 391 for (f = 0; f < Nf; ++f) tmpr[f] = prob->r[f]; 392 for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f]; 393 for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f]; 394 for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL; 395 for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL; 396 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL; 397 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgp[f] = NULL; 398 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgt[f] = NULL; 399 for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL; 400 for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL; 401 for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL; 402 ierr = PetscFree7(prob->obj, prob->f, prob->g, prob->gp, prob->gt, prob->r, prob->ctx);CHKERRQ(ierr); 403 ierr = PetscFree(prob->update);CHKERRQ(ierr); 404 prob->obj = tmpobj; 405 prob->f = tmpf; 406 prob->g = tmpg; 407 prob->gp = tmpgp; 408 prob->gt = tmpgt; 409 prob->r = tmpr; 410 prob->update = tmpup; 411 prob->ctx = tmpctx; 412 ierr = PetscCalloc3(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd, NfNew, &tmpexactSol);CHKERRQ(ierr); 413 for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f]; 414 for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f]; 415 for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f]; 416 for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL; 417 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL; 418 for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL; 419 ierr = PetscFree3(prob->fBd, prob->gBd, prob->exactSol);CHKERRQ(ierr); 420 prob->fBd = tmpfbd; 421 prob->gBd = tmpgbd; 422 prob->exactSol = tmpexactSol; 423 PetscFunctionReturn(0); 424 } 425 426 /*@ 427 PetscDSDestroy - Destroys a PetscDS object 428 429 Collective on PetscDS 430 431 Input Parameter: 432 . prob - the PetscDS object to destroy 433 434 Level: developer 435 436 .seealso PetscDSView() 437 @*/ 438 PetscErrorCode PetscDSDestroy(PetscDS *prob) 439 { 440 PetscInt f; 441 DSBoundary next; 442 PetscErrorCode ierr; 443 444 PetscFunctionBegin; 445 if (!*prob) PetscFunctionReturn(0); 446 PetscValidHeaderSpecific((*prob), PETSCDS_CLASSID, 1); 447 448 if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; PetscFunctionReturn(0);} 449 ((PetscObject) (*prob))->refct = 0; 450 ierr = PetscDSDestroyStructs_Static(*prob);CHKERRQ(ierr); 451 for (f = 0; f < (*prob)->Nf; ++f) { 452 ierr = PetscObjectDereference((*prob)->disc[f]);CHKERRQ(ierr); 453 } 454 ierr = PetscFree3((*prob)->disc, (*prob)->implicit, (*prob)->adjacency);CHKERRQ(ierr); 455 ierr = PetscFree7((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->gp,(*prob)->gt,(*prob)->r,(*prob)->ctx);CHKERRQ(ierr); 456 ierr = PetscFree((*prob)->update);CHKERRQ(ierr); 457 ierr = PetscFree3((*prob)->fBd,(*prob)->gBd,(*prob)->exactSol);CHKERRQ(ierr); 458 if ((*prob)->ops->destroy) {ierr = (*(*prob)->ops->destroy)(*prob);CHKERRQ(ierr);} 459 next = (*prob)->boundary; 460 while (next) { 461 DSBoundary b = next; 462 463 next = b->next; 464 ierr = PetscFree(b->comps);CHKERRQ(ierr); 465 ierr = PetscFree(b->ids);CHKERRQ(ierr); 466 ierr = PetscFree(b->name);CHKERRQ(ierr); 467 ierr = PetscFree(b->labelname);CHKERRQ(ierr); 468 ierr = PetscFree(b);CHKERRQ(ierr); 469 } 470 ierr = PetscFree((*prob)->constants);CHKERRQ(ierr); 471 ierr = PetscHeaderDestroy(prob);CHKERRQ(ierr); 472 PetscFunctionReturn(0); 473 } 474 475 /*@ 476 PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType(). 477 478 Collective on MPI_Comm 479 480 Input Parameter: 481 . comm - The communicator for the PetscDS object 482 483 Output Parameter: 484 . prob - The PetscDS object 485 486 Level: beginner 487 488 .seealso: PetscDSSetType(), PETSCDSBASIC 489 @*/ 490 PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob) 491 { 492 PetscDS p; 493 PetscErrorCode ierr; 494 495 PetscFunctionBegin; 496 PetscValidPointer(prob, 2); 497 *prob = NULL; 498 ierr = PetscDSInitializePackage();CHKERRQ(ierr); 499 500 ierr = PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);CHKERRQ(ierr); 501 502 p->Nf = 0; 503 p->setup = PETSC_FALSE; 504 p->numConstants = 0; 505 p->constants = NULL; 506 p->dimEmbed = -1; 507 508 *prob = p; 509 PetscFunctionReturn(0); 510 } 511 512 /*@ 513 PetscDSGetNumFields - Returns the number of fields in the DS 514 515 Not collective 516 517 Input Parameter: 518 . prob - The PetscDS object 519 520 Output Parameter: 521 . Nf - The number of fields 522 523 Level: beginner 524 525 .seealso: PetscDSGetSpatialDimension(), PetscDSCreate() 526 @*/ 527 PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf) 528 { 529 PetscFunctionBegin; 530 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 531 PetscValidPointer(Nf, 2); 532 *Nf = prob->Nf; 533 PetscFunctionReturn(0); 534 } 535 536 /*@ 537 PetscDSGetSpatialDimension - Returns the spatial dimension of the DS, meaning the topological dimension of the discretizations 538 539 Not collective 540 541 Input Parameter: 542 . prob - The PetscDS object 543 544 Output Parameter: 545 . dim - The spatial dimension 546 547 Level: beginner 548 549 .seealso: PetscDSGetCoordinateDimension(), PetscDSGetNumFields(), PetscDSCreate() 550 @*/ 551 PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim) 552 { 553 PetscErrorCode ierr; 554 555 PetscFunctionBegin; 556 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 557 PetscValidPointer(dim, 2); 558 *dim = 0; 559 if (prob->Nf) { 560 PetscObject obj; 561 PetscClassId id; 562 563 ierr = PetscDSGetDiscretization(prob, 0, &obj);CHKERRQ(ierr); 564 ierr = PetscObjectGetClassId(obj, &id);CHKERRQ(ierr); 565 if (id == PETSCFE_CLASSID) {ierr = PetscFEGetSpatialDimension((PetscFE) obj, dim);CHKERRQ(ierr);} 566 else if (id == PETSCFV_CLASSID) {ierr = PetscFVGetSpatialDimension((PetscFV) obj, dim);CHKERRQ(ierr);} 567 else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0); 568 } 569 PetscFunctionReturn(0); 570 } 571 572 /*@ 573 PetscDSGetCoordinateDimension - Returns the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded 574 575 Not collective 576 577 Input Parameter: 578 . prob - The PetscDS object 579 580 Output Parameter: 581 . dimEmbed - The coordinate dimension 582 583 Level: beginner 584 585 .seealso: PetscDSSetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate() 586 @*/ 587 PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed) 588 { 589 PetscFunctionBegin; 590 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 591 PetscValidPointer(dimEmbed, 2); 592 if (prob->dimEmbed < 0) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS"); 593 *dimEmbed = prob->dimEmbed; 594 PetscFunctionReturn(0); 595 } 596 597 /*@ 598 PetscDSSetCoordinateDimension - Set the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded 599 600 Not collective 601 602 Input Parameters: 603 + prob - The PetscDS object 604 - dimEmbed - The coordinate dimension 605 606 Level: beginner 607 608 .seealso: PetscDSGetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate() 609 @*/ 610 PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed) 611 { 612 PetscFunctionBegin; 613 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 614 prob->dimEmbed = dimEmbed; 615 PetscFunctionReturn(0); 616 } 617 618 /*@ 619 PetscDSGetTotalDimension - Returns the total size of the approximation space for this system 620 621 Not collective 622 623 Input Parameter: 624 . prob - The PetscDS object 625 626 Output Parameter: 627 . dim - The total problem dimension 628 629 Level: beginner 630 631 .seealso: PetscDSGetNumFields(), PetscDSCreate() 632 @*/ 633 PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim) 634 { 635 PetscErrorCode ierr; 636 637 PetscFunctionBegin; 638 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 639 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 640 PetscValidPointer(dim, 2); 641 *dim = prob->totDim; 642 PetscFunctionReturn(0); 643 } 644 645 /*@ 646 PetscDSGetTotalComponents - Returns the total number of components in this system 647 648 Not collective 649 650 Input Parameter: 651 . prob - The PetscDS object 652 653 Output Parameter: 654 . dim - The total number of components 655 656 Level: beginner 657 658 .seealso: PetscDSGetNumFields(), PetscDSCreate() 659 @*/ 660 PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc) 661 { 662 PetscErrorCode ierr; 663 664 PetscFunctionBegin; 665 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 666 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 667 PetscValidPointer(Nc, 2); 668 *Nc = prob->totComp; 669 PetscFunctionReturn(0); 670 } 671 672 /*@ 673 PetscDSGetDiscretization - Returns the discretization object for the given field 674 675 Not collective 676 677 Input Parameters: 678 + prob - The PetscDS object 679 - f - The field number 680 681 Output Parameter: 682 . disc - The discretization object 683 684 Level: beginner 685 686 .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 687 @*/ 688 PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc) 689 { 690 PetscFunctionBegin; 691 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 692 PetscValidPointer(disc, 3); 693 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); 694 *disc = prob->disc[f]; 695 PetscFunctionReturn(0); 696 } 697 698 /*@ 699 PetscDSSetDiscretization - Sets the discretization object for the given field 700 701 Not collective 702 703 Input Parameters: 704 + prob - The PetscDS object 705 . f - The field number 706 - disc - The discretization object 707 708 Level: beginner 709 710 .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 711 @*/ 712 PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc) 713 { 714 PetscErrorCode ierr; 715 716 PetscFunctionBegin; 717 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 718 PetscValidPointer(disc, 3); 719 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 720 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 721 if (prob->disc[f]) {ierr = PetscObjectDereference(prob->disc[f]);CHKERRQ(ierr);} 722 prob->disc[f] = disc; 723 ierr = PetscObjectReference(disc);CHKERRQ(ierr); 724 { 725 PetscClassId id; 726 727 ierr = PetscObjectGetClassId(disc, &id);CHKERRQ(ierr); 728 if (id == PETSCFV_CLASSID) { 729 prob->implicit[f] = PETSC_FALSE; 730 prob->adjacency[f*2+0] = PETSC_TRUE; 731 prob->adjacency[f*2+1] = PETSC_FALSE; 732 } 733 } 734 PetscFunctionReturn(0); 735 } 736 737 /*@ 738 PetscDSAddDiscretization - Adds a discretization object 739 740 Not collective 741 742 Input Parameters: 743 + prob - The PetscDS object 744 - disc - The boundary discretization object 745 746 Level: beginner 747 748 .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 749 @*/ 750 PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc) 751 { 752 PetscErrorCode ierr; 753 754 PetscFunctionBegin; 755 ierr = PetscDSSetDiscretization(prob, prob->Nf, disc);CHKERRQ(ierr); 756 PetscFunctionReturn(0); 757 } 758 759 /*@ 760 PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX 761 762 Not collective 763 764 Input Parameters: 765 + prob - The PetscDS object 766 - f - The field number 767 768 Output Parameter: 769 . implicit - The flag indicating what kind of solve to use for this field 770 771 Level: developer 772 773 .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 774 @*/ 775 PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit) 776 { 777 PetscFunctionBegin; 778 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 779 PetscValidPointer(implicit, 3); 780 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); 781 *implicit = prob->implicit[f]; 782 PetscFunctionReturn(0); 783 } 784 785 /*@ 786 PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for IMEX 787 788 Not collective 789 790 Input Parameters: 791 + prob - The PetscDS object 792 . f - The field number 793 - implicit - The flag indicating what kind of solve to use for this field 794 795 Level: developer 796 797 .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 798 @*/ 799 PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit) 800 { 801 PetscFunctionBegin; 802 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 803 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); 804 prob->implicit[f] = implicit; 805 PetscFunctionReturn(0); 806 } 807 808 /*@ 809 PetscDSGetAdjacency - Returns the flags for determining variable influence 810 811 Not collective 812 813 Input Parameters: 814 + prob - The PetscDS object 815 - f - The field number 816 817 Output Parameter: 818 + useCone - Flag for variable influence starting with the cone operation 819 - useClosure - Flag for variable influence using transitive closure 820 821 Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure() 822 823 Level: developer 824 825 .seealso: PetscDSSetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 826 @*/ 827 PetscErrorCode PetscDSGetAdjacency(PetscDS prob, PetscInt f, PetscBool *useCone, PetscBool *useClosure) 828 { 829 PetscFunctionBegin; 830 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 831 PetscValidPointer(useCone, 3); 832 PetscValidPointer(useClosure, 4); 833 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); 834 *useCone = prob->adjacency[f*2+0]; 835 *useClosure = prob->adjacency[f*2+1]; 836 PetscFunctionReturn(0); 837 } 838 839 /*@ 840 PetscDSSetAdjacency - Set the flags for determining variable influence 841 842 Not collective 843 844 Input Parameters: 845 + prob - The PetscDS object 846 . f - The field number 847 . useCone - Flag for variable influence starting with the cone operation 848 - useClosure - Flag for variable influence using transitive closure 849 850 Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure() 851 852 Level: developer 853 854 .seealso: PetscDSGetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 855 @*/ 856 PetscErrorCode PetscDSSetAdjacency(PetscDS prob, PetscInt f, PetscBool useCone, PetscBool useClosure) 857 { 858 PetscFunctionBegin; 859 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 860 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); 861 prob->adjacency[f*2+0] = useCone; 862 prob->adjacency[f*2+1] = useClosure; 863 PetscFunctionReturn(0); 864 } 865 866 PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f, 867 void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 868 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 869 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 870 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[])) 871 { 872 PetscFunctionBegin; 873 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 874 PetscValidPointer(obj, 2); 875 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); 876 *obj = prob->obj[f]; 877 PetscFunctionReturn(0); 878 } 879 880 PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f, 881 void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 882 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 883 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 884 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[])) 885 { 886 PetscErrorCode ierr; 887 888 PetscFunctionBegin; 889 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 890 if (obj) PetscValidFunction(obj, 2); 891 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 892 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 893 prob->obj[f] = obj; 894 PetscFunctionReturn(0); 895 } 896 897 /*@C 898 PetscDSGetResidual - Get the pointwise residual function for a given test field 899 900 Not collective 901 902 Input Parameters: 903 + prob - The PetscDS 904 - f - The test field number 905 906 Output Parameters: 907 + f0 - integrand for the test function term 908 - f1 - integrand for the test function gradient term 909 910 Note: We are using a first order FEM model for the weak form: 911 912 \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) 913 914 The calling sequence for the callbacks f0 and f1 is given by: 915 916 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 917 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 918 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 919 $ PetscReal t, const PetscReal x[], PetscScalar f0[]) 920 921 + dim - the spatial dimension 922 . Nf - the number of fields 923 . uOff - the offset into u[] and u_t[] for each field 924 . uOff_x - the offset into u_x[] for each field 925 . u - each field evaluated at the current point 926 . u_t - the time derivative of each field evaluated at the current point 927 . u_x - the gradient of each field evaluated at the current point 928 . aOff - the offset into a[] and a_t[] for each auxiliary field 929 . aOff_x - the offset into a_x[] for each auxiliary field 930 . a - each auxiliary field evaluated at the current point 931 . a_t - the time derivative of each auxiliary field evaluated at the current point 932 . a_x - the gradient of auxiliary each field evaluated at the current point 933 . t - current time 934 . x - coordinates of the current point 935 . numConstants - number of constant parameters 936 . constants - constant parameters 937 - f0 - output values at the current point 938 939 Level: intermediate 940 941 .seealso: PetscDSSetResidual() 942 @*/ 943 PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f, 944 void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 945 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 946 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 947 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), 948 void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 949 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 950 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 951 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])) 952 { 953 PetscFunctionBegin; 954 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 955 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); 956 if (f0) {PetscValidPointer(f0, 3); *f0 = prob->f[f*2+0];} 957 if (f1) {PetscValidPointer(f1, 4); *f1 = prob->f[f*2+1];} 958 PetscFunctionReturn(0); 959 } 960 961 /*@C 962 PetscDSSetResidual - Set the pointwise residual function for a given test field 963 964 Not collective 965 966 Input Parameters: 967 + prob - The PetscDS 968 . f - The test field number 969 . f0 - integrand for the test function term 970 - f1 - integrand for the test function gradient term 971 972 Note: We are using a first order FEM model for the weak form: 973 974 \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) 975 976 The calling sequence for the callbacks f0 and f1 is given by: 977 978 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 979 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 980 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 981 $ PetscReal t, const PetscReal x[], PetscScalar f0[]) 982 983 + dim - the spatial dimension 984 . Nf - the number of fields 985 . uOff - the offset into u[] and u_t[] for each field 986 . uOff_x - the offset into u_x[] for each field 987 . u - each field evaluated at the current point 988 . u_t - the time derivative of each field evaluated at the current point 989 . u_x - the gradient of each field evaluated at the current point 990 . aOff - the offset into a[] and a_t[] for each auxiliary field 991 . aOff_x - the offset into a_x[] for each auxiliary field 992 . a - each auxiliary field evaluated at the current point 993 . a_t - the time derivative of each auxiliary field evaluated at the current point 994 . a_x - the gradient of auxiliary each field evaluated at the current point 995 . t - current time 996 . x - coordinates of the current point 997 . numConstants - number of constant parameters 998 . constants - constant parameters 999 - f0 - output values at the current point 1000 1001 Level: intermediate 1002 1003 .seealso: PetscDSGetResidual() 1004 @*/ 1005 PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f, 1006 void (*f0)(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[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), 1010 void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1011 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1012 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1013 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])) 1014 { 1015 PetscErrorCode ierr; 1016 1017 PetscFunctionBegin; 1018 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1019 if (f0) PetscValidFunction(f0, 3); 1020 if (f1) PetscValidFunction(f1, 4); 1021 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1022 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1023 prob->f[f*2+0] = f0; 1024 prob->f[f*2+1] = f1; 1025 PetscFunctionReturn(0); 1026 } 1027 1028 /*@C 1029 PetscDSHasJacobian - Signals that Jacobian functions have been set 1030 1031 Not collective 1032 1033 Input Parameter: 1034 . prob - The PetscDS 1035 1036 Output Parameter: 1037 . hasJac - flag that pointwise function for the Jacobian has been set 1038 1039 Level: intermediate 1040 1041 .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian() 1042 @*/ 1043 PetscErrorCode PetscDSHasJacobian(PetscDS prob, PetscBool *hasJac) 1044 { 1045 PetscInt f, g, h; 1046 1047 PetscFunctionBegin; 1048 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1049 *hasJac = PETSC_FALSE; 1050 for (f = 0; f < prob->Nf; ++f) { 1051 for (g = 0; g < prob->Nf; ++g) { 1052 for (h = 0; h < 4; ++h) { 1053 if (prob->g[(f*prob->Nf + g)*4+h]) *hasJac = PETSC_TRUE; 1054 } 1055 } 1056 } 1057 PetscFunctionReturn(0); 1058 } 1059 1060 /*@C 1061 PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field 1062 1063 Not collective 1064 1065 Input Parameters: 1066 + prob - The PetscDS 1067 . f - The test field number 1068 - g - The field number 1069 1070 Output Parameters: 1071 + g0 - integrand for the test and basis function term 1072 . g1 - integrand for the test function and basis function gradient term 1073 . g2 - integrand for the test function gradient and basis function term 1074 - g3 - integrand for the test function gradient and basis function gradient term 1075 1076 Note: We are using a first order FEM model for the weak form: 1077 1078 \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 1079 1080 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1081 1082 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1083 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1084 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1085 $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]) 1086 1087 + dim - the spatial dimension 1088 . Nf - the number of fields 1089 . uOff - the offset into u[] and u_t[] for each field 1090 . uOff_x - the offset into u_x[] for each field 1091 . u - each field evaluated at the current point 1092 . u_t - the time derivative of each field evaluated at the current point 1093 . u_x - the gradient of each field evaluated at the current point 1094 . aOff - the offset into a[] and a_t[] for each auxiliary field 1095 . aOff_x - the offset into a_x[] for each auxiliary field 1096 . a - each auxiliary field evaluated at the current point 1097 . a_t - the time derivative of each auxiliary field evaluated at the current point 1098 . a_x - the gradient of auxiliary each field evaluated at the current point 1099 . t - current time 1100 . u_tShift - the multiplier a for dF/dU_t 1101 . x - coordinates of the current point 1102 . numConstants - number of constant parameters 1103 . constants - constant parameters 1104 - g0 - output values at the current point 1105 1106 Level: intermediate 1107 1108 .seealso: PetscDSSetJacobian() 1109 @*/ 1110 PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g, 1111 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1112 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1113 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1114 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1115 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1116 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1117 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1118 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1119 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1120 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1121 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1122 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1123 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1124 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1125 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1126 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1127 { 1128 PetscFunctionBegin; 1129 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1130 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); 1131 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); 1132 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->g[(f*prob->Nf + g)*4+0];} 1133 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->g[(f*prob->Nf + g)*4+1];} 1134 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->g[(f*prob->Nf + g)*4+2];} 1135 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->g[(f*prob->Nf + g)*4+3];} 1136 PetscFunctionReturn(0); 1137 } 1138 1139 /*@C 1140 PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields 1141 1142 Not collective 1143 1144 Input Parameters: 1145 + prob - The PetscDS 1146 . f - The test field number 1147 . g - The field number 1148 . g0 - integrand for the test and basis function term 1149 . g1 - integrand for the test function and basis function gradient term 1150 . g2 - integrand for the test function gradient and basis function term 1151 - g3 - integrand for the test function gradient and basis function gradient term 1152 1153 Note: We are using a first order FEM model for the weak form: 1154 1155 \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 1156 1157 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1158 1159 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1160 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1161 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1162 $ PetscReal t, const PetscReal x[], PetscScalar g0[]) 1163 1164 + dim - the spatial dimension 1165 . Nf - the number of fields 1166 . uOff - the offset into u[] and u_t[] for each field 1167 . uOff_x - the offset into u_x[] for each field 1168 . u - each field evaluated at the current point 1169 . u_t - the time derivative of each field evaluated at the current point 1170 . u_x - the gradient of each field evaluated at the current point 1171 . aOff - the offset into a[] and a_t[] for each auxiliary field 1172 . aOff_x - the offset into a_x[] for each auxiliary field 1173 . a - each auxiliary field evaluated at the current point 1174 . a_t - the time derivative of each auxiliary field evaluated at the current point 1175 . a_x - the gradient of auxiliary each field evaluated at the current point 1176 . t - current time 1177 . u_tShift - the multiplier a for dF/dU_t 1178 . x - coordinates of the current point 1179 . numConstants - number of constant parameters 1180 . constants - constant parameters 1181 - g0 - output values at the current point 1182 1183 Level: intermediate 1184 1185 .seealso: PetscDSGetJacobian() 1186 @*/ 1187 PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g, 1188 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1189 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1190 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1191 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1192 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1193 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1194 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1195 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1196 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1197 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1198 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1199 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1200 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1201 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1202 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1203 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1204 { 1205 PetscErrorCode ierr; 1206 1207 PetscFunctionBegin; 1208 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1209 if (g0) PetscValidFunction(g0, 4); 1210 if (g1) PetscValidFunction(g1, 5); 1211 if (g2) PetscValidFunction(g2, 6); 1212 if (g3) PetscValidFunction(g3, 7); 1213 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1214 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1215 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1216 prob->g[(f*prob->Nf + g)*4+0] = g0; 1217 prob->g[(f*prob->Nf + g)*4+1] = g1; 1218 prob->g[(f*prob->Nf + g)*4+2] = g2; 1219 prob->g[(f*prob->Nf + g)*4+3] = g3; 1220 PetscFunctionReturn(0); 1221 } 1222 1223 /*@C 1224 PetscDSHasJacobianPreconditioner - Signals that a Jacobian preconditioner matrix has been set 1225 1226 Not collective 1227 1228 Input Parameter: 1229 . prob - The PetscDS 1230 1231 Output Parameter: 1232 . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set 1233 1234 Level: intermediate 1235 1236 .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian() 1237 @*/ 1238 PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS prob, PetscBool *hasJacPre) 1239 { 1240 PetscInt f, g, h; 1241 1242 PetscFunctionBegin; 1243 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1244 *hasJacPre = PETSC_FALSE; 1245 for (f = 0; f < prob->Nf; ++f) { 1246 for (g = 0; g < prob->Nf; ++g) { 1247 for (h = 0; h < 4; ++h) { 1248 if (prob->gp[(f*prob->Nf + g)*4+h]) *hasJacPre = PETSC_TRUE; 1249 } 1250 } 1251 } 1252 PetscFunctionReturn(0); 1253 } 1254 1255 /*@C 1256 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. 1257 1258 Not collective 1259 1260 Input Parameters: 1261 + prob - The PetscDS 1262 . f - The test field number 1263 - g - The field number 1264 1265 Output Parameters: 1266 + g0 - integrand for the test and basis function term 1267 . g1 - integrand for the test function and basis function gradient term 1268 . g2 - integrand for the test function gradient and basis function term 1269 - g3 - integrand for the test function gradient and basis function gradient term 1270 1271 Note: We are using a first order FEM model for the weak form: 1272 1273 \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 1274 1275 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1276 1277 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1278 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1279 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1280 $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]) 1281 1282 + dim - the spatial dimension 1283 . Nf - the number of fields 1284 . uOff - the offset into u[] and u_t[] for each field 1285 . uOff_x - the offset into u_x[] for each field 1286 . u - each field evaluated at the current point 1287 . u_t - the time derivative of each field evaluated at the current point 1288 . u_x - the gradient of each field evaluated at the current point 1289 . aOff - the offset into a[] and a_t[] for each auxiliary field 1290 . aOff_x - the offset into a_x[] for each auxiliary field 1291 . a - each auxiliary field evaluated at the current point 1292 . a_t - the time derivative of each auxiliary field evaluated at the current point 1293 . a_x - the gradient of auxiliary each field evaluated at the current point 1294 . t - current time 1295 . u_tShift - the multiplier a for dF/dU_t 1296 . x - coordinates of the current point 1297 . numConstants - number of constant parameters 1298 . constants - constant parameters 1299 - g0 - output values at the current point 1300 1301 Level: intermediate 1302 1303 .seealso: PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian() 1304 @*/ 1305 PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g, 1306 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1307 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1308 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1309 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1310 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1311 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1312 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1313 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1314 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1315 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1316 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1317 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1318 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1319 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1320 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1321 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1322 { 1323 PetscFunctionBegin; 1324 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1325 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); 1326 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); 1327 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->gp[(f*prob->Nf + g)*4+0];} 1328 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->gp[(f*prob->Nf + g)*4+1];} 1329 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->gp[(f*prob->Nf + g)*4+2];} 1330 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->gp[(f*prob->Nf + g)*4+3];} 1331 PetscFunctionReturn(0); 1332 } 1333 1334 /*@C 1335 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. 1336 1337 Not collective 1338 1339 Input Parameters: 1340 + prob - The PetscDS 1341 . f - The test field number 1342 . g - The field number 1343 . g0 - integrand for the test and basis function term 1344 . g1 - integrand for the test function and basis function gradient term 1345 . g2 - integrand for the test function gradient and basis function term 1346 - g3 - integrand for the test function gradient and basis function gradient term 1347 1348 Note: We are using a first order FEM model for the weak form: 1349 1350 \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 1351 1352 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1353 1354 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1355 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1356 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1357 $ PetscReal t, const PetscReal x[], PetscScalar g0[]) 1358 1359 + dim - the spatial dimension 1360 . Nf - the number of fields 1361 . uOff - the offset into u[] and u_t[] for each field 1362 . uOff_x - the offset into u_x[] for each field 1363 . u - each field evaluated at the current point 1364 . u_t - the time derivative of each field evaluated at the current point 1365 . u_x - the gradient of each field evaluated at the current point 1366 . aOff - the offset into a[] and a_t[] for each auxiliary field 1367 . aOff_x - the offset into a_x[] for each auxiliary field 1368 . a - each auxiliary field evaluated at the current point 1369 . a_t - the time derivative of each auxiliary field evaluated at the current point 1370 . a_x - the gradient of auxiliary each field evaluated at the current point 1371 . t - current time 1372 . u_tShift - the multiplier a for dF/dU_t 1373 . x - coordinates of the current point 1374 . numConstants - number of constant parameters 1375 . constants - constant parameters 1376 - g0 - output values at the current point 1377 1378 Level: intermediate 1379 1380 .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobian() 1381 @*/ 1382 PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g, 1383 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1384 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1385 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1386 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1387 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1388 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1389 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1390 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1391 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1392 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1393 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1394 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1395 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1396 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1397 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1398 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1399 { 1400 PetscErrorCode ierr; 1401 1402 PetscFunctionBegin; 1403 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1404 if (g0) PetscValidFunction(g0, 4); 1405 if (g1) PetscValidFunction(g1, 5); 1406 if (g2) PetscValidFunction(g2, 6); 1407 if (g3) PetscValidFunction(g3, 7); 1408 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1409 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1410 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1411 prob->gp[(f*prob->Nf + g)*4+0] = g0; 1412 prob->gp[(f*prob->Nf + g)*4+1] = g1; 1413 prob->gp[(f*prob->Nf + g)*4+2] = g2; 1414 prob->gp[(f*prob->Nf + g)*4+3] = g3; 1415 PetscFunctionReturn(0); 1416 } 1417 1418 /*@C 1419 PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set 1420 1421 Not collective 1422 1423 Input Parameter: 1424 . prob - The PetscDS 1425 1426 Output Parameter: 1427 . hasDynJac - flag that pointwise function for dynamic Jacobian has been set 1428 1429 Level: intermediate 1430 1431 .seealso: PetscDSGetDynamicJacobian(), PetscDSSetDynamicJacobian(), PetscDSGetJacobian() 1432 @*/ 1433 PetscErrorCode PetscDSHasDynamicJacobian(PetscDS prob, PetscBool *hasDynJac) 1434 { 1435 PetscInt f, g, h; 1436 1437 PetscFunctionBegin; 1438 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1439 *hasDynJac = PETSC_FALSE; 1440 for (f = 0; f < prob->Nf; ++f) { 1441 for (g = 0; g < prob->Nf; ++g) { 1442 for (h = 0; h < 4; ++h) { 1443 if (prob->gt[(f*prob->Nf + g)*4+h]) *hasDynJac = PETSC_TRUE; 1444 } 1445 } 1446 } 1447 PetscFunctionReturn(0); 1448 } 1449 1450 /*@C 1451 PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field 1452 1453 Not collective 1454 1455 Input Parameters: 1456 + prob - The PetscDS 1457 . f - The test field number 1458 - g - The field number 1459 1460 Output Parameters: 1461 + g0 - integrand for the test and basis function term 1462 . g1 - integrand for the test function and basis function gradient term 1463 . g2 - integrand for the test function gradient and basis function term 1464 - g3 - integrand for the test function gradient and basis function gradient term 1465 1466 Note: We are using a first order FEM model for the weak form: 1467 1468 \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 1469 1470 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1471 1472 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1473 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1474 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1475 $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]) 1476 1477 + dim - the spatial dimension 1478 . Nf - the number of fields 1479 . uOff - the offset into u[] and u_t[] for each field 1480 . uOff_x - the offset into u_x[] for each field 1481 . u - each field evaluated at the current point 1482 . u_t - the time derivative of each field evaluated at the current point 1483 . u_x - the gradient of each field evaluated at the current point 1484 . aOff - the offset into a[] and a_t[] for each auxiliary field 1485 . aOff_x - the offset into a_x[] for each auxiliary field 1486 . a - each auxiliary field evaluated at the current point 1487 . a_t - the time derivative of each auxiliary field evaluated at the current point 1488 . a_x - the gradient of auxiliary each field evaluated at the current point 1489 . t - current time 1490 . u_tShift - the multiplier a for dF/dU_t 1491 . x - coordinates of the current point 1492 . numConstants - number of constant parameters 1493 . constants - constant parameters 1494 - g0 - output values at the current point 1495 1496 Level: intermediate 1497 1498 .seealso: PetscDSSetJacobian() 1499 @*/ 1500 PetscErrorCode PetscDSGetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g, 1501 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1502 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1503 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1504 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1505 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1506 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1507 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1508 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1509 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1510 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1511 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1512 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1513 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1514 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1515 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1516 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1517 { 1518 PetscFunctionBegin; 1519 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1520 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); 1521 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); 1522 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->gt[(f*prob->Nf + g)*4+0];} 1523 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->gt[(f*prob->Nf + g)*4+1];} 1524 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->gt[(f*prob->Nf + g)*4+2];} 1525 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->gt[(f*prob->Nf + g)*4+3];} 1526 PetscFunctionReturn(0); 1527 } 1528 1529 /*@C 1530 PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields 1531 1532 Not collective 1533 1534 Input Parameters: 1535 + prob - The PetscDS 1536 . f - The test field number 1537 . g - The field number 1538 . g0 - integrand for the test and basis function term 1539 . g1 - integrand for the test function and basis function gradient term 1540 . g2 - integrand for the test function gradient and basis function term 1541 - g3 - integrand for the test function gradient and basis function gradient term 1542 1543 Note: We are using a first order FEM model for the weak form: 1544 1545 \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 1546 1547 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1548 1549 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1550 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1551 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1552 $ PetscReal t, const PetscReal x[], PetscScalar g0[]) 1553 1554 + dim - the spatial dimension 1555 . Nf - the number of fields 1556 . uOff - the offset into u[] and u_t[] for each field 1557 . uOff_x - the offset into u_x[] for each field 1558 . u - each field evaluated at the current point 1559 . u_t - the time derivative of each field evaluated at the current point 1560 . u_x - the gradient of each field evaluated at the current point 1561 . aOff - the offset into a[] and a_t[] for each auxiliary field 1562 . aOff_x - the offset into a_x[] for each auxiliary field 1563 . a - each auxiliary field evaluated at the current point 1564 . a_t - the time derivative of each auxiliary field evaluated at the current point 1565 . a_x - the gradient of auxiliary each field evaluated at the current point 1566 . t - current time 1567 . u_tShift - the multiplier a for dF/dU_t 1568 . x - coordinates of the current point 1569 . numConstants - number of constant parameters 1570 . constants - constant parameters 1571 - g0 - output values at the current point 1572 1573 Level: intermediate 1574 1575 .seealso: PetscDSGetJacobian() 1576 @*/ 1577 PetscErrorCode PetscDSSetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g, 1578 void (*g0)(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, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1582 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1583 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1584 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1585 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1586 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1587 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1588 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1589 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1590 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1591 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1592 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1593 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1594 { 1595 PetscErrorCode ierr; 1596 1597 PetscFunctionBegin; 1598 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1599 if (g0) PetscValidFunction(g0, 4); 1600 if (g1) PetscValidFunction(g1, 5); 1601 if (g2) PetscValidFunction(g2, 6); 1602 if (g3) PetscValidFunction(g3, 7); 1603 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1604 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1605 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1606 prob->gt[(f*prob->Nf + g)*4+0] = g0; 1607 prob->gt[(f*prob->Nf + g)*4+1] = g1; 1608 prob->gt[(f*prob->Nf + g)*4+2] = g2; 1609 prob->gt[(f*prob->Nf + g)*4+3] = g3; 1610 PetscFunctionReturn(0); 1611 } 1612 1613 /*@C 1614 PetscDSGetRiemannSolver - Returns the Riemann solver for the given field 1615 1616 Not collective 1617 1618 Input Arguments: 1619 + prob - The PetscDS object 1620 - f - The field number 1621 1622 Output Argument: 1623 . r - Riemann solver 1624 1625 Calling sequence for r: 1626 1627 $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx) 1628 1629 + dim - The spatial dimension 1630 . Nf - The number of fields 1631 . x - The coordinates at a point on the interface 1632 . n - The normal vector to the interface 1633 . uL - The state vector to the left of the interface 1634 . uR - The state vector to the right of the interface 1635 . flux - output array of flux through the interface 1636 . numConstants - number of constant parameters 1637 . constants - constant parameters 1638 - ctx - optional user context 1639 1640 Level: intermediate 1641 1642 .seealso: PetscDSSetRiemannSolver() 1643 @*/ 1644 PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f, 1645 void (**r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx)) 1646 { 1647 PetscFunctionBegin; 1648 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1649 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); 1650 PetscValidPointer(r, 3); 1651 *r = prob->r[f]; 1652 PetscFunctionReturn(0); 1653 } 1654 1655 /*@C 1656 PetscDSSetRiemannSolver - Sets the Riemann solver for the given field 1657 1658 Not collective 1659 1660 Input Arguments: 1661 + prob - The PetscDS object 1662 . f - The field number 1663 - r - Riemann solver 1664 1665 Calling sequence for r: 1666 1667 $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx) 1668 1669 + dim - The spatial dimension 1670 . Nf - The number of fields 1671 . x - The coordinates at a point on the interface 1672 . n - The normal vector to the interface 1673 . uL - The state vector to the left of the interface 1674 . uR - The state vector to the right of the interface 1675 . flux - output array of flux through the interface 1676 . numConstants - number of constant parameters 1677 . constants - constant parameters 1678 - ctx - optional user context 1679 1680 Level: intermediate 1681 1682 .seealso: PetscDSGetRiemannSolver() 1683 @*/ 1684 PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f, 1685 void (*r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx)) 1686 { 1687 PetscErrorCode ierr; 1688 1689 PetscFunctionBegin; 1690 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1691 if (r) PetscValidFunction(r, 3); 1692 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1693 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1694 prob->r[f] = r; 1695 PetscFunctionReturn(0); 1696 } 1697 1698 /*@C 1699 PetscDSGetUpdate - Get the pointwise update function for a given field 1700 1701 Not collective 1702 1703 Input Parameters: 1704 + prob - The PetscDS 1705 - f - The field number 1706 1707 Output Parameters: 1708 + update - update function 1709 1710 Note: The calling sequence for the callback update is given by: 1711 1712 $ update(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[], PetscScalar uNew[]) 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 . x - coordinates of the current point 1731 - uNew - new value for field at the current point 1732 1733 Level: intermediate 1734 1735 .seealso: PetscDSSetUpdate(), PetscDSSetResidual() 1736 @*/ 1737 PetscErrorCode PetscDSGetUpdate(PetscDS prob, PetscInt f, 1738 void (**update)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1739 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1740 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1741 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[])) 1742 { 1743 PetscFunctionBegin; 1744 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1745 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); 1746 if (update) {PetscValidPointer(update, 3); *update = prob->update[f];} 1747 PetscFunctionReturn(0); 1748 } 1749 1750 /*@C 1751 PetscDSSetUpdate - Set the pointwise update function for a given field 1752 1753 Not collective 1754 1755 Input Parameters: 1756 + prob - The PetscDS 1757 . f - The field number 1758 - update - update function 1759 1760 Note: The calling sequence for the callback update is given by: 1761 1762 $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1763 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1764 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1765 $ PetscReal t, const PetscReal x[], PetscScalar uNew[]) 1766 1767 + dim - the spatial dimension 1768 . Nf - the number of fields 1769 . uOff - the offset into u[] and u_t[] for each field 1770 . uOff_x - the offset into u_x[] for each field 1771 . u - each field evaluated at the current point 1772 . u_t - the time derivative of each field evaluated at the current point 1773 . u_x - the gradient of each field evaluated at the current point 1774 . aOff - the offset into a[] and a_t[] for each auxiliary field 1775 . aOff_x - the offset into a_x[] for each auxiliary field 1776 . a - each auxiliary field evaluated at the current point 1777 . a_t - the time derivative of each auxiliary field evaluated at the current point 1778 . a_x - the gradient of auxiliary each field evaluated at the current point 1779 . t - current time 1780 . x - coordinates of the current point 1781 - uNew - new field values at the current point 1782 1783 Level: intermediate 1784 1785 .seealso: PetscDSGetResidual() 1786 @*/ 1787 PetscErrorCode PetscDSSetUpdate(PetscDS prob, PetscInt f, 1788 void (*update)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1789 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1790 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1791 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[])) 1792 { 1793 PetscErrorCode ierr; 1794 1795 PetscFunctionBegin; 1796 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1797 if (update) PetscValidFunction(update, 3); 1798 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1799 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1800 prob->update[f] = update; 1801 PetscFunctionReturn(0); 1802 } 1803 1804 PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx) 1805 { 1806 PetscFunctionBegin; 1807 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1808 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); 1809 PetscValidPointer(ctx, 3); 1810 *ctx = prob->ctx[f]; 1811 PetscFunctionReturn(0); 1812 } 1813 1814 PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx) 1815 { 1816 PetscErrorCode ierr; 1817 1818 PetscFunctionBegin; 1819 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1820 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1821 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1822 prob->ctx[f] = ctx; 1823 PetscFunctionReturn(0); 1824 } 1825 1826 /*@C 1827 PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field 1828 1829 Not collective 1830 1831 Input Parameters: 1832 + prob - The PetscDS 1833 - f - The test field number 1834 1835 Output Parameters: 1836 + f0 - boundary integrand for the test function term 1837 - f1 - boundary integrand for the test function gradient term 1838 1839 Note: We are using a first order FEM model for the weak form: 1840 1841 \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 1842 1843 The calling sequence for the callbacks f0 and f1 is given by: 1844 1845 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1846 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1847 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1848 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]) 1849 1850 + dim - the spatial dimension 1851 . Nf - the number of fields 1852 . uOff - the offset into u[] and u_t[] for each field 1853 . uOff_x - the offset into u_x[] for each field 1854 . u - each field evaluated at the current point 1855 . u_t - the time derivative of each field evaluated at the current point 1856 . u_x - the gradient of each field evaluated at the current point 1857 . aOff - the offset into a[] and a_t[] for each auxiliary field 1858 . aOff_x - the offset into a_x[] for each auxiliary field 1859 . a - each auxiliary field evaluated at the current point 1860 . a_t - the time derivative of each auxiliary field evaluated at the current point 1861 . a_x - the gradient of auxiliary each field evaluated at the current point 1862 . t - current time 1863 . x - coordinates of the current point 1864 . n - unit normal at the current point 1865 . numConstants - number of constant parameters 1866 . constants - constant parameters 1867 - f0 - output values at the current point 1868 1869 Level: intermediate 1870 1871 .seealso: PetscDSSetBdResidual() 1872 @*/ 1873 PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f, 1874 void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1875 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1876 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1877 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), 1878 void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1879 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1880 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1881 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])) 1882 { 1883 PetscFunctionBegin; 1884 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1885 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); 1886 if (f0) {PetscValidPointer(f0, 3); *f0 = prob->fBd[f*2+0];} 1887 if (f1) {PetscValidPointer(f1, 4); *f1 = prob->fBd[f*2+1];} 1888 PetscFunctionReturn(0); 1889 } 1890 1891 /*@C 1892 PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field 1893 1894 Not collective 1895 1896 Input Parameters: 1897 + prob - The PetscDS 1898 . f - The test field number 1899 . f0 - boundary integrand for the test function term 1900 - f1 - boundary integrand for the test function gradient term 1901 1902 Note: We are using a first order FEM model for the weak form: 1903 1904 \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 1905 1906 The calling sequence for the callbacks f0 and f1 is given by: 1907 1908 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1909 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1910 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1911 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]) 1912 1913 + dim - the spatial dimension 1914 . Nf - the number of fields 1915 . uOff - the offset into u[] and u_t[] for each field 1916 . uOff_x - the offset into u_x[] for each field 1917 . u - each field evaluated at the current point 1918 . u_t - the time derivative of each field evaluated at the current point 1919 . u_x - the gradient of each field evaluated at the current point 1920 . aOff - the offset into a[] and a_t[] for each auxiliary field 1921 . aOff_x - the offset into a_x[] for each auxiliary field 1922 . a - each auxiliary field evaluated at the current point 1923 . a_t - the time derivative of each auxiliary field evaluated at the current point 1924 . a_x - the gradient of auxiliary each field evaluated at the current point 1925 . t - current time 1926 . x - coordinates of the current point 1927 . n - unit normal at the current point 1928 . numConstants - number of constant parameters 1929 . constants - constant parameters 1930 - f0 - output values at the current point 1931 1932 Level: intermediate 1933 1934 .seealso: PetscDSGetBdResidual() 1935 @*/ 1936 PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f, 1937 void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1938 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1939 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1940 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), 1941 void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1942 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1943 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1944 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])) 1945 { 1946 PetscErrorCode ierr; 1947 1948 PetscFunctionBegin; 1949 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1950 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1951 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1952 if (f0) {PetscValidFunction(f0, 3); prob->fBd[f*2+0] = f0;} 1953 if (f1) {PetscValidFunction(f1, 4); prob->fBd[f*2+1] = f1;} 1954 PetscFunctionReturn(0); 1955 } 1956 1957 /*@C 1958 PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field 1959 1960 Not collective 1961 1962 Input Parameters: 1963 + prob - The PetscDS 1964 . f - The test field number 1965 - g - The field number 1966 1967 Output Parameters: 1968 + g0 - integrand for the test and basis function term 1969 . g1 - integrand for the test function and basis function gradient term 1970 . g2 - integrand for the test function gradient and basis function term 1971 - g3 - integrand for the test function gradient and basis function gradient term 1972 1973 Note: We are using a first order FEM model for the weak form: 1974 1975 \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 1976 1977 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1978 1979 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1980 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1981 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1982 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[]) 1983 1984 + dim - the spatial dimension 1985 . Nf - the number of fields 1986 . uOff - the offset into u[] and u_t[] for each field 1987 . uOff_x - the offset into u_x[] for each field 1988 . u - each field evaluated at the current point 1989 . u_t - the time derivative of each field evaluated at the current point 1990 . u_x - the gradient of each field evaluated at the current point 1991 . aOff - the offset into a[] and a_t[] for each auxiliary field 1992 . aOff_x - the offset into a_x[] for each auxiliary field 1993 . a - each auxiliary field evaluated at the current point 1994 . a_t - the time derivative of each auxiliary field evaluated at the current point 1995 . a_x - the gradient of auxiliary each field evaluated at the current point 1996 . t - current time 1997 . u_tShift - the multiplier a for dF/dU_t 1998 . x - coordinates of the current point 1999 . n - normal at the current point 2000 . numConstants - number of constant parameters 2001 . constants - constant parameters 2002 - g0 - output values at the current point 2003 2004 Level: intermediate 2005 2006 .seealso: PetscDSSetBdJacobian() 2007 @*/ 2008 PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g, 2009 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 2010 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 2011 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 2012 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 2013 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 2014 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 2015 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 2016 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 2017 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 2018 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 2019 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 2020 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 2021 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 2022 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 2023 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 2024 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 2025 { 2026 PetscFunctionBegin; 2027 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2028 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); 2029 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); 2030 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->gBd[(f*prob->Nf + g)*4+0];} 2031 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->gBd[(f*prob->Nf + g)*4+1];} 2032 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->gBd[(f*prob->Nf + g)*4+2];} 2033 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->gBd[(f*prob->Nf + g)*4+3];} 2034 PetscFunctionReturn(0); 2035 } 2036 2037 /*@C 2038 PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field 2039 2040 Not collective 2041 2042 Input Parameters: 2043 + prob - The PetscDS 2044 . f - The test field number 2045 . g - The field number 2046 . g0 - integrand for the test and basis function term 2047 . g1 - integrand for the test function and basis function gradient term 2048 . g2 - integrand for the test function gradient and basis function term 2049 - g3 - integrand for the test function gradient and basis function gradient term 2050 2051 Note: We are using a first order FEM model for the weak form: 2052 2053 \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 2054 2055 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 2056 2057 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 2058 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 2059 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 2060 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[]) 2061 2062 + dim - the spatial dimension 2063 . Nf - the number of fields 2064 . uOff - the offset into u[] and u_t[] for each field 2065 . uOff_x - the offset into u_x[] for each field 2066 . u - each field evaluated at the current point 2067 . u_t - the time derivative of each field evaluated at the current point 2068 . u_x - the gradient of each field evaluated at the current point 2069 . aOff - the offset into a[] and a_t[] for each auxiliary field 2070 . aOff_x - the offset into a_x[] for each auxiliary field 2071 . a - each auxiliary field evaluated at the current point 2072 . a_t - the time derivative of each auxiliary field evaluated at the current point 2073 . a_x - the gradient of auxiliary each field evaluated at the current point 2074 . t - current time 2075 . u_tShift - the multiplier a for dF/dU_t 2076 . x - coordinates of the current point 2077 . n - normal at the current point 2078 . numConstants - number of constant parameters 2079 . constants - constant parameters 2080 - g0 - output values at the current point 2081 2082 Level: intermediate 2083 2084 .seealso: PetscDSGetBdJacobian() 2085 @*/ 2086 PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g, 2087 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 2088 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 2089 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 2090 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 2091 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 2092 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 2093 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 2094 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 2095 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 2096 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 2097 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 2098 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 2099 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 2100 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 2101 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 2102 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 2103 { 2104 PetscErrorCode ierr; 2105 2106 PetscFunctionBegin; 2107 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2108 if (g0) PetscValidFunction(g0, 4); 2109 if (g1) PetscValidFunction(g1, 5); 2110 if (g2) PetscValidFunction(g2, 6); 2111 if (g3) PetscValidFunction(g3, 7); 2112 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 2113 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 2114 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 2115 prob->gBd[(f*prob->Nf + g)*4+0] = g0; 2116 prob->gBd[(f*prob->Nf + g)*4+1] = g1; 2117 prob->gBd[(f*prob->Nf + g)*4+2] = g2; 2118 prob->gBd[(f*prob->Nf + g)*4+3] = g3; 2119 PetscFunctionReturn(0); 2120 } 2121 2122 /*@C 2123 PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field 2124 2125 Not collective 2126 2127 Input Parameters: 2128 + prob - The PetscDS 2129 - f - The test field number 2130 2131 Output Parameter: 2132 . exactSol - exact solution for the test field 2133 2134 Note: The calling sequence for the solution functions is given by: 2135 2136 $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) 2137 2138 + dim - the spatial dimension 2139 . t - current time 2140 . x - coordinates of the current point 2141 . Nc - the number of field components 2142 . u - the solution field evaluated at the current point 2143 - ctx - a user context 2144 2145 Level: intermediate 2146 2147 .seealso: PetscDSSetExactSolution() 2148 @*/ 2149 PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)) 2150 { 2151 PetscFunctionBegin; 2152 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2153 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); 2154 if (sol) {PetscValidPointer(sol, 3); *sol = prob->exactSol[f];} 2155 PetscFunctionReturn(0); 2156 } 2157 2158 /*@C 2159 PetscDSSetExactSolution - Get the pointwise exact solution function for a given test field 2160 2161 Not collective 2162 2163 Input Parameters: 2164 + prob - The PetscDS 2165 . f - The test field number 2166 - sol - solution function for the test fields 2167 2168 Note: The calling sequence for solution functions is given by: 2169 2170 $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) 2171 2172 + dim - the spatial dimension 2173 . t - current time 2174 . x - coordinates of the current point 2175 . Nc - the number of field components 2176 . u - the solution field evaluated at the current point 2177 - ctx - a user context 2178 2179 Level: intermediate 2180 2181 .seealso: PetscDSGetExactSolution() 2182 @*/ 2183 PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)) 2184 { 2185 PetscErrorCode ierr; 2186 2187 PetscFunctionBegin; 2188 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2189 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 2190 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 2191 if (sol) {PetscValidFunction(sol, 3); prob->exactSol[f] = sol;} 2192 PetscFunctionReturn(0); 2193 } 2194 2195 /*@C 2196 PetscDSGetConstants - Returns the array of constants passed to point functions 2197 2198 Not collective 2199 2200 Input Parameter: 2201 . prob - The PetscDS object 2202 2203 Output Parameters: 2204 + numConstants - The number of constants 2205 - constants - The array of constants, NULL if there are none 2206 2207 Level: intermediate 2208 2209 .seealso: PetscDSSetConstants(), PetscDSCreate() 2210 @*/ 2211 PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[]) 2212 { 2213 PetscFunctionBegin; 2214 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2215 if (numConstants) {PetscValidPointer(numConstants, 2); *numConstants = prob->numConstants;} 2216 if (constants) {PetscValidPointer(constants, 3); *constants = prob->constants;} 2217 PetscFunctionReturn(0); 2218 } 2219 2220 /*@C 2221 PetscDSSetConstants - Set the array of constants passed to point functions 2222 2223 Not collective 2224 2225 Input Parameters: 2226 + prob - The PetscDS object 2227 . numConstants - The number of constants 2228 - constants - The array of constants, NULL if there are none 2229 2230 Level: intermediate 2231 2232 .seealso: PetscDSGetConstants(), PetscDSCreate() 2233 @*/ 2234 PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[]) 2235 { 2236 PetscErrorCode ierr; 2237 2238 PetscFunctionBegin; 2239 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2240 if (numConstants != prob->numConstants) { 2241 ierr = PetscFree(prob->constants);CHKERRQ(ierr); 2242 prob->numConstants = numConstants; 2243 if (prob->numConstants) { 2244 ierr = PetscMalloc1(prob->numConstants, &prob->constants);CHKERRQ(ierr); 2245 } else { 2246 prob->constants = NULL; 2247 } 2248 } 2249 if (prob->numConstants) { 2250 PetscValidPointer(constants, 3); 2251 ierr = PetscMemcpy(prob->constants, constants, prob->numConstants * sizeof(PetscScalar));CHKERRQ(ierr); 2252 } 2253 PetscFunctionReturn(0); 2254 } 2255 2256 /*@ 2257 PetscDSGetFieldIndex - Returns the index of the given field 2258 2259 Not collective 2260 2261 Input Parameters: 2262 + prob - The PetscDS object 2263 - disc - The discretization object 2264 2265 Output Parameter: 2266 . f - The field number 2267 2268 Level: beginner 2269 2270 .seealso: PetscGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 2271 @*/ 2272 PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f) 2273 { 2274 PetscInt g; 2275 2276 PetscFunctionBegin; 2277 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2278 PetscValidPointer(f, 3); 2279 *f = -1; 2280 for (g = 0; g < prob->Nf; ++g) {if (disc == prob->disc[g]) break;} 2281 if (g == prob->Nf) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS."); 2282 *f = g; 2283 PetscFunctionReturn(0); 2284 } 2285 2286 /*@ 2287 PetscDSGetFieldSize - Returns the size of the given field in the full space basis 2288 2289 Not collective 2290 2291 Input Parameters: 2292 + prob - The PetscDS object 2293 - f - The field number 2294 2295 Output Parameter: 2296 . size - The size 2297 2298 Level: beginner 2299 2300 .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 2301 @*/ 2302 PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size) 2303 { 2304 PetscErrorCode ierr; 2305 2306 PetscFunctionBegin; 2307 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2308 PetscValidPointer(size, 3); 2309 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); 2310 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2311 *size = prob->Nb[f]; 2312 PetscFunctionReturn(0); 2313 } 2314 2315 /*@ 2316 PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis 2317 2318 Not collective 2319 2320 Input Parameters: 2321 + prob - The PetscDS object 2322 - f - The field number 2323 2324 Output Parameter: 2325 . off - The offset 2326 2327 Level: beginner 2328 2329 .seealso: PetscDSGetFieldSize(), PetscDSGetNumFields(), PetscDSCreate() 2330 @*/ 2331 PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off) 2332 { 2333 PetscInt size, g; 2334 PetscErrorCode ierr; 2335 2336 PetscFunctionBegin; 2337 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2338 PetscValidPointer(off, 3); 2339 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); 2340 *off = 0; 2341 for (g = 0; g < f; ++g) { 2342 ierr = PetscDSGetFieldSize(prob, g, &size);CHKERRQ(ierr); 2343 *off += size; 2344 } 2345 PetscFunctionReturn(0); 2346 } 2347 2348 /*@ 2349 PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point 2350 2351 Not collective 2352 2353 Input Parameter: 2354 . prob - The PetscDS object 2355 2356 Output Parameter: 2357 . dimensions - The number of dimensions 2358 2359 Level: beginner 2360 2361 .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate() 2362 @*/ 2363 PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[]) 2364 { 2365 PetscErrorCode ierr; 2366 2367 PetscFunctionBegin; 2368 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2369 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2370 PetscValidPointer(dimensions, 2); 2371 *dimensions = prob->Nb; 2372 PetscFunctionReturn(0); 2373 } 2374 2375 /*@ 2376 PetscDSGetComponents - Returns the number of components for each field on an evaluation point 2377 2378 Not collective 2379 2380 Input Parameter: 2381 . prob - The PetscDS object 2382 2383 Output Parameter: 2384 . components - The number of components 2385 2386 Level: beginner 2387 2388 .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate() 2389 @*/ 2390 PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[]) 2391 { 2392 PetscErrorCode ierr; 2393 2394 PetscFunctionBegin; 2395 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2396 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2397 PetscValidPointer(components, 2); 2398 *components = prob->Nc; 2399 PetscFunctionReturn(0); 2400 } 2401 2402 /*@ 2403 PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point 2404 2405 Not collective 2406 2407 Input Parameters: 2408 + prob - The PetscDS object 2409 - f - The field number 2410 2411 Output Parameter: 2412 . off - The offset 2413 2414 Level: beginner 2415 2416 .seealso: PetscDSGetNumFields(), PetscDSCreate() 2417 @*/ 2418 PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off) 2419 { 2420 PetscFunctionBegin; 2421 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2422 PetscValidPointer(off, 3); 2423 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); 2424 *off = prob->off[f]; 2425 PetscFunctionReturn(0); 2426 } 2427 2428 /*@ 2429 PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point 2430 2431 Not collective 2432 2433 Input Parameter: 2434 . prob - The PetscDS object 2435 2436 Output Parameter: 2437 . offsets - The offsets 2438 2439 Level: beginner 2440 2441 .seealso: PetscDSGetNumFields(), PetscDSCreate() 2442 @*/ 2443 PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[]) 2444 { 2445 PetscFunctionBegin; 2446 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2447 PetscValidPointer(offsets, 2); 2448 *offsets = prob->off; 2449 PetscFunctionReturn(0); 2450 } 2451 2452 /*@ 2453 PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point 2454 2455 Not collective 2456 2457 Input Parameter: 2458 . prob - The PetscDS object 2459 2460 Output Parameter: 2461 . offsets - The offsets 2462 2463 Level: beginner 2464 2465 .seealso: PetscDSGetNumFields(), PetscDSCreate() 2466 @*/ 2467 PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[]) 2468 { 2469 PetscFunctionBegin; 2470 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2471 PetscValidPointer(offsets, 2); 2472 *offsets = prob->offDer; 2473 PetscFunctionReturn(0); 2474 } 2475 2476 /*@C 2477 PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization 2478 2479 Not collective 2480 2481 Input Parameter: 2482 . prob - The PetscDS object 2483 2484 Output Parameters: 2485 + basis - The basis function tabulation at quadrature points 2486 - basisDer - The basis function derivative tabulation at quadrature points 2487 2488 Level: intermediate 2489 2490 .seealso: PetscDSCreate() 2491 @*/ 2492 PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer) 2493 { 2494 PetscErrorCode ierr; 2495 2496 PetscFunctionBegin; 2497 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2498 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2499 if (basis) {PetscValidPointer(basis, 2); *basis = prob->basis;} 2500 if (basisDer) {PetscValidPointer(basisDer, 3); *basisDer = prob->basisDer;} 2501 PetscFunctionReturn(0); 2502 } 2503 2504 /*@C 2505 PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces 2506 2507 Not collective 2508 2509 Input Parameter: 2510 . prob - The PetscDS object 2511 2512 Output Parameters: 2513 + basisFace - The basis function tabulation at quadrature points 2514 - basisDerFace - The basis function derivative tabulation at quadrature points 2515 2516 Level: intermediate 2517 2518 .seealso: PetscDSGetTabulation(), PetscDSCreate() 2519 @*/ 2520 PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer) 2521 { 2522 PetscErrorCode ierr; 2523 2524 PetscFunctionBegin; 2525 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2526 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2527 if (basis) {PetscValidPointer(basis, 2); *basis = prob->basisFace;} 2528 if (basisDer) {PetscValidPointer(basisDer, 3); *basisDer = prob->basisDerFace;} 2529 PetscFunctionReturn(0); 2530 } 2531 2532 PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x) 2533 { 2534 PetscErrorCode ierr; 2535 2536 PetscFunctionBegin; 2537 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2538 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2539 if (u) {PetscValidPointer(u, 2); *u = prob->u;} 2540 if (u_t) {PetscValidPointer(u_t, 3); *u_t = prob->u_t;} 2541 if (u_x) {PetscValidPointer(u_x, 4); *u_x = prob->u_x;} 2542 PetscFunctionReturn(0); 2543 } 2544 2545 PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3) 2546 { 2547 PetscErrorCode ierr; 2548 2549 PetscFunctionBegin; 2550 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2551 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2552 if (f0) {PetscValidPointer(f0, 2); *f0 = prob->f0;} 2553 if (f1) {PetscValidPointer(f1, 3); *f1 = prob->f1;} 2554 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->g0;} 2555 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->g1;} 2556 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->g2;} 2557 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->g3;} 2558 PetscFunctionReturn(0); 2559 } 2560 2561 PetscErrorCode PetscDSGetRefCoordArrays(PetscDS prob, PetscReal **x, PetscScalar **refSpaceDer) 2562 { 2563 PetscErrorCode ierr; 2564 2565 PetscFunctionBegin; 2566 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2567 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2568 if (x) {PetscValidPointer(x, 2); *x = prob->x;} 2569 if (refSpaceDer) {PetscValidPointer(refSpaceDer, 3); *refSpaceDer = prob->refSpaceDer;} 2570 PetscFunctionReturn(0); 2571 } 2572 2573 /*@C 2574 PetscDSAddBoundary - Add a boundary condition to the model 2575 2576 Input Parameters: 2577 + ds - The PetscDS object 2578 . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann) 2579 . name - The BC name 2580 . labelname - The label defining constrained points 2581 . field - The field to constrain 2582 . numcomps - The number of constrained field components 2583 . comps - An array of constrained component numbers 2584 . bcFunc - A pointwise function giving boundary values 2585 . numids - The number of DMLabel ids for constrained points 2586 . ids - An array of ids for constrained points 2587 - ctx - An optional user context for bcFunc 2588 2589 Options Database Keys: 2590 + -bc_<boundary name> <num> - Overrides the boundary ids 2591 - -bc_<boundary name>_comp <num> - Overrides the boundary components 2592 2593 Level: developer 2594 2595 .seealso: PetscDSGetBoundary() 2596 @*/ 2597 PetscErrorCode PetscDSAddBoundary(PetscDS ds, DMBoundaryConditionType type, const char name[], const char labelname[], PetscInt field, PetscInt numcomps, const PetscInt *comps, void (*bcFunc)(void), PetscInt numids, const PetscInt *ids, void *ctx) 2598 { 2599 DSBoundary b; 2600 PetscErrorCode ierr; 2601 2602 PetscFunctionBegin; 2603 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2604 ierr = PetscNew(&b);CHKERRQ(ierr); 2605 ierr = PetscStrallocpy(name, (char **) &b->name);CHKERRQ(ierr); 2606 ierr = PetscStrallocpy(labelname, (char **) &b->labelname);CHKERRQ(ierr); 2607 ierr = PetscMalloc1(numcomps, &b->comps);CHKERRQ(ierr); 2608 if (numcomps) {ierr = PetscMemcpy(b->comps, comps, numcomps*sizeof(PetscInt));CHKERRQ(ierr);} 2609 ierr = PetscMalloc1(numids, &b->ids);CHKERRQ(ierr); 2610 if (numids) {ierr = PetscMemcpy(b->ids, ids, numids*sizeof(PetscInt));CHKERRQ(ierr);} 2611 b->type = type; 2612 b->field = field; 2613 b->numcomps = numcomps; 2614 b->func = bcFunc; 2615 b->numids = numids; 2616 b->ctx = ctx; 2617 b->next = ds->boundary; 2618 ds->boundary = b; 2619 PetscFunctionReturn(0); 2620 } 2621 2622 /*@ 2623 PetscDSGetNumBoundary - Get the number of registered BC 2624 2625 Input Parameters: 2626 . ds - The PetscDS object 2627 2628 Output Parameters: 2629 . numBd - The number of BC 2630 2631 Level: intermediate 2632 2633 .seealso: PetscDSAddBoundary(), PetscDSGetBoundary() 2634 @*/ 2635 PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd) 2636 { 2637 DSBoundary b = ds->boundary; 2638 2639 PetscFunctionBegin; 2640 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2641 PetscValidPointer(numBd, 2); 2642 *numBd = 0; 2643 while (b) {++(*numBd); b = b->next;} 2644 PetscFunctionReturn(0); 2645 } 2646 2647 /*@C 2648 PetscDSGetBoundary - Add a boundary condition to the model 2649 2650 Input Parameters: 2651 + ds - The PetscDS object 2652 - bd - The BC number 2653 2654 Output Parameters: 2655 + type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann) 2656 . name - The BC name 2657 . labelname - The label defining constrained points 2658 . field - The field to constrain 2659 . numcomps - The number of constrained field components 2660 . comps - An array of constrained component numbers 2661 . bcFunc - A pointwise function giving boundary values 2662 . numids - The number of DMLabel ids for constrained points 2663 . ids - An array of ids for constrained points 2664 - ctx - An optional user context for bcFunc 2665 2666 Options Database Keys: 2667 + -bc_<boundary name> <num> - Overrides the boundary ids 2668 - -bc_<boundary name>_comp <num> - Overrides the boundary components 2669 2670 Level: developer 2671 2672 .seealso: PetscDSAddBoundary() 2673 @*/ 2674 PetscErrorCode PetscDSGetBoundary(PetscDS ds, PetscInt bd, DMBoundaryConditionType *type, const char **name, const char **labelname, PetscInt *field, PetscInt *numcomps, const PetscInt **comps, void (**func)(void), PetscInt *numids, const PetscInt **ids, void **ctx) 2675 { 2676 DSBoundary b = ds->boundary; 2677 PetscInt n = 0; 2678 2679 PetscFunctionBegin; 2680 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2681 while (b) { 2682 if (n == bd) break; 2683 b = b->next; 2684 ++n; 2685 } 2686 if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n); 2687 if (type) { 2688 PetscValidPointer(type, 3); 2689 *type = b->type; 2690 } 2691 if (name) { 2692 PetscValidPointer(name, 4); 2693 *name = b->name; 2694 } 2695 if (labelname) { 2696 PetscValidPointer(labelname, 5); 2697 *labelname = b->labelname; 2698 } 2699 if (field) { 2700 PetscValidPointer(field, 6); 2701 *field = b->field; 2702 } 2703 if (numcomps) { 2704 PetscValidPointer(numcomps, 7); 2705 *numcomps = b->numcomps; 2706 } 2707 if (comps) { 2708 PetscValidPointer(comps, 8); 2709 *comps = b->comps; 2710 } 2711 if (func) { 2712 PetscValidPointer(func, 9); 2713 *func = b->func; 2714 } 2715 if (numids) { 2716 PetscValidPointer(numids, 10); 2717 *numids = b->numids; 2718 } 2719 if (ids) { 2720 PetscValidPointer(ids, 11); 2721 *ids = b->ids; 2722 } 2723 if (ctx) { 2724 PetscValidPointer(ctx, 12); 2725 *ctx = b->ctx; 2726 } 2727 PetscFunctionReturn(0); 2728 } 2729 2730 PetscErrorCode PetscDSCopyBoundary(PetscDS probA, PetscDS probB) 2731 { 2732 DSBoundary b, next, *lastnext; 2733 PetscErrorCode ierr; 2734 2735 PetscFunctionBegin; 2736 PetscValidHeaderSpecific(probA, PETSCDS_CLASSID, 1); 2737 PetscValidHeaderSpecific(probB, PETSCDS_CLASSID, 2); 2738 if (probA == probB) PetscFunctionReturn(0); 2739 next = probB->boundary; 2740 while (next) { 2741 DSBoundary b = next; 2742 2743 next = b->next; 2744 ierr = PetscFree(b->comps);CHKERRQ(ierr); 2745 ierr = PetscFree(b->ids);CHKERRQ(ierr); 2746 ierr = PetscFree(b->name);CHKERRQ(ierr); 2747 ierr = PetscFree(b->labelname);CHKERRQ(ierr); 2748 ierr = PetscFree(b);CHKERRQ(ierr); 2749 } 2750 lastnext = &(probB->boundary); 2751 for (b = probA->boundary; b; b = b->next) { 2752 DSBoundary bNew; 2753 2754 ierr = PetscNew(&bNew);CHKERRQ(ierr); 2755 bNew->numcomps = b->numcomps; 2756 ierr = PetscMalloc1(bNew->numcomps, &bNew->comps);CHKERRQ(ierr); 2757 ierr = PetscMemcpy(bNew->comps, b->comps, bNew->numcomps*sizeof(PetscInt));CHKERRQ(ierr); 2758 bNew->numids = b->numids; 2759 ierr = PetscMalloc1(bNew->numids, &bNew->ids);CHKERRQ(ierr); 2760 ierr = PetscMemcpy(bNew->ids, b->ids, bNew->numids*sizeof(PetscInt));CHKERRQ(ierr); 2761 ierr = PetscStrallocpy(b->labelname,(char **) &(bNew->labelname));CHKERRQ(ierr); 2762 ierr = PetscStrallocpy(b->name,(char **) &(bNew->name));CHKERRQ(ierr); 2763 bNew->ctx = b->ctx; 2764 bNew->type = b->type; 2765 bNew->field = b->field; 2766 bNew->func = b->func; 2767 2768 *lastnext = bNew; 2769 lastnext = &(bNew->next); 2770 } 2771 PetscFunctionReturn(0); 2772 } 2773 2774 /*@ 2775 PetscDSCopyEquations - Copy all pointwise function pointers to the new problem 2776 2777 Not collective 2778 2779 Input Parameter: 2780 . prob - The PetscDS object 2781 2782 Output Parameter: 2783 . newprob - The PetscDS copy 2784 2785 Level: intermediate 2786 2787 .seealso: PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate() 2788 @*/ 2789 PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob) 2790 { 2791 PetscInt Nf, Ng, f, g; 2792 PetscErrorCode ierr; 2793 2794 PetscFunctionBegin; 2795 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2796 PetscValidHeaderSpecific(newprob, PETSCDS_CLASSID, 2); 2797 ierr = PetscDSGetNumFields(prob, &Nf);CHKERRQ(ierr); 2798 ierr = PetscDSGetNumFields(newprob, &Ng);CHKERRQ(ierr); 2799 if (Nf != Ng) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %D != %D", Nf, Ng); 2800 for (f = 0; f < Nf; ++f) { 2801 PetscPointFunc obj; 2802 PetscPointFunc f0, f1; 2803 PetscPointJac g0, g1, g2, g3; 2804 PetscBdPointFunc f0Bd, f1Bd; 2805 PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd; 2806 PetscRiemannFunc r; 2807 2808 ierr = PetscDSGetObjective(prob, f, &obj);CHKERRQ(ierr); 2809 ierr = PetscDSGetResidual(prob, f, &f0, &f1);CHKERRQ(ierr); 2810 ierr = PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);CHKERRQ(ierr); 2811 ierr = PetscDSGetRiemannSolver(prob, f, &r);CHKERRQ(ierr); 2812 ierr = PetscDSSetObjective(newprob, f, obj);CHKERRQ(ierr); 2813 ierr = PetscDSSetResidual(newprob, f, f0, f1);CHKERRQ(ierr); 2814 ierr = PetscDSSetBdResidual(newprob, f, f0Bd, f1Bd);CHKERRQ(ierr); 2815 ierr = PetscDSSetRiemannSolver(newprob, f, r);CHKERRQ(ierr); 2816 for (g = 0; g < Nf; ++g) { 2817 ierr = PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);CHKERRQ(ierr); 2818 ierr = PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);CHKERRQ(ierr); 2819 ierr = PetscDSSetJacobian(newprob, f, g, g0, g1, g2, g3);CHKERRQ(ierr); 2820 ierr = PetscDSSetBdJacobian(newprob, f, g, g0Bd, g1Bd, g2Bd, g3Bd);CHKERRQ(ierr); 2821 } 2822 } 2823 PetscFunctionReturn(0); 2824 } 2825 2826 static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob) 2827 { 2828 PetscErrorCode ierr; 2829 2830 PetscFunctionBegin; 2831 ierr = PetscFree(prob->data);CHKERRQ(ierr); 2832 PetscFunctionReturn(0); 2833 } 2834 2835 static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob) 2836 { 2837 PetscFunctionBegin; 2838 prob->ops->setfromoptions = NULL; 2839 prob->ops->setup = NULL; 2840 prob->ops->view = NULL; 2841 prob->ops->destroy = PetscDSDestroy_Basic; 2842 PetscFunctionReturn(0); 2843 } 2844 2845 /*MC 2846 PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions 2847 2848 Level: intermediate 2849 2850 .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType() 2851 M*/ 2852 2853 PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob) 2854 { 2855 PetscDS_Basic *b; 2856 PetscErrorCode ierr; 2857 2858 PetscFunctionBegin; 2859 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2860 ierr = PetscNewLog(prob, &b);CHKERRQ(ierr); 2861 prob->data = b; 2862 2863 ierr = PetscDSInitialize_Basic(prob);CHKERRQ(ierr); 2864 PetscFunctionReturn(0); 2865 } 2866