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: #if defined(PETSC_USE_LOG)
  8: PETSC_INTERN PetscObject *PetscObjects;
  9: PETSC_INTERN PetscInt     PetscObjectsCounts;
 10: PETSC_INTERN PetscInt     PetscObjectsMaxCounts;
 11: PETSC_INTERN PetscBool    PetscObjectsLog;
 12: #endif

 14: #if defined(PETSC_USE_LOG)
 15: PetscObject *PetscObjects       = NULL;
 16: PetscInt     PetscObjectsCounts = 0, PetscObjectsMaxCounts = 0;
 17: PetscBool    PetscObjectsLog = PETSC_FALSE;
 18: #endif

 20: PETSC_EXTERN PetscErrorCode PetscObjectCompose_Petsc(PetscObject, const char[], PetscObject);
 21: PETSC_EXTERN PetscErrorCode PetscObjectQuery_Petsc(PetscObject, const char[], PetscObject *);
 22: PETSC_EXTERN PetscErrorCode PetscObjectComposeFunction_Petsc(PetscObject, const char[], void (*)(void));
 23: PETSC_EXTERN PetscErrorCode PetscObjectQueryFunction_Petsc(PetscObject, const char[], void (**)(void));

 25: PetscObjectId PetscObjectNewId_Internal(void)
 26: {
 27:   static PetscObjectId idcnt = 1;
 28:   return idcnt++;
 29: }

 31: /*
 32:    PetscHeaderCreate_Private - Creates a base PETSc object header and fills
 33:    in the default values.  Called by the macro PetscHeaderCreate().
 34: */
 35: PetscErrorCode PetscHeaderCreate_Private(PetscObject h, PetscClassId classid, const char class_name[], const char descr[], const char mansec[], MPI_Comm comm, PetscObjectDestroyFunction destroy, PetscObjectViewFunction view)
 36: {
 37:   void       *get_tmp;
 38:   PetscInt64 *cidx;
 39:   PetscMPIInt flg;

 41:   PetscFunctionBegin;
 42:   h->classid               = classid;
 43:   h->class_name            = (char *)class_name;
 44:   h->description           = (char *)descr;
 45:   h->mansec                = (char *)mansec;
 46:   h->refct                 = 1;
 47:   h->non_cyclic_references = NULL;
 48:   h->id                    = PetscObjectNewId_Internal();
 49:   h->bops->destroy         = destroy;
 50:   h->bops->view            = view;
 51:   h->bops->compose         = PetscObjectCompose_Petsc;
 52:   h->bops->query           = PetscObjectQuery_Petsc;
 53:   h->bops->composefunction = PetscObjectComposeFunction_Petsc;
 54:   h->bops->queryfunction   = PetscObjectQueryFunction_Petsc;

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

 58:   /* Increment and store current object creation index */
 59:   PetscCallMPI(MPI_Comm_get_attr(h->comm, Petsc_CreationIdx_keyval, &get_tmp, &flg));
 60:   PetscCheck(flg, h->comm, PETSC_ERR_ARG_CORRUPT, "MPI_Comm does not have an object creation index");
 61:   cidx    = (PetscInt64 *)get_tmp;
 62:   h->cidx = (*cidx)++;
 63:   PetscCallMPI(MPI_Comm_set_attr(h->comm, Petsc_CreationIdx_keyval, cidx));

 65: #if defined(PETSC_USE_LOG)
 66:   /* Keep a record of object created */
 67:   if (PetscObjectsLog) {
 68:     PetscObject *newPetscObjects;
 69:     PetscInt     newPetscObjectsMaxCounts;

 71:     PetscObjectsCounts++;
 72:     for (PetscInt i = 0; i < PetscObjectsMaxCounts; ++i) {
 73:       if (!PetscObjects[i]) {
 74:         PetscObjects[i] = h;
 75:         PetscFunctionReturn(PETSC_SUCCESS);
 76:       }
 77:     }
 78:     /* Need to increase the space for storing PETSc objects */
 79:     if (!PetscObjectsMaxCounts) newPetscObjectsMaxCounts = 100;
 80:     else newPetscObjectsMaxCounts = 2 * PetscObjectsMaxCounts;
 81:     PetscCall(PetscCalloc1(newPetscObjectsMaxCounts, &newPetscObjects));
 82:     PetscCall(PetscArraycpy(newPetscObjects, PetscObjects, PetscObjectsMaxCounts));
 83:     PetscCall(PetscFree(PetscObjects));

 85:     PetscObjects                        = newPetscObjects;
 86:     PetscObjects[PetscObjectsMaxCounts] = h;
 87:     PetscObjectsMaxCounts               = newPetscObjectsMaxCounts;
 88:   }
 89: #endif
 90:   PetscFunctionReturn(PETSC_SUCCESS);
 91: }

 93: PETSC_INTERN PetscBool      PetscMemoryCollectMaximumUsage;
 94: PETSC_INTERN PetscLogDouble PetscMemoryMaximumUsage;

 96: /*
 97:     PetscHeaderDestroy_Private - Destroys a base PETSc object header. Called by
 98:     the macro PetscHeaderDestroy().
 99: */
100: PetscErrorCode PetscHeaderDestroy_Private(PetscObject obj, PetscBool clear_for_reuse)
101: {
102:   PetscFunctionBegin;
104:   PetscCall(PetscLogObjectDestroy(obj));
105:   PetscCall(PetscComposedQuantitiesDestroy(obj));
106:   if (PetscMemoryCollectMaximumUsage) {
107:     PetscLogDouble usage;

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

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

124:   /* destroy allocated quantities */
125:   if (PetscPrintFunctionList) PetscCall(PetscFunctionListPrintNonEmpty(obj->qlist));
126:   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);
127:   PetscCall(PetscFree(obj->name));
128:   PetscCall(PetscFree(obj->prefix));
129:   PetscCall(PetscFree(obj->type_name));

131:   if (clear_for_reuse) {
132:     /* we will assume that obj->bops->view and destroy are safe to leave as-is */
133:     obj->bops->compose         = PetscObjectCompose_Petsc;
134:     obj->bops->query           = PetscObjectQuery_Petsc;
135:     obj->bops->composefunction = PetscObjectComposeFunction_Petsc;
136:     obj->bops->queryfunction   = PetscObjectQueryFunction_Petsc;

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

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

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

187:   - classid
188:   - bops->view
189:   - bops->destroy
190:   - comm
191:   - tag
192:   - class_name
193:   - description
194:   - mansec
195:   - cpp

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

208: /*@C
209:    PetscObjectCopyFortranFunctionPointers - Copy function pointers to another object

211:    Logically Collective

213:    Input Parameters:
214: +  src - source object
215: -  dest - destination object

217:    Level: developer

219:    Note:
220:    Both objects must have the same class.

222:    This is used to help manage user callback functions that were provided in Fortran
223: @*/
224: PetscErrorCode PetscObjectCopyFortranFunctionPointers(PetscObject src, PetscObject dest)
225: {
226:   PetscFortranCallbackId cbtype, numcb[PETSC_FORTRAN_CALLBACK_MAXTYPE];

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

233:   PetscCall(PetscFree(dest->fortran_func_pointers));
234:   PetscCall(PetscMalloc(src->num_fortran_func_pointers * sizeof(void (*)(void)), &dest->fortran_func_pointers));
235:   PetscCall(PetscMemcpy(dest->fortran_func_pointers, src->fortran_func_pointers, src->num_fortran_func_pointers * sizeof(void (*)(void))));

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

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

249: /*@C
250:    PetscObjectSetFortranCallback - set fortran callback function pointer and context

252:    Logically Collective

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

261:    Level: developer

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

266: .seealso: `PetscObjectGetFortranCallback()`
267: @*/
268: PetscErrorCode PetscObjectSetFortranCallback(PetscObject obj, PetscFortranCallbackType cbtype, PetscFortranCallbackId *cid, void (*func)(void), void *ctx)
269: {
270:   const char *subtype = NULL;

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

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

292: /*@C
293:    PetscObjectGetFortranCallback - get fortran callback function pointer and context

295:    Logically Collective

297:    Input Parameters:
298: +  obj - object on which to get callback
299: .  cbtype - callback type
300: -  cid - address of callback Id

302:    Output Parameters:
303: +  func - Fortran function (or NULL if not needed)
304: -  ctx - Fortran context (or NULL if not needed)

306:    Level: developer

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

311: .seealso: `PetscObjectSetFortranCallback()`
312: @*/
313: PetscErrorCode PetscObjectGetFortranCallback(PetscObject obj, PetscFortranCallbackType cbtype, PetscFortranCallbackId cid, void (**func)(void), void **ctx)
314: {
315:   PetscFortranCallback *cb;

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

327: #if defined(PETSC_USE_LOG)
328: /*@C
329:    PetscObjectsDump - Prints all the currently existing objects.

331:    On rank 0 of `PETSC_COMM_WORLD` prints the values

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

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

340:    Level: advanced

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

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

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

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

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

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

390:    Logically Collective

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

395:    Level: advanced

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
422: -  classname - the name of the class

424:    Level: advanced

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

433:   PetscFunctionBegin;
436:   *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:         *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 Parameters:
456: .  obj  - the `PetscObject`

458:    Level: developer

460:    Developer Note:
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()`
465: @*/
466: PetscErrorCode PetscObjectSetPrintedOptions(PetscObject obj)
467: {
468:   PetscFunctionBegin;
470:   obj->optionsprinted = PETSC_TRUE;
471:   PetscFunctionReturn(PETSC_SUCCESS);
472: }

474: /*@
475:    PetscObjectInheritPrintedOptions - If the child object is not on the 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()`
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 context if provided
513: -   ctx - optional context for check function

515:     Level: developer

517: .seealso: `KSPSetFromOptions()`, `PCSetFromOptions()`, `SNESSetFromOptions()`, `PetscObjectProcessOptionsHandlers()`, `PetscObjectDestroyOptionsHandlers()`
518: @*/
519: PetscErrorCode PetscObjectAddOptionsHandler(PetscObject obj, PetscErrorCode (*handle)(PetscObject, PetscOptionItems *, void *), PetscErrorCode (*destroy)(PetscObject, void *), void *ctx)
520: {
521:   PetscFunctionBegin;
523:   PetscCheck(obj->noptionhandler < PETSC_MAX_OPTIONS_HANDLER, obj->comm, PETSC_ERR_ARG_OUTOFRANGE, "To many options handlers added");
524:   obj->optionhandler[obj->noptionhandler] = handle;
525:   obj->optiondestroy[obj->noptionhandler] = destroy;
526:   obj->optionctx[obj->noptionhandler++]   = ctx;
527:   PetscFunctionReturn(PETSC_SUCCESS);
528: }

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

533:     Not Collective

535:     Input Parameters:
536: +   obj - the PETSc object
537: -   PetscOptionsObject - the options context

539:     Level: developer

541: .seealso: `KSPSetFromOptions()`, `PCSetFromOptions()`, `SNESSetFromOptions()`, `PetscObjectAddOptionsHandler()`, `PetscObjectDestroyOptionsHandlers()`
542: @*/
543: PetscErrorCode PetscObjectProcessOptionsHandlers(PetscObject obj, PetscOptionItems *PetscOptionsObject)
544: {
545:   PetscFunctionBegin;
547:   for (PetscInt i = 0; i < obj->noptionhandler; i++) PetscCall((*obj->optionhandler[i])(obj, PetscOptionsObject, obj->optionctx[i]));
548:   PetscFunctionReturn(PETSC_SUCCESS);
549: }

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

554:     Not Collective

556:     Input Parameter:
557: .   obj - the PETSc object

559:     Level: developer

561: .seealso: `KSPSetFromOptions()`, `PCSetFromOptions()`, `SNESSetFromOptions()`, `PetscObjectAddOptionsHandler()`, `PetscObjectProcessOptionsHandlers()`
562: @*/
563: PetscErrorCode PetscObjectDestroyOptionsHandlers(PetscObject obj)
564: {
565:   PetscFunctionBegin;
567:   for (PetscInt i = 0; i < obj->noptionhandler; i++) {
568:     if (obj->optiondestroy[i]) PetscCall((*obj->optiondestroy[i])(obj, obj->optionctx[i]));
569:   }
570:   obj->noptionhandler = 0;
571:   PetscFunctionReturn(PETSC_SUCCESS);
572: }

574: /*@C
575:    PetscObjectReference - Indicates to any `PetscObject` that it is being
576:    referenced by another `PetscObject`. This increases the reference
577:    count for that object by one.

579:    Logically Collective

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

585:    Level: advanced

587: .seealso: `PetscObjectCompose()`, `PetscObjectDereference()`
588: @*/
589: PetscErrorCode PetscObjectReference(PetscObject obj)
590: {
591:   PetscFunctionBegin;
592:   if (!obj) PetscFunctionReturn(PETSC_SUCCESS);
594:   obj->refct++;
595:   PetscFunctionReturn(PETSC_SUCCESS);
596: }

598: /*@C
599:    PetscObjectGetReference - Gets the current reference count for
600:    any PETSc object.

602:    Not Collective

604:    Input Parameter:
605: .  obj - the PETSc object; this must be cast with (`PetscObject`), for example,
606:          `PetscObjectGetReference`((`PetscObject`)mat,&cnt);

608:    Output Parameter:
609: .  cnt - the reference count

611:    Level: advanced

613: .seealso: `PetscObjectCompose()`, `PetscObjectDereference()`, `PetscObjectReference()`
614: @*/
615: PetscErrorCode PetscObjectGetReference(PetscObject obj, PetscInt *cnt)
616: {
617:   PetscFunctionBegin;
620:   *cnt = obj->refct;
621:   PetscFunctionReturn(PETSC_SUCCESS);
622: }

624: /*@C
625:    PetscObjectDereference - Indicates to any `PetscObject` that it is being
626:    referenced by one less `PetscObject`. This decreases the reference
627:    count for that object by one.

629:    Collective on obj if reference reaches 0 otherwise Logically Collective

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

635:    Note:
636:     `PetscObjectDestroy()` sets the obj pointer to null after the call, this routine does not.

638:    Level: advanced

640: .seealso: `PetscObjectCompose()`, `PetscObjectReference()`, `PetscObjectDestroy()`
641: @*/
642: PetscErrorCode PetscObjectDereference(PetscObject obj)
643: {
644:   PetscFunctionBegin;
645:   if (!obj) PetscFunctionReturn(PETSC_SUCCESS);
647:   if (obj->bops->destroy) PetscCall((*obj->bops->destroy)(&obj));
648:   else PetscCheck(--(obj->refct), PETSC_COMM_SELF, PETSC_ERR_SUP, "This PETSc object does not have a generic destroy routine");
649:   PetscFunctionReturn(PETSC_SUCCESS);
650: }

652: /* ----------------------------------------------------------------------- */
653: /*
654:      The following routines are the versions private to the PETSc object
655:      data structures.
656: */
657: PetscErrorCode PetscObjectRemoveReference(PetscObject obj, const char name[])
658: {
659:   PetscFunctionBegin;
661:   PetscCall(PetscObjectListRemoveReference(&obj->olist, name));
662:   PetscFunctionReturn(PETSC_SUCCESS);
663: }

665: PetscErrorCode PetscObjectCompose_Petsc(PetscObject obj, const char name[], PetscObject ptr)
666: {
667:   PetscFunctionBegin;
668:   if (ptr) {
669:     char     *tname;
670:     PetscBool skipreference;

672:     PetscCall(PetscObjectListReverseFind(ptr->olist, obj, &tname, &skipreference));
673:     if (tname) PetscCheck(skipreference, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "An object cannot be composed with an object that was composed with it");
674:   }
675:   PetscCall(PetscObjectListAdd(&obj->olist, name, ptr));
676:   PetscFunctionReturn(PETSC_SUCCESS);
677: }

679: PetscErrorCode PetscObjectQuery_Petsc(PetscObject obj, const char name[], PetscObject *ptr)
680: {
681:   PetscFunctionBegin;
683:   PetscCall(PetscObjectListFind(obj->olist, name, ptr));
684:   PetscFunctionReturn(PETSC_SUCCESS);
685: }

687: PetscErrorCode PetscObjectComposeFunction_Petsc(PetscObject obj, const char name[], void (*ptr)(void))
688: {
689:   PetscFunctionBegin;
691:   PetscCall(PetscFunctionListAdd(&obj->qlist, name, ptr));
692:   PetscFunctionReturn(PETSC_SUCCESS);
693: }

695: PetscErrorCode PetscObjectQueryFunction_Petsc(PetscObject obj, const char name[], void (**ptr)(void))
696: {
697:   PetscFunctionBegin;
699:   PetscCall(PetscFunctionListFind(obj->qlist, name, ptr));
700:   PetscFunctionReturn(PETSC_SUCCESS);
701: }

703: /*@C
704:    PetscObjectCompose - Associates another PETSc object with a given PETSc object.

706:    Not Collective

708:    Input Parameters:
709: +  obj - the PETSc object; this must be cast with (`PetscObject`), for example,
710:          `PetscObjectCompose`((`PetscObject`)mat,...);
711: .  name - name associated with the child object
712: -  ptr - the other PETSc object to associate with the PETSc object; this must also be
713:          cast with (`PetscObject`)

715:    Level: advanced

717:    Notes:
718:    The second objects reference count is automatically increased by one when it is
719:    composed.

721:    Replaces any previous object that had the same name.

723:    If ptr is null and name has previously been composed using an object, then that
724:    entry is removed from the obj.

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

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

732: .seealso: `PetscObjectQuery()`, `PetscContainerCreate()`, `PetscObjectComposeFunction()`, `PetscObjectQueryFunction()`, `PetscContainer`,
733:           `PetscContainerSetPointer()`
734: @*/
735: PetscErrorCode PetscObjectCompose(PetscObject obj, const char name[], PetscObject ptr)
736: {
737:   PetscFunctionBegin;
741:   PetscCheck(obj != ptr, PetscObjectComm((PetscObject)obj), PETSC_ERR_SUP, "Cannot compose object with itself");
742:   PetscCall((*obj->bops->compose)(obj, name, ptr));
743:   PetscFunctionReturn(PETSC_SUCCESS);
744: }

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

749:    Not Collective

751:    Input Parameters:
752: +  obj - the PETSc object
753:          Thus must be cast with a (`PetscObject`), for example,
754:          `PetscObjectCompose`((`PetscObject`)mat,...);
755: .  name - name associated with child object
756: -  ptr - the other PETSc object associated with the PETSc object, this must be
757:          cast with (`PetscObject`*)

759:    Level: advanced

761:    Note:
762:    The reference count of neither object is increased in this call

764: .seealso: `PetscObjectCompose()`, `PetscObjectComposeFunction()`, `PetscObjectQueryFunction()`, `PetscContainer`
765:           `PetscContainerGetPointer()`
766: @*/
767: PetscErrorCode PetscObjectQuery(PetscObject obj, const char name[], PetscObject *ptr)
768: {
769:   PetscFunctionBegin;
773:   PetscCall((*obj->bops->query)(obj, name, ptr));
774:   PetscFunctionReturn(PETSC_SUCCESS);
775: }

777: /*MC
778:    PetscObjectComposeFunction - Associates a function with a given PETSc object.

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

784:    Logically Collective

786:    Input Parameters:
787: +  obj - the PETSc object; this must be cast with a (`PetscObject`), for example,
788:          `PetscObjectCompose`((`PetscObject`)mat,...);
789: .  name - name associated with the child function
790: -  fptr - function pointer

792:    Level: advanced

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

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

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

803:    `PetscCallMethod()` is used to call a function that is stored in the objects `obj->ops` table.

805: .seealso: `PetscObjectQueryFunction()`, `PetscContainerCreate()` `PetscObjectCompose()`, `PetscObjectQuery()`, `PetscTryMethod()`, `PetscUseMethod()`,
806:           `PetscCallMethod()`
807: M*/

809: PetscErrorCode PetscObjectComposeFunction_Private(PetscObject obj, const char name[], void (*fptr)(void))
810: {
811:   PetscFunctionBegin;
814:   PetscCall((*obj->bops->composefunction)(obj, name, fptr));
815:   PetscFunctionReturn(PETSC_SUCCESS);
816: }

818: /*MC
819:    PetscObjectQueryFunction - Gets a function associated with a given object.

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

825:    Logically Collective

827:    Input Parameters:
828: +  obj - the PETSc object; this must be cast with (`PetscObject`), for example,
829:          `PetscObjectQueryFunction`((`PetscObject`)ksp,...);
830: -  name - name associated with the child function

832:    Output Parameter:
833: .  fptr - function pointer

835:    Level: advanced

837: .seealso: `PetscObjectComposeFunction()`, `PetscFunctionListFind()`, `PetscObjectCompose()`, `PetscObjectQuery()`
838: M*/
839: PETSC_EXTERN PetscErrorCode PetscObjectQueryFunction_Private(PetscObject obj, const char name[], void (**ptr)(void))
840: {
841:   PetscFunctionBegin;
844:   PetscCall((*obj->bops->queryfunction)(obj, name, ptr));
845:   PetscFunctionReturn(PETSC_SUCCESS);
846: }

848: struct _p_PetscContainer {
849:   PETSCHEADER(int);
850:   void *ptr;
851:   PetscErrorCode (*userdestroy)(void *);
852: };

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

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

860:    Input Parameter:
861: .  ctx - pointer to user-provided data

863:    Level: advanced

865: .seealso: `PetscContainerDestroy()`, `PetscContainerSetUserDestroy()`
866: @*/
867: PetscErrorCode PetscContainerUserDestroyDefault(void *ctx)
868: {
869:   PetscFunctionBegin;
870:   PetscCall(PetscFree(ctx));
871:   PetscFunctionReturn(PETSC_SUCCESS);
872: }

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

877:    Not Collective

879:    Input Parameter:
880: .  obj - the object created with `PetscContainerCreate()`

882:    Output Parameter:
883: .  ptr - the pointer value

885:    Level: advanced

887: .seealso: `PetscContainerCreate()`, `PetscContainerDestroy()`,
888:           `PetscContainerSetPointer()`
889: @*/
890: PetscErrorCode PetscContainerGetPointer(PetscContainer obj, void **ptr)
891: {
892:   PetscFunctionBegin;
895:   *ptr = obj->ptr;
896:   PetscFunctionReturn(PETSC_SUCCESS);
897: }

899: /*@C
900:    PetscContainerSetPointer - Sets the pointer value contained in the container.

902:    Logically Collective

904:    Input Parameters:
905: +  obj - the object created with `PetscContainerCreate()`
906: -  ptr - the pointer value

908:    Level: advanced

910: .seealso: `PetscContainerCreate()`, `PetscContainerDestroy()`, `PetscObjectCompose()`, `PetscObjectQuery()`,
911:           `PetscContainerGetPointer()`
912: @*/
913: PetscErrorCode PetscContainerSetPointer(PetscContainer obj, void *ptr)
914: {
915:   PetscFunctionBegin;
918:   obj->ptr = ptr;
919:   PetscFunctionReturn(PETSC_SUCCESS);
920: }

922: /*@C
923:    PetscContainerDestroy - Destroys a PETSc container object.

925:    Collective

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

930:    Level: advanced

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

936: .seealso: `PetscContainerCreate()`, `PetscContainerSetUserDestroy()`
937: @*/
938: PetscErrorCode PetscContainerDestroy(PetscContainer *obj)
939: {
940:   PetscFunctionBegin;
941:   if (!*obj) PetscFunctionReturn(PETSC_SUCCESS);
943:   if (--((PetscObject)(*obj))->refct > 0) {
944:     *obj = NULL;
945:     PetscFunctionReturn(PETSC_SUCCESS);
946:   }
947:   if ((*obj)->userdestroy) PetscCall((*(*obj)->userdestroy)((*obj)->ptr));
948:   PetscCall(PetscHeaderDestroy(obj));
949:   PetscFunctionReturn(PETSC_SUCCESS);
950: }

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

955:    Logically Collective

957:    Input Parameters:
958: +  obj - an object that was created with `PetscContainerCreate()`
959: -  des - name of the user destroy function

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

964:    Level: advanced

966: .seealso: `PetscContainerDestroy()`, `PetscContainerUserDestroyDefault()`, `PetscMalloc()`, `PetscMalloc1()`, `PetscCalloc()`, `PetscCalloc1()`
967: @*/
968: PetscErrorCode PetscContainerSetUserDestroy(PetscContainer obj, PetscErrorCode (*des)(void *))
969: {
970:   PetscFunctionBegin;
972:   obj->userdestroy = des;
973:   PetscFunctionReturn(PETSC_SUCCESS);
974: }

976: PetscClassId PETSC_CONTAINER_CLASSID;

978: /*@C
979:    PetscContainerCreate - Creates a PETSc object that has room to hold
980:    a single pointer. This allows one to attach any type of data (accessible
981:    through a pointer) with the `PetscObjectCompose()` function to a `PetscObject`.
982:    The data item itself is attached by a call to `PetscContainerSetPointer()`.

984:    Collective

986:    Input Parameters:
987: .  comm - MPI communicator that shares the object

989:    Output Parameters:
990: .  container - the container created

992:    Level: advanced

994: .seealso: `PetscContainerDestroy()`, `PetscContainerSetPointer()`, `PetscContainerGetPointer()`, `PetscObjectCompose()`, `PetscObjectQuery()`,
995:           `PetscContainerSetUserDestroy()`
996: @*/
997: PetscErrorCode PetscContainerCreate(MPI_Comm comm, PetscContainer *container)
998: {
999:   PetscFunctionBegin;
1001:   PetscCall(PetscSysInitializePackage());
1002:   PetscCall(PetscHeaderCreate(*container, PETSC_CONTAINER_CLASSID, "PetscContainer", "Container", "Sys", comm, PetscContainerDestroy, NULL));
1003:   PetscFunctionReturn(PETSC_SUCCESS);
1004: }

1006: /*@
1007:    PetscObjectSetFromOptions - Sets generic parameters from user options.

1009:    Collective

1011:    Input Parameter:
1012: .  obj - the `PetscObject`

1014:    Note:
1015:    We have no generic options at present, so this does nothing

1017:    Level: beginner

1019: .seealso: `PetscObjectSetOptionsPrefix()`, `PetscObjectGetOptionsPrefix()`
1020: @*/
1021: PetscErrorCode PetscObjectSetFromOptions(PetscObject obj)
1022: {
1023:   PetscFunctionBegin;
1025:   PetscFunctionReturn(PETSC_SUCCESS);
1026: }

1028: /*@
1029:    PetscObjectSetUp - Sets up the internal data structures for the later use.

1031:    Collective

1033:    Input Parameters:
1034: .  obj - the `PetscObject`

1036:    Note:
1037:    This does nothing at present.

1039:    Level: advanced

1041: .seealso: `PetscObjectDestroy()`
1042: @*/
1043: PetscErrorCode PetscObjectSetUp(PetscObject obj)
1044: {
1045:   PetscFunctionBegin;
1047:   PetscFunctionReturn(PETSC_SUCCESS);
1048: }