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