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((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; 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 = PetscCalloc5(NfNew, &tmpobj, NfNew*2, &tmpf, NfNew*NfNew*4, &tmpg, 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; ++f) tmpr[f] = prob->r[f]; 378 for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f]; 379 for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL; 380 for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL; 381 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL; 382 for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL; 383 for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL; 384 ierr = PetscFree5(prob->obj, prob->f, prob->g, prob->r, prob->ctx);CHKERRQ(ierr); 385 prob->obj = tmpobj; 386 prob->f = tmpf; 387 prob->g = tmpg; 388 prob->r = tmpr; 389 prob->ctx = tmpctx; 390 ierr = PetscCalloc2(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd);CHKERRQ(ierr); 391 for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f]; 392 for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f]; 393 for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL; 394 for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL; 395 ierr = PetscFree2(prob->fBd, prob->gBd);CHKERRQ(ierr); 396 prob->fBd = tmpfbd; 397 prob->gBd = tmpgbd; 398 PetscFunctionReturn(0); 399 } 400 401 #undef __FUNCT__ 402 #define __FUNCT__ "PetscDSDestroy" 403 /*@ 404 PetscDSDestroy - Destroys a PetscDS object 405 406 Collective on PetscDS 407 408 Input Parameter: 409 . prob - the PetscDS object to destroy 410 411 Level: developer 412 413 .seealso PetscDSView() 414 @*/ 415 PetscErrorCode PetscDSDestroy(PetscDS *prob) 416 { 417 PetscInt f; 418 PetscErrorCode ierr; 419 420 PetscFunctionBegin; 421 if (!*prob) PetscFunctionReturn(0); 422 PetscValidHeaderSpecific((*prob), PETSCDS_CLASSID, 1); 423 424 if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; PetscFunctionReturn(0);} 425 ((PetscObject) (*prob))->refct = 0; 426 ierr = PetscDSDestroyStructs_Static(*prob);CHKERRQ(ierr); 427 for (f = 0; f < (*prob)->Nf; ++f) { 428 ierr = PetscObjectDereference((*prob)->disc[f]);CHKERRQ(ierr); 429 ierr = PetscObjectDereference((*prob)->discBd[f]);CHKERRQ(ierr); 430 } 431 ierr = PetscFree4((*prob)->disc, (*prob)->discBd, (*prob)->implicit, (*prob)->adjacency);CHKERRQ(ierr); 432 ierr = PetscFree5((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->r,(*prob)->ctx);CHKERRQ(ierr); 433 ierr = PetscFree2((*prob)->fBd,(*prob)->gBd);CHKERRQ(ierr); 434 if ((*prob)->ops->destroy) {ierr = (*(*prob)->ops->destroy)(*prob);CHKERRQ(ierr);} 435 ierr = PetscHeaderDestroy(prob);CHKERRQ(ierr); 436 PetscFunctionReturn(0); 437 } 438 439 #undef __FUNCT__ 440 #define __FUNCT__ "PetscDSCreate" 441 /*@ 442 PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType(). 443 444 Collective on MPI_Comm 445 446 Input Parameter: 447 . comm - The communicator for the PetscDS object 448 449 Output Parameter: 450 . prob - The PetscDS object 451 452 Level: beginner 453 454 .seealso: PetscDSSetType(), PETSCDSBASIC 455 @*/ 456 PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob) 457 { 458 PetscDS p; 459 PetscErrorCode ierr; 460 461 PetscFunctionBegin; 462 PetscValidPointer(prob, 2); 463 *prob = NULL; 464 ierr = PetscDSInitializePackage();CHKERRQ(ierr); 465 466 ierr = PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);CHKERRQ(ierr); 467 468 p->Nf = 0; 469 p->setup = PETSC_FALSE; 470 471 *prob = p; 472 PetscFunctionReturn(0); 473 } 474 475 #undef __FUNCT__ 476 #define __FUNCT__ "PetscDSGetNumFields" 477 /*@ 478 PetscDSGetNumFields - Returns the number of fields in the DS 479 480 Not collective 481 482 Input Parameter: 483 . prob - The PetscDS object 484 485 Output Parameter: 486 . Nf - The number of fields 487 488 Level: beginner 489 490 .seealso: PetscDSGetSpatialDimension(), PetscDSCreate() 491 @*/ 492 PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf) 493 { 494 PetscFunctionBegin; 495 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 496 PetscValidPointer(Nf, 2); 497 *Nf = prob->Nf; 498 PetscFunctionReturn(0); 499 } 500 501 #undef __FUNCT__ 502 #define __FUNCT__ "PetscDSGetSpatialDimension" 503 /*@ 504 PetscDSGetSpatialDimension - Returns the spatial dimension of the DS 505 506 Not collective 507 508 Input Parameter: 509 . prob - The PetscDS object 510 511 Output Parameter: 512 . dim - The spatial dimension 513 514 Level: beginner 515 516 .seealso: PetscDSGetNumFields(), PetscDSCreate() 517 @*/ 518 PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim) 519 { 520 PetscErrorCode ierr; 521 522 PetscFunctionBegin; 523 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 524 PetscValidPointer(dim, 2); 525 *dim = 0; 526 if (prob->Nf) { 527 PetscObject obj; 528 PetscClassId id; 529 530 ierr = PetscDSGetDiscretization(prob, 0, &obj);CHKERRQ(ierr); 531 ierr = PetscObjectGetClassId(obj, &id);CHKERRQ(ierr); 532 if (id == PETSCFE_CLASSID) {ierr = PetscFEGetSpatialDimension((PetscFE) obj, dim);CHKERRQ(ierr);} 533 else if (id == PETSCFV_CLASSID) {ierr = PetscFVGetSpatialDimension((PetscFV) obj, dim);CHKERRQ(ierr);} 534 else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0); 535 } 536 PetscFunctionReturn(0); 537 } 538 539 #undef __FUNCT__ 540 #define __FUNCT__ "PetscDSGetTotalDimension" 541 /*@ 542 PetscDSGetTotalDimension - Returns the total size of the approximation space for this system 543 544 Not collective 545 546 Input Parameter: 547 . prob - The PetscDS object 548 549 Output Parameter: 550 . dim - The total problem dimension 551 552 Level: beginner 553 554 .seealso: PetscDSGetNumFields(), PetscDSCreate() 555 @*/ 556 PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim) 557 { 558 PetscErrorCode ierr; 559 560 PetscFunctionBegin; 561 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 562 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 563 PetscValidPointer(dim, 2); 564 *dim = prob->totDim; 565 PetscFunctionReturn(0); 566 } 567 568 #undef __FUNCT__ 569 #define __FUNCT__ "PetscDSGetTotalBdDimension" 570 /*@ 571 PetscDSGetTotalBdDimension - Returns the total size of the boundary approximation space for this system 572 573 Not collective 574 575 Input Parameter: 576 . prob - The PetscDS object 577 578 Output Parameter: 579 . dim - The total boundary problem dimension 580 581 Level: beginner 582 583 .seealso: PetscDSGetNumFields(), PetscDSCreate() 584 @*/ 585 PetscErrorCode PetscDSGetTotalBdDimension(PetscDS prob, PetscInt *dim) 586 { 587 PetscErrorCode ierr; 588 589 PetscFunctionBegin; 590 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 591 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 592 PetscValidPointer(dim, 2); 593 *dim = prob->totDimBd; 594 PetscFunctionReturn(0); 595 } 596 597 #undef __FUNCT__ 598 #define __FUNCT__ "PetscDSGetTotalComponents" 599 /*@ 600 PetscDSGetTotalComponents - Returns the total number of components in this system 601 602 Not collective 603 604 Input Parameter: 605 . prob - The PetscDS object 606 607 Output Parameter: 608 . dim - The total number of components 609 610 Level: beginner 611 612 .seealso: PetscDSGetNumFields(), PetscDSCreate() 613 @*/ 614 PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc) 615 { 616 PetscErrorCode ierr; 617 618 PetscFunctionBegin; 619 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 620 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 621 PetscValidPointer(Nc, 2); 622 *Nc = prob->totComp; 623 PetscFunctionReturn(0); 624 } 625 626 #undef __FUNCT__ 627 #define __FUNCT__ "PetscDSGetDiscretization" 628 /*@ 629 PetscDSGetDiscretization - Returns the discretization object for the given field 630 631 Not collective 632 633 Input Parameters: 634 + prob - The PetscDS object 635 - f - The field number 636 637 Output Parameter: 638 . disc - The discretization object 639 640 Level: beginner 641 642 .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 643 @*/ 644 PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc) 645 { 646 PetscFunctionBegin; 647 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 648 PetscValidPointer(disc, 3); 649 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); 650 *disc = prob->disc[f]; 651 PetscFunctionReturn(0); 652 } 653 654 #undef __FUNCT__ 655 #define __FUNCT__ "PetscDSGetBdDiscretization" 656 /*@ 657 PetscDSGetBdDiscretization - Returns the boundary discretization object for the given field 658 659 Not collective 660 661 Input Parameters: 662 + prob - The PetscDS object 663 - f - The field number 664 665 Output Parameter: 666 . disc - The boundary discretization object 667 668 Level: beginner 669 670 .seealso: PetscDSSetBdDiscretization(), PetscDSAddBdDiscretization(), PetscDSGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 671 @*/ 672 PetscErrorCode PetscDSGetBdDiscretization(PetscDS prob, PetscInt f, PetscObject *disc) 673 { 674 PetscFunctionBegin; 675 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 676 PetscValidPointer(disc, 3); 677 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); 678 *disc = prob->discBd[f]; 679 PetscFunctionReturn(0); 680 } 681 682 #undef __FUNCT__ 683 #define __FUNCT__ "PetscDSSetDiscretization" 684 /*@ 685 PetscDSSetDiscretization - Sets the discretization object for the given field 686 687 Not collective 688 689 Input Parameters: 690 + prob - The PetscDS object 691 . f - The field number 692 - disc - The discretization object 693 694 Level: beginner 695 696 .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 697 @*/ 698 PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc) 699 { 700 PetscErrorCode ierr; 701 702 PetscFunctionBegin; 703 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 704 PetscValidPointer(disc, 3); 705 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 706 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 707 if (prob->disc[f]) {ierr = PetscObjectDereference(prob->disc[f]);CHKERRQ(ierr);} 708 prob->disc[f] = disc; 709 ierr = PetscObjectReference(disc);CHKERRQ(ierr); 710 { 711 PetscClassId id; 712 713 ierr = PetscObjectGetClassId(disc, &id);CHKERRQ(ierr); 714 if (id == PETSCFV_CLASSID) { 715 prob->implicit[f] = PETSC_FALSE; 716 prob->adjacency[f*2+0] = PETSC_TRUE; 717 prob->adjacency[f*2+1] = PETSC_FALSE; 718 } 719 } 720 PetscFunctionReturn(0); 721 } 722 723 #undef __FUNCT__ 724 #define __FUNCT__ "PetscDSSetBdDiscretization" 725 /*@ 726 PetscDSSetBdDiscretization - Sets the boundary discretization object for the given field 727 728 Not collective 729 730 Input Parameters: 731 + prob - The PetscDS object 732 . f - The field number 733 - disc - The boundary discretization object 734 735 Level: beginner 736 737 .seealso: PetscDSGetBdDiscretization(), PetscDSAddBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 738 @*/ 739 PetscErrorCode PetscDSSetBdDiscretization(PetscDS prob, PetscInt f, PetscObject disc) 740 { 741 PetscErrorCode ierr; 742 743 PetscFunctionBegin; 744 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 745 if (disc) PetscValidPointer(disc, 3); 746 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 747 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 748 if (prob->discBd[f]) {ierr = PetscObjectDereference(prob->discBd[f]);CHKERRQ(ierr);} 749 prob->discBd[f] = disc; 750 ierr = PetscObjectReference(disc);CHKERRQ(ierr); 751 PetscFunctionReturn(0); 752 } 753 754 #undef __FUNCT__ 755 #define __FUNCT__ "PetscDSAddDiscretization" 756 /*@ 757 PetscDSAddDiscretization - Adds a discretization object 758 759 Not collective 760 761 Input Parameters: 762 + prob - The PetscDS object 763 - disc - The boundary discretization object 764 765 Level: beginner 766 767 .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 768 @*/ 769 PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc) 770 { 771 PetscErrorCode ierr; 772 773 PetscFunctionBegin; 774 ierr = PetscDSSetDiscretization(prob, prob->Nf, disc);CHKERRQ(ierr); 775 PetscFunctionReturn(0); 776 } 777 778 #undef __FUNCT__ 779 #define __FUNCT__ "PetscDSAddBdDiscretization" 780 /*@ 781 PetscDSAddBdDiscretization - Adds a boundary discretization object 782 783 Not collective 784 785 Input Parameters: 786 + prob - The PetscDS object 787 - disc - The boundary discretization object 788 789 Level: beginner 790 791 .seealso: PetscDSGetBdDiscretization(), PetscDSSetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 792 @*/ 793 PetscErrorCode PetscDSAddBdDiscretization(PetscDS prob, PetscObject disc) 794 { 795 PetscErrorCode ierr; 796 797 PetscFunctionBegin; 798 ierr = PetscDSSetBdDiscretization(prob, prob->Nf, disc);CHKERRQ(ierr); 799 PetscFunctionReturn(0); 800 } 801 802 #undef __FUNCT__ 803 #define __FUNCT__ "PetscDSGetImplicit" 804 /*@ 805 PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX 806 807 Not collective 808 809 Input Parameters: 810 + prob - The PetscDS object 811 - f - The field number 812 813 Output Parameter: 814 . implicit - The flag indicating what kind of solve to use for this field 815 816 Level: developer 817 818 .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 819 @*/ 820 PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit) 821 { 822 PetscFunctionBegin; 823 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 824 PetscValidPointer(implicit, 3); 825 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); 826 *implicit = prob->implicit[f]; 827 PetscFunctionReturn(0); 828 } 829 830 #undef __FUNCT__ 831 #define __FUNCT__ "PetscDSSetImplicit" 832 /*@ 833 PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for IMEX 834 835 Not collective 836 837 Input Parameters: 838 + prob - The PetscDS object 839 . f - The field number 840 - implicit - The flag indicating what kind of solve to use for this field 841 842 Level: developer 843 844 .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 845 @*/ 846 PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit) 847 { 848 PetscFunctionBegin; 849 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 850 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); 851 prob->implicit[f] = implicit; 852 PetscFunctionReturn(0); 853 } 854 855 #undef __FUNCT__ 856 #define __FUNCT__ "PetscDSGetAdjacency" 857 /*@ 858 PetscDSGetAdjacency - Returns the flags for determining variable influence 859 860 Not collective 861 862 Input Parameters: 863 + prob - The PetscDS object 864 - f - The field number 865 866 Output Parameter: 867 + useCone - Flag for variable influence starting with the cone operation 868 - useClosure - Flag for variable influence using transitive closure 869 870 Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure() 871 872 Level: developer 873 874 .seealso: PetscDSSetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 875 @*/ 876 PetscErrorCode PetscDSGetAdjacency(PetscDS prob, PetscInt f, PetscBool *useCone, PetscBool *useClosure) 877 { 878 PetscFunctionBegin; 879 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 880 PetscValidPointer(useCone, 3); 881 PetscValidPointer(useClosure, 4); 882 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); 883 *useCone = prob->adjacency[f*2+0]; 884 *useClosure = prob->adjacency[f*2+1]; 885 PetscFunctionReturn(0); 886 } 887 888 #undef __FUNCT__ 889 #define __FUNCT__ "PetscDSSetAdjacency" 890 /*@ 891 PetscDSSetAdjacency - Set the flags for determining variable influence 892 893 Not collective 894 895 Input Parameters: 896 + prob - The PetscDS object 897 . f - The field number 898 . useCone - Flag for variable influence starting with the cone operation 899 - useClosure - Flag for variable influence using transitive closure 900 901 Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure() 902 903 Level: developer 904 905 .seealso: PetscDSGetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetBdDiscretization(), PetscDSGetNumFields(), PetscDSCreate() 906 @*/ 907 PetscErrorCode PetscDSSetAdjacency(PetscDS prob, PetscInt f, PetscBool useCone, PetscBool useClosure) 908 { 909 PetscFunctionBegin; 910 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 911 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); 912 prob->adjacency[f*2+0] = useCone; 913 prob->adjacency[f*2+1] = useClosure; 914 PetscFunctionReturn(0); 915 } 916 917 #undef __FUNCT__ 918 #define __FUNCT__ "PetscDSGetObjective" 919 PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f, 920 void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 921 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 922 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 923 PetscReal t, const PetscReal x[], PetscScalar obj[])) 924 { 925 PetscFunctionBegin; 926 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 927 PetscValidPointer(obj, 2); 928 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); 929 *obj = prob->obj[f]; 930 PetscFunctionReturn(0); 931 } 932 933 #undef __FUNCT__ 934 #define __FUNCT__ "PetscDSSetObjective" 935 PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f, 936 void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 937 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 938 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 939 PetscReal t, const PetscReal x[], PetscScalar obj[])) 940 { 941 PetscErrorCode ierr; 942 943 PetscFunctionBegin; 944 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 945 PetscValidFunction(obj, 2); 946 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 947 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 948 prob->obj[f] = obj; 949 PetscFunctionReturn(0); 950 } 951 952 #undef __FUNCT__ 953 #define __FUNCT__ "PetscDSGetResidual" 954 /*@C 955 PetscDSGetResidual - Get the pointwise residual function for a given test field 956 957 Not collective 958 959 Input Parameters: 960 + prob - The PetscDS 961 - f - The test field number 962 963 Output Parameters: 964 + f0 - integrand for the test function term 965 - f1 - integrand for the test function gradient term 966 967 Note: We are using a first order FEM model for the weak form: 968 969 \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) 970 971 The calling sequence for the callbacks f0 and f1 is given by: 972 973 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 974 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 975 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 976 $ PetscReal t, const PetscReal x[], PetscScalar f0[]) 977 978 + dim - the spatial dimension 979 . Nf - the number of fields 980 . uOff - the offset into u[] and u_t[] for each field 981 . uOff_x - the offset into u_x[] for each field 982 . u - each field evaluated at the current point 983 . u_t - the time derivative of each field evaluated at the current point 984 . u_x - the gradient of each field evaluated at the current point 985 . aOff - the offset into a[] and a_t[] for each auxiliary field 986 . aOff_x - the offset into a_x[] for each auxiliary field 987 . a - each auxiliary field evaluated at the current point 988 . a_t - the time derivative of each auxiliary field evaluated at the current point 989 . a_x - the gradient of auxiliary each field evaluated at the current point 990 . t - current time 991 . x - coordinates of the current point 992 - f0 - output values at the current point 993 994 Level: intermediate 995 996 .seealso: PetscDSSetResidual() 997 @*/ 998 PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f, 999 void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1000 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1001 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1002 PetscReal t, const PetscReal x[], PetscScalar f0[]), 1003 void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1004 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1005 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1006 PetscReal t, const PetscReal x[], PetscScalar f1[])) 1007 { 1008 PetscFunctionBegin; 1009 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1010 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); 1011 if (f0) {PetscValidPointer(f0, 3); *f0 = prob->f[f*2+0];} 1012 if (f1) {PetscValidPointer(f1, 4); *f1 = prob->f[f*2+1];} 1013 PetscFunctionReturn(0); 1014 } 1015 1016 #undef __FUNCT__ 1017 #define __FUNCT__ "PetscDSSetResidual" 1018 /*@C 1019 PetscDSSetResidual - Set the pointwise residual function for a given test field 1020 1021 Not collective 1022 1023 Input Parameters: 1024 + prob - The PetscDS 1025 . f - The test field number 1026 . f0 - integrand for the test function term 1027 - f1 - integrand for the test function gradient term 1028 1029 Note: We are using a first order FEM model for the weak form: 1030 1031 \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) 1032 1033 The calling sequence for the callbacks f0 and f1 is given by: 1034 1035 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1036 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1037 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1038 $ PetscReal t, const PetscReal x[], PetscScalar f0[]) 1039 1040 + dim - the spatial dimension 1041 . Nf - the number of fields 1042 . uOff - the offset into u[] and u_t[] for each field 1043 . uOff_x - the offset into u_x[] for each field 1044 . u - each field evaluated at the current point 1045 . u_t - the time derivative of each field evaluated at the current point 1046 . u_x - the gradient of each field evaluated at the current point 1047 . aOff - the offset into a[] and a_t[] for each auxiliary field 1048 . aOff_x - the offset into a_x[] for each auxiliary field 1049 . a - each auxiliary field evaluated at the current point 1050 . a_t - the time derivative of each auxiliary field evaluated at the current point 1051 . a_x - the gradient of auxiliary each field evaluated at the current point 1052 . t - current time 1053 . x - coordinates of the current point 1054 - f0 - output values at the current point 1055 1056 Level: intermediate 1057 1058 .seealso: PetscDSGetResidual() 1059 @*/ 1060 PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f, 1061 void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1062 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1063 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1064 PetscReal t, const PetscReal x[], PetscScalar f0[]), 1065 void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1066 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1067 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1068 PetscReal t, const PetscReal x[], PetscScalar f1[])) 1069 { 1070 PetscErrorCode ierr; 1071 1072 PetscFunctionBegin; 1073 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1074 if (f0) PetscValidFunction(f0, 3); 1075 if (f1) PetscValidFunction(f1, 4); 1076 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1077 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1078 prob->f[f*2+0] = f0; 1079 prob->f[f*2+1] = f1; 1080 PetscFunctionReturn(0); 1081 } 1082 1083 #undef __FUNCT__ 1084 #define __FUNCT__ "PetscDSGetJacobian" 1085 /*@C 1086 PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field 1087 1088 Not collective 1089 1090 Input Parameters: 1091 + prob - The PetscDS 1092 . f - The test field number 1093 - g - The field number 1094 1095 Output Parameters: 1096 + g0 - integrand for the test and basis function term 1097 . g1 - integrand for the test function and basis function gradient term 1098 . g2 - integrand for the test function gradient and basis function term 1099 - g3 - integrand for the test function gradient and basis function gradient term 1100 1101 Note: We are using a first order FEM model for the weak form: 1102 1103 \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 1104 1105 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1106 1107 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1108 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1109 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1110 $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]) 1111 1112 + dim - the spatial dimension 1113 . Nf - the number of fields 1114 . uOff - the offset into u[] and u_t[] for each field 1115 . uOff_x - the offset into u_x[] for each field 1116 . u - each field evaluated at the current point 1117 . u_t - the time derivative of each field evaluated at the current point 1118 . u_x - the gradient of each field evaluated at the current point 1119 . aOff - the offset into a[] and a_t[] for each auxiliary field 1120 . aOff_x - the offset into a_x[] for each auxiliary field 1121 . a - each auxiliary field evaluated at the current point 1122 . a_t - the time derivative of each auxiliary field evaluated at the current point 1123 . a_x - the gradient of auxiliary each field evaluated at the current point 1124 . t - current time 1125 . u_tShift - the multiplier a for dF/dU_t 1126 . x - coordinates of the current point 1127 - g0 - output values at the current point 1128 1129 Level: intermediate 1130 1131 .seealso: PetscDSSetJacobian() 1132 @*/ 1133 PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g, 1134 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1135 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1136 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1137 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]), 1138 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1139 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1140 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1141 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g1[]), 1142 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1143 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1144 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1145 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g2[]), 1146 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1147 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1148 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1149 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g3[])) 1150 { 1151 PetscFunctionBegin; 1152 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1153 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); 1154 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); 1155 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->g[(f*prob->Nf + g)*4+0];} 1156 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->g[(f*prob->Nf + g)*4+1];} 1157 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->g[(f*prob->Nf + g)*4+2];} 1158 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->g[(f*prob->Nf + g)*4+3];} 1159 PetscFunctionReturn(0); 1160 } 1161 1162 #undef __FUNCT__ 1163 #define __FUNCT__ "PetscDSSetJacobian" 1164 /*@C 1165 PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields 1166 1167 Not collective 1168 1169 Input Parameters: 1170 + prob - The PetscDS 1171 . f - The test field number 1172 . g - The field number 1173 . g0 - integrand for the test and basis function term 1174 . g1 - integrand for the test function and basis function gradient term 1175 . g2 - integrand for the test function gradient and basis function term 1176 - g3 - integrand for the test function gradient and basis function gradient term 1177 1178 Note: We are using a first order FEM model for the weak form: 1179 1180 \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 1181 1182 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1183 1184 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1185 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1186 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1187 $ PetscReal t, const PetscReal x[], PetscScalar g0[]) 1188 1189 + dim - the spatial dimension 1190 . Nf - the number of fields 1191 . uOff - the offset into u[] and u_t[] for each field 1192 . uOff_x - the offset into u_x[] for each field 1193 . u - each field evaluated at the current point 1194 . u_t - the time derivative of each field evaluated at the current point 1195 . u_x - the gradient of each field evaluated at the current point 1196 . aOff - the offset into a[] and a_t[] for each auxiliary field 1197 . aOff_x - the offset into a_x[] for each auxiliary field 1198 . a - each auxiliary field evaluated at the current point 1199 . a_t - the time derivative of each auxiliary field evaluated at the current point 1200 . a_x - the gradient of auxiliary each field evaluated at the current point 1201 . t - current time 1202 . u_tShift - the multiplier a for dF/dU_t 1203 . x - coordinates of the current point 1204 - g0 - output values at the current point 1205 1206 Level: intermediate 1207 1208 .seealso: PetscDSGetJacobian() 1209 @*/ 1210 PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g, 1211 void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1212 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1213 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1214 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g0[]), 1215 void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1216 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1217 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1218 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g1[]), 1219 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1220 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1221 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1222 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g2[]), 1223 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1224 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1225 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1226 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscScalar g3[])) 1227 { 1228 PetscErrorCode ierr; 1229 1230 PetscFunctionBegin; 1231 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1232 if (g0) PetscValidFunction(g0, 4); 1233 if (g1) PetscValidFunction(g1, 5); 1234 if (g2) PetscValidFunction(g2, 6); 1235 if (g3) PetscValidFunction(g3, 7); 1236 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1237 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1238 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1239 prob->g[(f*prob->Nf + g)*4+0] = g0; 1240 prob->g[(f*prob->Nf + g)*4+1] = g1; 1241 prob->g[(f*prob->Nf + g)*4+2] = g2; 1242 prob->g[(f*prob->Nf + g)*4+3] = g3; 1243 PetscFunctionReturn(0); 1244 } 1245 1246 #undef __FUNCT__ 1247 #define __FUNCT__ "PetscDSGetRiemannSolver" 1248 /*@C 1249 PetscDSGetRiemannSolver - Returns the Riemann solver for the given field 1250 1251 Not collective 1252 1253 Input Arguments: 1254 + prob - The PetscDS object 1255 - f - The field number 1256 1257 Output Argument: 1258 . r - Riemann solver 1259 1260 Calling sequence for r: 1261 1262 $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx) 1263 1264 + dim - The spatial dimension 1265 . Nf - The number of fields 1266 . x - The coordinates at a point on the interface 1267 . n - The normal vector to the interface 1268 . uL - The state vector to the left of the interface 1269 . uR - The state vector to the right of the interface 1270 . flux - output array of flux through the interface 1271 - ctx - optional user context 1272 1273 Level: intermediate 1274 1275 .seealso: PetscDSSetRiemannSolver() 1276 @*/ 1277 PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f, 1278 void (**r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)) 1279 { 1280 PetscFunctionBegin; 1281 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1282 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); 1283 PetscValidPointer(r, 3); 1284 *r = prob->r[f]; 1285 PetscFunctionReturn(0); 1286 } 1287 1288 #undef __FUNCT__ 1289 #define __FUNCT__ "PetscDSSetRiemannSolver" 1290 /*@C 1291 PetscDSSetRiemannSolver - Sets the Riemann solver for the given field 1292 1293 Not collective 1294 1295 Input Arguments: 1296 + prob - The PetscDS object 1297 . f - The field number 1298 - r - Riemann solver 1299 1300 Calling sequence for r: 1301 1302 $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx) 1303 1304 + dim - The spatial dimension 1305 . Nf - The number of fields 1306 . x - The coordinates at a point on the interface 1307 . n - The normal vector to the interface 1308 . uL - The state vector to the left of the interface 1309 . uR - The state vector to the right of the interface 1310 . flux - output array of flux through the interface 1311 - ctx - optional user context 1312 1313 Level: intermediate 1314 1315 .seealso: PetscDSGetRiemannSolver() 1316 @*/ 1317 PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f, 1318 void (*r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)) 1319 { 1320 PetscErrorCode ierr; 1321 1322 PetscFunctionBegin; 1323 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1324 PetscValidFunction(r, 3); 1325 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1326 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1327 prob->r[f] = r; 1328 PetscFunctionReturn(0); 1329 } 1330 1331 #undef __FUNCT__ 1332 #define __FUNCT__ "PetscDSGetContext" 1333 PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx) 1334 { 1335 PetscFunctionBegin; 1336 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1337 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); 1338 PetscValidPointer(ctx, 3); 1339 *ctx = prob->ctx[f]; 1340 PetscFunctionReturn(0); 1341 } 1342 1343 #undef __FUNCT__ 1344 #define __FUNCT__ "PetscDSSetContext" 1345 PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx) 1346 { 1347 PetscErrorCode ierr; 1348 1349 PetscFunctionBegin; 1350 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1351 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1352 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1353 prob->ctx[f] = ctx; 1354 PetscFunctionReturn(0); 1355 } 1356 1357 #undef __FUNCT__ 1358 #define __FUNCT__ "PetscDSGetBdResidual" 1359 /*@C 1360 PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field 1361 1362 Not collective 1363 1364 Input Parameters: 1365 + prob - The PetscDS 1366 - f - The test field number 1367 1368 Output Parameters: 1369 + f0 - boundary integrand for the test function term 1370 - f1 - boundary integrand for the test function gradient term 1371 1372 Note: We are using a first order FEM model for the weak form: 1373 1374 \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 1375 1376 The calling sequence for the callbacks f0 and f1 is given by: 1377 1378 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1379 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1380 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1381 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]) 1382 1383 + dim - the spatial dimension 1384 . Nf - the number of fields 1385 . uOff - the offset into u[] and u_t[] for each field 1386 . uOff_x - the offset into u_x[] for each field 1387 . u - each field evaluated at the current point 1388 . u_t - the time derivative of each field evaluated at the current point 1389 . u_x - the gradient of each field evaluated at the current point 1390 . aOff - the offset into a[] and a_t[] for each auxiliary field 1391 . aOff_x - the offset into a_x[] for each auxiliary field 1392 . a - each auxiliary field evaluated at the current point 1393 . a_t - the time derivative of each auxiliary field evaluated at the current point 1394 . a_x - the gradient of auxiliary each field evaluated at the current point 1395 . t - current time 1396 . x - coordinates of the current point 1397 . n - unit normal at the current point 1398 - f0 - output values at the current point 1399 1400 Level: intermediate 1401 1402 .seealso: PetscDSSetBdResidual() 1403 @*/ 1404 PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f, 1405 void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1406 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1407 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1408 PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]), 1409 void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1410 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1411 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1412 PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f1[])) 1413 { 1414 PetscFunctionBegin; 1415 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1416 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); 1417 if (f0) {PetscValidPointer(f0, 3); *f0 = prob->fBd[f*2+0];} 1418 if (f1) {PetscValidPointer(f1, 4); *f1 = prob->fBd[f*2+1];} 1419 PetscFunctionReturn(0); 1420 } 1421 1422 #undef __FUNCT__ 1423 #define __FUNCT__ "PetscDSSetBdResidual" 1424 /*@C 1425 PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field 1426 1427 Not collective 1428 1429 Input Parameters: 1430 + prob - The PetscDS 1431 . f - The test field number 1432 . f0 - boundary integrand for the test function term 1433 - f1 - boundary integrand for the test function gradient term 1434 1435 Note: We are using a first order FEM model for the weak form: 1436 1437 \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 1438 1439 The calling sequence for the callbacks f0 and f1 is given by: 1440 1441 $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1442 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1443 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1444 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]) 1445 1446 + dim - the spatial dimension 1447 . Nf - the number of fields 1448 . uOff - the offset into u[] and u_t[] for each field 1449 . uOff_x - the offset into u_x[] for each field 1450 . u - each field evaluated at the current point 1451 . u_t - the time derivative of each field evaluated at the current point 1452 . u_x - the gradient of each field evaluated at the current point 1453 . aOff - the offset into a[] and a_t[] for each auxiliary field 1454 . aOff_x - the offset into a_x[] for each auxiliary field 1455 . a - each auxiliary field evaluated at the current point 1456 . a_t - the time derivative of each auxiliary field evaluated at the current point 1457 . a_x - the gradient of auxiliary each field evaluated at the current point 1458 . t - current time 1459 . x - coordinates of the current point 1460 . n - unit normal at the current point 1461 - f0 - output values at the current point 1462 1463 Level: intermediate 1464 1465 .seealso: PetscDSGetBdResidual() 1466 @*/ 1467 PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f, 1468 void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1469 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1470 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1471 PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[]), 1472 void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1473 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1474 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1475 PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f1[])) 1476 { 1477 PetscErrorCode ierr; 1478 1479 PetscFunctionBegin; 1480 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1481 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1482 ierr = PetscDSEnlarge_Static(prob, f+1);CHKERRQ(ierr); 1483 if (f0) {PetscValidFunction(f0, 3); prob->fBd[f*2+0] = f0;} 1484 if (f1) {PetscValidFunction(f1, 4); prob->fBd[f*2+1] = f1;} 1485 PetscFunctionReturn(0); 1486 } 1487 1488 #undef __FUNCT__ 1489 #define __FUNCT__ "PetscDSGetBdJacobian" 1490 /*@C 1491 PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field 1492 1493 Not collective 1494 1495 Input Parameters: 1496 + prob - The PetscDS 1497 . f - The test field number 1498 - g - The field number 1499 1500 Output Parameters: 1501 + g0 - integrand for the test and basis function term 1502 . g1 - integrand for the test function and basis function gradient term 1503 . g2 - integrand for the test function gradient and basis function term 1504 - g3 - integrand for the test function gradient and basis function gradient term 1505 1506 Note: We are using a first order FEM model for the weak form: 1507 1508 \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 1509 1510 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1511 1512 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1513 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1514 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1515 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[]) 1516 1517 + dim - the spatial dimension 1518 . Nf - the number of fields 1519 . uOff - the offset into u[] and u_t[] for each field 1520 . uOff_x - the offset into u_x[] for each field 1521 . u - each field evaluated at the current point 1522 . u_t - the time derivative of each field evaluated at the current point 1523 . u_x - the gradient of each field evaluated at the current point 1524 . aOff - the offset into a[] and a_t[] for each auxiliary field 1525 . aOff_x - the offset into a_x[] for each auxiliary field 1526 . a - each auxiliary field evaluated at the current point 1527 . a_t - the time derivative of each auxiliary field evaluated at the current point 1528 . a_x - the gradient of auxiliary each field evaluated at the current point 1529 . t - current time 1530 . u_tShift - the multiplier a for dF/dU_t 1531 . x - coordinates of the current point 1532 . n - normal at the current point 1533 - g0 - output values at the current point 1534 1535 Level: intermediate 1536 1537 .seealso: PetscDSSetBdJacobian() 1538 @*/ 1539 PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g, 1540 void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1541 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1542 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1543 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g0[]), 1544 void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1545 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1546 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1547 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g1[]), 1548 void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1549 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1550 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1551 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g2[]), 1552 void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1553 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1554 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1555 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g3[])) 1556 { 1557 PetscFunctionBegin; 1558 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1559 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); 1560 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); 1561 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->gBd[(f*prob->Nf + g)*4+0];} 1562 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->gBd[(f*prob->Nf + g)*4+1];} 1563 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->gBd[(f*prob->Nf + g)*4+2];} 1564 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->gBd[(f*prob->Nf + g)*4+3];} 1565 PetscFunctionReturn(0); 1566 } 1567 1568 #undef __FUNCT__ 1569 #define __FUNCT__ "PetscDSSetBdJacobian" 1570 /*@C 1571 PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field 1572 1573 Not collective 1574 1575 Input Parameters: 1576 + prob - The PetscDS 1577 . f - The test field number 1578 . g - The field number 1579 . g0 - integrand for the test and basis function term 1580 . g1 - integrand for the test function and basis function gradient term 1581 . g2 - integrand for the test function gradient and basis function term 1582 - g3 - integrand for the test function gradient and basis function gradient term 1583 1584 Note: We are using a first order FEM model for the weak form: 1585 1586 \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 1587 1588 The calling sequence for the callbacks g0, g1, g2 and g3 is given by: 1589 1590 $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1591 $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1592 $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1593 $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[]) 1594 1595 + dim - the spatial dimension 1596 . Nf - the number of fields 1597 . uOff - the offset into u[] and u_t[] for each field 1598 . uOff_x - the offset into u_x[] for each field 1599 . u - each field evaluated at the current point 1600 . u_t - the time derivative of each field evaluated at the current point 1601 . u_x - the gradient of each field evaluated at the current point 1602 . aOff - the offset into a[] and a_t[] for each auxiliary field 1603 . aOff_x - the offset into a_x[] for each auxiliary field 1604 . a - each auxiliary field evaluated at the current point 1605 . a_t - the time derivative of each auxiliary field evaluated at the current point 1606 . a_x - the gradient of auxiliary each field evaluated at the current point 1607 . t - current time 1608 . u_tShift - the multiplier a for dF/dU_t 1609 . x - coordinates of the current point 1610 . n - normal at the current point 1611 - g0 - output values at the current point 1612 1613 Level: intermediate 1614 1615 .seealso: PetscDSGetBdJacobian() 1616 @*/ 1617 PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g, 1618 void (*g0)(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[], const PetscReal n[], PetscScalar g0[]), 1622 void (*g1)(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[], const PetscReal n[], PetscScalar g1[]), 1626 void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1627 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1628 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1629 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g2[]), 1630 void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, 1631 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 1632 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 1633 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscScalar g3[])) 1634 { 1635 PetscErrorCode ierr; 1636 1637 PetscFunctionBegin; 1638 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1639 if (g0) PetscValidFunction(g0, 4); 1640 if (g1) PetscValidFunction(g1, 5); 1641 if (g2) PetscValidFunction(g2, 6); 1642 if (g3) PetscValidFunction(g3, 7); 1643 if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f); 1644 if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g); 1645 ierr = PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);CHKERRQ(ierr); 1646 prob->gBd[(f*prob->Nf + g)*4+0] = g0; 1647 prob->gBd[(f*prob->Nf + g)*4+1] = g1; 1648 prob->gBd[(f*prob->Nf + g)*4+2] = g2; 1649 prob->gBd[(f*prob->Nf + g)*4+3] = g3; 1650 PetscFunctionReturn(0); 1651 } 1652 1653 #undef __FUNCT__ 1654 #define __FUNCT__ "PetscDSGetFieldOffset" 1655 /*@ 1656 PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis 1657 1658 Not collective 1659 1660 Input Parameters: 1661 + prob - The PetscDS object 1662 - f - The field number 1663 1664 Output Parameter: 1665 . off - The offset 1666 1667 Level: beginner 1668 1669 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1670 @*/ 1671 PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off) 1672 { 1673 PetscInt g; 1674 PetscErrorCode ierr; 1675 1676 PetscFunctionBegin; 1677 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1678 PetscValidPointer(off, 3); 1679 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); 1680 *off = 0; 1681 for (g = 0; g < f; ++g) { 1682 PetscFE fe = (PetscFE) prob->disc[g]; 1683 PetscInt Nb, Nc; 1684 1685 ierr = PetscFEGetDimension(fe, &Nb);CHKERRQ(ierr); 1686 ierr = PetscFEGetNumComponents(fe, &Nc);CHKERRQ(ierr); 1687 *off += Nb*Nc; 1688 } 1689 PetscFunctionReturn(0); 1690 } 1691 1692 #undef __FUNCT__ 1693 #define __FUNCT__ "PetscDSGetBdFieldOffset" 1694 /*@ 1695 PetscDSGetBdFieldOffset - Returns the offset of the given field in the full space boundary basis 1696 1697 Not collective 1698 1699 Input Parameters: 1700 + prob - The PetscDS object 1701 - f - The field number 1702 1703 Output Parameter: 1704 . off - The boundary offset 1705 1706 Level: beginner 1707 1708 .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1709 @*/ 1710 PetscErrorCode PetscDSGetBdFieldOffset(PetscDS prob, PetscInt f, PetscInt *off) 1711 { 1712 PetscInt g; 1713 PetscErrorCode ierr; 1714 1715 PetscFunctionBegin; 1716 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1717 PetscValidPointer(off, 3); 1718 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); 1719 *off = 0; 1720 for (g = 0; g < f; ++g) { 1721 PetscFE fe = (PetscFE) prob->discBd[g]; 1722 PetscInt Nb, Nc; 1723 1724 ierr = PetscFEGetDimension(fe, &Nb);CHKERRQ(ierr); 1725 ierr = PetscFEGetNumComponents(fe, &Nc);CHKERRQ(ierr); 1726 *off += Nb*Nc; 1727 } 1728 PetscFunctionReturn(0); 1729 } 1730 1731 #undef __FUNCT__ 1732 #define __FUNCT__ "PetscDSGetComponentOffset" 1733 /*@ 1734 PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point 1735 1736 Not collective 1737 1738 Input Parameters: 1739 + prob - The PetscDS object 1740 - f - The field number 1741 1742 Output Parameter: 1743 . off - The offset 1744 1745 Level: beginner 1746 1747 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1748 @*/ 1749 PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off) 1750 { 1751 PetscInt g; 1752 PetscErrorCode ierr; 1753 1754 PetscFunctionBegin; 1755 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1756 PetscValidPointer(off, 3); 1757 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); 1758 *off = 0; 1759 for (g = 0; g < f; ++g) { 1760 PetscFE fe = (PetscFE) prob->disc[g]; 1761 PetscInt Nc; 1762 1763 ierr = PetscFEGetNumComponents(fe, &Nc);CHKERRQ(ierr); 1764 *off += Nc; 1765 } 1766 PetscFunctionReturn(0); 1767 } 1768 1769 #undef __FUNCT__ 1770 #define __FUNCT__ "PetscDSGetComponentOffsets" 1771 /*@ 1772 PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point 1773 1774 Not collective 1775 1776 Input Parameter: 1777 . prob - The PetscDS object 1778 1779 Output Parameter: 1780 . offsets - The offsets 1781 1782 Level: beginner 1783 1784 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1785 @*/ 1786 PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[]) 1787 { 1788 PetscFunctionBegin; 1789 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1790 PetscValidPointer(offsets, 2); 1791 *offsets = prob->off; 1792 PetscFunctionReturn(0); 1793 } 1794 1795 #undef __FUNCT__ 1796 #define __FUNCT__ "PetscDSGetComponentDerivativeOffsets" 1797 /*@ 1798 PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point 1799 1800 Not collective 1801 1802 Input Parameter: 1803 . prob - The PetscDS object 1804 1805 Output Parameter: 1806 . offsets - The offsets 1807 1808 Level: beginner 1809 1810 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1811 @*/ 1812 PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[]) 1813 { 1814 PetscFunctionBegin; 1815 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1816 PetscValidPointer(offsets, 2); 1817 *offsets = prob->offDer; 1818 PetscFunctionReturn(0); 1819 } 1820 1821 #undef __FUNCT__ 1822 #define __FUNCT__ "PetscDSGetComponentBdOffsets" 1823 /*@ 1824 PetscDSGetComponentBdOffsets - Returns the offset of each field on a boundary evaluation point 1825 1826 Not collective 1827 1828 Input Parameter: 1829 . prob - The PetscDS object 1830 1831 Output Parameter: 1832 . offsets - The offsets 1833 1834 Level: beginner 1835 1836 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1837 @*/ 1838 PetscErrorCode PetscDSGetComponentBdOffsets(PetscDS prob, PetscInt *offsets[]) 1839 { 1840 PetscFunctionBegin; 1841 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1842 PetscValidPointer(offsets, 2); 1843 *offsets = prob->offBd; 1844 PetscFunctionReturn(0); 1845 } 1846 1847 #undef __FUNCT__ 1848 #define __FUNCT__ "PetscDSGetComponentBdDerivativeOffsets" 1849 /*@ 1850 PetscDSGetComponentBdDerivativeOffsets - Returns the offset of each field derivative on a boundary evaluation point 1851 1852 Not collective 1853 1854 Input Parameter: 1855 . prob - The PetscDS object 1856 1857 Output Parameter: 1858 . offsets - The offsets 1859 1860 Level: beginner 1861 1862 .seealso: PetscDSGetBdFieldOffset(), PetscDSGetNumFields(), PetscDSCreate() 1863 @*/ 1864 PetscErrorCode PetscDSGetComponentBdDerivativeOffsets(PetscDS prob, PetscInt *offsets[]) 1865 { 1866 PetscFunctionBegin; 1867 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1868 PetscValidPointer(offsets, 2); 1869 *offsets = prob->offDerBd; 1870 PetscFunctionReturn(0); 1871 } 1872 1873 #undef __FUNCT__ 1874 #define __FUNCT__ "PetscDSGetTabulation" 1875 /*@C 1876 PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization 1877 1878 Not collective 1879 1880 Input Parameter: 1881 . prob - The PetscDS object 1882 1883 Output Parameters: 1884 + basis - The basis function tabulation at quadrature points 1885 - basisDer - The basis function derivative tabulation at quadrature points 1886 1887 Level: intermediate 1888 1889 .seealso: PetscDSGetBdTabulation(), PetscDSCreate() 1890 @*/ 1891 PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer) 1892 { 1893 PetscErrorCode ierr; 1894 1895 PetscFunctionBegin; 1896 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1897 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 1898 if (basis) {PetscValidPointer(basis, 2); *basis = prob->basis;} 1899 if (basisDer) {PetscValidPointer(basisDer, 3); *basisDer = prob->basisDer;} 1900 PetscFunctionReturn(0); 1901 } 1902 1903 #undef __FUNCT__ 1904 #define __FUNCT__ "PetscDSGetBdTabulation" 1905 /*@C 1906 PetscDSGetBdTabulation - Return the basis tabulation at quadrature points for the boundary discretization 1907 1908 Not collective 1909 1910 Input Parameter: 1911 . prob - The PetscDS object 1912 1913 Output Parameters: 1914 + basis - The basis function tabulation at quadrature points 1915 - basisDer - The basis function derivative tabulation at quadrature points 1916 1917 Level: intermediate 1918 1919 .seealso: PetscDSGetTabulation(), PetscDSCreate() 1920 @*/ 1921 PetscErrorCode PetscDSGetBdTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer) 1922 { 1923 PetscErrorCode ierr; 1924 1925 PetscFunctionBegin; 1926 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1927 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 1928 if (basis) {PetscValidPointer(basis, 2); *basis = prob->basisBd;} 1929 if (basisDer) {PetscValidPointer(basisDer, 3); *basisDer = prob->basisDerBd;} 1930 PetscFunctionReturn(0); 1931 } 1932 1933 #undef __FUNCT__ 1934 #define __FUNCT__ "PetscDSGetEvaluationArrays" 1935 PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x) 1936 { 1937 PetscErrorCode ierr; 1938 1939 PetscFunctionBegin; 1940 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1941 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 1942 if (u) {PetscValidPointer(u, 2); *u = prob->u;} 1943 if (u_t) {PetscValidPointer(u_t, 3); *u_t = prob->u_t;} 1944 if (u_x) {PetscValidPointer(u_x, 4); *u_x = prob->u_x;} 1945 PetscFunctionReturn(0); 1946 } 1947 1948 #undef __FUNCT__ 1949 #define __FUNCT__ "PetscDSGetWeakFormArrays" 1950 PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3) 1951 { 1952 PetscErrorCode ierr; 1953 1954 PetscFunctionBegin; 1955 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1956 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 1957 if (f0) {PetscValidPointer(f0, 2); *f0 = prob->f0;} 1958 if (f1) {PetscValidPointer(f1, 3); *f1 = prob->f1;} 1959 if (g0) {PetscValidPointer(g0, 4); *g0 = prob->g0;} 1960 if (g1) {PetscValidPointer(g1, 5); *g1 = prob->g1;} 1961 if (g2) {PetscValidPointer(g2, 6); *g2 = prob->g2;} 1962 if (g3) {PetscValidPointer(g3, 7); *g3 = prob->g3;} 1963 PetscFunctionReturn(0); 1964 } 1965 1966 #undef __FUNCT__ 1967 #define __FUNCT__ "PetscDSGetRefCoordArrays" 1968 PetscErrorCode PetscDSGetRefCoordArrays(PetscDS prob, PetscReal **x, PetscScalar **refSpaceDer) 1969 { 1970 PetscErrorCode ierr; 1971 1972 PetscFunctionBegin; 1973 PetscValidHeaderSpecific(prob, PETSCDS_CLASSID, 1); 1974 ierr = PetscDSSetUp(prob);CHKERRQ(ierr); 1975 if (x) {PetscValidPointer(x, 2); *x = prob->x;} 1976 if (refSpaceDer) {PetscValidPointer(refSpaceDer, 3); *refSpaceDer = prob->refSpaceDer;} 1977 PetscFunctionReturn(0); 1978 } 1979 1980 #undef __FUNCT__ 1981 #define __FUNCT__ "PetscDSDestroy_Basic" 1982 static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob) 1983 { 1984 PetscFunctionBegin; 1985 PetscFunctionReturn(0); 1986 } 1987 1988 #undef __FUNCT__ 1989 #define __FUNCT__ "PetscDSInitialize_Basic" 1990 static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob) 1991 { 1992 PetscFunctionBegin; 1993 prob->ops->setfromoptions = NULL; 1994 prob->ops->setup = NULL; 1995 prob->ops->view = NULL; 1996 prob->ops->destroy = PetscDSDestroy_Basic; 1997 PetscFunctionReturn(0); 1998 } 1999 2000 /*MC 2001 PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions 2002 2003 Level: intermediate 2004 2005 .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType() 2006 M*/ 2007 2008 #undef __FUNCT__ 2009 #define __FUNCT__ "PetscDSCreate_Basic" 2010 PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob) 2011 { 2012 PetscDS_Basic *b; 2013 PetscErrorCode ierr; 2014 2015 PetscFunctionBegin; 2016 PetscValidHeaderSpecific(prob, PETSCSPACE_CLASSID, 1); 2017 ierr = PetscNewLog(prob, &b);CHKERRQ(ierr); 2018 prob->data = b; 2019 2020 ierr = PetscDSInitialize_Basic(prob);CHKERRQ(ierr); 2021 PetscFunctionReturn(0); 2022 } 2023