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", 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; 364 PetscPointJac *tmpg, *tmpgp, *tmpgt; 365 PetscBdPointFunc *tmpfbd; 366 PetscBdPointJac *tmpgbd; 367 PetscRiemannFunc *tmpr; 368 void **tmpctx; 369 PetscInt Nf = prob->Nf, f, i; 370 PetscErrorCode ierr; 371 372 PetscFunctionBegin; 373 if (Nf >= NfNew) PetscFunctionReturn(0); 374 prob->setup = PETSC_FALSE; 375 ierr = PetscDSDestroyStructs_Static(prob);CHKERRQ(ierr); 376 ierr = PetscMalloc3(NfNew, &tmpd, NfNew, &tmpi, NfNew*2, &tmpa);CHKERRQ(ierr); 377 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];} 378 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;} 379 ierr = PetscFree3(prob->disc, prob->implicit, prob->adjacency);CHKERRQ(ierr); 380 prob->Nf = NfNew; 381 prob->disc = tmpd; 382 prob->implicit = tmpi; 383 prob->adjacency = tmpa; 384 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); 385 for (f = 0; f < Nf; ++f) tmpobj[f] = prob->obj[f]; 386 for (f = 0; f < Nf*2; ++f) tmpf[f] = prob->f[f]; 387 for (f = 0; f < Nf*Nf*4; ++f) tmpg[f] = prob->g[f]; 388 for (f = 0; f < Nf*Nf*4; ++f) tmpgp[f] = prob->gp[f]; 389 for (f = 0; f < Nf; ++f) tmpr[f] = prob->r[f]; 390 for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f]; 391 for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL; 392 for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL; 393 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL; 394 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgp[f] = NULL; 395 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgt[f] = NULL; 396 for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL; 397 for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL; 398 ierr = PetscFree7(prob->obj, prob->f, prob->g, prob->gp, prob->gt, prob->r, prob->ctx);CHKERRQ(ierr); 399 prob->obj = tmpobj; 400 prob->f = tmpf; 401 prob->g = tmpg; 402 prob->gp = tmpgp; 403 prob->gt = tmpgt; 404 prob->r = tmpr; 405 prob->ctx = tmpctx; 406 ierr = PetscCalloc2(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd);CHKERRQ(ierr); 407 for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f]; 408 for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f]; 409 for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL; 410 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL; 411 ierr = PetscFree2(prob->fBd, prob->gBd);CHKERRQ(ierr); 412 prob->fBd = tmpfbd; 413 prob->gBd = tmpgbd; 414 PetscFunctionReturn(0); 415 } 416 417 /*@ 418 PetscDSDestroy - Destroys a PetscDS object 419 420 Collective on PetscDS 421 422 Input Parameter: 423 . prob - the PetscDS object to destroy 424 425 Level: developer 426 427 .seealso PetscDSView() 428 @*/ 429 PetscErrorCode PetscDSDestroy(PetscDS *prob) 430 { 431 PetscInt f; 432 DSBoundary next; 433 PetscErrorCode ierr; 434 435 PetscFunctionBegin; 436 if (!*prob) PetscFunctionReturn(0); 437 PetscValidHeaderSpecific((*prob), PETSCDS_CLASSID, 1); 438 439 if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; PetscFunctionReturn(0);} 440 ((PetscObject) (*prob))->refct = 0; 441 ierr = PetscDSDestroyStructs_Static(*prob);CHKERRQ(ierr); 442 for (f = 0; f < (*prob)->Nf; ++f) { 443 ierr = PetscObjectDereference((*prob)->disc[f]);CHKERRQ(ierr); 444 } 445 ierr = PetscFree3((*prob)->disc, (*prob)->implicit, (*prob)->adjacency);CHKERRQ(ierr); 446 ierr = PetscFree7((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->gp,(*prob)->gt,(*prob)->r,(*prob)->ctx);CHKERRQ(ierr); 447 ierr = PetscFree2((*prob)->fBd,(*prob)->gBd);CHKERRQ(ierr); 448 if ((*prob)->ops->destroy) {ierr = (*(*prob)->ops->destroy)(*prob);CHKERRQ(ierr);} 449 next = (*prob)->boundary; 450 while (next) { 451 DSBoundary b = next; 452 453 next = b->next; 454 ierr = PetscFree(b->comps);CHKERRQ(ierr); 455 ierr = PetscFree(b->ids);CHKERRQ(ierr); 456 ierr = PetscFree(b->name);CHKERRQ(ierr); 457 ierr = PetscFree(b->labelname);CHKERRQ(ierr); 458 ierr = PetscFree(b);CHKERRQ(ierr); 459 } 460 ierr = PetscFree((*prob)->constants);CHKERRQ(ierr); 461 ierr = PetscHeaderDestroy(prob);CHKERRQ(ierr); 462 PetscFunctionReturn(0); 463 } 464 465 /*@ 466 PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType(). 467 468 Collective on MPI_Comm 469 470 Input Parameter: 471 . comm - The communicator for the PetscDS object 472 473 Output Parameter: 474 . prob - The PetscDS object 475 476 Level: beginner 477 478 .seealso: PetscDSSetType(), PETSCDSBASIC 479 @*/ 480 PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob) 481 { 482 PetscDS p; 483 PetscErrorCode ierr; 484 485 PetscFunctionBegin; 486 PetscValidPointer(prob, 2); 487 *prob = NULL; 488 ierr = PetscDSInitializePackage();CHKERRQ(ierr); 489 490 ierr = PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);CHKERRQ(ierr); 491 492 p->Nf = 0; 493 p->setup = PETSC_FALSE; 494 p->numConstants = 0; 495 p->constants = NULL; 496 497 *prob = p; 498 PetscFunctionReturn(0); 499 } 500 501 /*@ 502 PetscDSGetNumFields - Returns the number of fields in the DS 503 504 Not collective 505 506 Input Parameter: 507 . prob - The PetscDS object 508 509 Output Parameter: 510 . Nf - The number of fields 511 512 Level: beginner 513 514 .seealso: PetscDSGetSpatialDimension(), PetscDSCreate() 515 @*/ 516 PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf) 517 { 518 PetscFunctionBegin; 519 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 520 PetscValidPointer(Nf, 2); 521 *Nf = prob->Nf; 522 PetscFunctionReturn(0); 523 } 524 525 /*@ 526 PetscDSGetSpatialDimension - Returns the spatial dimension of the DS 527 528 Not collective 529 530 Input Parameter: 531 . prob - The PetscDS object 532 533 Output Parameter: 534 . dim - The spatial dimension 535 536 Level: beginner 537 538 .seealso: PetscDSGetNumFields(), PetscDSCreate() 539 @*/ 540 PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim) 541 { 542 PetscErrorCode ierr; 543 544 PetscFunctionBegin; 545 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 546 PetscValidPointer(dim, 2); 547 *dim = 0; 548 if (prob->Nf) { 549 PetscObject obj; 550 PetscClassId id; 551 552 ierr = PetscDSGetDiscretization(prob, 0, &obj);CHKERRQ(ierr); 553 ierr = PetscObjectGetClassId(obj, &id);CHKERRQ(ierr); 554 if (id == PETSCFE_CLASSID) {ierr = PetscFEGetSpatialDimension((PetscFE) obj, dim);CHKERRQ(ierr);} 555 else if (id == PETSCFV_CLASSID) {ierr = PetscFVGetSpatialDimension((PetscFV) obj, dim);CHKERRQ(ierr);} 556 else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0); 557 } 558 PetscFunctionReturn(0); 559 } 560 561 /*@ 562 PetscDSGetTotalDimension - Returns the total size of the approximation space for this system 563 564 Not collective 565 566 Input Parameter: 567 . prob - The PetscDS object 568 569 Output Parameter: 570 . dim - The total problem dimension 571 572 Level: beginner 573 574 .seealso: PetscDSGetNumFields(), PetscDSCreate() 575 @*/ 576 PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim) 577 { 578 PetscErrorCode ierr; 579 580 PetscFunctionBegin; 581 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 582 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 583 PetscValidPointer(dim, 2); 584 *dim = prob->totDim; 585 PetscFunctionReturn(0); 586 } 587 588 /*@ 589 PetscDSGetTotalComponents - Returns the total number of components in this system 590 591 Not collective 592 593 Input Parameter: 594 . prob - The PetscDS object 595 596 Output Parameter: 597 . dim - The total number of components 598 599 Level: beginner 600 601 .seealso: PetscDSGetNumFields(), PetscDSCreate() 602 @*/ 603 PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc) 604 { 605 PetscErrorCode ierr; 606 607 PetscFunctionBegin; 608 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 609 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 610 PetscValidPointer(Nc, 2); 611 *Nc = prob->totComp; 612 PetscFunctionReturn(0); 613 } 614 615 /*@ 616 PetscDSGetDiscretization - Returns the discretization object for the given field 617 618 Not collective 619 620 Input Parameters: 621 + prob - The PetscDS object 622 - f - The field number 623 624 Output Parameter: 625 . disc - The discretization object 626 627 Level: beginner 628 629 .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 630 @*/ 631 PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc) 632 { 633 PetscFunctionBegin; 634 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 635 PetscValidPointer(disc, 3); 636 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); 637 *disc = prob->disc[f]; 638 PetscFunctionReturn(0); 639 } 640 641 /*@ 642 PetscDSSetDiscretization - Sets the discretization object for the given field 643 644 Not collective 645 646 Input Parameters: 647 + prob - The PetscDS object 648 . f - The field number 649 - disc - The discretization object 650 651 Level: beginner 652 653 .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 654 @*/ 655 PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc) 656 { 657 PetscErrorCode ierr; 658 659 PetscFunctionBegin; 660 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 661 PetscValidPointer(disc, 3); 662 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 663 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 664 if (prob->disc[f]) {ierr = PetscObjectDereference(prob->disc[f]);CHKERRQ(ierr);} 665 prob->disc[f] = disc; 666 ierr = PetscObjectReference(disc);CHKERRQ(ierr); 667 { 668 PetscClassId id; 669 670 ierr = PetscObjectGetClassId(disc, &id);CHKERRQ(ierr); 671 if (id == PETSCFV_CLASSID) { 672 prob->implicit[f] = PETSC_FALSE; 673 prob->adjacency[f*2+0] = PETSC_TRUE; 674 prob->adjacency[f*2+1] = PETSC_FALSE; 675 } 676 } 677 PetscFunctionReturn(0); 678 } 679 680 /*@ 681 PetscDSAddDiscretization - Adds a discretization object 682 683 Not collective 684 685 Input Parameters: 686 + prob - The PetscDS object 687 - disc - The boundary discretization object 688 689 Level: beginner 690 691 .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 692 @*/ 693 PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc) 694 { 695 PetscErrorCode ierr; 696 697 PetscFunctionBegin; 698 ierr = PetscDSSetDiscretization(prob, prob->Nf, disc);CHKERRQ(ierr); 699 PetscFunctionReturn(0); 700 } 701 702 /*@ 703 PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX 704 705 Not collective 706 707 Input Parameters: 708 + prob - The PetscDS object 709 - f - The field number 710 711 Output Parameter: 712 . implicit - The flag indicating what kind of solve to use for this field 713 714 Level: developer 715 716 .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 717 @*/ 718 PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit) 719 { 720 PetscFunctionBegin; 721 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 722 PetscValidPointer(implicit, 3); 723 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); 724 *implicit = prob->implicit[f]; 725 PetscFunctionReturn(0); 726 } 727 728 /*@ 729 PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for IMEX 730 731 Not collective 732 733 Input Parameters: 734 + prob - The PetscDS object 735 . f - The field number 736 - implicit - The flag indicating what kind of solve to use for this field 737 738 Level: developer 739 740 .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 741 @*/ 742 PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit) 743 { 744 PetscFunctionBegin; 745 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 746 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); 747 prob->implicit[f] = implicit; 748 PetscFunctionReturn(0); 749 } 750 751 /*@ 752 PetscDSGetAdjacency - Returns the flags for determining variable influence 753 754 Not collective 755 756 Input Parameters: 757 + prob - The PetscDS object 758 - f - The field number 759 760 Output Parameter: 761 + useCone - Flag for variable influence starting with the cone operation 762 - useClosure - Flag for variable influence using transitive closure 763 764 Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure() 765 766 Level: developer 767 768 .seealso: PetscDSSetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 769 @*/ 770 PetscErrorCode PetscDSGetAdjacency(PetscDS prob, PetscInt f, PetscBool *useCone, PetscBool *useClosure) 771 { 772 PetscFunctionBegin; 773 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 774 PetscValidPointer(useCone, 3); 775 PetscValidPointer(useClosure, 4); 776 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); 777 *useCone = prob->adjacency[f*2+0]; 778 *useClosure = prob->adjacency[f*2+1]; 779 PetscFunctionReturn(0); 780 } 781 782 /*@ 783 PetscDSSetAdjacency - Set the flags for determining variable influence 784 785 Not collective 786 787 Input Parameters: 788 + prob - The PetscDS object 789 . f - The field number 790 . useCone - Flag for variable influence starting with the cone operation 791 - useClosure - Flag for variable influence using transitive closure 792 793 Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure() 794 795 Level: developer 796 797 .seealso: PetscDSGetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 798 @*/ 799 PetscErrorCode PetscDSSetAdjacency(PetscDS prob, PetscInt f, PetscBool useCone, PetscBool useClosure) 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->adjacency[f*2+0] = useCone; 805 prob->adjacency[f*2+1] = useClosure; 806 PetscFunctionReturn(0); 807 } 808 809 PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f, 810 void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 811 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 812 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 813 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[])) 814 { 815 PetscFunctionBegin; 816 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 817 PetscValidPointer(obj, 2); 818 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); 819 *obj = prob->obj[f]; 820 PetscFunctionReturn(0); 821 } 822 823 PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f, 824 void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 825 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 826 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 827 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[])) 828 { 829 PetscErrorCode ierr; 830 831 PetscFunctionBegin; 832 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 833 if (obj) PetscValidFunction(obj, 2); 834 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 835 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 836 prob->obj[f] = obj; 837 PetscFunctionReturn(0); 838 } 839 840 /*@C 841 PetscDSGetResidual - Get the pointwise residual function for a given test field 842 843 Not collective 844 845 Input Parameters: 846 + prob - The PetscDS 847 - f - The test field number 848 849 Output Parameters: 850 + f0 - integrand for the test function term 851 - f1 - integrand for the test function gradient term 852 853 Note: We are using a first order FEM model for the weak form: 854 855 \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) 856 857 The calling sequence for the callbacks f0 and f1 is given by: 858 859 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 860 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 861 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 862 $ PetscReal t, const PetscReal x[], PetscScalar f0[]) 863 864 + dim - the spatial dimension 865 . Nf - the number of fields 866 . uOff - the offset into u[] and u_t[] for each field 867 . uOff_x - the offset into u_x[] for each field 868 . u - each field evaluated at the current point 869 . u_t - the time derivative of each field evaluated at the current point 870 . u_x - the gradient of each field evaluated at the current point 871 . aOff - the offset into a[] and a_t[] for each auxiliary field 872 . aOff_x - the offset into a_x[] for each auxiliary field 873 . a - each auxiliary field evaluated at the current point 874 . a_t - the time derivative of each auxiliary field evaluated at the current point 875 . a_x - the gradient of auxiliary each field evaluated at the current point 876 . t - current time 877 . x - coordinates of the current point 878 . numConstants - number of constant parameters 879 . constants - constant parameters 880 - f0 - output values at the current point 881 882 Level: intermediate 883 884 .seealso: PetscDSSetResidual() 885 @*/ 886 PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f, 887 void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 888 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 889 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 890 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), 891 void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 892 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 893 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 894 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])) 895 { 896 PetscFunctionBegin; 897 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 898 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); 899 if (f0) {PetscValidPointer(f0, 3); *f0 = prob->f[f*2+0];} 900 if (f1) {PetscValidPointer(f1, 4); *f1 = prob->f[f*2+1];} 901 PetscFunctionReturn(0); 902 } 903 904 /*@C 905 PetscDSSetResidual - Set the pointwise residual function for a given test field 906 907 Not collective 908 909 Input Parameters: 910 + prob - The PetscDS 911 . f - The test field number 912 . f0 - integrand for the test function term 913 - f1 - integrand for the test function gradient term 914 915 Note: We are using a first order FEM model for the weak form: 916 917 \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) 918 919 The calling sequence for the callbacks f0 and f1 is given by: 920 921 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 922 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 923 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 924 $ PetscReal t, const PetscReal x[], PetscScalar f0[]) 925 926 + dim - the spatial dimension 927 . Nf - the number of fields 928 . uOff - the offset into u[] and u_t[] for each field 929 . uOff_x - the offset into u_x[] for each field 930 . u - each field evaluated at the current point 931 . u_t - the time derivative of each field evaluated at the current point 932 . u_x - the gradient of each field evaluated at the current point 933 . aOff - the offset into a[] and a_t[] for each auxiliary field 934 . aOff_x - the offset into a_x[] for each auxiliary field 935 . a - each auxiliary field evaluated at the current point 936 . a_t - the time derivative of each auxiliary field evaluated at the current point 937 . a_x - the gradient of auxiliary each field evaluated at the current point 938 . t - current time 939 . x - coordinates of the current point 940 . numConstants - number of constant parameters 941 . constants - constant parameters 942 - f0 - output values at the current point 943 944 Level: intermediate 945 946 .seealso: PetscDSGetResidual() 947 @*/ 948 PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f, 949 void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 950 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 951 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 952 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), 953 void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 954 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 955 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 956 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])) 957 { 958 PetscErrorCode ierr; 959 960 PetscFunctionBegin; 961 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 962 if (f0) PetscValidFunction(f0, 3); 963 if (f1) PetscValidFunction(f1, 4); 964 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 965 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 966 prob->f[f*2+0] = f0; 967 prob->f[f*2+1] = f1; 968 PetscFunctionReturn(0); 969 } 970 971 /*@C 972 PetscDSHasJacobian - Signals that Jacobian functions have been set 973 974 Not collective 975 976 Input Parameter: 977 . prob - The PetscDS 978 979 Output Parameter: 980 . hasJac - flag that pointwise function for the Jacobian has been set 981 982 Level: intermediate 983 984 .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian() 985 @*/ 986 PetscErrorCode PetscDSHasJacobian(PetscDS prob, PetscBool *hasJac) 987 { 988 PetscInt f, g, h; 989 990 PetscFunctionBegin; 991 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 992 *hasJac = PETSC_FALSE; 993 for (f = 0; f < prob->Nf; ++f) { 994 for (g = 0; g < prob->Nf; ++g) { 995 for (h = 0; h < 4; ++h) { 996 if (prob->g[(f*prob->Nf + g)*4+h]) *hasJac = PETSC_TRUE; 997 } 998 } 999 } 1000 PetscFunctionReturn(0); 1001 } 1002 1003 /*@C 1004 PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field 1005 1006 Not collective 1007 1008 Input Parameters: 1009 + prob - The PetscDS 1010 . f - The test field number 1011 - g - The field number 1012 1013 Output Parameters: 1014 + g0 - integrand for the test and basis function term 1015 . g1 - integrand for the test function and basis function gradient term 1016 . g2 - integrand for the test function gradient and basis function term 1017 - g3 - integrand for the test function gradient and basis function gradient term 1018 1019 Note: We are using a first order FEM model for the weak form: 1020 1021 \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 1022 1023 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1024 1025 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1026 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1027 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1028 $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]) 1029 1030 + dim - the spatial dimension 1031 . Nf - the number of fields 1032 . uOff - the offset into u[] and u_t[] for each field 1033 . uOff_x - the offset into u_x[] for each field 1034 . u - each field evaluated at the current point 1035 . u_t - the time derivative of each field evaluated at the current point 1036 . u_x - the gradient of each field evaluated at the current point 1037 . aOff - the offset into a[] and a_t[] for each auxiliary field 1038 . aOff_x - the offset into a_x[] for each auxiliary field 1039 . a - each auxiliary field evaluated at the current point 1040 . a_t - the time derivative of each auxiliary field evaluated at the current point 1041 . a_x - the gradient of auxiliary each field evaluated at the current point 1042 . t - current time 1043 . u_tShift - the multiplier a for dF/dU_t 1044 . x - coordinates of the current point 1045 . numConstants - number of constant parameters 1046 . constants - constant parameters 1047 - g0 - output values at the current point 1048 1049 Level: intermediate 1050 1051 .seealso: PetscDSSetJacobian() 1052 @*/ 1053 PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g, 1054 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1055 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1056 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1057 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1058 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1059 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1060 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1061 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1062 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1063 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1064 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1065 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1066 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1067 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1068 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1069 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1070 { 1071 PetscFunctionBegin; 1072 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1073 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); 1074 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); 1075 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->g[(f*prob->Nf + g)*4+0];} 1076 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->g[(f*prob->Nf + g)*4+1];} 1077 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->g[(f*prob->Nf + g)*4+2];} 1078 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->g[(f*prob->Nf + g)*4+3];} 1079 PetscFunctionReturn(0); 1080 } 1081 1082 /*@C 1083 PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields 1084 1085 Not collective 1086 1087 Input Parameters: 1088 + prob - The PetscDS 1089 . f - The test field number 1090 . g - The field number 1091 . g0 - integrand for the test and basis function term 1092 . g1 - integrand for the test function and basis function gradient term 1093 . g2 - integrand for the test function gradient and basis function term 1094 - g3 - integrand for the test function gradient and basis function gradient term 1095 1096 Note: We are using a first order FEM model for the weak form: 1097 1098 \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 1099 1100 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1101 1102 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1103 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1104 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1105 $ PetscReal t, const PetscReal x[], PetscScalar g0[]) 1106 1107 + dim - the spatial dimension 1108 . Nf - the number of fields 1109 . uOff - the offset into u[] and u_t[] for each field 1110 . uOff_x - the offset into u_x[] for each field 1111 . u - each field evaluated at the current point 1112 . u_t - the time derivative of each field evaluated at the current point 1113 . u_x - the gradient of each field evaluated at the current point 1114 . aOff - the offset into a[] and a_t[] for each auxiliary field 1115 . aOff_x - the offset into a_x[] for each auxiliary field 1116 . a - each auxiliary field evaluated at the current point 1117 . a_t - the time derivative of each auxiliary field evaluated at the current point 1118 . a_x - the gradient of auxiliary each field evaluated at the current point 1119 . t - current time 1120 . u_tShift - the multiplier a for dF/dU_t 1121 . x - coordinates of the current point 1122 . numConstants - number of constant parameters 1123 . constants - constant parameters 1124 - g0 - output values at the current point 1125 1126 Level: intermediate 1127 1128 .seealso: PetscDSGetJacobian() 1129 @*/ 1130 PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g, 1131 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1132 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1133 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1134 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1135 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1136 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1137 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1138 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1139 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1140 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1141 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1142 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1143 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1144 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1145 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1146 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1147 { 1148 PetscErrorCode ierr; 1149 1150 PetscFunctionBegin; 1151 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1152 if (g0) PetscValidFunction(g0, 4); 1153 if (g1) PetscValidFunction(g1, 5); 1154 if (g2) PetscValidFunction(g2, 6); 1155 if (g3) PetscValidFunction(g3, 7); 1156 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1157 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1158 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1159 prob->g[(f*prob->Nf + g)*4+0] = g0; 1160 prob->g[(f*prob->Nf + g)*4+1] = g1; 1161 prob->g[(f*prob->Nf + g)*4+2] = g2; 1162 prob->g[(f*prob->Nf + g)*4+3] = g3; 1163 PetscFunctionReturn(0); 1164 } 1165 1166 /*@C 1167 PetscDSHasJacobianPreconditioner - Signals that a Jacobian preconditioner matrix has been set 1168 1169 Not collective 1170 1171 Input Parameter: 1172 . prob - The PetscDS 1173 1174 Output Parameter: 1175 . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set 1176 1177 Level: intermediate 1178 1179 .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian() 1180 @*/ 1181 PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS prob, PetscBool *hasJacPre) 1182 { 1183 PetscInt f, g, h; 1184 1185 PetscFunctionBegin; 1186 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1187 *hasJacPre = PETSC_FALSE; 1188 for (f = 0; f < prob->Nf; ++f) { 1189 for (g = 0; g < prob->Nf; ++g) { 1190 for (h = 0; h < 4; ++h) { 1191 if (prob->gp[(f*prob->Nf + g)*4+h]) *hasJacPre = PETSC_TRUE; 1192 } 1193 } 1194 } 1195 PetscFunctionReturn(0); 1196 } 1197 1198 /*@C 1199 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. 1200 1201 Not collective 1202 1203 Input Parameters: 1204 + prob - The PetscDS 1205 . f - The test field number 1206 - g - The field number 1207 1208 Output Parameters: 1209 + g0 - integrand for the test and basis function term 1210 . g1 - integrand for the test function and basis function gradient term 1211 . g2 - integrand for the test function gradient and basis function term 1212 - g3 - integrand for the test function gradient and basis function gradient term 1213 1214 Note: We are using a first order FEM model for the weak form: 1215 1216 \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 1217 1218 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1219 1220 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1221 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1222 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1223 $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]) 1224 1225 + dim - the spatial dimension 1226 . Nf - the number of fields 1227 . uOff - the offset into u[] and u_t[] for each field 1228 . uOff_x - the offset into u_x[] for each field 1229 . u - each field evaluated at the current point 1230 . u_t - the time derivative of each field evaluated at the current point 1231 . u_x - the gradient of each field evaluated at the current point 1232 . aOff - the offset into a[] and a_t[] for each auxiliary field 1233 . aOff_x - the offset into a_x[] for each auxiliary field 1234 . a - each auxiliary field evaluated at the current point 1235 . a_t - the time derivative of each auxiliary field evaluated at the current point 1236 . a_x - the gradient of auxiliary each field evaluated at the current point 1237 . t - current time 1238 . u_tShift - the multiplier a for dF/dU_t 1239 . x - coordinates of the current point 1240 . numConstants - number of constant parameters 1241 . constants - constant parameters 1242 - g0 - output values at the current point 1243 1244 Level: intermediate 1245 1246 .seealso: PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian() 1247 @*/ 1248 PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g, 1249 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1250 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1251 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1252 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1253 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1254 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1255 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1256 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1257 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1258 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1259 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1260 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1261 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1262 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1263 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1264 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1265 { 1266 PetscFunctionBegin; 1267 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1268 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); 1269 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); 1270 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->gp[(f*prob->Nf + g)*4+0];} 1271 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->gp[(f*prob->Nf + g)*4+1];} 1272 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->gp[(f*prob->Nf + g)*4+2];} 1273 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->gp[(f*prob->Nf + g)*4+3];} 1274 PetscFunctionReturn(0); 1275 } 1276 1277 /*@C 1278 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. 1279 1280 Not collective 1281 1282 Input Parameters: 1283 + prob - The PetscDS 1284 . f - The test field number 1285 . g - The field number 1286 . g0 - integrand for the test and basis function term 1287 . g1 - integrand for the test function and basis function gradient term 1288 . g2 - integrand for the test function gradient and basis function term 1289 - g3 - integrand for the test function gradient and basis function gradient term 1290 1291 Note: We are using a first order FEM model for the weak form: 1292 1293 \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 1294 1295 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1296 1297 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1298 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1299 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1300 $ PetscReal t, const PetscReal x[], PetscScalar g0[]) 1301 1302 + dim - the spatial dimension 1303 . Nf - the number of fields 1304 . uOff - the offset into u[] and u_t[] for each field 1305 . uOff_x - the offset into u_x[] for each field 1306 . u - each field evaluated at the current point 1307 . u_t - the time derivative of each field evaluated at the current point 1308 . u_x - the gradient of each field evaluated at the current point 1309 . aOff - the offset into a[] and a_t[] for each auxiliary field 1310 . aOff_x - the offset into a_x[] for each auxiliary field 1311 . a - each auxiliary field evaluated at the current point 1312 . a_t - the time derivative of each auxiliary field evaluated at the current point 1313 . a_x - the gradient of auxiliary each field evaluated at the current point 1314 . t - current time 1315 . u_tShift - the multiplier a for dF/dU_t 1316 . x - coordinates of the current point 1317 . numConstants - number of constant parameters 1318 . constants - constant parameters 1319 - g0 - output values at the current point 1320 1321 Level: intermediate 1322 1323 .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobian() 1324 @*/ 1325 PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g, 1326 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1327 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1328 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1329 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1330 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1331 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1332 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1333 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1334 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1335 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1336 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1337 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1338 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1339 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1340 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1341 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1342 { 1343 PetscErrorCode ierr; 1344 1345 PetscFunctionBegin; 1346 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1347 if (g0) PetscValidFunction(g0, 4); 1348 if (g1) PetscValidFunction(g1, 5); 1349 if (g2) PetscValidFunction(g2, 6); 1350 if (g3) PetscValidFunction(g3, 7); 1351 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1352 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1353 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1354 prob->gp[(f*prob->Nf + g)*4+0] = g0; 1355 prob->gp[(f*prob->Nf + g)*4+1] = g1; 1356 prob->gp[(f*prob->Nf + g)*4+2] = g2; 1357 prob->gp[(f*prob->Nf + g)*4+3] = g3; 1358 PetscFunctionReturn(0); 1359 } 1360 1361 /*@C 1362 PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set 1363 1364 Not collective 1365 1366 Input Parameter: 1367 . prob - The PetscDS 1368 1369 Output Parameter: 1370 . hasDynJac - flag that pointwise function for dynamic Jacobian has been set 1371 1372 Level: intermediate 1373 1374 .seealso: PetscDSGetDynamicJacobian(), PetscDSSetDynamicJacobian(), PetscDSGetJacobian() 1375 @*/ 1376 PetscErrorCode PetscDSHasDynamicJacobian(PetscDS prob, PetscBool *hasDynJac) 1377 { 1378 PetscInt f, g, h; 1379 1380 PetscFunctionBegin; 1381 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1382 *hasDynJac = PETSC_FALSE; 1383 for (f = 0; f < prob->Nf; ++f) { 1384 for (g = 0; g < prob->Nf; ++g) { 1385 for (h = 0; h < 4; ++h) { 1386 if (prob->gt[(f*prob->Nf + g)*4+h]) *hasDynJac = PETSC_TRUE; 1387 } 1388 } 1389 } 1390 PetscFunctionReturn(0); 1391 } 1392 1393 /*@C 1394 PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field 1395 1396 Not collective 1397 1398 Input Parameters: 1399 + prob - The PetscDS 1400 . f - The test field number 1401 - g - The field number 1402 1403 Output Parameters: 1404 + g0 - integrand for the test and basis function term 1405 . g1 - integrand for the test function and basis function gradient term 1406 . g2 - integrand for the test function gradient and basis function term 1407 - g3 - integrand for the test function gradient and basis function gradient term 1408 1409 Note: We are using a first order FEM model for the weak form: 1410 1411 \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 1412 1413 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1414 1415 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1416 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1417 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1418 $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]) 1419 1420 + dim - the spatial dimension 1421 . Nf - the number of fields 1422 . uOff - the offset into u[] and u_t[] for each field 1423 . uOff_x - the offset into u_x[] for each field 1424 . u - each field evaluated at the current point 1425 . u_t - the time derivative of each field evaluated at the current point 1426 . u_x - the gradient of each field evaluated at the current point 1427 . aOff - the offset into a[] and a_t[] for each auxiliary field 1428 . aOff_x - the offset into a_x[] for each auxiliary field 1429 . a - each auxiliary field evaluated at the current point 1430 . a_t - the time derivative of each auxiliary field evaluated at the current point 1431 . a_x - the gradient of auxiliary each field evaluated at the current point 1432 . t - current time 1433 . u_tShift - the multiplier a for dF/dU_t 1434 . x - coordinates of the current point 1435 . numConstants - number of constant parameters 1436 . constants - constant parameters 1437 - g0 - output values at the current point 1438 1439 Level: intermediate 1440 1441 .seealso: PetscDSSetJacobian() 1442 @*/ 1443 PetscErrorCode PetscDSGetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g, 1444 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1445 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1446 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1447 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1448 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1449 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1450 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1451 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1452 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1453 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1454 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1455 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1456 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1457 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1458 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1459 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1460 { 1461 PetscFunctionBegin; 1462 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1463 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); 1464 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); 1465 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->gt[(f*prob->Nf + g)*4+0];} 1466 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->gt[(f*prob->Nf + g)*4+1];} 1467 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->gt[(f*prob->Nf + g)*4+2];} 1468 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->gt[(f*prob->Nf + g)*4+3];} 1469 PetscFunctionReturn(0); 1470 } 1471 1472 /*@C 1473 PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields 1474 1475 Not collective 1476 1477 Input Parameters: 1478 + prob - The PetscDS 1479 . f - The test field number 1480 . g - The field number 1481 . g0 - integrand for the test and basis function term 1482 . g1 - integrand for the test function and basis function gradient term 1483 . g2 - integrand for the test function gradient and basis function term 1484 - g3 - integrand for the test function gradient and basis function gradient term 1485 1486 Note: We are using a first order FEM model for the weak form: 1487 1488 \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 1489 1490 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1491 1492 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1493 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1494 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1495 $ PetscReal t, const PetscReal x[], PetscScalar g0[]) 1496 1497 + dim - the spatial dimension 1498 . Nf - the number of fields 1499 . uOff - the offset into u[] and u_t[] for each field 1500 . uOff_x - the offset into u_x[] for each field 1501 . u - each field evaluated at the current point 1502 . u_t - the time derivative of each field evaluated at the current point 1503 . u_x - the gradient of each field evaluated at the current point 1504 . aOff - the offset into a[] and a_t[] for each auxiliary field 1505 . aOff_x - the offset into a_x[] for each auxiliary field 1506 . a - each auxiliary field evaluated at the current point 1507 . a_t - the time derivative of each auxiliary field evaluated at the current point 1508 . a_x - the gradient of auxiliary each field evaluated at the current point 1509 . t - current time 1510 . u_tShift - the multiplier a for dF/dU_t 1511 . x - coordinates of the current point 1512 . numConstants - number of constant parameters 1513 . constants - constant parameters 1514 - g0 - output values at the current point 1515 1516 Level: intermediate 1517 1518 .seealso: PetscDSGetJacobian() 1519 @*/ 1520 PetscErrorCode PetscDSSetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g, 1521 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1522 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1523 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1524 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1525 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1526 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1527 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1528 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1529 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1530 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1531 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1532 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1533 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1534 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1535 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1536 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1537 { 1538 PetscErrorCode ierr; 1539 1540 PetscFunctionBegin; 1541 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1542 if (g0) PetscValidFunction(g0, 4); 1543 if (g1) PetscValidFunction(g1, 5); 1544 if (g2) PetscValidFunction(g2, 6); 1545 if (g3) PetscValidFunction(g3, 7); 1546 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1547 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1548 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1549 prob->gt[(f*prob->Nf + g)*4+0] = g0; 1550 prob->gt[(f*prob->Nf + g)*4+1] = g1; 1551 prob->gt[(f*prob->Nf + g)*4+2] = g2; 1552 prob->gt[(f*prob->Nf + g)*4+3] = g3; 1553 PetscFunctionReturn(0); 1554 } 1555 1556 /*@C 1557 PetscDSGetRiemannSolver - Returns the Riemann solver for the given field 1558 1559 Not collective 1560 1561 Input Arguments: 1562 + prob - The PetscDS object 1563 - f - The field number 1564 1565 Output Argument: 1566 . r - Riemann solver 1567 1568 Calling sequence for r: 1569 1570 $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx) 1571 1572 + dim - The spatial dimension 1573 . Nf - The number of fields 1574 . x - The coordinates at a point on the interface 1575 . n - The normal vector to the interface 1576 . uL - The state vector to the left of the interface 1577 . uR - The state vector to the right of the interface 1578 . flux - output array of flux through the interface 1579 . numConstants - number of constant parameters 1580 . constants - constant parameters 1581 - ctx - optional user context 1582 1583 Level: intermediate 1584 1585 .seealso: PetscDSSetRiemannSolver() 1586 @*/ 1587 PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f, 1588 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)) 1589 { 1590 PetscFunctionBegin; 1591 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1592 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); 1593 PetscValidPointer(r, 3); 1594 *r = prob->r[f]; 1595 PetscFunctionReturn(0); 1596 } 1597 1598 /*@C 1599 PetscDSSetRiemannSolver - Sets the Riemann solver for the given field 1600 1601 Not collective 1602 1603 Input Arguments: 1604 + prob - The PetscDS object 1605 . f - The field number 1606 - r - Riemann solver 1607 1608 Calling sequence for r: 1609 1610 $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx) 1611 1612 + dim - The spatial dimension 1613 . Nf - The number of fields 1614 . x - The coordinates at a point on the interface 1615 . n - The normal vector to the interface 1616 . uL - The state vector to the left of the interface 1617 . uR - The state vector to the right of the interface 1618 . flux - output array of flux through the interface 1619 . numConstants - number of constant parameters 1620 . constants - constant parameters 1621 - ctx - optional user context 1622 1623 Level: intermediate 1624 1625 .seealso: PetscDSGetRiemannSolver() 1626 @*/ 1627 PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f, 1628 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)) 1629 { 1630 PetscErrorCode ierr; 1631 1632 PetscFunctionBegin; 1633 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1634 if (r) PetscValidFunction(r, 3); 1635 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1636 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1637 prob->r[f] = r; 1638 PetscFunctionReturn(0); 1639 } 1640 1641 PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx) 1642 { 1643 PetscFunctionBegin; 1644 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1645 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); 1646 PetscValidPointer(ctx, 3); 1647 *ctx = prob->ctx[f]; 1648 PetscFunctionReturn(0); 1649 } 1650 1651 PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx) 1652 { 1653 PetscErrorCode ierr; 1654 1655 PetscFunctionBegin; 1656 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1657 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1658 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1659 prob->ctx[f] = ctx; 1660 PetscFunctionReturn(0); 1661 } 1662 1663 /*@C 1664 PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field 1665 1666 Not collective 1667 1668 Input Parameters: 1669 + prob - The PetscDS 1670 - f - The test field number 1671 1672 Output Parameters: 1673 + f0 - boundary integrand for the test function term 1674 - f1 - boundary integrand for the test function gradient term 1675 1676 Note: We are using a first order FEM model for the weak form: 1677 1678 \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 1679 1680 The calling sequence for the callbacks f0 and f1 is given by: 1681 1682 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1683 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1684 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1685 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]) 1686 1687 + dim - the spatial dimension 1688 . Nf - the number of fields 1689 . uOff - the offset into u[] and u_t[] for each field 1690 . uOff_x - the offset into u_x[] for each field 1691 . u - each field evaluated at the current point 1692 . u_t - the time derivative of each field evaluated at the current point 1693 . u_x - the gradient of each field evaluated at the current point 1694 . aOff - the offset into a[] and a_t[] for each auxiliary field 1695 . aOff_x - the offset into a_x[] for each auxiliary field 1696 . a - each auxiliary field evaluated at the current point 1697 . a_t - the time derivative of each auxiliary field evaluated at the current point 1698 . a_x - the gradient of auxiliary each field evaluated at the current point 1699 . t - current time 1700 . x - coordinates of the current point 1701 . n - unit normal at the current point 1702 . numConstants - number of constant parameters 1703 . constants - constant parameters 1704 - f0 - output values at the current point 1705 1706 Level: intermediate 1707 1708 .seealso: PetscDSSetBdResidual() 1709 @*/ 1710 PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f, 1711 void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1712 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1713 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1714 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), 1715 void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1716 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1717 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1718 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])) 1719 { 1720 PetscFunctionBegin; 1721 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1722 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); 1723 if (f0) {PetscValidPointer(f0, 3); *f0 = prob->fBd[f*2+0];} 1724 if (f1) {PetscValidPointer(f1, 4); *f1 = prob->fBd[f*2+1];} 1725 PetscFunctionReturn(0); 1726 } 1727 1728 /*@C 1729 PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field 1730 1731 Not collective 1732 1733 Input Parameters: 1734 + prob - The PetscDS 1735 . f - The test field number 1736 . f0 - boundary integrand for the test function term 1737 - f1 - boundary integrand for the test function gradient term 1738 1739 Note: We are using a first order FEM model for the weak form: 1740 1741 \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 1742 1743 The calling sequence for the callbacks f0 and f1 is given by: 1744 1745 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1746 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1747 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1748 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]) 1749 1750 + dim - the spatial dimension 1751 . Nf - the number of fields 1752 . uOff - the offset into u[] and u_t[] for each field 1753 . uOff_x - the offset into u_x[] for each field 1754 . u - each field evaluated at the current point 1755 . u_t - the time derivative of each field evaluated at the current point 1756 . u_x - the gradient of each field evaluated at the current point 1757 . aOff - the offset into a[] and a_t[] for each auxiliary field 1758 . aOff_x - the offset into a_x[] for each auxiliary field 1759 . a - each auxiliary field evaluated at the current point 1760 . a_t - the time derivative of each auxiliary field evaluated at the current point 1761 . a_x - the gradient of auxiliary each field evaluated at the current point 1762 . t - current time 1763 . x - coordinates of the current point 1764 . n - unit normal at the current point 1765 . numConstants - number of constant parameters 1766 . constants - constant parameters 1767 - f0 - output values at the current point 1768 1769 Level: intermediate 1770 1771 .seealso: PetscDSGetBdResidual() 1772 @*/ 1773 PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f, 1774 void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1775 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1776 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1777 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), 1778 void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1779 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1780 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1781 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])) 1782 { 1783 PetscErrorCode ierr; 1784 1785 PetscFunctionBegin; 1786 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1787 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1788 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1789 if (f0) {PetscValidFunction(f0, 3); prob->fBd[f*2+0] = f0;} 1790 if (f1) {PetscValidFunction(f1, 4); prob->fBd[f*2+1] = f1;} 1791 PetscFunctionReturn(0); 1792 } 1793 1794 /*@C 1795 PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field 1796 1797 Not collective 1798 1799 Input Parameters: 1800 + prob - The PetscDS 1801 . f - The test field number 1802 - g - The field number 1803 1804 Output Parameters: 1805 + g0 - integrand for the test and basis function term 1806 . g1 - integrand for the test function and basis function gradient term 1807 . g2 - integrand for the test function gradient and basis function term 1808 - g3 - integrand for the test function gradient and basis function gradient term 1809 1810 Note: We are using a first order FEM model for the weak form: 1811 1812 \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 1813 1814 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1815 1816 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1817 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1818 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1819 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[]) 1820 1821 + dim - the spatial dimension 1822 . Nf - the number of fields 1823 . uOff - the offset into u[] and u_t[] for each field 1824 . uOff_x - the offset into u_x[] for each field 1825 . u - each field evaluated at the current point 1826 . u_t - the time derivative of each field evaluated at the current point 1827 . u_x - the gradient of each field evaluated at the current point 1828 . aOff - the offset into a[] and a_t[] for each auxiliary field 1829 . aOff_x - the offset into a_x[] for each auxiliary field 1830 . a - each auxiliary field evaluated at the current point 1831 . a_t - the time derivative of each auxiliary field evaluated at the current point 1832 . a_x - the gradient of auxiliary each field evaluated at the current point 1833 . t - current time 1834 . u_tShift - the multiplier a for dF/dU_t 1835 . x - coordinates of the current point 1836 . n - normal at the current point 1837 . numConstants - number of constant parameters 1838 . constants - constant parameters 1839 - g0 - output values at the current point 1840 1841 Level: intermediate 1842 1843 .seealso: PetscDSSetBdJacobian() 1844 @*/ 1845 PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g, 1846 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1847 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1848 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1849 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1850 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1851 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1852 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1853 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1854 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1855 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1856 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1857 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1858 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1859 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1860 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1861 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1862 { 1863 PetscFunctionBegin; 1864 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1865 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); 1866 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); 1867 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->gBd[(f*prob->Nf + g)*4+0];} 1868 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->gBd[(f*prob->Nf + g)*4+1];} 1869 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->gBd[(f*prob->Nf + g)*4+2];} 1870 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->gBd[(f*prob->Nf + g)*4+3];} 1871 PetscFunctionReturn(0); 1872 } 1873 1874 /*@C 1875 PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field 1876 1877 Not collective 1878 1879 Input Parameters: 1880 + prob - The PetscDS 1881 . f - The test field number 1882 . g - The field number 1883 . g0 - integrand for the test and basis function term 1884 . g1 - integrand for the test function and basis function gradient term 1885 . g2 - integrand for the test function gradient and basis function term 1886 - g3 - integrand for the test function gradient and basis function gradient term 1887 1888 Note: We are using a first order FEM model for the weak form: 1889 1890 \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 1891 1892 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1893 1894 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1895 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1896 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1897 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[]) 1898 1899 + dim - the spatial dimension 1900 . Nf - the number of fields 1901 . uOff - the offset into u[] and u_t[] for each field 1902 . uOff_x - the offset into u_x[] for each field 1903 . u - each field evaluated at the current point 1904 . u_t - the time derivative of each field evaluated at the current point 1905 . u_x - the gradient of each field evaluated at the current point 1906 . aOff - the offset into a[] and a_t[] for each auxiliary field 1907 . aOff_x - the offset into a_x[] for each auxiliary field 1908 . a - each auxiliary field evaluated at the current point 1909 . a_t - the time derivative of each auxiliary field evaluated at the current point 1910 . a_x - the gradient of auxiliary each field evaluated at the current point 1911 . t - current time 1912 . u_tShift - the multiplier a for dF/dU_t 1913 . x - coordinates of the current point 1914 . n - normal at the current point 1915 . numConstants - number of constant parameters 1916 . constants - constant parameters 1917 - g0 - output values at the current point 1918 1919 Level: intermediate 1920 1921 .seealso: PetscDSGetBdJacobian() 1922 @*/ 1923 PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g, 1924 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1925 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1926 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1927 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), 1928 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1929 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1930 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1931 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), 1932 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1933 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1934 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1935 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), 1936 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1937 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1938 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1939 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])) 1940 { 1941 PetscErrorCode ierr; 1942 1943 PetscFunctionBegin; 1944 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1945 if (g0) PetscValidFunction(g0, 4); 1946 if (g1) PetscValidFunction(g1, 5); 1947 if (g2) PetscValidFunction(g2, 6); 1948 if (g3) PetscValidFunction(g3, 7); 1949 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1950 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1951 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1952 prob->gBd[(f*prob->Nf + g)*4+0] = g0; 1953 prob->gBd[(f*prob->Nf + g)*4+1] = g1; 1954 prob->gBd[(f*prob->Nf + g)*4+2] = g2; 1955 prob->gBd[(f*prob->Nf + g)*4+3] = g3; 1956 PetscFunctionReturn(0); 1957 } 1958 1959 1960 /*@ 1961 PetscDSGetConstants - Returns the array of constants passed to point functions 1962 1963 Not collective 1964 1965 Input Parameter: 1966 . prob - The PetscDS object 1967 1968 Output Parameters: 1969 + numConstants - The number of constants 1970 - constants - The array of constants, NULL if there are none 1971 1972 Level: intermediate 1973 1974 .seealso: PetscDSSetConstants(), PetscDSCreate() 1975 @*/ 1976 PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[]) 1977 { 1978 PetscFunctionBegin; 1979 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1980 if (numConstants) {PetscValidPointer(numConstants, 2); *numConstants = prob->numConstants;} 1981 if (constants) {PetscValidPointer(constants, 3); *constants = prob->constants;} 1982 PetscFunctionReturn(0); 1983 } 1984 1985 /*@ 1986 PetscDSSetConstants - Set the array of constants passed to point functions 1987 1988 Not collective 1989 1990 Input Parameters: 1991 + prob - The PetscDS object 1992 . numConstants - The number of constants 1993 - constants - The array of constants, NULL if there are none 1994 1995 Level: intermediate 1996 1997 .seealso: PetscDSGetConstants(), PetscDSCreate() 1998 @*/ 1999 PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[]) 2000 { 2001 PetscErrorCode ierr; 2002 2003 PetscFunctionBegin; 2004 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2005 if (numConstants != prob->numConstants) { 2006 ierr = PetscFree(prob->constants);CHKERRQ(ierr); 2007 prob->constants = NULL; 2008 } 2009 prob->numConstants = numConstants; 2010 if (prob->numConstants) { 2011 PetscValidPointer(constants, 3); 2012 ierr = PetscMalloc1(prob->numConstants, &prob->constants);CHKERRQ(ierr); 2013 ierr = PetscMemcpy(prob->constants, constants, prob->numConstants * sizeof(PetscScalar));CHKERRQ(ierr); 2014 } 2015 PetscFunctionReturn(0); 2016 } 2017 2018 /*@ 2019 PetscDSGetFieldIndex - Returns the index of the given field 2020 2021 Not collective 2022 2023 Input Parameters: 2024 + prob - The PetscDS object 2025 - disc - The discretization object 2026 2027 Output Parameter: 2028 . f - The field number 2029 2030 Level: beginner 2031 2032 .seealso: PetscGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 2033 @*/ 2034 PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f) 2035 { 2036 PetscInt g; 2037 2038 PetscFunctionBegin; 2039 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2040 PetscValidPointer(f, 3); 2041 *f = -1; 2042 for (g = 0; g < prob->Nf; ++g) {if (disc == prob->disc[g]) break;} 2043 if (g == prob->Nf) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS."); 2044 *f = g; 2045 PetscFunctionReturn(0); 2046 } 2047 2048 /*@ 2049 PetscDSGetFieldSize - Returns the size of the given field in the full space basis 2050 2051 Not collective 2052 2053 Input Parameters: 2054 + prob - The PetscDS object 2055 - f - The field number 2056 2057 Output Parameter: 2058 . size - The size 2059 2060 Level: beginner 2061 2062 .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 2063 @*/ 2064 PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size) 2065 { 2066 PetscErrorCode ierr; 2067 2068 PetscFunctionBegin; 2069 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2070 PetscValidPointer(size, 3); 2071 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); 2072 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2073 *size = prob->Nb[f]; 2074 PetscFunctionReturn(0); 2075 } 2076 2077 /*@ 2078 PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis 2079 2080 Not collective 2081 2082 Input Parameters: 2083 + prob - The PetscDS object 2084 - f - The field number 2085 2086 Output Parameter: 2087 . off - The offset 2088 2089 Level: beginner 2090 2091 .seealso: PetscDSGetFieldSize(), PetscDSGetNumFields(), PetscDSCreate() 2092 @*/ 2093 PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off) 2094 { 2095 PetscInt size, g; 2096 PetscErrorCode ierr; 2097 2098 PetscFunctionBegin; 2099 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2100 PetscValidPointer(off, 3); 2101 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); 2102 *off = 0; 2103 for (g = 0; g < f; ++g) { 2104 ierr = PetscDSGetFieldSize(prob, g, &size);CHKERRQ(ierr); 2105 *off += size; 2106 } 2107 PetscFunctionReturn(0); 2108 } 2109 2110 /*@ 2111 PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point 2112 2113 Not collective 2114 2115 Input Parameter: 2116 . prob - The PetscDS object 2117 2118 Output Parameter: 2119 . dimensions - The number of dimensions 2120 2121 Level: beginner 2122 2123 .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate() 2124 @*/ 2125 PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[]) 2126 { 2127 PetscErrorCode ierr; 2128 2129 PetscFunctionBegin; 2130 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2131 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2132 PetscValidPointer(dimensions, 2); 2133 *dimensions = prob->Nb; 2134 PetscFunctionReturn(0); 2135 } 2136 2137 /*@ 2138 PetscDSGetComponents - Returns the number of components for each field on an evaluation point 2139 2140 Not collective 2141 2142 Input Parameter: 2143 . prob - The PetscDS object 2144 2145 Output Parameter: 2146 . components - The number of components 2147 2148 Level: beginner 2149 2150 .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate() 2151 @*/ 2152 PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[]) 2153 { 2154 PetscErrorCode ierr; 2155 2156 PetscFunctionBegin; 2157 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2158 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2159 PetscValidPointer(components, 2); 2160 *components = prob->Nc; 2161 PetscFunctionReturn(0); 2162 } 2163 2164 /*@ 2165 PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point 2166 2167 Not collective 2168 2169 Input Parameters: 2170 + prob - The PetscDS object 2171 - f - The field number 2172 2173 Output Parameter: 2174 . off - The offset 2175 2176 Level: beginner 2177 2178 .seealso: PetscDSGetNumFields(), PetscDSCreate() 2179 @*/ 2180 PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off) 2181 { 2182 PetscFunctionBegin; 2183 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2184 PetscValidPointer(off, 3); 2185 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); 2186 *off = prob->off[f]; 2187 PetscFunctionReturn(0); 2188 } 2189 2190 /*@ 2191 PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point 2192 2193 Not collective 2194 2195 Input Parameter: 2196 . prob - The PetscDS object 2197 2198 Output Parameter: 2199 . offsets - The offsets 2200 2201 Level: beginner 2202 2203 .seealso: PetscDSGetNumFields(), PetscDSCreate() 2204 @*/ 2205 PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[]) 2206 { 2207 PetscFunctionBegin; 2208 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2209 PetscValidPointer(offsets, 2); 2210 *offsets = prob->off; 2211 PetscFunctionReturn(0); 2212 } 2213 2214 /*@ 2215 PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point 2216 2217 Not collective 2218 2219 Input Parameter: 2220 . prob - The PetscDS object 2221 2222 Output Parameter: 2223 . offsets - The offsets 2224 2225 Level: beginner 2226 2227 .seealso: PetscDSGetNumFields(), PetscDSCreate() 2228 @*/ 2229 PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[]) 2230 { 2231 PetscFunctionBegin; 2232 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2233 PetscValidPointer(offsets, 2); 2234 *offsets = prob->offDer; 2235 PetscFunctionReturn(0); 2236 } 2237 2238 /*@C 2239 PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization 2240 2241 Not collective 2242 2243 Input Parameter: 2244 . prob - The PetscDS object 2245 2246 Output Parameters: 2247 + basis - The basis function tabulation at quadrature points 2248 - basisDer - The basis function derivative tabulation at quadrature points 2249 2250 Level: intermediate 2251 2252 .seealso: PetscDSCreate() 2253 @*/ 2254 PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer) 2255 { 2256 PetscErrorCode ierr; 2257 2258 PetscFunctionBegin; 2259 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2260 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2261 if (basis) {PetscValidPointer(basis, 2); *basis = prob->basis;} 2262 if (basisDer) {PetscValidPointer(basisDer, 3); *basisDer = prob->basisDer;} 2263 PetscFunctionReturn(0); 2264 } 2265 2266 /*@C 2267 PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces 2268 2269 Not collective 2270 2271 Input Parameter: 2272 . prob - The PetscDS object 2273 2274 Output Parameters: 2275 + basisFace - The basis function tabulation at quadrature points 2276 - basisDerFace - The basis function derivative tabulation at quadrature points 2277 2278 Level: intermediate 2279 2280 .seealso: PetscDSGetTabulation(), PetscDSCreate() 2281 @*/ 2282 PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer) 2283 { 2284 PetscErrorCode ierr; 2285 2286 PetscFunctionBegin; 2287 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2288 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2289 if (basis) {PetscValidPointer(basis, 2); *basis = prob->basisFace;} 2290 if (basisDer) {PetscValidPointer(basisDer, 3); *basisDer = prob->basisDerFace;} 2291 PetscFunctionReturn(0); 2292 } 2293 2294 PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x) 2295 { 2296 PetscErrorCode ierr; 2297 2298 PetscFunctionBegin; 2299 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2300 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2301 if (u) {PetscValidPointer(u, 2); *u = prob->u;} 2302 if (u_t) {PetscValidPointer(u_t, 3); *u_t = prob->u_t;} 2303 if (u_x) {PetscValidPointer(u_x, 4); *u_x = prob->u_x;} 2304 PetscFunctionReturn(0); 2305 } 2306 2307 PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3) 2308 { 2309 PetscErrorCode ierr; 2310 2311 PetscFunctionBegin; 2312 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2313 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2314 if (f0) {PetscValidPointer(f0, 2); *f0 = prob->f0;} 2315 if (f1) {PetscValidPointer(f1, 3); *f1 = prob->f1;} 2316 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->g0;} 2317 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->g1;} 2318 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->g2;} 2319 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->g3;} 2320 PetscFunctionReturn(0); 2321 } 2322 2323 PetscErrorCode PetscDSGetRefCoordArrays(PetscDS prob, PetscReal **x, PetscScalar **refSpaceDer) 2324 { 2325 PetscErrorCode ierr; 2326 2327 PetscFunctionBegin; 2328 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2329 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 2330 if (x) {PetscValidPointer(x, 2); *x = prob->x;} 2331 if (refSpaceDer) {PetscValidPointer(refSpaceDer, 3); *refSpaceDer = prob->refSpaceDer;} 2332 PetscFunctionReturn(0); 2333 } 2334 2335 /*@C 2336 PetscDSAddBoundary - Add a boundary condition to the model 2337 2338 Input Parameters: 2339 + ds - The PetscDS object 2340 . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann) 2341 . name - The BC name 2342 . labelname - The label defining constrained points 2343 . field - The field to constrain 2344 . numcomps - The number of constrained field components 2345 . comps - An array of constrained component numbers 2346 . bcFunc - A pointwise function giving boundary values 2347 . numids - The number of DMLabel ids for constrained points 2348 . ids - An array of ids for constrained points 2349 - ctx - An optional user context for bcFunc 2350 2351 Options Database Keys: 2352 + -bc_<boundary name> <num> - Overrides the boundary ids 2353 - -bc_<boundary name>_comp <num> - Overrides the boundary components 2354 2355 Level: developer 2356 2357 .seealso: PetscDSGetBoundary() 2358 @*/ 2359 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) 2360 { 2361 DSBoundary b; 2362 PetscErrorCode ierr; 2363 2364 PetscFunctionBegin; 2365 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2366 ierr = PetscNew(&b);CHKERRQ(ierr); 2367 ierr = PetscStrallocpy(name, (char **) &b->name);CHKERRQ(ierr); 2368 ierr = PetscStrallocpy(labelname, (char **) &b->labelname);CHKERRQ(ierr); 2369 ierr = PetscMalloc1(numcomps, &b->comps);CHKERRQ(ierr); 2370 if (numcomps) {ierr = PetscMemcpy(b->comps, comps, numcomps*sizeof(PetscInt));CHKERRQ(ierr);} 2371 ierr = PetscMalloc1(numids, &b->ids);CHKERRQ(ierr); 2372 if (numids) {ierr = PetscMemcpy(b->ids, ids, numids*sizeof(PetscInt));CHKERRQ(ierr);} 2373 b->type = type; 2374 b->field = field; 2375 b->numcomps = numcomps; 2376 b->func = bcFunc; 2377 b->numids = numids; 2378 b->ctx = ctx; 2379 b->next = ds->boundary; 2380 ds->boundary = b; 2381 PetscFunctionReturn(0); 2382 } 2383 2384 /*@ 2385 PetscDSGetNumBoundary - Get the number of registered BC 2386 2387 Input Parameters: 2388 . ds - The PetscDS object 2389 2390 Output Parameters: 2391 . numBd - The number of BC 2392 2393 Level: intermediate 2394 2395 .seealso: PetscDSAddBoundary(), PetscDSGetBoundary() 2396 @*/ 2397 PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd) 2398 { 2399 DSBoundary b = ds->boundary; 2400 2401 PetscFunctionBegin; 2402 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2403 PetscValidPointer(numBd, 2); 2404 *numBd = 0; 2405 while (b) {++(*numBd); b = b->next;} 2406 PetscFunctionReturn(0); 2407 } 2408 2409 /*@C 2410 PetscDSGetBoundary - Add a boundary condition to the model 2411 2412 Input Parameters: 2413 + ds - The PetscDS object 2414 - bd - The BC number 2415 2416 Output Parameters: 2417 + type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann) 2418 . name - The BC name 2419 . labelname - The label defining constrained points 2420 . field - The field to constrain 2421 . numcomps - The number of constrained field components 2422 . comps - An array of constrained component numbers 2423 . bcFunc - A pointwise function giving boundary values 2424 . numids - The number of DMLabel ids for constrained points 2425 . ids - An array of ids for constrained points 2426 - ctx - An optional user context for bcFunc 2427 2428 Options Database Keys: 2429 + -bc_<boundary name> <num> - Overrides the boundary ids 2430 - -bc_<boundary name>_comp <num> - Overrides the boundary components 2431 2432 Level: developer 2433 2434 .seealso: PetscDSAddBoundary() 2435 @*/ 2436 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) 2437 { 2438 DSBoundary b = ds->boundary; 2439 PetscInt n = 0; 2440 2441 PetscFunctionBegin; 2442 PetscValidHeaderSpecific(ds, PETSCDS_CLASSID, 1); 2443 while (b) { 2444 if (n == bd) break; 2445 b = b->next; 2446 ++n; 2447 } 2448 if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n); 2449 if (type) { 2450 PetscValidPointer(type, 3); 2451 *type = b->type; 2452 } 2453 if (name) { 2454 PetscValidPointer(name, 4); 2455 *name = b->name; 2456 } 2457 if (labelname) { 2458 PetscValidPointer(labelname, 5); 2459 *labelname = b->labelname; 2460 } 2461 if (field) { 2462 PetscValidPointer(field, 6); 2463 *field = b->field; 2464 } 2465 if (numcomps) { 2466 PetscValidPointer(numcomps, 7); 2467 *numcomps = b->numcomps; 2468 } 2469 if (comps) { 2470 PetscValidPointer(comps, 8); 2471 *comps = b->comps; 2472 } 2473 if (func) { 2474 PetscValidPointer(func, 9); 2475 *func = b->func; 2476 } 2477 if (numids) { 2478 PetscValidPointer(numids, 10); 2479 *numids = b->numids; 2480 } 2481 if (ids) { 2482 PetscValidPointer(ids, 11); 2483 *ids = b->ids; 2484 } 2485 if (ctx) { 2486 PetscValidPointer(ctx, 12); 2487 *ctx = b->ctx; 2488 } 2489 PetscFunctionReturn(0); 2490 } 2491 2492 PetscErrorCode PetscDSCopyBoundary(PetscDS probA, PetscDS probB) 2493 { 2494 DSBoundary b, next, *lastnext; 2495 PetscErrorCode ierr; 2496 2497 PetscFunctionBegin; 2498 PetscValidHeaderSpecific(probA, PETSCDS_CLASSID, 1); 2499 PetscValidHeaderSpecific(probB, PETSCDS_CLASSID, 2); 2500 if (probA == probB) PetscFunctionReturn(0); 2501 next = probB->boundary; 2502 while (next) { 2503 DSBoundary b = next; 2504 2505 next = b->next; 2506 ierr = PetscFree(b->comps);CHKERRQ(ierr); 2507 ierr = PetscFree(b->ids);CHKERRQ(ierr); 2508 ierr = PetscFree(b->name);CHKERRQ(ierr); 2509 ierr = PetscFree(b->labelname);CHKERRQ(ierr); 2510 ierr = PetscFree(b);CHKERRQ(ierr); 2511 } 2512 lastnext = &(probB->boundary); 2513 for (b = probA->boundary; b; b = b->next) { 2514 DSBoundary bNew; 2515 2516 ierr = PetscNew(&bNew);CHKERRQ(ierr); 2517 bNew->numcomps = b->numcomps; 2518 ierr = PetscMalloc1(bNew->numcomps, &bNew->comps);CHKERRQ(ierr); 2519 ierr = PetscMemcpy(bNew->comps, b->comps, bNew->numcomps*sizeof(PetscInt));CHKERRQ(ierr); 2520 bNew->numids = b->numids; 2521 ierr = PetscMalloc1(bNew->numids, &bNew->ids);CHKERRQ(ierr); 2522 ierr = PetscMemcpy(bNew->ids, b->ids, bNew->numids*sizeof(PetscInt));CHKERRQ(ierr); 2523 ierr = PetscStrallocpy(b->labelname,(char **) &(bNew->labelname));CHKERRQ(ierr); 2524 ierr = PetscStrallocpy(b->name,(char **) &(bNew->name));CHKERRQ(ierr); 2525 bNew->ctx = b->ctx; 2526 bNew->type = b->type; 2527 bNew->field = b->field; 2528 bNew->func = b->func; 2529 2530 *lastnext = bNew; 2531 lastnext = &(bNew->next); 2532 } 2533 PetscFunctionReturn(0); 2534 } 2535 2536 /*@ 2537 PetscDSCopyEquations - Copy all pointwise function pointers to the new problem 2538 2539 Not collective 2540 2541 Input Parameter: 2542 . prob - The PetscDS object 2543 2544 Output Parameter: 2545 . newprob - The PetscDS copy 2546 2547 Level: intermediate 2548 2549 .seealso: PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate() 2550 @*/ 2551 PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob) 2552 { 2553 PetscInt Nf, Ng, f, g; 2554 PetscErrorCode ierr; 2555 2556 PetscFunctionBegin; 2557 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2558 PetscValidHeaderSpecific(newprob, PETSCDS_CLASSID, 2); 2559 ierr = PetscDSGetNumFields(prob, &Nf);CHKERRQ(ierr); 2560 ierr = PetscDSGetNumFields(newprob, &Ng);CHKERRQ(ierr); 2561 if (Nf != Ng) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %D != %D", Nf, Ng);CHKERRQ(ierr); 2562 for (f = 0; f < Nf; ++f) { 2563 PetscPointFunc obj; 2564 PetscPointFunc f0, f1; 2565 PetscPointJac g0, g1, g2, g3; 2566 PetscBdPointFunc f0Bd, f1Bd; 2567 PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd; 2568 PetscRiemannFunc r; 2569 2570 ierr = PetscDSGetObjective(prob, f, &obj);CHKERRQ(ierr); 2571 ierr = PetscDSGetResidual(prob, f, &f0, &f1);CHKERRQ(ierr); 2572 ierr = PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);CHKERRQ(ierr); 2573 ierr = PetscDSGetRiemannSolver(prob, f, &r);CHKERRQ(ierr); 2574 ierr = PetscDSSetObjective(newprob, f, obj);CHKERRQ(ierr); 2575 ierr = PetscDSSetResidual(newprob, f, f0, f1);CHKERRQ(ierr); 2576 ierr = PetscDSSetBdResidual(newprob, f, f0Bd, f1Bd);CHKERRQ(ierr); 2577 ierr = PetscDSSetRiemannSolver(newprob, f, r);CHKERRQ(ierr); 2578 for (g = 0; g < Nf; ++g) { 2579 ierr = PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);CHKERRQ(ierr); 2580 ierr = PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);CHKERRQ(ierr); 2581 ierr = PetscDSSetJacobian(newprob, f, g, g0, g1, g2, g3);CHKERRQ(ierr); 2582 ierr = PetscDSSetBdJacobian(newprob, f, g, g0Bd, g1Bd, g2Bd, g3Bd);CHKERRQ(ierr); 2583 } 2584 } 2585 PetscFunctionReturn(0); 2586 } 2587 2588 static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob) 2589 { 2590 PetscErrorCode ierr; 2591 2592 PetscFunctionBegin; 2593 ierr = PetscFree(prob->data);CHKERRQ(ierr); 2594 PetscFunctionReturn(0); 2595 } 2596 2597 static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob) 2598 { 2599 PetscFunctionBegin; 2600 prob->ops->setfromoptions = NULL; 2601 prob->ops->setup = NULL; 2602 prob->ops->view = NULL; 2603 prob->ops->destroy = PetscDSDestroy_Basic; 2604 PetscFunctionReturn(0); 2605 } 2606 2607 /*MC 2608 PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions 2609 2610 Level: intermediate 2611 2612 .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType() 2613 M*/ 2614 2615 PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob) 2616 { 2617 PetscDS_Basic *b; 2618 PetscErrorCode ierr; 2619 2620 PetscFunctionBegin; 2621 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 2622 ierr = PetscNewLog(prob, &b);CHKERRQ(ierr); 2623 prob->data = b; 2624 2625 ierr = PetscDSInitialize_Basic(prob);CHKERRQ(ierr); 2626 PetscFunctionReturn(0); 2627 } 2628