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