Actual source code: inherit.c

  1: /*
  2:      Provides utility routines for manipulating any type of PETSc object.
  3: */
  4: #include <petsc/private/petscimpl.h>
  5: #include <petscviewer.h>

  7: PETSC_INTERN PetscObject *PetscObjects;
  8: PETSC_INTERN PetscInt     PetscObjectsCounts;
  9: PETSC_INTERN PetscInt     PetscObjectsMaxCounts;
 10: PETSC_INTERN PetscBool    PetscObjectsLog;

 12: PetscObject *PetscObjects       = NULL;
 13: PetscInt     PetscObjectsCounts = 0, PetscObjectsMaxCounts = 0;
 14: PetscBool    PetscObjectsLog = PETSC_FALSE;

 16: PetscObjectId PetscObjectNewId_Internal(void)
 17: {
 18:   static PetscObjectId idcnt = 1;
 19:   return idcnt++;
 20: }

 22: PetscErrorCode PetscHeaderCreate_Function(PetscErrorCode ierr, PetscObject *h, PetscClassId classid, const char class_name[], const char descr[], const char mansec[], MPI_Comm comm, PetscObjectDestroyFn *destroy, PetscObjectViewFn *view)
 23: {
 24:   PetscFunctionBegin;
 25:   if (ierr) PetscFunctionReturn(ierr);
 26:   PetscCall(PetscHeaderCreate_Private(*h, classid, class_name, descr, mansec, comm, destroy, view));
 27:   PetscCall(PetscLogObjectCreate(*h));
 28:   PetscFunctionReturn(PETSC_SUCCESS);
 29: }

 31: /*
 32:    PetscHeaderCreate_Private - Fills in the default values.
 33: */
 34: PetscErrorCode PetscHeaderCreate_Private(PetscObject h, PetscClassId classid, const char class_name[], const char descr[], const char mansec[], MPI_Comm comm, PetscObjectDestroyFn *destroy, PetscObjectViewFn *view)
 35: {
 36:   void       *get_tmp;
 37:   PetscInt64 *cidx;
 38:   PetscMPIInt flg;

 40:   PetscFunctionBegin;
 41:   h->classid               = classid;
 42:   h->class_name            = (char *)class_name;
 43:   h->description           = (char *)descr;
 44:   h->mansec                = (char *)mansec;
 45:   h->refct                 = 1;
 46:   h->non_cyclic_references = NULL;
 47:   h->id                    = PetscObjectNewId_Internal();
 48:   h->bops->destroy         = destroy;
 49:   h->bops->view            = view;

 51:   PetscCall(PetscCommDuplicate(comm, &h->comm, &h->tag));

 53:   /* Increment and store current object creation index */
 54:   PetscCallMPI(MPI_Comm_get_attr(h->comm, Petsc_CreationIdx_keyval, &get_tmp, &flg));
 55:   PetscCheck(flg, h->comm, PETSC_ERR_ARG_CORRUPT, "MPI_Comm does not have an object creation index");
 56:   cidx    = (PetscInt64 *)get_tmp;
 57:   h->cidx = (*cidx)++;

 59:   /* Keep a record of object created */
 60:   if (PetscDefined(USE_LOG) && PetscObjectsLog) {
 61:     PetscObject *newPetscObjects;
 62:     PetscInt     newPetscObjectsMaxCounts;

 64:     PetscObjectsCounts++;
 65:     for (PetscInt i = 0; i < PetscObjectsMaxCounts; ++i) {
 66:       if (!PetscObjects[i]) {
 67:         PetscObjects[i] = h;
 68:         PetscFunctionReturn(PETSC_SUCCESS);
 69:       }
 70:     }
 71:     /* Need to increase the space for storing PETSc objects */
 72:     if (!PetscObjectsMaxCounts) newPetscObjectsMaxCounts = 100;
 73:     else newPetscObjectsMaxCounts = 2 * PetscObjectsMaxCounts;
 74:     PetscCall(PetscCalloc1(newPetscObjectsMaxCounts, &newPetscObjects));
 75:     PetscCall(PetscArraycpy(newPetscObjects, PetscObjects, PetscObjectsMaxCounts));
 76:     PetscCall(PetscFree(PetscObjects));

 78:     PetscObjects                        = newPetscObjects;
 79:     PetscObjects[PetscObjectsMaxCounts] = h;
 80:     PetscObjectsMaxCounts               = newPetscObjectsMaxCounts;
 81:   }
 82:   PetscFunctionReturn(PETSC_SUCCESS);
 83: }

 85: PETSC_INTERN PetscBool      PetscMemoryCollectMaximumUsage;
 86: PETSC_INTERN PetscLogDouble PetscMemoryMaximumUsage;

 88: PetscErrorCode PetscHeaderDestroy_Function(PetscObject *h)
 89: {
 90:   PetscFunctionBegin;
 91:   PetscCall(PetscLogObjectDestroy(*h));
 92:   PetscCall(PetscHeaderDestroy_Private(*h, PETSC_FALSE));
 93:   PetscCall(PetscFree(*h));
 94:   PetscFunctionReturn(PETSC_SUCCESS);
 95: }

 97: /*
 98:     PetscHeaderDestroy_Private - Destroys a base PETSc object header. Called by
 99:     the macro PetscHeaderDestroy().
100: */
101: PetscErrorCode PetscHeaderDestroy_Private(PetscObject obj, PetscBool clear_for_reuse)
102: {
103:   PetscFunctionBegin;
105:   PetscCheck(!obj->persistent, PetscObjectComm((PetscObject)obj), PETSC_ERR_ARG_WRONGSTATE, "Cannot destroy this object, it is destroyed automatically in PetscFinalize()");
106:   PetscCall(PetscComposedQuantitiesDestroy(obj));
107:   if (PetscMemoryCollectMaximumUsage) {
108:     PetscLogDouble usage;

110:     PetscCall(PetscMemoryGetCurrentUsage(&usage));
111:     if (usage > PetscMemoryMaximumUsage) PetscMemoryMaximumUsage = usage;
112:   }
113:   /* first destroy things that could execute arbitrary code */
114:   if (obj->python_destroy) {
115:     void *python_context                     = obj->python_context;
116:     PetscErrorCode (*python_destroy)(void *) = obj->python_destroy;

118:     obj->python_context = NULL;
119:     obj->python_destroy = NULL;
120:     PetscCall((*python_destroy)(python_context));
121:   }
122:   PetscCall(PetscObjectDestroyOptionsHandlers(obj));
123:   PetscCall(PetscObjectListDestroy(&obj->olist));

125:   /* destroy allocated quantities */
126:   if (PetscPrintFunctionList) PetscCall(PetscFunctionListPrintNonEmpty(obj->qlist));
127:   PetscCheck(--(obj->refct) <= 0, obj->comm, PETSC_ERR_PLIB, "Destroying a PetscObject (%s) with reference count %" PetscInt_FMT " >= 1", obj->name ? obj->name : "unnamed", obj->refct);
128:   PetscCall(PetscFree(obj->name));
129:   PetscCall(PetscFree(obj->prefix));
130:   PetscCall(PetscFree(obj->type_name));

132:   if (clear_for_reuse) {
133:     /* we will assume that obj->bops->view and destroy are safe to leave as-is */

135:     /* reset quantities, in order of appearance in _p_PetscObject */
136:     obj->id       = PetscObjectNewId_Internal();
137:     obj->refct    = 1;
138:     obj->tablevel = 0;
139:     obj->state    = 0;
140:     /* don't deallocate, zero these out instead */
141:     PetscCall(PetscFunctionListClear(obj->qlist));
142:     PetscCall(PetscArrayzero(obj->fortran_func_pointers, obj->num_fortran_func_pointers));
143:     PetscCall(PetscArrayzero(obj->fortrancallback[PETSC_FORTRAN_CALLBACK_CLASS], obj->num_fortrancallback[PETSC_FORTRAN_CALLBACK_CLASS]));
144:     PetscCall(PetscArrayzero(obj->fortrancallback[PETSC_FORTRAN_CALLBACK_SUBTYPE], obj->num_fortrancallback[PETSC_FORTRAN_CALLBACK_SUBTYPE]));
145:     obj->optionsprinted = PETSC_FALSE;
146: #if PetscDefined(HAVE_SAWS)
147:     obj->amsmem          = PETSC_FALSE;
148:     obj->amspublishblock = PETSC_FALSE;
149: #endif
150:     obj->options                                  = NULL;
151:     obj->donotPetscObjectPrintClassNamePrefixType = PETSC_FALSE;
152:   } else {
153:     PetscCall(PetscFunctionListDestroy(&obj->qlist));
154:     PetscCall(PetscFree(obj->fortran_func_pointers));
155:     PetscCall(PetscFree(obj->fortrancallback[PETSC_FORTRAN_CALLBACK_CLASS]));
156:     PetscCall(PetscFree(obj->fortrancallback[PETSC_FORTRAN_CALLBACK_SUBTYPE]));
157:     PetscCall(PetscCommDestroy(&obj->comm));
158:     obj->classid = PETSCFREEDHEADER;

160:     if (PetscDefined(USE_LOG) && PetscObjectsLog) {
161:       /* Record object removal from list of all objects */
162:       for (PetscInt i = 0; i < PetscObjectsMaxCounts; ++i) {
163:         if (PetscObjects[i] == obj) {
164:           PetscObjects[i] = NULL;
165:           --PetscObjectsCounts;
166:           break;
167:         }
168:       }
169:       if (!PetscObjectsCounts) {
170:         PetscCall(PetscFree(PetscObjects));
171:         PetscObjectsMaxCounts = 0;
172:       }
173:     }
174:   }
175:   PetscFunctionReturn(PETSC_SUCCESS);
176: }

178: /*
179:   PetscHeaderReset_Internal - "Reset" a PetscObject header. This is tantamount to destroying
180:   the object but does not free all resources. The object retains its:

182:   - classid
183:   - bops->view
184:   - bops->destroy
185:   - comm
186:   - tag
187:   - class_name
188:   - description
189:   - mansec
190:   - cpp

192:   Note that while subclass information is lost, superclass info remains. Thus this function is
193:   intended to be used to reuse a PetscObject within the same class to avoid reallocating its
194:   resources.
195: */
196: PetscErrorCode PetscHeaderReset_Internal(PetscObject obj)
197: {
198:   PetscFunctionBegin;
199:   PetscCall(PetscHeaderDestroy_Private(obj, PETSC_TRUE));
200:   PetscFunctionReturn(PETSC_SUCCESS);
201: }

203: /*@C
204:   PetscObjectCopyFortranFunctionPointers - Copy function pointers to another object

206:   Logically Collective

208:   Input Parameters:
209: + src  - source object
210: - dest - destination object

212:   Level: developer

214:   Note:
215:   Both objects must have the same class.

217:   This is used to help manage user callback functions that were provided in Fortran

219: .seealso: `PetscFortranCallbackRegister()`, `PetscFortranCallbackGetSizes()`
220: @*/
221: PetscErrorCode PetscObjectCopyFortranFunctionPointers(PetscObject src, PetscObject dest)
222: {
223:   PetscFortranCallbackId cbtype, numcb[PETSC_FORTRAN_CALLBACK_MAXTYPE];

225:   PetscFunctionBegin;
228:   PetscCheck(src->classid == dest->classid, src->comm, PETSC_ERR_ARG_INCOMP, "Objects must be of the same class");

230:   PetscCall(PetscFree(dest->fortran_func_pointers));
231:   PetscCall(PetscMalloc(src->num_fortran_func_pointers * sizeof(void (*)(void)), &dest->fortran_func_pointers));
232:   PetscCall(PetscMemcpy(dest->fortran_func_pointers, src->fortran_func_pointers, src->num_fortran_func_pointers * sizeof(void (*)(void))));

234:   dest->num_fortran_func_pointers = src->num_fortran_func_pointers;

236:   PetscCall(PetscFortranCallbackGetSizes(src->classid, &numcb[PETSC_FORTRAN_CALLBACK_CLASS], &numcb[PETSC_FORTRAN_CALLBACK_SUBTYPE]));
237:   for (cbtype = PETSC_FORTRAN_CALLBACK_CLASS; cbtype < PETSC_FORTRAN_CALLBACK_MAXTYPE; cbtype++) {
238:     PetscCall(PetscFree(dest->fortrancallback[cbtype]));
239:     PetscCall(PetscCalloc1(numcb[cbtype], &dest->fortrancallback[cbtype]));
240:     PetscCall(PetscMemcpy(dest->fortrancallback[cbtype], src->fortrancallback[cbtype], src->num_fortrancallback[cbtype] * sizeof(PetscFortranCallback)));
241:     dest->num_fortrancallback[cbtype] = src->num_fortrancallback[cbtype];
242:   }
243:   PetscFunctionReturn(PETSC_SUCCESS);
244: }

246: /*@C
247:   PetscObjectSetFortranCallback - set Fortran callback function pointer and context

249:   Logically Collective

251:   Input Parameters:
252: + obj    - object on which to set callback
253: . cbtype - callback type (class or subtype)
254: . cid    - address of callback Id, updated if not yet initialized (zero)
255: . func   - Fortran function
256: - ctx    - Fortran context

258:   Level: developer

260:   Note:
261:   This is used to help manage user callback functions that were provided in Fortran

263: .seealso: `PetscObjectGetFortranCallback()`, `PetscFortranCallbackRegister()`, `PetscFortranCallbackGetSizes()`
264: @*/
265: PetscErrorCode PetscObjectSetFortranCallback(PetscObject obj, PetscFortranCallbackType cbtype, PetscFortranCallbackId *cid, void (*func)(void), void *ctx)
266: {
267:   const char *subtype = NULL;

269:   PetscFunctionBegin;
271:   if (cbtype == PETSC_FORTRAN_CALLBACK_SUBTYPE) subtype = obj->type_name;
272:   if (!*cid) PetscCall(PetscFortranCallbackRegister(obj->classid, subtype, cid));
273:   if (*cid >= PETSC_SMALLEST_FORTRAN_CALLBACK + obj->num_fortrancallback[cbtype]) {
274:     PetscFortranCallbackId oldnum = obj->num_fortrancallback[cbtype];
275:     PetscFortranCallbackId newnum = PetscMax(*cid - PETSC_SMALLEST_FORTRAN_CALLBACK + 1, 2 * oldnum);
276:     PetscFortranCallback  *callback;
277:     PetscCall(PetscMalloc1(newnum, &callback));
278:     PetscCall(PetscMemcpy(callback, obj->fortrancallback[cbtype], oldnum * sizeof(*obj->fortrancallback[cbtype])));
279:     PetscCall(PetscFree(obj->fortrancallback[cbtype]));

281:     obj->fortrancallback[cbtype]     = callback;
282:     obj->num_fortrancallback[cbtype] = newnum;
283:   }
284:   obj->fortrancallback[cbtype][*cid - PETSC_SMALLEST_FORTRAN_CALLBACK].func = func;
285:   obj->fortrancallback[cbtype][*cid - PETSC_SMALLEST_FORTRAN_CALLBACK].ctx  = ctx;
286:   PetscFunctionReturn(PETSC_SUCCESS);
287: }

289: /*@C
290:   PetscObjectGetFortranCallback - get Fortran callback function pointer and context

292:   Logically Collective

294:   Input Parameters:
295: + obj    - object on which to get callback
296: . cbtype - callback type
297: - cid    - address of callback Id

299:   Output Parameters:
300: + func - Fortran function (or `NULL` if not needed)
301: - ctx  - Fortran context (or `NULL` if not needed)

303:   Level: developer

305:   Note:
306:   This is used to help manage user callback functions that were provided in Fortran

308: .seealso: `PetscObjectSetFortranCallback()`, `PetscFortranCallbackRegister()`, `PetscFortranCallbackGetSizes()`
309: @*/
310: PetscErrorCode PetscObjectGetFortranCallback(PetscObject obj, PetscFortranCallbackType cbtype, PetscFortranCallbackId cid, void (**func)(void), void **ctx)
311: {
312:   PetscFortranCallback *cb;

314:   PetscFunctionBegin;
316:   PetscCheck(cid >= PETSC_SMALLEST_FORTRAN_CALLBACK, obj->comm, PETSC_ERR_ARG_CORRUPT, "Fortran callback Id invalid");
317:   PetscCheck(cid < PETSC_SMALLEST_FORTRAN_CALLBACK + obj->num_fortrancallback[cbtype], obj->comm, PETSC_ERR_ARG_CORRUPT, "Fortran callback not set on this object");
318:   cb = &obj->fortrancallback[cbtype][cid - PETSC_SMALLEST_FORTRAN_CALLBACK];
319:   if (func) *func = cb->func;
320:   if (ctx) *ctx = cb->ctx;
321:   PetscFunctionReturn(PETSC_SUCCESS);
322: }

324: #if defined(PETSC_USE_LOG)
325: /*@C
326:   PetscObjectsDump - Prints all the currently existing objects.

328:   Input Parameters:
329: + fd  - file pointer
330: - all - by default only tries to display objects created explicitly by the user, if all is `PETSC_TRUE` then lists all outstanding objects

332:   Options Database Key:
333: . -objects_dump <all> - print information about all the objects that exist at the end of the programs run

335:   Level: advanced

337:   Note:
338:   Only MPI rank 0 of `PETSC_COMM_WORLD` prints the values

340: .seealso: `PetscObject`
341: @*/
342: PetscErrorCode PetscObjectsDump(FILE *fd, PetscBool all)
343: {
344:   PetscInt    i, j, k = 0;
345:   PetscObject h;

347:   PetscFunctionBegin;
348:   if (PetscObjectsCounts) {
349:     PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "The following objects were never freed\n"));
350:     PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "-----------------------------------------\n"));
351:     for (i = 0; i < PetscObjectsMaxCounts; i++) {
352:       if ((h = PetscObjects[i])) {
353:         PetscCall(PetscObjectName(h));
354:         {
355:           PetscStack *stack  = NULL;
356:           char       *create = NULL, *rclass = NULL;

358:           /* if the PETSc function the user calls is not a create then this object was NOT directly created by them */
359:           PetscCall(PetscMallocGetStack(h, &stack));
360:           if (stack) {
361:             k = stack->currentsize - 2;
362:             if (!all) {
363:               k = 0;
364:               while (!stack->petscroutine[k]) k++;
365:               PetscCall(PetscStrstr(stack->function[k], "Create", &create));
366:               if (!create) PetscCall(PetscStrstr(stack->function[k], "Get", &create));
367:               PetscCall(PetscStrstr(stack->function[k], h->class_name, &rclass));
368:               if (!create) continue;
369:               if (!rclass) continue;
370:             }
371:           }

373:           PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "[%d] %s %s %s\n", PetscGlobalRank, h->class_name, h->type_name, h->name));

375:           PetscCall(PetscMallocGetStack(h, &stack));
376:           if (stack) {
377:             for (j = k; j >= 0; j--) fprintf(fd, "      [%d]  %s() in %s\n", PetscGlobalRank, stack->function[j], stack->file[j]);
378:           }
379:         }
380:       }
381:     }
382:   }
383:   PetscFunctionReturn(PETSC_SUCCESS);
384: }

386: /*@C
387:   PetscObjectsView - Prints the currently existing objects.

389:   Logically Collective

391:   Input Parameter:
392: . viewer - must be an `PETSCVIEWERASCII` viewer

394:   Level: advanced

396: .seealso: `PetscObject`
397: @*/
398: PetscErrorCode PetscObjectsView(PetscViewer viewer)
399: {
400:   PetscBool isascii;
401:   FILE     *fd;

403:   PetscFunctionBegin;
404:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_WORLD;
405:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
406:   PetscCheck(isascii, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Only supports ASCII viewer");
407:   PetscCall(PetscViewerASCIIGetPointer(viewer, &fd));
408:   PetscCall(PetscObjectsDump(fd, PETSC_TRUE));
409:   PetscFunctionReturn(PETSC_SUCCESS);
410: }

412: /*@C
413:   PetscObjectsGetObject - Get a pointer to a named object

415:   Not Collective

417:   Input Parameter:
418: . name - the name of an object

420:   Output Parameters:
421: + obj       - the object or `NULL` if there is no object, optional, pass in `NULL` if not needed
422: - classname - the name of the class of the object, optional, pass in `NULL` if not needed

424:   Level: advanced

426: .seealso: `PetscObject`
427: @*/
428: PetscErrorCode PetscObjectsGetObject(const char *name, PetscObject *obj, char **classname)
429: {
430:   PetscInt    i;
431:   PetscObject h;
432:   PetscBool   flg;

434:   PetscFunctionBegin;
435:   PetscAssertPointer(name, 1);
436:   if (obj) *obj = NULL;
437:   for (i = 0; i < PetscObjectsMaxCounts; i++) {
438:     if ((h = PetscObjects[i])) {
439:       PetscCall(PetscObjectName(h));
440:       PetscCall(PetscStrcmp(h->name, name, &flg));
441:       if (flg) {
442:         if (obj) *obj = h;
443:         if (classname) *classname = h->class_name;
444:         PetscFunctionReturn(PETSC_SUCCESS);
445:       }
446:     }
447:   }
448:   PetscFunctionReturn(PETSC_SUCCESS);
449: }
450: #endif

452: /*@
453:   PetscObjectSetPrintedOptions - indicate to an object that it should behave as if it has already printed the help for its options so it will not display the help message

455:   Input Parameter:
456: . obj - the `PetscObject`

458:   Level: developer

460:   Developer Notes:
461:   This is used, for example to prevent sequential objects that are created from a parallel object; such as the `KSP` created by
462:   `PCBJACOBI` from all printing the same help messages to the screen

464: .seealso: `PetscOptionsInsert()`, `PetscObject`
465: @*/
466: PetscErrorCode PetscObjectSetPrintedOptions(PetscObject obj)
467: {
468:   PetscFunctionBegin;
469:   PetscAssertPointer(obj, 1);
470:   obj->optionsprinted = PETSC_TRUE;
471:   PetscFunctionReturn(PETSC_SUCCESS);
472: }

474: /*@
475:   PetscObjectInheritPrintedOptions - If the child object is not on the MPI rank 0 process of the parent object and the child is sequential then the child gets it set.

477:   Input Parameters:
478: + pobj - the parent object
479: - obj  - the `PetscObject`

481:   Level: developer

483:   Developer Notes:
484:   This is used, for example to prevent sequential objects that are created from a parallel object; such as the `KSP` created by
485:   `PCBJACOBI` from all printing the same help messages to the screen

487:   This will not handle more complicated situations like with `PCGASM` where children may live on any subset of the parent's processes and overlap

489: .seealso: `PetscOptionsInsert()`, `PetscObjectSetPrintedOptions()`, `PetscObject`
490: @*/
491: PetscErrorCode PetscObjectInheritPrintedOptions(PetscObject pobj, PetscObject obj)
492: {
493:   PetscMPIInt prank, size;

495:   PetscFunctionBegin;
498:   PetscCallMPI(MPI_Comm_rank(pobj->comm, &prank));
499:   PetscCallMPI(MPI_Comm_size(obj->comm, &size));
500:   if (size == 1 && prank > 0) obj->optionsprinted = PETSC_TRUE;
501:   PetscFunctionReturn(PETSC_SUCCESS);
502: }

504: /*@C
505:   PetscObjectAddOptionsHandler - Adds an additional function to check for options when `XXXSetFromOptions()` is called.

507:   Not Collective

509:   Input Parameters:
510: + obj     - the PETSc object
511: . handle  - function that checks for options
512: . destroy - function to destroy `ctx` if provided
513: - ctx     - optional context for check function

515:   Calling sequence of `handle`:
516: + obj                - the PETSc object
517: . PetscOptionsObject - the `PetscOptionItems` object
518: - ctx                - optional context for `handle`

520:   Calling sequence of `destroy`:
521: + obj - the PETSc object
522: - ctx - optional context for `handle`

524:   Level: developer

526: .seealso: `KSPSetFromOptions()`, `PCSetFromOptions()`, `SNESSetFromOptions()`, `PetscObjectProcessOptionsHandlers()`, `PetscObjectDestroyOptionsHandlers()`,
527:           `PetscObject`
528: @*/
529: PetscErrorCode PetscObjectAddOptionsHandler(PetscObject obj, PetscErrorCode (*handle)(PetscObject obj, PetscOptionItems *PetscOptionsObject, void *ctx), PetscErrorCode (*destroy)(PetscObject obj, void *ctx), void *ctx)
530: {
531:   PetscFunctionBegin;
533:   for (PetscInt i = 0; i < obj->noptionhandler; i++) {
534:     PetscBool identical = (PetscBool)(obj->optionhandler[i] == handle && obj->optiondestroy[i] == destroy && obj->optionctx[i] == ctx);
535:     if (identical) PetscFunctionReturn(PETSC_SUCCESS);
536:   }
537:   PetscCheck(obj->noptionhandler < PETSC_MAX_OPTIONS_HANDLER, obj->comm, PETSC_ERR_ARG_OUTOFRANGE, "Too many options handlers added");
538:   obj->optionhandler[obj->noptionhandler] = handle;
539:   obj->optiondestroy[obj->noptionhandler] = destroy;
540:   obj->optionctx[obj->noptionhandler++]   = ctx;
541:   PetscFunctionReturn(PETSC_SUCCESS);
542: }

544: /*@C
545:   PetscObjectProcessOptionsHandlers - Calls all the options handlers attached to an object

547:   Not Collective

549:   Input Parameters:
550: + obj                - the PETSc object
551: - PetscOptionsObject - the options context

553:   Level: developer

555: .seealso: `KSPSetFromOptions()`, `PCSetFromOptions()`, `SNESSetFromOptions()`, `PetscObjectAddOptionsHandler()`, `PetscObjectDestroyOptionsHandlers()`,
556:           `PetscObject`
557: @*/
558: PetscErrorCode PetscObjectProcessOptionsHandlers(PetscObject obj, PetscOptionItems *PetscOptionsObject)
559: {
560:   PetscFunctionBegin;
562:   for (PetscInt i = 0; i < obj->noptionhandler; i++) PetscCall((*obj->optionhandler[i])(obj, PetscOptionsObject, obj->optionctx[i]));
563:   PetscFunctionReturn(PETSC_SUCCESS);
564: }

566: /*@C
567:   PetscObjectDestroyOptionsHandlers - Destroys all the option handlers attached to an object

569:   Not Collective

571:   Input Parameter:
572: . obj - the PETSc object

574:   Level: developer

576: .seealso: `KSPSetFromOptions()`, `PCSetFromOptions()`, `SNESSetFromOptions()`, `PetscObjectAddOptionsHandler()`, `PetscObjectProcessOptionsHandlers()`,
577:           `PetscObject`
578: @*/
579: PetscErrorCode PetscObjectDestroyOptionsHandlers(PetscObject obj)
580: {
581:   PetscFunctionBegin;
583:   for (PetscInt i = 0; i < obj->noptionhandler; i++) {
584:     if (obj->optiondestroy[i]) PetscCall((*obj->optiondestroy[i])(obj, obj->optionctx[i]));
585:   }
586:   obj->noptionhandler = 0;
587:   PetscFunctionReturn(PETSC_SUCCESS);
588: }

590: /*@C
591:   PetscObjectReference - Indicates to a `PetscObject` that it is being
592:   referenced by another `PetscObject`. This increases the reference
593:   count for that object by one.

595:   Logically Collective

597:   Input Parameter:
598: . obj - the PETSc object. This must be cast with (`PetscObject`), for example, `PetscObjectReference`((`PetscObject`)mat);

600:   Level: advanced

602:   Note:
603:   If `obj` is `NULL` this function returns without doing anything.

605: .seealso: `PetscObjectCompose()`, `PetscObjectDereference()`, `PetscObject`
606: @*/
607: PetscErrorCode PetscObjectReference(PetscObject obj)
608: {
609:   PetscFunctionBegin;
610:   if (!obj) PetscFunctionReturn(PETSC_SUCCESS);
612:   obj->refct++;
613:   PetscFunctionReturn(PETSC_SUCCESS);
614: }

616: /*@C
617:   PetscObjectGetReference - Gets the current reference count for a PETSc object.

619:   Not Collective

621:   Input Parameter:
622: . obj - the PETSc object; this must be cast with (`PetscObject`), for example,
623:         `PetscObjectGetReference`((`PetscObject`)mat,&cnt); `obj` cannot be `NULL`

625:   Output Parameter:
626: . cnt - the reference count

628:   Level: advanced

630: .seealso: `PetscObjectCompose()`, `PetscObjectDereference()`, `PetscObjectReference()`, `PetscObject`
631: @*/
632: PetscErrorCode PetscObjectGetReference(PetscObject obj, PetscInt *cnt)
633: {
634:   PetscFunctionBegin;
636:   PetscAssertPointer(cnt, 2);
637:   *cnt = obj->refct;
638:   PetscFunctionReturn(PETSC_SUCCESS);
639: }

641: /*@C
642:   PetscObjectDereference - Indicates to any `PetscObject` that it is being
643:   referenced by one less `PetscObject`. This decreases the reference
644:   count for that object by one.

646:   Collective on `obj` if reference reaches 0 otherwise Logically Collective

648:   Input Parameter:
649: . obj - the PETSc object; this must be cast with (`PetscObject`), for example,
650:         `PetscObjectDereference`((`PetscObject`)mat);

652:   Level: advanced

654:   Notes:
655:   `PetscObjectDestroy()` sets the `obj` pointer to `NULL` after the call, this routine does not.

657:   If `obj` is `NULL` this function returns without doing anything.

659: .seealso: `PetscObjectCompose()`, `PetscObjectReference()`, `PetscObjectDestroy()`, `PetscObject`
660: @*/
661: PetscErrorCode PetscObjectDereference(PetscObject obj)
662: {
663:   PetscFunctionBegin;
664:   if (!obj) PetscFunctionReturn(PETSC_SUCCESS);
666:   if (obj->bops->destroy) PetscCall((*obj->bops->destroy)(&obj));
667:   else PetscCheck(--(obj->refct), PETSC_COMM_SELF, PETSC_ERR_SUP, "This PETSc object does not have a generic destroy routine");
668:   PetscFunctionReturn(PETSC_SUCCESS);
669: }

671: /*
672:      The following routines are the versions private to the PETSc object
673:      data structures.
674: */
675: PetscErrorCode PetscObjectRemoveReference(PetscObject obj, const char name[])
676: {
677:   PetscFunctionBegin;
679:   PetscCall(PetscObjectListRemoveReference(&obj->olist, name));
680:   PetscFunctionReturn(PETSC_SUCCESS);
681: }

683: /*@C
684:   PetscObjectCompose - Associates another PETSc object with a given PETSc object.

686:   Not Collective

688:   Input Parameters:
689: + obj  - the PETSc object; this must be cast with (`PetscObject`), for example,
690:          `PetscObjectCompose`((`PetscObject`)mat,...);
691: . name - name associated with the child object
692: - ptr  - the other PETSc object to associate with the PETSc object; this must also be
693:          cast with (`PetscObject`)

695:   Level: advanced

697:   Notes:
698:   The second objects reference count is automatically increased by one when it is
699:   composed.

701:   Replaces any previous object that had been composed with the same name.

703:   If `ptr` is `NULL` and `name` has previously been composed using an object, then that
704:   entry is removed from `obj`.

706:   `PetscObjectCompose()` can be used with any PETSc object (such as
707:   `Mat`, `Vec`, `KSP`, `SNES`, etc.) or any user-provided object.

709:   `PetscContainerCreate()` can be used to create an object from a
710:   user-provided pointer that may then be composed with PETSc objects using `PetscObjectCompose()`

712: .seealso: `PetscObjectQuery()`, `PetscContainerCreate()`, `PetscObjectComposeFunction()`, `PetscObjectQueryFunction()`, `PetscContainer`,
713:           `PetscContainerSetPointer()`, `PetscObject`
714: @*/
715: PetscErrorCode PetscObjectCompose(PetscObject obj, const char name[], PetscObject ptr)
716: {
717:   PetscFunctionBegin;
719:   PetscAssertPointer(name, 2);
721:   PetscCheck(obj != ptr, PetscObjectComm((PetscObject)obj), PETSC_ERR_SUP, "Cannot compose object with itself");
722:   if (ptr) {
723:     char     *tname;
724:     PetscBool skipreference;

726:     PetscCall(PetscObjectListReverseFind(ptr->olist, obj, &tname, &skipreference));
727:     if (tname) PetscCheck(skipreference, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "An object cannot be composed with an object that was composed with it");
728:   }
729:   PetscCall(PetscObjectListAdd(&obj->olist, name, ptr));
730:   PetscFunctionReturn(PETSC_SUCCESS);
731: }

733: /*@C
734:   PetscObjectQuery  - Gets a PETSc object associated with a given object that was composed with `PetscObjectCompose()`

736:   Not Collective

738:   Input Parameters:
739: + obj  - the PETSc object. It must be cast with a (`PetscObject`), for example,
740:          `PetscObjectCompose`((`PetscObject`)mat,...);
741: . name - name associated with child object
742: - ptr  - the other PETSc object associated with the PETSc object, this must be
743:          cast with (`PetscObject`*)

745:   Level: advanced

747:   Note:
748:   The reference count of neither object is increased in this call

750: .seealso: `PetscObjectCompose()`, `PetscObjectComposeFunction()`, `PetscObjectQueryFunction()`, `PetscContainer`
751:           `PetscContainerGetPointer()`, `PetscObject`
752: @*/
753: PetscErrorCode PetscObjectQuery(PetscObject obj, const char name[], PetscObject *ptr)
754: {
755:   PetscFunctionBegin;
757:   PetscAssertPointer(name, 2);
758:   PetscAssertPointer(ptr, 3);
759:   PetscCall(PetscObjectListFind(obj->olist, name, ptr));
760:   PetscFunctionReturn(PETSC_SUCCESS);
761: }

763: /*MC
764:   PetscObjectComposeFunction - Associates a function with a given PETSc object.

766:   Synopsis:
767: #include <petscsys.h>
768:   PetscErrorCode PetscObjectComposeFunction(PetscObject obj, const char name[], void (*fptr)(void))

770:   Logically Collective

772:   Input Parameters:
773: + obj  - the PETSc object; this must be cast with a (`PetscObject`), for example,
774:          `PetscObjectCompose`((`PetscObject`)mat,...);
775: . name - name associated with the child function
776: - fptr - function pointer

778:   Level: advanced

780:   Notes:
781:   When the first argument of `fptr` is (or is derived from) a `PetscObject` then `PetscTryMethod()` and `PetscUseMethod()`
782:   can be used to call the function directly with error checking.

784:   To remove a registered routine, pass in `NULL` for `fptr`.

786:   `PetscObjectComposeFunction()` can be used with any PETSc object (such as
787:   `Mat`, `Vec`, `KSP`, `SNES`, etc.) or any user-provided object.

789:   `PetscUseTypeMethod()` and `PetscTryTypeMethod()` are used to call a function that is stored in the objects `obj->ops` table.

791: .seealso: `PetscObjectQueryFunction()`, `PetscContainerCreate()` `PetscObjectCompose()`, `PetscObjectQuery()`, `PetscTryMethod()`, `PetscUseMethod()`,
792:           `PetscUseTypeMethod()`, `PetscTryTypeMethod()`, `PetscObject`
793: M*/
794: PetscErrorCode PetscObjectComposeFunction_Private(PetscObject obj, const char name[], void (*fptr)(void))
795: {
796:   PetscFunctionBegin;
798:   PetscAssertPointer(name, 2);
799:   PetscCall(PetscFunctionListAdd(&obj->qlist, name, fptr));
800:   PetscFunctionReturn(PETSC_SUCCESS);
801: }

803: /*MC
804:   PetscObjectQueryFunction - Gets a function associated with a given object.

806:   Synopsis:
807: #include <petscsys.h>
808:   PetscErrorCode PetscObjectQueryFunction(PetscObject obj, const char name[], void (**fptr)(void))

810:   Logically Collective

812:   Input Parameters:
813: + obj  - the PETSc object; this must be cast with (`PetscObject`), for example,
814:          `PetscObjectQueryFunction`((`PetscObject`)ksp,...);
815: - name - name associated with the child function

817:   Output Parameter:
818: . fptr - function pointer

820:   Level: advanced

822: .seealso: `PetscObjectComposeFunction()`, `PetscFunctionListFind()`, `PetscObjectCompose()`, `PetscObjectQuery()`, `PetscObject`
823: M*/
824: PETSC_EXTERN PetscErrorCode PetscObjectQueryFunction_Private(PetscObject obj, const char name[], void (**fptr)(void))
825: {
826:   PetscFunctionBegin;
828:   PetscAssertPointer(name, 2);
829:   PetscCall(PetscFunctionListFind(obj->qlist, name, fptr));
830:   PetscFunctionReturn(PETSC_SUCCESS);
831: }

833: struct _p_PetscContainer {
834:   PETSCHEADER(int);
835:   void *ptr;
836:   PetscErrorCode (*userdestroy)(void *);
837: };

839: /*@C
840:   PetscContainerUserDestroyDefault - Default destroy routine for user-provided data that simply calls `PetscFree()` in the data
841:   provided with `PetscContainerSetPointer()`

843:   Logically Collective on the `PetscContainer` containing the user data

845:   Input Parameter:
846: . ctx - pointer to user-provided data

848:   Level: advanced

850: .seealso: `PetscContainerDestroy()`, `PetscContainerSetUserDestroy()`, `PetscObject`
851: @*/
852: PetscErrorCode PetscContainerUserDestroyDefault(void *ctx)
853: {
854:   PetscFunctionBegin;
855:   PetscCall(PetscFree(ctx));
856:   PetscFunctionReturn(PETSC_SUCCESS);
857: }

859: /*@C
860:   PetscContainerGetPointer - Gets the pointer value contained in the container that was provided with `PetscContainerSetPointer()`

862:   Not Collective

864:   Input Parameter:
865: . obj - the object created with `PetscContainerCreate()`

867:   Output Parameter:
868: . ptr - the pointer value

870:   Level: advanced

872: .seealso: `PetscContainerCreate()`, `PetscContainerDestroy()`, `PetscObject`,
873:           `PetscContainerSetPointer()`
874: @*/
875: PetscErrorCode PetscContainerGetPointer(PetscContainer obj, void **ptr)
876: {
877:   PetscFunctionBegin;
879:   PetscAssertPointer(ptr, 2);
880:   *ptr = obj->ptr;
881:   PetscFunctionReturn(PETSC_SUCCESS);
882: }

884: /*@C
885:   PetscContainerSetPointer - Sets the pointer value contained in the container.

887:   Logically Collective

889:   Input Parameters:
890: + obj - the object created with `PetscContainerCreate()`
891: - ptr - the pointer value

893:   Level: advanced

895: .seealso: `PetscContainerCreate()`, `PetscContainerDestroy()`, `PetscObjectCompose()`, `PetscObjectQuery()`, `PetscObject`,
896:           `PetscContainerGetPointer()`
897: @*/
898: PetscErrorCode PetscContainerSetPointer(PetscContainer obj, void *ptr)
899: {
900:   PetscFunctionBegin;
902:   if (ptr) PetscAssertPointer(ptr, 2);
903:   obj->ptr = ptr;
904:   PetscFunctionReturn(PETSC_SUCCESS);
905: }

907: /*@C
908:   PetscContainerDestroy - Destroys a PETSc container object.

910:   Collective

912:   Input Parameter:
913: . obj - an object that was created with `PetscContainerCreate()`

915:   Level: advanced

917:   Note:
918:   If `PetscContainerSetUserDestroy()` was used to provide a user destroy object for the data provided with `PetscContainerSetPointer()`
919:   then that function is called to destroy the data.

921: .seealso: `PetscContainerCreate()`, `PetscContainerSetUserDestroy()`, `PetscObject`
922: @*/
923: PetscErrorCode PetscContainerDestroy(PetscContainer *obj)
924: {
925:   PetscFunctionBegin;
926:   if (!*obj) PetscFunctionReturn(PETSC_SUCCESS);
928:   if (--((PetscObject)*obj)->refct > 0) {
929:     *obj = NULL;
930:     PetscFunctionReturn(PETSC_SUCCESS);
931:   }
932:   if ((*obj)->userdestroy) PetscCall((*(*obj)->userdestroy)((*obj)->ptr));
933:   PetscCall(PetscHeaderDestroy(obj));
934:   PetscFunctionReturn(PETSC_SUCCESS);
935: }

937: /*@C
938:   PetscContainerSetUserDestroy - Sets name of the user destroy function for the data provided to the `PetscContainer` with `PetscContainerSetPointer()`

940:   Logically Collective

942:   Input Parameters:
943: + obj - an object that was created with `PetscContainerCreate()`
944: - des - name of the user destroy function

946:   Level: advanced

948:   Note:
949:   Use `PetscContainerUserDestroyDefault()` if the memory was obtained by calling `PetscMalloc()` or one of its variants for single memory allocation.

951: .seealso: `PetscContainerDestroy()`, `PetscContainerUserDestroyDefault()`, `PetscMalloc()`, `PetscMalloc1()`, `PetscCalloc()`, `PetscCalloc1()`, `PetscObject`
952: @*/
953: PetscErrorCode PetscContainerSetUserDestroy(PetscContainer obj, PetscErrorCode (*des)(void *))
954: {
955:   PetscFunctionBegin;
957:   obj->userdestroy = des;
958:   PetscFunctionReturn(PETSC_SUCCESS);
959: }

961: PetscClassId PETSC_CONTAINER_CLASSID;

963: /*@C
964:   PetscContainerCreate - Creates a PETSc object that has room to hold a single pointer.

966:   Collective

968:   Input Parameter:
969: . comm - MPI communicator that shares the object

971:   Output Parameter:
972: . container - the container created

974:   Level: advanced

976:   Notes:
977:   This allows one to attach any type of data (accessible through a pointer) with the
978:   `PetscObjectCompose()` function to a `PetscObject`. The data item itself is attached by a
979:   call to `PetscContainerSetPointer()`.

981: .seealso: `PetscContainerDestroy()`, `PetscContainerSetPointer()`, `PetscContainerGetPointer()`, `PetscObjectCompose()`, `PetscObjectQuery()`,
982:           `PetscContainerSetUserDestroy()`, `PetscObject`
983: @*/
984: PetscErrorCode PetscContainerCreate(MPI_Comm comm, PetscContainer *container)
985: {
986:   PetscFunctionBegin;
987:   PetscAssertPointer(container, 2);
988:   PetscCall(PetscSysInitializePackage());
989:   PetscCall(PetscHeaderCreate(*container, PETSC_CONTAINER_CLASSID, "PetscContainer", "Container", "Sys", comm, PetscContainerDestroy, NULL));
990:   PetscFunctionReturn(PETSC_SUCCESS);
991: }

993: /*@
994:   PetscObjectSetFromOptions - Sets generic parameters from user options.

996:   Collective

998:   Input Parameter:
999: . obj - the `PetscObject`

1001:   Level: beginner

1003:   Note:
1004:   We have no generic options at present, so this does nothing

1006: .seealso: `PetscObjectSetOptionsPrefix()`, `PetscObjectGetOptionsPrefix()`, `PetscObject`
1007: @*/
1008: PetscErrorCode PetscObjectSetFromOptions(PetscObject obj)
1009: {
1010:   PetscFunctionBegin;
1012:   PetscFunctionReturn(PETSC_SUCCESS);
1013: }

1015: /*@
1016:   PetscObjectSetUp - Sets up the internal data structures for later use of the object

1018:   Collective

1020:   Input Parameter:
1021: . obj - the `PetscObject`

1023:   Level: advanced

1025:   Note:
1026:   This does nothing at present.

1028: .seealso: `PetscObjectDestroy()`, `PetscObject`
1029: @*/
1030: PetscErrorCode PetscObjectSetUp(PetscObject obj)
1031: {
1032:   PetscFunctionBegin;
1034:   PetscFunctionReturn(PETSC_SUCCESS);
1035: }