Actual source code: rgbasic.c
slepc-3.21.1 2024-04-26
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 = NULL;
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: PetscFunctionBegin;
32: PetscCall(PetscFunctionListDestroy(&RGList));
33: RGPackageInitialized = PETSC_FALSE;
34: RGRegisterAllCalled = PETSC_FALSE;
35: PetscFunctionReturn(PETSC_SUCCESS);
36: }
38: /*@C
39: RGInitializePackage - This function initializes everything in the RG package.
40: It is called from PetscDLLibraryRegister() when using dynamic libraries, and
41: on the first call to RGCreate() when using static libraries.
43: Level: developer
45: .seealso: SlepcInitialize()
46: @*/
47: PetscErrorCode RGInitializePackage(void)
48: {
49: char logList[256];
50: PetscBool opt,pkg;
51: PetscClassId classids[1];
53: PetscFunctionBegin;
54: if (RGPackageInitialized) PetscFunctionReturn(PETSC_SUCCESS);
55: RGPackageInitialized = PETSC_TRUE;
56: /* Register Classes */
57: PetscCall(PetscClassIdRegister("Region",&RG_CLASSID));
58: /* Register Constructors */
59: PetscCall(RGRegisterAll());
60: /* Process Info */
61: classids[0] = RG_CLASSID;
62: PetscCall(PetscInfoProcessClass("rg",1,&classids[0]));
63: /* Process summary exclusions */
64: PetscCall(PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt));
65: if (opt) {
66: PetscCall(PetscStrInList("rg",logList,',',&pkg));
67: if (pkg) PetscCall(PetscLogEventDeactivateClass(RG_CLASSID));
68: }
69: /* Register package finalizer */
70: PetscCall(PetscRegisterFinalize(RGFinalizePackage));
71: PetscFunctionReturn(PETSC_SUCCESS);
72: }
74: /*@
75: RGCreate - Creates an RG context.
77: Collective
79: Input Parameter:
80: . comm - MPI communicator
82: Output Parameter:
83: . newrg - location to put the RG context
85: Level: beginner
87: .seealso: RGDestroy(), RG
88: @*/
89: PetscErrorCode RGCreate(MPI_Comm comm,RG *newrg)
90: {
91: RG rg;
93: PetscFunctionBegin;
94: PetscAssertPointer(newrg,2);
95: *newrg = NULL;
96: PetscCall(RGInitializePackage());
97: PetscCall(SlepcHeaderCreate(rg,RG_CLASSID,"RG","Region","RG",comm,RGDestroy,RGView));
98: rg->complement = PETSC_FALSE;
99: rg->sfactor = 1.0;
100: rg->osfactor = 0.0;
101: rg->data = NULL;
103: *newrg = rg;
104: PetscFunctionReturn(PETSC_SUCCESS);
105: }
107: /*@C
108: RGSetOptionsPrefix - Sets the prefix used for searching for all
109: RG options in the database.
111: Logically Collective
113: Input Parameters:
114: + rg - the region context
115: - prefix - the prefix string to prepend to all RG option requests
117: Notes:
118: A hyphen (-) must NOT be given at the beginning of the prefix name.
119: The first character of all runtime options is AUTOMATICALLY the
120: hyphen.
122: Level: advanced
124: .seealso: RGAppendOptionsPrefix()
125: @*/
126: PetscErrorCode RGSetOptionsPrefix(RG rg,const char *prefix)
127: {
128: PetscFunctionBegin;
130: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)rg,prefix));
131: PetscFunctionReturn(PETSC_SUCCESS);
132: }
134: /*@C
135: RGAppendOptionsPrefix - Appends to the prefix used for searching for all
136: RG options in the database.
138: Logically Collective
140: Input Parameters:
141: + rg - the region context
142: - prefix - the prefix string to prepend to all RG option requests
144: Notes:
145: A hyphen (-) must NOT be given at the beginning of the prefix name.
146: The first character of all runtime options is AUTOMATICALLY the hyphen.
148: Level: advanced
150: .seealso: RGSetOptionsPrefix()
151: @*/
152: PetscErrorCode RGAppendOptionsPrefix(RG rg,const char *prefix)
153: {
154: PetscFunctionBegin;
156: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)rg,prefix));
157: PetscFunctionReturn(PETSC_SUCCESS);
158: }
160: /*@C
161: RGGetOptionsPrefix - Gets the prefix used for searching for all
162: RG options in the database.
164: Not Collective
166: Input Parameters:
167: . rg - the region context
169: Output Parameters:
170: . prefix - pointer to the prefix string used is returned
172: Note:
173: On the Fortran side, the user should pass in a string 'prefix' of
174: sufficient length to hold the prefix.
176: Level: advanced
178: .seealso: RGSetOptionsPrefix(), RGAppendOptionsPrefix()
179: @*/
180: PetscErrorCode RGGetOptionsPrefix(RG rg,const char *prefix[])
181: {
182: PetscFunctionBegin;
184: PetscAssertPointer(prefix,2);
185: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)rg,prefix));
186: PetscFunctionReturn(PETSC_SUCCESS);
187: }
189: /*@C
190: RGSetType - Selects the type for the RG object.
192: Logically Collective
194: Input Parameters:
195: + rg - the region context
196: - type - a known type
198: Level: intermediate
200: .seealso: RGGetType()
201: @*/
202: PetscErrorCode RGSetType(RG rg,RGType type)
203: {
204: PetscErrorCode (*r)(RG);
205: PetscBool match;
207: PetscFunctionBegin;
209: PetscAssertPointer(type,2);
211: PetscCall(PetscObjectTypeCompare((PetscObject)rg,type,&match));
212: if (match) PetscFunctionReturn(PETSC_SUCCESS);
214: PetscCall(PetscFunctionListFind(RGList,type,&r));
215: PetscCheck(r,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested RG type %s",type);
217: PetscTryTypeMethod(rg,destroy);
218: PetscCall(PetscMemzero(rg->ops,sizeof(struct _RGOps)));
220: PetscCall(PetscObjectChangeTypeName((PetscObject)rg,type));
221: PetscCall((*r)(rg));
222: PetscFunctionReturn(PETSC_SUCCESS);
223: }
225: /*@C
226: RGGetType - Gets the RG type name (as a string) from the RG context.
228: Not Collective
230: Input Parameter:
231: . rg - the region context
233: Output Parameter:
234: . type - name of the region
236: Level: intermediate
238: .seealso: RGSetType()
239: @*/
240: PetscErrorCode RGGetType(RG rg,RGType *type)
241: {
242: PetscFunctionBegin;
244: PetscAssertPointer(type,2);
245: *type = ((PetscObject)rg)->type_name;
246: PetscFunctionReturn(PETSC_SUCCESS);
247: }
249: /*@
250: RGSetFromOptions - Sets RG options from the options database.
252: Collective
254: Input Parameters:
255: . rg - the region context
257: Notes:
258: To see all options, run your program with the -help option.
260: Level: beginner
262: .seealso: RGSetOptionsPrefix()
263: @*/
264: PetscErrorCode RGSetFromOptions(RG rg)
265: {
266: char type[256];
267: PetscBool flg;
268: PetscReal sfactor;
270: PetscFunctionBegin;
272: PetscCall(RGRegisterAll());
273: PetscObjectOptionsBegin((PetscObject)rg);
274: PetscCall(PetscOptionsFList("-rg_type","Region type","RGSetType",RGList,(char*)(((PetscObject)rg)->type_name?((PetscObject)rg)->type_name:RGINTERVAL),type,sizeof(type),&flg));
275: if (flg) PetscCall(RGSetType(rg,type));
276: else if (!((PetscObject)rg)->type_name) PetscCall(RGSetType(rg,RGINTERVAL));
278: PetscCall(PetscOptionsBool("-rg_complement","Whether region is complemented or not","RGSetComplement",rg->complement,&rg->complement,NULL));
280: PetscCall(PetscOptionsReal("-rg_scale","Scaling factor","RGSetScale",1.0,&sfactor,&flg));
281: if (flg) PetscCall(RGSetScale(rg,sfactor));
283: PetscTryTypeMethod(rg,setfromoptions,PetscOptionsObject);
284: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)rg,PetscOptionsObject));
285: PetscOptionsEnd();
286: PetscFunctionReturn(PETSC_SUCCESS);
287: }
289: /*@C
290: RGView - Prints the RG data structure.
292: Collective
294: Input Parameters:
295: + rg - the region context
296: - viewer - optional visualization context
298: Note:
299: The available visualization contexts include
300: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
301: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
302: output where only the first processor opens
303: the file. All other processors send their
304: data to the first processor to print.
306: The user can open an alternative visualization context with
307: PetscViewerASCIIOpen() - output to a specified file.
309: Level: beginner
311: .seealso: RGCreate()
312: @*/
313: PetscErrorCode RGView(RG rg,PetscViewer viewer)
314: {
315: PetscBool isdraw,isascii;
317: PetscFunctionBegin;
319: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)rg),&viewer));
321: PetscCheckSameComm(rg,1,viewer,2);
322: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
323: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
324: if (isascii) {
325: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)rg,viewer));
326: PetscCall(PetscViewerASCIIPushTab(viewer));
327: PetscTryTypeMethod(rg,view,viewer);
328: PetscCall(PetscViewerASCIIPopTab(viewer));
329: if (rg->complement) PetscCall(PetscViewerASCIIPrintf(viewer," selected region is the complement of the specified one\n"));
330: if (rg->sfactor!=1.0) PetscCall(PetscViewerASCIIPrintf(viewer," scaling factor = %g\n",(double)rg->sfactor));
331: } else if (isdraw) PetscTryTypeMethod(rg,view,viewer);
332: PetscFunctionReturn(PETSC_SUCCESS);
333: }
335: /*@C
336: RGViewFromOptions - View from options
338: Collective
340: Input Parameters:
341: + rg - the region context
342: . obj - optional object
343: - name - command line option
345: Level: intermediate
347: .seealso: RGView(), RGCreate()
348: @*/
349: PetscErrorCode RGViewFromOptions(RG rg,PetscObject obj,const char name[])
350: {
351: PetscFunctionBegin;
353: PetscCall(PetscObjectViewFromOptions((PetscObject)rg,obj,name));
354: PetscFunctionReturn(PETSC_SUCCESS);
355: }
357: /*@
358: RGIsTrivial - Whether it is the trivial region (whole complex plane).
360: Not Collective
362: Input Parameter:
363: . rg - the region context
365: Output Parameter:
366: . trivial - true if the region is equal to the whole complex plane, e.g.,
367: an interval region with all four endpoints unbounded or an
368: ellipse with infinite radius.
370: Level: beginner
372: .seealso: RGCheckInside()
373: @*/
374: PetscErrorCode RGIsTrivial(RG rg,PetscBool *trivial)
375: {
376: PetscFunctionBegin;
379: PetscAssertPointer(trivial,2);
380: *trivial = PETSC_FALSE;
381: PetscTryTypeMethod(rg,istrivial,trivial);
382: PetscFunctionReturn(PETSC_SUCCESS);
383: }
385: /*@
386: RGCheckInside - Determines if a set of given points are inside the region or not.
388: Not Collective
390: Input Parameters:
391: + rg - the region context
392: . n - number of points to check
393: . ar - array of real parts
394: - ai - array of imaginary parts
396: Output Parameter:
397: . inside - array of results (1=inside, 0=on the contour, -1=outside)
399: Note:
400: The point a is expressed as a couple of PetscScalar variables ar,ai.
401: If built with complex scalars, the point is supposed to be stored in ar,
402: otherwise ar,ai contain the real and imaginary parts, respectively.
404: If a scaling factor was set, the points are scaled before checking.
406: Level: intermediate
408: .seealso: RGSetScale(), RGSetComplement()
409: @*/
410: PetscErrorCode RGCheckInside(RG rg,PetscInt n,PetscScalar *ar,PetscScalar *ai,PetscInt *inside)
411: {
412: PetscReal px,py;
413: PetscInt i;
415: PetscFunctionBegin;
418: PetscAssertPointer(ar,3);
419: #if !defined(PETSC_USE_COMPLEX)
420: PetscAssertPointer(ai,4);
421: #endif
422: PetscAssertPointer(inside,5);
424: for (i=0;i<n;i++) {
425: #if defined(PETSC_USE_COMPLEX)
426: px = PetscRealPart(ar[i]);
427: py = PetscImaginaryPart(ar[i]);
428: #else
429: px = ar[i];
430: py = ai[i];
431: #endif
432: if (PetscUnlikely(rg->sfactor != 1.0)) {
433: px /= rg->sfactor;
434: py /= rg->sfactor;
435: }
436: PetscUseTypeMethod(rg,checkinside,px,py,inside+i);
437: if (PetscUnlikely(rg->complement)) inside[i] = -inside[i];
438: }
439: PetscFunctionReturn(PETSC_SUCCESS);
440: }
442: /*@
443: RGIsAxisymmetric - Determines if the region is symmetric with respect
444: to the real or imaginary axis.
446: Not Collective
448: Input Parameters:
449: + rg - the region context
450: - vertical - true if symmetry must be checked against the vertical axis
452: Output Parameter:
453: . symm - true if the region is axisymmetric
455: Note:
456: If the vertical argument is true, symmetry is checked with respect to
457: the vertical axis, otherwise with respect to the horizontal axis.
459: Level: intermediate
461: .seealso: RGCanUseConjugates()
462: @*/
463: PetscErrorCode RGIsAxisymmetric(RG rg,PetscBool vertical,PetscBool *symm)
464: {
465: PetscFunctionBegin;
468: PetscAssertPointer(symm,3);
469: *symm = PETSC_FALSE;
470: PetscTryTypeMethod(rg,isaxisymmetric,vertical,symm);
471: PetscFunctionReturn(PETSC_SUCCESS);
472: }
474: /*@
475: RGCanUseConjugates - Used in contour integral methods to determine whether
476: half of integration points can be avoided (use their conjugates).
478: Not Collective
480: Input Parameters:
481: + rg - the region context
482: - realmats - true if the problem matrices are real
484: Output Parameter:
485: . useconj - whether it is possible to use conjugates
487: Notes:
488: If some integration points are the conjugates of other points, then the
489: associated computational cost can be saved. This depends on the problem
490: matrices being real and also the region being symmetric with respect to
491: the horizontal axis. The result is false if using real arithmetic or
492: in the case of a flat region (height equal to zero).
494: Level: developer
496: .seealso: RGIsAxisymmetric()
497: @*/
498: PetscErrorCode RGCanUseConjugates(RG rg,PetscBool realmats,PetscBool *useconj)
499: {
500: #if defined(PETSC_USE_COMPLEX)
501: PetscReal c,d;
502: PetscBool isaxisymm;
503: #endif
505: PetscFunctionBegin;
508: PetscAssertPointer(useconj,3);
509: *useconj = PETSC_FALSE;
510: #if defined(PETSC_USE_COMPLEX)
511: if (realmats) {
512: PetscCall(RGIsAxisymmetric(rg,PETSC_FALSE,&isaxisymm));
513: if (isaxisymm) {
514: PetscCall(RGComputeBoundingBox(rg,NULL,NULL,&c,&d));
515: if (c!=d) *useconj = PETSC_TRUE;
516: }
517: }
518: #endif
519: PetscFunctionReturn(PETSC_SUCCESS);
520: }
522: /*@
523: RGComputeContour - Computes the coordinates of several points lying on the
524: contour of the region.
526: Not Collective
528: Input Parameters:
529: + rg - the region context
530: - n - number of points to compute
532: Output Parameters:
533: + cr - location to store real parts
534: - ci - location to store imaginary parts
536: Notes:
537: In real scalars, either cr or ci can be NULL (but not both). In complex
538: scalars, the coordinates are stored in cr, which cannot be NULL (ci is
539: not referenced).
541: Level: intermediate
543: .seealso: RGComputeBoundingBox()
544: @*/
545: PetscErrorCode RGComputeContour(RG rg,PetscInt n,PetscScalar cr[],PetscScalar ci[])
546: {
547: PetscInt i;
549: PetscFunctionBegin;
552: #if defined(PETSC_USE_COMPLEX)
553: PetscAssertPointer(cr,3);
554: #else
555: PetscCheck(cr || ci,PetscObjectComm((PetscObject)rg),PETSC_ERR_SUP,"cr and ci cannot be NULL simultaneously");
556: #endif
557: PetscCheck(!rg->complement,PetscObjectComm((PetscObject)rg),PETSC_ERR_SUP,"Cannot compute contour of region with complement flag set");
558: PetscUseTypeMethod(rg,computecontour,n,cr,ci);
559: for (i=0;i<n;i++) {
560: if (cr) cr[i] *= rg->sfactor;
561: if (ci) ci[i] *= rg->sfactor;
562: }
563: PetscFunctionReturn(PETSC_SUCCESS);
564: }
566: /*@
567: RGComputeBoundingBox - Determines the endpoints of a rectangle in the complex plane that
568: contains the region.
570: Not Collective
572: Input Parameter:
573: . rg - the region context
575: Output Parameters:
576: + a - left endpoint of the bounding box in the real axis
577: . b - right endpoint of the bounding box in the real axis
578: . c - bottom endpoint of the bounding box in the imaginary axis
579: - d - top endpoint of the bounding box in the imaginary axis
581: Notes:
582: The bounding box is defined as [a,b]x[c,d]. In regions that are not bounded (e.g. an
583: open interval) or with the complement flag set, it makes no sense to compute a bounding
584: box, so the return values are infinite.
586: Level: intermediate
588: .seealso: RGComputeContour()
589: @*/
590: PetscErrorCode RGComputeBoundingBox(RG rg,PetscReal *a,PetscReal *b,PetscReal *c,PetscReal *d)
591: {
592: PetscFunctionBegin;
596: if (rg->complement) { /* cannot compute bounding box */
597: if (a) *a = -PETSC_MAX_REAL;
598: if (b) *b = PETSC_MAX_REAL;
599: if (c) *c = -PETSC_MAX_REAL;
600: if (d) *d = PETSC_MAX_REAL;
601: } else {
602: PetscUseTypeMethod(rg,computebbox,a,b,c,d);
603: if (a && *a!=-PETSC_MAX_REAL) *a *= rg->sfactor;
604: if (b && *b!= PETSC_MAX_REAL) *b *= rg->sfactor;
605: if (c && *c!=-PETSC_MAX_REAL) *c *= rg->sfactor;
606: if (d && *d!= PETSC_MAX_REAL) *d *= rg->sfactor;
607: }
608: PetscFunctionReturn(PETSC_SUCCESS);
609: }
611: /*@
612: RGComputeQuadrature - Computes the values of the parameters used in a
613: quadrature rule for a contour integral around the boundary of the region.
615: Not Collective
617: Input Parameters:
618: + rg - the region context
619: . quad - the type of quadrature
620: - n - number of quadrature points to compute
622: Output Parameters:
623: + z - quadrature points
624: . zn - normalized quadrature points
625: - w - quadrature weights
627: Notes:
628: In complex scalars, the values returned in z are often the same as those
629: computed by RGComputeContour(), but this is not the case in real scalars
630: where all output arguments are real.
632: The computed values change for different quadrature rules (trapezoidal
633: or Chebyshev).
635: Level: intermediate
637: .seealso: RGComputeContour()
638: @*/
639: PetscErrorCode RGComputeQuadrature(RG rg,RGQuadRule quad,PetscInt n,PetscScalar z[],PetscScalar zn[],PetscScalar w[])
640: {
641: PetscFunctionBegin;
644: PetscAssertPointer(z,4);
645: PetscAssertPointer(zn,5);
646: PetscAssertPointer(w,6);
648: PetscCall(RGComputeContour(rg,n,z,NULL));
649: PetscUseTypeMethod(rg,computequadrature,quad,n,z,zn,w);
650: PetscFunctionReturn(PETSC_SUCCESS);
651: }
653: /*@
654: RGSetComplement - Sets a flag to indicate that the region is the complement
655: of the specified one.
657: Logically Collective
659: Input Parameters:
660: + rg - the region context
661: - flg - the boolean flag
663: Options Database Key:
664: . -rg_complement <bool> - Activate/deactivate the complementation of the region
666: Level: intermediate
668: .seealso: RGGetComplement()
669: @*/
670: PetscErrorCode RGSetComplement(RG rg,PetscBool flg)
671: {
672: PetscFunctionBegin;
675: rg->complement = flg;
676: PetscFunctionReturn(PETSC_SUCCESS);
677: }
679: /*@
680: RGGetComplement - Gets a flag that indicates whether the region
681: is complemented or not.
683: Not Collective
685: Input Parameter:
686: . rg - the region context
688: Output Parameter:
689: . flg - the flag
691: Level: intermediate
693: .seealso: RGSetComplement()
694: @*/
695: PetscErrorCode RGGetComplement(RG rg,PetscBool *flg)
696: {
697: PetscFunctionBegin;
699: PetscAssertPointer(flg,2);
700: *flg = rg->complement;
701: PetscFunctionReturn(PETSC_SUCCESS);
702: }
704: /*@
705: RGSetScale - Sets the scaling factor to be used when checking that a
706: point is inside the region and when computing the contour.
708: Logically Collective
710: Input Parameters:
711: + rg - the region context
712: - sfactor - the scaling factor
714: Options Database Key:
715: . -rg_scale <real> - Sets the scaling factor
717: Level: advanced
719: .seealso: RGGetScale(), RGCheckInside()
720: @*/
721: PetscErrorCode RGSetScale(RG rg,PetscReal sfactor)
722: {
723: PetscFunctionBegin;
726: if (sfactor == (PetscReal)PETSC_DEFAULT || sfactor == (PetscReal)PETSC_DECIDE) sfactor = 1.0;
727: PetscCheck(sfactor>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of scaling factor. Must be > 0");
728: rg->sfactor = sfactor;
729: PetscFunctionReturn(PETSC_SUCCESS);
730: }
732: /*@
733: RGGetScale - Gets the scaling factor.
735: Not Collective
737: Input Parameter:
738: . rg - the region context
740: Output Parameter:
741: . sfactor - the scaling factor
743: Level: advanced
745: .seealso: RGSetScale()
746: @*/
747: PetscErrorCode RGGetScale(RG rg,PetscReal *sfactor)
748: {
749: PetscFunctionBegin;
751: PetscAssertPointer(sfactor,2);
752: *sfactor = rg->sfactor;
753: PetscFunctionReturn(PETSC_SUCCESS);
754: }
756: /*@
757: RGPushScale - Sets an additional scaling factor, that will multiply the
758: user-defined scaling factor.
760: Logically Collective
762: Input Parameters:
763: + rg - the region context
764: - sfactor - the scaling factor
766: Notes:
767: The current implementation does not allow pushing several scaling factors.
769: This is intended for internal use, for instance in polynomial eigensolvers
770: that use parameter scaling.
772: Level: developer
774: .seealso: RGPopScale(), RGSetScale()
775: @*/
776: PetscErrorCode RGPushScale(RG rg,PetscReal sfactor)
777: {
778: PetscFunctionBegin;
781: PetscCheck(sfactor>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of scaling factor. Must be > 0");
782: PetscCheck(!rg->osfactor,PetscObjectComm((PetscObject)rg),PETSC_ERR_SUP,"Current implementation does not allow pushing several scaling factors");
783: rg->osfactor = rg->sfactor;
784: rg->sfactor *= sfactor;
785: PetscFunctionReturn(PETSC_SUCCESS);
786: }
788: /*@
789: RGPopScale - Pops the scaling factor set with RGPushScale().
791: Logically Collective
793: Input Parameter:
794: . rg - the region context
796: Level: developer
798: .seealso: RGPushScale()
799: @*/
800: PetscErrorCode RGPopScale(RG rg)
801: {
802: PetscFunctionBegin;
804: PetscCheck(rg->osfactor,PetscObjectComm((PetscObject)rg),PETSC_ERR_ORDER,"Must call RGPushScale first");
805: rg->sfactor = rg->osfactor;
806: rg->osfactor = 0.0;
807: PetscFunctionReturn(PETSC_SUCCESS);
808: }
810: /*@C
811: RGDestroy - Destroys RG context that was created with RGCreate().
813: Collective
815: Input Parameter:
816: . rg - the region context
818: Level: beginner
820: .seealso: RGCreate()
821: @*/
822: PetscErrorCode RGDestroy(RG *rg)
823: {
824: PetscFunctionBegin;
825: if (!*rg) PetscFunctionReturn(PETSC_SUCCESS);
827: if (--((PetscObject)*rg)->refct > 0) { *rg = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
828: PetscTryTypeMethod(*rg,destroy);
829: PetscCall(PetscHeaderDestroy(rg));
830: PetscFunctionReturn(PETSC_SUCCESS);
831: }
833: /*@C
834: RGRegister - Adds a region to the RG package.
836: Not Collective
838: Input Parameters:
839: + name - name of a new user-defined RG
840: - function - routine to create context
842: Notes:
843: RGRegister() may be called multiple times to add several user-defined regions.
845: Level: advanced
847: .seealso: RGRegisterAll()
848: @*/
849: PetscErrorCode RGRegister(const char *name,PetscErrorCode (*function)(RG))
850: {
851: PetscFunctionBegin;
852: PetscCall(RGInitializePackage());
853: PetscCall(PetscFunctionListAdd(&RGList,name,function));
854: PetscFunctionReturn(PETSC_SUCCESS);
855: }