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