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