Actual source code: rgbasic.c

slepc-3.17.2 2022-08-09
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    Basic RG routines
 12: */

 14: #include <slepc/private/rgimpl.h>

 16: PetscFunctionList RGList = 0;
 17: PetscBool         RGRegisterAllCalled = PETSC_FALSE;
 18: PetscClassId      RG_CLASSID = 0;
 19: static PetscBool  RGPackageInitialized = PETSC_FALSE;

 21: /*@C
 22:    RGFinalizePackage - This function destroys everything in the Slepc interface
 23:    to the RG package. It is called from SlepcFinalize().

 25:    Level: developer

 27: .seealso: SlepcFinalize()
 28: @*/
 29: PetscErrorCode RGFinalizePackage(void)
 30: {
 31:   PetscFunctionListDestroy(&RGList);
 32:   RGPackageInitialized = PETSC_FALSE;
 33:   RGRegisterAllCalled  = PETSC_FALSE;
 34:   PetscFunctionReturn(0);
 35: }

 37: /*@C
 38:   RGInitializePackage - This function initializes everything in the RG package.
 39:   It is called from PetscDLLibraryRegister() when using dynamic libraries, and
 40:   on the first call to RGCreate() when using static libraries.

 42:   Level: developer

 44: .seealso: SlepcInitialize()
 45: @*/
 46: PetscErrorCode RGInitializePackage(void)
 47: {
 48:   char           logList[256];
 49:   PetscBool      opt,pkg;
 50:   PetscClassId   classids[1];

 52:   if (RGPackageInitialized) PetscFunctionReturn(0);
 53:   RGPackageInitialized = PETSC_TRUE;
 54:   /* Register Classes */
 55:   PetscClassIdRegister("Region",&RG_CLASSID);
 56:   /* Register Constructors */
 57:   RGRegisterAll();
 58:   /* Process Info */
 59:   classids[0] = RG_CLASSID;
 60:   PetscInfoProcessClass("rg",1,&classids[0]);
 61:   /* Process summary exclusions */
 62:   PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt);
 63:   if (opt) {
 64:     PetscStrInList("rg",logList,',',&pkg);
 65:     if (pkg) PetscLogEventDeactivateClass(RG_CLASSID);
 66:   }
 67:   /* Register package finalizer */
 68:   PetscRegisterFinalize(RGFinalizePackage);
 69:   PetscFunctionReturn(0);
 70: }

 72: /*@
 73:    RGCreate - Creates an RG context.

 75:    Collective

 77:    Input Parameter:
 78: .  comm - MPI communicator

 80:    Output Parameter:
 81: .  newrg - location to put the RG context

 83:    Level: beginner

 85: .seealso: RGDestroy(), RG
 86: @*/
 87: PetscErrorCode RGCreate(MPI_Comm comm,RG *newrg)
 88: {
 89:   RG             rg;

 92:   *newrg = 0;
 93:   RGInitializePackage();
 94:   SlepcHeaderCreate(rg,RG_CLASSID,"RG","Region","RG",comm,RGDestroy,RGView);
 95:   rg->complement = PETSC_FALSE;
 96:   rg->sfactor    = 1.0;
 97:   rg->osfactor   = 0.0;
 98:   rg->data       = NULL;

100:   *newrg = rg;
101:   PetscFunctionReturn(0);
102: }

104: /*@C
105:    RGSetOptionsPrefix - Sets the prefix used for searching for all
106:    RG options in the database.

108:    Logically Collective on rg

110:    Input Parameters:
111: +  rg     - the region context
112: -  prefix - the prefix string to prepend to all RG option requests

114:    Notes:
115:    A hyphen (-) must NOT be given at the beginning of the prefix name.
116:    The first character of all runtime options is AUTOMATICALLY the
117:    hyphen.

119:    Level: advanced

121: .seealso: RGAppendOptionsPrefix()
122: @*/
123: PetscErrorCode RGSetOptionsPrefix(RG rg,const char *prefix)
124: {
126:   PetscObjectSetOptionsPrefix((PetscObject)rg,prefix);
127:   PetscFunctionReturn(0);
128: }

130: /*@C
131:    RGAppendOptionsPrefix - Appends to the prefix used for searching for all
132:    RG options in the database.

134:    Logically Collective on rg

136:    Input Parameters:
137: +  rg     - the region context
138: -  prefix - the prefix string to prepend to all RG option requests

140:    Notes:
141:    A hyphen (-) must NOT be given at the beginning of the prefix name.
142:    The first character of all runtime options is AUTOMATICALLY the hyphen.

144:    Level: advanced

146: .seealso: RGSetOptionsPrefix()
147: @*/
148: PetscErrorCode RGAppendOptionsPrefix(RG rg,const char *prefix)
149: {
151:   PetscObjectAppendOptionsPrefix((PetscObject)rg,prefix);
152:   PetscFunctionReturn(0);
153: }

155: /*@C
156:    RGGetOptionsPrefix - Gets the prefix used for searching for all
157:    RG options in the database.

159:    Not Collective

161:    Input Parameters:
162: .  rg - the region context

164:    Output Parameters:
165: .  prefix - pointer to the prefix string used is returned

167:    Note:
168:    On the Fortran side, the user should pass in a string 'prefix' of
169:    sufficient length to hold the prefix.

171:    Level: advanced

173: .seealso: RGSetOptionsPrefix(), RGAppendOptionsPrefix()
174: @*/
175: PetscErrorCode RGGetOptionsPrefix(RG rg,const char *prefix[])
176: {
179:   PetscObjectGetOptionsPrefix((PetscObject)rg,prefix);
180:   PetscFunctionReturn(0);
181: }

183: /*@C
184:    RGSetType - Selects the type for the RG object.

186:    Logically Collective on rg

188:    Input Parameters:
189: +  rg   - the region context
190: -  type - a known type

192:    Level: intermediate

194: .seealso: RGGetType()
195: @*/
196: PetscErrorCode RGSetType(RG rg,RGType type)
197: {
198:   PetscErrorCode (*r)(RG);
199:   PetscBool      match;


204:   PetscObjectTypeCompare((PetscObject)rg,type,&match);
205:   if (match) PetscFunctionReturn(0);

207:   PetscFunctionListFind(RGList,type,&r);

210:   if (rg->ops->destroy) (*rg->ops->destroy)(rg);
211:   PetscMemzero(rg->ops,sizeof(struct _RGOps));

213:   PetscObjectChangeTypeName((PetscObject)rg,type);
214:   (*r)(rg);
215:   PetscFunctionReturn(0);
216: }

218: /*@C
219:    RGGetType - Gets the RG type name (as a string) from the RG context.

221:    Not Collective

223:    Input Parameter:
224: .  rg - the region context

226:    Output Parameter:
227: .  type - name of the region

229:    Level: intermediate

231: .seealso: RGSetType()
232: @*/
233: PetscErrorCode RGGetType(RG rg,RGType *type)
234: {
237:   *type = ((PetscObject)rg)->type_name;
238:   PetscFunctionReturn(0);
239: }

241: /*@
242:    RGSetFromOptions - Sets RG options from the options database.

244:    Collective on rg

246:    Input Parameters:
247: .  rg - the region context

249:    Notes:
250:    To see all options, run your program with the -help option.

252:    Level: beginner

254: .seealso: RGSetOptionsPrefix()
255: @*/
256: PetscErrorCode RGSetFromOptions(RG rg)
257: {
259:   char           type[256];
260:   PetscBool      flg;
261:   PetscReal      sfactor;

264:   RGRegisterAll();
265:   ierr = PetscObjectOptionsBegin((PetscObject)rg);
266:     PetscOptionsFList("-rg_type","Region type","RGSetType",RGList,(char*)(((PetscObject)rg)->type_name?((PetscObject)rg)->type_name:RGINTERVAL),type,sizeof(type),&flg);
267:     if (flg) RGSetType(rg,type);
268:     else if (!((PetscObject)rg)->type_name) RGSetType(rg,RGINTERVAL);

270:     PetscOptionsBool("-rg_complement","Whether region is complemented or not","RGSetComplement",rg->complement,&rg->complement,NULL);

272:     PetscOptionsReal("-rg_scale","Scaling factor","RGSetScale",1.0,&sfactor,&flg);
273:     if (flg) RGSetScale(rg,sfactor);

275:     if (rg->ops->setfromoptions) (*rg->ops->setfromoptions)(PetscOptionsObject,rg);
276:     PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)rg);
277:   ierr = PetscOptionsEnd();
278:   PetscFunctionReturn(0);
279: }

281: /*@C
282:    RGView - Prints the RG data structure.

284:    Collective on rg

286:    Input Parameters:
287: +  rg - the region context
288: -  viewer - optional visualization context

290:    Note:
291:    The available visualization contexts include
292: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
293: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
294:          output where only the first processor opens
295:          the file.  All other processors send their
296:          data to the first processor to print.

298:    The user can open an alternative visualization context with
299:    PetscViewerASCIIOpen() - output to a specified file.

301:    Level: beginner

303: .seealso: RGCreate()
304: @*/
305: PetscErrorCode RGView(RG rg,PetscViewer viewer)
306: {
307:   PetscBool      isdraw,isascii;

310:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)rg),&viewer);
313:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
314:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
315:   if (isascii) {
316:     PetscObjectPrintClassNamePrefixType((PetscObject)rg,viewer);
317:     if (rg->ops->view) {
318:       PetscViewerASCIIPushTab(viewer);
319:       (*rg->ops->view)(rg,viewer);
320:       PetscViewerASCIIPopTab(viewer);
321:     }
322:     if (rg->complement) PetscViewerASCIIPrintf(viewer,"  selected region is the complement of the specified one\n");
323:     if (rg->sfactor!=1.0) PetscViewerASCIIPrintf(viewer,"  scaling factor = %g\n",(double)rg->sfactor);
324:   } else if (isdraw) {
325:     if (rg->ops->view) (*rg->ops->view)(rg,viewer);
326:   }
327:   PetscFunctionReturn(0);
328: }

330: /*@C
331:    RGViewFromOptions - View from options

333:    Collective on RG

335:    Input Parameters:
336: +  rg   - the region context
337: .  obj  - optional object
338: -  name - command line option

340:    Level: intermediate

342: .seealso: RGView(), RGCreate()
343: @*/
344: PetscErrorCode RGViewFromOptions(RG rg,PetscObject obj,const char name[])
345: {
347:   PetscObjectViewFromOptions((PetscObject)rg,obj,name);
348:   PetscFunctionReturn(0);
349: }

351: /*@
352:    RGIsTrivial - Whether it is the trivial region (whole complex plane).

354:    Not Collective

356:    Input Parameter:
357: .  rg - the region context

359:    Output Parameter:
360: .  trivial - true if the region is equal to the whole complex plane, e.g.,
361:              an interval region with all four endpoints unbounded or an
362:              ellipse with infinite radius.

364:    Level: beginner

366: .seealso: RGCheckInside()
367: @*/
368: PetscErrorCode RGIsTrivial(RG rg,PetscBool *trivial)
369: {
373:   if (rg->ops->istrivial) (*rg->ops->istrivial)(rg,trivial);
374:   else *trivial = PETSC_FALSE;
375:   PetscFunctionReturn(0);
376: }

378: /*@
379:    RGCheckInside - Determines if a set of given points are inside the region or not.

381:    Not Collective

383:    Input Parameters:
384: +  rg - the region context
385: .  n  - number of points to check
386: .  ar - array of real parts
387: -  ai - array of imaginary parts

389:    Output Parameter:
390: .  inside - array of results (1=inside, 0=on the contour, -1=outside)

392:    Note:
393:    The point a is expressed as a couple of PetscScalar variables ar,ai.
394:    If built with complex scalars, the point is supposed to be stored in ar,
395:    otherwise ar,ai contain the real and imaginary parts, respectively.

397:    If a scaling factor was set, the points are scaled before checking.

399:    Level: intermediate

401: .seealso: RGSetScale(), RGSetComplement()
402: @*/
403: PetscErrorCode RGCheckInside(RG rg,PetscInt n,PetscScalar *ar,PetscScalar *ai,PetscInt *inside)
404: {
405:   PetscReal      px,py;
406:   PetscInt       i;

411: #if !defined(PETSC_USE_COMPLEX)
413: #endif

416:   for (i=0;i<n;i++) {
417: #if defined(PETSC_USE_COMPLEX)
418:     px = PetscRealPart(ar[i]);
419:     py = PetscImaginaryPart(ar[i]);
420: #else
421:     px = ar[i];
422:     py = ai[i];
423: #endif
424:     if (PetscUnlikely(rg->sfactor != 1.0)) {
425:       px /= rg->sfactor;
426:       py /= rg->sfactor;
427:     }
428:     (*rg->ops->checkinside)(rg,px,py,inside+i);
429:     if (PetscUnlikely(rg->complement)) inside[i] = -inside[i];
430:   }
431:   PetscFunctionReturn(0);
432: }

434: /*@
435:    RGIsAxisymmetric - Determines if the region is symmetric with respect
436:    to the real or imaginary axis.

438:    Not Collective

440:    Input Parameters:
441: +  rg       - the region context
442: -  vertical - true if symmetry must be checked against the vertical axis

444:    Output Parameter:
445: .  symm - true if the region is axisymmetric

447:    Note:
448:    If the vertical argument is true, symmetry is checked with respect to
449:    the vertical axis, otherwise with respect to the horizontal axis.

451:    Level: intermediate

453: .seealso: RGCanUseConjugates()
454: @*/
455: PetscErrorCode RGIsAxisymmetric(RG rg,PetscBool vertical,PetscBool *symm)
456: {

461:   if (rg->ops->isaxisymmetric) (*rg->ops->isaxisymmetric)(rg,vertical,symm);
462:   else *symm = PETSC_FALSE;
463:   PetscFunctionReturn(0);
464: }

466: /*@
467:    RGCanUseConjugates - Used in contour integral methods to determine whether
468:    half of integration points can be avoided (use their conjugates).

470:    Not Collective

472:    Input Parameters:
473: +  rg       - the region context
474: -  realmats - true if the problem matrices are real

476:    Output Parameter:
477: .  useconj  - whether it is possible to use conjugates

479:    Notes:
480:    If some integration points are the conjugates of other points, then the
481:    associated computational cost can be saved. This depends on the problem
482:    matrices being real and also the region being symmetric with respect to
483:    the horizontal axis. The result is false if using real arithmetic or
484:    in the case of a flat region (height equal to zero).

486:    Level: developer

488: .seealso: RGIsAxisymmetric()
489: @*/
490: PetscErrorCode RGCanUseConjugates(RG rg,PetscBool realmats,PetscBool *useconj)
491: {
492: #if defined(PETSC_USE_COMPLEX)
493:   PetscReal      c,d;
494:   PetscBool      isaxisymm;
495: #endif

500:   *useconj = PETSC_FALSE;
501: #if defined(PETSC_USE_COMPLEX)
502:   if (realmats) {
503:     RGIsAxisymmetric(rg,PETSC_FALSE,&isaxisymm);
504:     if (isaxisymm) {
505:       RGComputeBoundingBox(rg,NULL,NULL,&c,&d);
506:       if (c!=d) *useconj = PETSC_TRUE;
507:     }
508:   }
509: #endif
510:   PetscFunctionReturn(0);
511: }

513: /*@
514:    RGComputeContour - Computes the coordinates of several points lying on the
515:    contour of the region.

517:    Not Collective

519:    Input Parameters:
520: +  rg - the region context
521: -  n  - number of points to compute

523:    Output Parameters:
524: +  cr - location to store real parts
525: -  ci - location to store imaginary parts

527:    Notes:
528:    In real scalars, either cr or ci can be NULL (but not both). In complex
529:    scalars, the coordinates are stored in cr, which cannot be NULL (ci is
530:    not referenced).

532:    Level: intermediate

534: .seealso: RGComputeBoundingBox()
535: @*/
536: PetscErrorCode RGComputeContour(RG rg,PetscInt n,PetscScalar cr[],PetscScalar ci[])
537: {
538:   PetscInt       i;

542: #if defined(PETSC_USE_COMPLEX)
544: #else
546: #endif
548:   (*rg->ops->computecontour)(rg,n,cr,ci);
549:   for (i=0;i<n;i++) {
550:     if (cr) cr[i] *= rg->sfactor;
551:     if (ci) ci[i] *= rg->sfactor;
552:   }
553:   PetscFunctionReturn(0);
554: }

556: /*@
557:    RGComputeBoundingBox - Determines the endpoints of a rectangle in the complex plane that
558:    contains the region.

560:    Not Collective

562:    Input Parameter:
563: .  rg - the region context

565:    Output Parameters:
566: +  a - left endpoint of the bounding box in the real axis
567: .  b - right endpoint of the bounding box in the real axis
568: .  c - bottom endpoint of the bounding box in the imaginary axis
569: -  d - top endpoint of the bounding box in the imaginary axis

571:    Notes:
572:    The bounding box is defined as [a,b]x[c,d]. In regions that are not bounded (e.g. an
573:    open interval) or with the complement flag set, it makes no sense to compute a bounding
574:    box, so the return values are infinite.

576:    Level: intermediate

578: .seealso: RGComputeContour()
579: @*/
580: PetscErrorCode RGComputeBoundingBox(RG rg,PetscReal *a,PetscReal *b,PetscReal *c,PetscReal *d)
581: {

585:   if (rg->complement) {  /* cannot compute bounding box */
586:     if (a) *a = -PETSC_MAX_REAL;
587:     if (b) *b =  PETSC_MAX_REAL;
588:     if (c) *c = -PETSC_MAX_REAL;
589:     if (d) *d =  PETSC_MAX_REAL;
590:   } else {
591:     (*rg->ops->computebbox)(rg,a,b,c,d);
592:     if (a && *a!=-PETSC_MAX_REAL) *a *= rg->sfactor;
593:     if (b && *b!= PETSC_MAX_REAL) *b *= rg->sfactor;
594:     if (c && *c!=-PETSC_MAX_REAL) *c *= rg->sfactor;
595:     if (d && *d!= PETSC_MAX_REAL) *d *= rg->sfactor;
596:   }
597:   PetscFunctionReturn(0);
598: }

600: /*@
601:    RGComputeQuadrature - Computes the values of the parameters used in a
602:    quadrature rule for a contour integral around the boundary of the region.

604:    Not Collective

606:    Input Parameters:
607: +  rg   - the region context
608: .  quad - the type of quadrature
609: -  n    - number of quadrature points to compute

611:    Output Parameters:
612: +  z  - quadrature points
613: .  zn - normalized quadrature points
614: -  w  - quadrature weights

616:    Notes:
617:    In complex scalars, the values returned in z are often the same as those
618:    computed by RGComputeContour(), but this is not the case in real scalars
619:    where all output arguments are real.

621:    The computed values change for different quadrature rules (trapezoidal
622:    or Chebyshev).

624:    Level: intermediate

626: .seealso: RGComputeContour()
627: @*/
628: PetscErrorCode RGComputeQuadrature(RG rg,RGQuadRule quad,PetscInt n,PetscScalar z[],PetscScalar zn[],PetscScalar w[])
629: {

636:   RGComputeContour(rg,n,z,NULL);
638:   (*rg->ops->computequadrature)(rg,quad,n,z,zn,w);
639:   PetscFunctionReturn(0);
640: }

642: /*@
643:    RGSetComplement - Sets a flag to indicate that the region is the complement
644:    of the specified one.

646:    Logically Collective on rg

648:    Input Parameters:
649: +  rg  - the region context
650: -  flg - the boolean flag

652:    Options Database Key:
653: .  -rg_complement <bool> - Activate/deactivate the complementation of the region

655:    Level: intermediate

657: .seealso: RGGetComplement()
658: @*/
659: PetscErrorCode RGSetComplement(RG rg,PetscBool flg)
660: {
663:   rg->complement = flg;
664:   PetscFunctionReturn(0);
665: }

667: /*@
668:    RGGetComplement - Gets a flag that that indicates whether the region
669:    is complemented or not.

671:    Not Collective

673:    Input Parameter:
674: .  rg - the region context

676:    Output Parameter:
677: .  flg - the flag

679:    Level: intermediate

681: .seealso: RGSetComplement()
682: @*/
683: PetscErrorCode RGGetComplement(RG rg,PetscBool *flg)
684: {
687:   *flg = rg->complement;
688:   PetscFunctionReturn(0);
689: }

691: /*@
692:    RGSetScale - Sets the scaling factor to be used when checking that a
693:    point is inside the region and when computing the contour.

695:    Logically Collective on rg

697:    Input Parameters:
698: +  rg      - the region context
699: -  sfactor - the scaling factor

701:    Options Database Key:
702: .  -rg_scale <real> - Sets the scaling factor

704:    Level: advanced

706: .seealso: RGGetScale(), RGCheckInside()
707: @*/
708: PetscErrorCode RGSetScale(RG rg,PetscReal sfactor)
709: {
712:   if (sfactor == PETSC_DEFAULT || sfactor == PETSC_DECIDE) rg->sfactor = 1.0;
713:   else {
715:     rg->sfactor = sfactor;
716:   }
717:   PetscFunctionReturn(0);
718: }

720: /*@
721:    RGGetScale - Gets the scaling factor.

723:    Not Collective

725:    Input Parameter:
726: .  rg - the region context

728:    Output Parameter:
729: .  sfactor - the scaling factor

731:    Level: advanced

733: .seealso: RGSetScale()
734: @*/
735: PetscErrorCode RGGetScale(RG rg,PetscReal *sfactor)
736: {
739:   *sfactor = rg->sfactor;
740:   PetscFunctionReturn(0);
741: }

743: /*@
744:    RGPushScale - Sets an additional scaling factor, that will multiply the
745:    user-defined scaling factor.

747:    Logically Collective on rg

749:    Input Parameters:
750: +  rg      - the region context
751: -  sfactor - the scaling factor

753:    Notes:
754:    The current implementation does not allow pushing several scaling factors.

756:    This is intended for internal use, for instance in polynomial eigensolvers
757:    that use parameter scaling.

759:    Level: developer

761: .seealso: RGPopScale(), RGSetScale()
762: @*/
763: PetscErrorCode RGPushScale(RG rg,PetscReal sfactor)
764: {
769:   rg->osfactor = rg->sfactor;
770:   rg->sfactor *= sfactor;
771:   PetscFunctionReturn(0);
772: }

774: /*@
775:    RGPopScale - Pops the scaling factor set with RGPushScale().

777:    Not Collective

779:    Input Parameter:
780: .  rg - the region context

782:    Level: developer

784: .seealso: RGPushScale()
785: @*/
786: PetscErrorCode RGPopScale(RG rg)
787: {
790:   rg->sfactor  = rg->osfactor;
791:   rg->osfactor = 0.0;
792:   PetscFunctionReturn(0);
793: }

795: /*@C
796:    RGDestroy - Destroys RG context that was created with RGCreate().

798:    Collective on rg

800:    Input Parameter:
801: .  rg - the region context

803:    Level: beginner

805: .seealso: RGCreate()
806: @*/
807: PetscErrorCode RGDestroy(RG *rg)
808: {
809:   if (!*rg) PetscFunctionReturn(0);
811:   if (--((PetscObject)(*rg))->refct > 0) { *rg = 0; PetscFunctionReturn(0); }
812:   if ((*rg)->ops->destroy) (*(*rg)->ops->destroy)(*rg);
813:   PetscHeaderDestroy(rg);
814:   PetscFunctionReturn(0);
815: }

817: /*@C
818:    RGRegister - Adds a region to the RG package.

820:    Not collective

822:    Input Parameters:
823: +  name - name of a new user-defined RG
824: -  function - routine to create context

826:    Notes:
827:    RGRegister() may be called multiple times to add several user-defined regions.

829:    Level: advanced

831: .seealso: RGRegisterAll()
832: @*/
833: PetscErrorCode RGRegister(const char *name,PetscErrorCode (*function)(RG))
834: {
835:   RGInitializePackage();
836:   PetscFunctionListAdd(&RGList,name,function);
837:   PetscFunctionReturn(0);
838: }