Actual source code: rgbasic.c
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: [](sec:rg), `SlepcFinalize()`, `BVInitializePackage()`
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_slepc()` when using dynamic libraries, and
41: on the first call to `RGCreate()` when using shared or static libraries.
43: Note:
44: This function never needs to be called by SLEPc users.
46: Level: developer
48: .seealso: [](sec:rg), `RG`, `SlepcInitialize()`, `RGFinalizePackage()`
49: @*/
50: PetscErrorCode RGInitializePackage(void)
51: {
52: char logList[256];
53: PetscBool opt,pkg;
54: PetscClassId classids[1];
56: PetscFunctionBegin;
57: if (RGPackageInitialized) PetscFunctionReturn(PETSC_SUCCESS);
58: RGPackageInitialized = PETSC_TRUE;
59: /* Register Classes */
60: PetscCall(PetscClassIdRegister("Region",&RG_CLASSID));
61: /* Register Constructors */
62: PetscCall(RGRegisterAll());
63: /* Process Info */
64: classids[0] = RG_CLASSID;
65: PetscCall(PetscInfoProcessClass("rg",1,&classids[0]));
66: /* Process summary exclusions */
67: PetscCall(PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt));
68: if (opt) {
69: PetscCall(PetscStrInList("rg",logList,',',&pkg));
70: if (pkg) PetscCall(PetscLogEventDeactivateClass(RG_CLASSID));
71: }
72: /* Register package finalizer */
73: PetscCall(PetscRegisterFinalize(RGFinalizePackage));
74: PetscFunctionReturn(PETSC_SUCCESS);
75: }
77: /*@
78: RGCreate - Creates an `RG` context.
80: Collective
82: Input Parameter:
83: . comm - MPI communicator
85: Output Parameter:
86: . newrg - location to put the `RG` context
88: Level: beginner
90: .seealso: [](sec:rg), `RG`, `RGDestroy()`
91: @*/
92: PetscErrorCode RGCreate(MPI_Comm comm,RG *newrg)
93: {
94: RG rg;
96: PetscFunctionBegin;
97: PetscAssertPointer(newrg,2);
98: PetscCall(RGInitializePackage());
99: PetscCall(SlepcHeaderCreate(rg,RG_CLASSID,"RG","Region","RG",comm,RGDestroy,RGView));
100: rg->complement = PETSC_FALSE;
101: rg->sfactor = 1.0;
102: rg->osfactor = 0.0;
103: rg->data = NULL;
105: *newrg = rg;
106: PetscFunctionReturn(PETSC_SUCCESS);
107: }
109: /*@
110: RGSetOptionsPrefix - Sets the prefix used for searching for all
111: `RG` options in the database.
113: Logically Collective
115: Input Parameters:
116: + rg - the region context
117: - prefix - the prefix string to prepend to all `RG` option requests
119: Notes:
120: A hyphen (-) must NOT be given at the beginning of the prefix name.
121: The first character of all runtime options is AUTOMATICALLY the
122: hyphen.
124: Level: advanced
126: .seealso: [](sec:rg), `RGAppendOptionsPrefix()`
127: @*/
128: PetscErrorCode RGSetOptionsPrefix(RG rg,const char prefix[])
129: {
130: PetscFunctionBegin;
132: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)rg,prefix));
133: PetscFunctionReturn(PETSC_SUCCESS);
134: }
136: /*@
137: RGAppendOptionsPrefix - Appends to the prefix used for searching for all
138: `RG` options in the database.
140: Logically Collective
142: Input Parameters:
143: + rg - the region context
144: - prefix - the prefix string to prepend to all `RG` option requests
146: Notes:
147: A hyphen (-) must NOT be given at the beginning of the prefix name.
148: The first character of all runtime options is AUTOMATICALLY the hyphen.
150: Level: advanced
152: .seealso: [](sec:rg), `RGSetOptionsPrefix()`
153: @*/
154: PetscErrorCode RGAppendOptionsPrefix(RG rg,const char prefix[])
155: {
156: PetscFunctionBegin;
158: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)rg,prefix));
159: PetscFunctionReturn(PETSC_SUCCESS);
160: }
162: /*@
163: RGGetOptionsPrefix - Gets the prefix used for searching for all
164: `RG` options in the database.
166: Not Collective
168: Input Parameter:
169: . rg - the region context
171: Output Parameter:
172: . prefix - pointer to the prefix string used is returned
174: Level: advanced
176: .seealso: [](sec:rg), `RGSetOptionsPrefix()`, `RGAppendOptionsPrefix()`
177: @*/
178: PetscErrorCode RGGetOptionsPrefix(RG rg,const char *prefix[])
179: {
180: PetscFunctionBegin;
182: PetscAssertPointer(prefix,2);
183: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)rg,prefix));
184: PetscFunctionReturn(PETSC_SUCCESS);
185: }
187: /*@
188: RGSetType - Selects the type for the `RG` object.
190: Logically Collective
192: Input Parameters:
193: + rg - the region context
194: - type - a known type
196: Options Database Key:
197: . -rg_type <type> - Sets `RG` type
199: Level: beginner
201: .seealso: [](sec:rg), `RGGetType()`
202: @*/
203: PetscErrorCode RGSetType(RG rg,RGType type)
204: {
205: PetscErrorCode (*r)(RG);
206: PetscBool match;
208: PetscFunctionBegin;
210: PetscAssertPointer(type,2);
212: PetscCall(PetscObjectTypeCompare((PetscObject)rg,type,&match));
213: if (match) PetscFunctionReturn(PETSC_SUCCESS);
215: PetscCall(PetscFunctionListFind(RGList,type,&r));
216: PetscCheck(r,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested RG type %s",type);
218: PetscTryTypeMethod(rg,destroy);
219: PetscCall(PetscMemzero(rg->ops,sizeof(struct _RGOps)));
221: PetscCall(PetscObjectChangeTypeName((PetscObject)rg,type));
222: PetscCall((*r)(rg));
223: PetscFunctionReturn(PETSC_SUCCESS);
224: }
226: /*@
227: RGGetType - Gets the `RG` type name (as a string) from the `RG` context.
229: Not Collective
231: Input Parameter:
232: . rg - the region context
234: Output Parameter:
235: . type - name of the region
237: Level: beginner
239: .seealso: [](sec:rg), `RGSetType()`
240: @*/
241: PetscErrorCode RGGetType(RG rg,RGType *type)
242: {
243: PetscFunctionBegin;
245: PetscAssertPointer(type,2);
246: *type = ((PetscObject)rg)->type_name;
247: PetscFunctionReturn(PETSC_SUCCESS);
248: }
250: /*@
251: RGSetFromOptions - Sets `RG` options from the options database.
253: Collective
255: Input Parameter:
256: . rg - the region context
258: Note:
259: To see all options, run your program with the `-help` option.
261: Level: beginner
263: .seealso: [](sec:rg), `RGSetOptionsPrefix()`
264: @*/
265: PetscErrorCode RGSetFromOptions(RG rg)
266: {
267: char type[256];
268: PetscBool flg;
269: PetscReal sfactor;
271: PetscFunctionBegin;
273: PetscCall(RGRegisterAll());
274: PetscObjectOptionsBegin((PetscObject)rg);
275: PetscCall(PetscOptionsFList("-rg_type","Region type","RGSetType",RGList,(char*)(((PetscObject)rg)->type_name?((PetscObject)rg)->type_name:RGINTERVAL),type,sizeof(type),&flg));
276: if (flg) PetscCall(RGSetType(rg,type));
277: else if (!((PetscObject)rg)->type_name) PetscCall(RGSetType(rg,RGINTERVAL));
279: PetscCall(PetscOptionsBool("-rg_complement","Whether region is complemented or not","RGSetComplement",rg->complement,&rg->complement,NULL));
281: PetscCall(PetscOptionsReal("-rg_scale","Scaling factor","RGSetScale",1.0,&sfactor,&flg));
282: if (flg) PetscCall(RGSetScale(rg,sfactor));
284: PetscTryTypeMethod(rg,setfromoptions,PetscOptionsObject);
285: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)rg,PetscOptionsObject));
286: PetscOptionsEnd();
287: PetscFunctionReturn(PETSC_SUCCESS);
288: }
290: /*@
291: RGView - Prints the `RG` data structure.
293: Collective
295: Input Parameters:
296: + rg - the region context
297: - viewer - optional visualization context
299: Notes:
300: The available visualization contexts include
301: + `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
302: - `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard output where only the
303: first process opens the file; all other processes send their data to the
304: first one to print
306: The user can open an alternative visualization context with `PetscViewerASCIIOpen()`
307: to output to a specified file.
309: Use `RGViewFromOptions()` to allow the user to select many different `PetscViewerType`
310: and formats from the options database.
312: Level: beginner
314: .seealso: [](sec:rg), `RGCreate()`, `RGViewFromOptions()`
315: @*/
316: PetscErrorCode RGView(RG rg,PetscViewer viewer)
317: {
318: PetscBool isdraw,isascii;
320: PetscFunctionBegin;
322: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)rg),&viewer));
324: PetscCheckSameComm(rg,1,viewer,2);
325: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
326: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
327: if (isascii) {
328: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)rg,viewer));
329: PetscCall(PetscViewerASCIIPushTab(viewer));
330: PetscTryTypeMethod(rg,view,viewer);
331: PetscCall(PetscViewerASCIIPopTab(viewer));
332: if (rg->complement) PetscCall(PetscViewerASCIIPrintf(viewer," selected region is the complement of the specified one\n"));
333: if (rg->sfactor!=1.0) PetscCall(PetscViewerASCIIPrintf(viewer," scaling factor = %g\n",(double)rg->sfactor));
334: } else if (isdraw) PetscTryTypeMethod(rg,view,viewer);
335: PetscFunctionReturn(PETSC_SUCCESS);
336: }
338: /*@
339: RGViewFromOptions - View (print) an `RG` object based on values in the options database.
341: Collective
343: Input Parameters:
344: + rg - the region context
345: . obj - optional object that provides the options prefix used to query the options database
346: - name - command line option
348: Level: intermediate
350: .seealso: [](sec:rg), `RGView()`, `RGCreate()`, `PetscObjectViewFromOptions()`
351: @*/
352: PetscErrorCode RGViewFromOptions(RG rg,PetscObject obj,const char name[])
353: {
354: PetscFunctionBegin;
356: PetscCall(PetscObjectViewFromOptions((PetscObject)rg,obj,name));
357: PetscFunctionReturn(PETSC_SUCCESS);
358: }
360: /*@
361: RGIsTrivial - Whether it is the trivial region (whole complex plane).
363: Not Collective
365: Input Parameter:
366: . rg - the region context
368: Output Parameter:
369: . trivial - true if the region is equal to the whole complex plane, e.g.,
370: an interval region with all four endpoints unbounded or an
371: ellipse with infinite radius.
373: Level: beginner
375: .seealso: [](sec:rg), `RGCheckInside()`
376: @*/
377: PetscErrorCode RGIsTrivial(RG rg,PetscBool *trivial)
378: {
379: PetscFunctionBegin;
382: PetscAssertPointer(trivial,2);
383: *trivial = PETSC_FALSE;
384: PetscTryTypeMethod(rg,istrivial,trivial);
385: PetscFunctionReturn(PETSC_SUCCESS);
386: }
388: /*@
389: RGCheckInside - Determines if a set of given points are inside the region or not.
391: Not Collective
393: Input Parameters:
394: + rg - the region context
395: . n - number of points to check
396: . ar - array of real parts
397: - ai - array of imaginary parts
399: Output Parameter:
400: . inside - array of results (1=inside, 0=on the contour, -1=outside)
402: Notes:
403: The point `a` is expressed as a couple of `PetscScalar` variables `ar`, `ai`.
404: If built with complex scalars, the point is supposed to be stored in `ar`,
405: otherwise `ar`, `ai` contain the real and imaginary parts, respectively.
407: If a scaling factor was set, the points are scaled before checking.
409: Level: intermediate
411: .seealso: [](sec:rg), `RGSetScale()`, `RGSetComplement()`
412: @*/
413: PetscErrorCode RGCheckInside(RG rg,PetscInt n,PetscScalar ar[],PetscScalar ai[],PetscInt inside[])
414: {
415: PetscReal px,py;
416: PetscInt i;
418: PetscFunctionBegin;
421: PetscAssertPointer(ar,3);
422: #if !defined(PETSC_USE_COMPLEX)
423: PetscAssertPointer(ai,4);
424: #endif
425: PetscAssertPointer(inside,5);
427: for (i=0;i<n;i++) {
428: #if defined(PETSC_USE_COMPLEX)
429: px = PetscRealPart(ar[i]);
430: py = PetscImaginaryPart(ar[i]);
431: #else
432: px = ar[i];
433: py = ai[i];
434: #endif
435: if (PetscUnlikely(rg->sfactor != 1.0)) {
436: px /= rg->sfactor;
437: py /= rg->sfactor;
438: }
439: PetscUseTypeMethod(rg,checkinside,px,py,inside+i);
440: if (PetscUnlikely(rg->complement)) inside[i] = -inside[i];
441: }
442: PetscFunctionReturn(PETSC_SUCCESS);
443: }
445: /*@
446: RGIsAxisymmetric - Determines if the region is symmetric with respect
447: to the real or imaginary axis.
449: Not Collective
451: Input Parameters:
452: + rg - the region context
453: - vertical - `PETSC_TRUE` if symmetry must be checked against the vertical axis
455: Output Parameter:
456: . symm - `PETSC_TRUE` if the region is axisymmetric
458: Note:
459: If the vertical argument is true, symmetry is checked with respect to
460: the vertical axis, otherwise with respect to the horizontal axis.
462: Level: intermediate
464: .seealso: [](sec:rg), `RGCanUseConjugates()`
465: @*/
466: PetscErrorCode RGIsAxisymmetric(RG rg,PetscBool vertical,PetscBool *symm)
467: {
468: PetscFunctionBegin;
471: PetscAssertPointer(symm,3);
472: *symm = PETSC_FALSE;
473: PetscTryTypeMethod(rg,isaxisymmetric,vertical,symm);
474: PetscFunctionReturn(PETSC_SUCCESS);
475: }
477: /*@
478: RGCanUseConjugates - Used in contour integral methods to determine whether
479: half of integration points can be avoided (use their conjugates).
481: Not Collective
483: Input Parameters:
484: + rg - the region context
485: - realmats - `PETSC_TRUE` if the problem matrices are real
487: Output Parameter:
488: . useconj - whether it is possible to use conjugates
490: Notes:
491: If some integration points are the conjugates of other points, then the
492: associated computational cost can be saved. This depends on the problem
493: matrices being real and also the region being symmetric with respect to
494: the horizontal axis. The result is false if using real arithmetic or
495: in the case of a flat region (height equal to zero).
497: Level: developer
499: .seealso: [](sec:rg), `RGIsAxisymmetric()`
500: @*/
501: PetscErrorCode RGCanUseConjugates(RG rg,PetscBool realmats,PetscBool *useconj)
502: {
503: #if defined(PETSC_USE_COMPLEX)
504: PetscReal c,d;
505: PetscBool isaxisymm;
506: #endif
508: PetscFunctionBegin;
511: PetscAssertPointer(useconj,3);
512: *useconj = PETSC_FALSE;
513: #if defined(PETSC_USE_COMPLEX)
514: if (realmats) {
515: PetscCall(RGIsAxisymmetric(rg,PETSC_FALSE,&isaxisymm));
516: if (isaxisymm) {
517: PetscCall(RGComputeBoundingBox(rg,NULL,NULL,&c,&d));
518: if (c!=d) *useconj = PETSC_TRUE;
519: }
520: }
521: #endif
522: PetscFunctionReturn(PETSC_SUCCESS);
523: }
525: /*@
526: RGComputeContour - Computes the coordinates of several points lying on the
527: contour of the region.
529: Not Collective
531: Input Parameters:
532: + rg - the region context
533: - n - number of points to compute
535: Output Parameters:
536: + cr - location to store real parts
537: - ci - location to store imaginary parts
539: Note:
540: In real scalars, either `cr` or `ci` can be `NULL` (but not both). In complex
541: scalars, the coordinates are stored in `cr`, which cannot be `NULL` (`ci` is
542: not referenced).
544: Level: intermediate
546: .seealso: [](sec:rg), `RGComputeBoundingBox()`, `RGSetScale()`
547: @*/
548: PetscErrorCode RGComputeContour(RG rg,PetscInt n,PetscScalar cr[],PetscScalar ci[])
549: {
550: PetscInt i;
552: PetscFunctionBegin;
555: #if defined(PETSC_USE_COMPLEX)
556: PetscAssertPointer(cr,3);
557: #else
558: PetscCheck(cr || ci,PetscObjectComm((PetscObject)rg),PETSC_ERR_SUP,"cr and ci cannot be NULL simultaneously");
559: #endif
560: PetscCheck(!rg->complement,PetscObjectComm((PetscObject)rg),PETSC_ERR_SUP,"Cannot compute contour of region with complement flag set");
561: PetscUseTypeMethod(rg,computecontour,n,cr,ci);
562: for (i=0;i<n;i++) {
563: if (cr) cr[i] *= rg->sfactor;
564: if (ci) ci[i] *= rg->sfactor;
565: }
566: PetscFunctionReturn(PETSC_SUCCESS);
567: }
569: /*@
570: RGComputeBoundingBox - Determines the endpoints of a rectangle in the complex plane that
571: contains the region.
573: Not Collective
575: Input Parameter:
576: . rg - the region context
578: Output Parameters:
579: + a - left endpoint of the bounding box in the real axis
580: . b - right endpoint of the bounding box in the real axis
581: . c - bottom endpoint of the bounding box in the imaginary axis
582: - d - top endpoint of the bounding box in the imaginary axis
584: Note:
585: The bounding box is defined as $[a,b]\times[c,d]$. In regions that are not bounded (e.g., an
586: open interval) or with the complement flag set, it makes no sense to compute a bounding
587: box, so the return values are infinite.
589: Level: intermediate
591: .seealso: [](sec:rg), `RGComputeContour()`, `RGSetScale()`, `RGSetComplement()`
592: @*/
593: PetscErrorCode RGComputeBoundingBox(RG rg,PetscReal *a,PetscReal *b,PetscReal *c,PetscReal *d)
594: {
595: PetscFunctionBegin;
599: if (rg->complement) { /* cannot compute bounding box */
600: if (a) *a = -PETSC_MAX_REAL;
601: if (b) *b = PETSC_MAX_REAL;
602: if (c) *c = -PETSC_MAX_REAL;
603: if (d) *d = PETSC_MAX_REAL;
604: } else {
605: PetscUseTypeMethod(rg,computebbox,a,b,c,d);
606: if (a && *a!=-PETSC_MAX_REAL) *a *= rg->sfactor;
607: if (b && *b!= PETSC_MAX_REAL) *b *= rg->sfactor;
608: if (c && *c!=-PETSC_MAX_REAL) *c *= rg->sfactor;
609: if (d && *d!= PETSC_MAX_REAL) *d *= rg->sfactor;
610: }
611: PetscFunctionReturn(PETSC_SUCCESS);
612: }
614: /*@
615: RGComputeQuadrature - Computes the values of the parameters used in a
616: quadrature rule for a contour integral around the boundary of the region.
618: Not Collective
620: Input Parameters:
621: + rg - the region context
622: . quad - the type of quadrature
623: - n - number of quadrature points to compute
625: Output Parameters:
626: + z - quadrature points
627: . zn - normalized quadrature points
628: - w - quadrature weights
630: Notes:
631: In complex scalars, the values returned in `z` are often the same as those
632: computed by `RGComputeContour()`, but this is not the case in real scalars
633: where all output arguments are real.
635: The computed values change for different quadrature rules.
637: Level: advanced
639: .seealso: [](sec:rg), `RGComputeContour()`
640: @*/
641: PetscErrorCode RGComputeQuadrature(RG rg,RGQuadRule quad,PetscInt n,PetscScalar z[],PetscScalar zn[],PetscScalar w[])
642: {
643: PetscFunctionBegin;
646: PetscAssertPointer(z,4);
647: PetscAssertPointer(zn,5);
648: PetscAssertPointer(w,6);
650: PetscCall(RGComputeContour(rg,n,z,NULL));
651: PetscUseTypeMethod(rg,computequadrature,quad,n,z,zn,w);
652: PetscFunctionReturn(PETSC_SUCCESS);
653: }
655: /*@
656: RGSetComplement - Sets a flag to indicate that the region is the complement
657: of the specified one.
659: Logically Collective
661: Input Parameters:
662: + rg - the region context
663: - flg - the boolean flag
665: Options Database Key:
666: . -rg_complement <bool> - Activate/deactivate complementing of the region
668: Level: intermediate
670: .seealso: [](sec:rg), `RGGetComplement()`
671: @*/
672: PetscErrorCode RGSetComplement(RG rg,PetscBool flg)
673: {
674: PetscFunctionBegin;
677: rg->complement = flg;
678: PetscFunctionReturn(PETSC_SUCCESS);
679: }
681: /*@
682: RGGetComplement - Gets a flag that indicates whether the region
683: is complemented or not.
685: Not Collective
687: Input Parameter:
688: . rg - the region context
690: Output Parameter:
691: . flg - the flag
693: Level: intermediate
695: .seealso: [](sec:rg), `RGSetComplement()`
696: @*/
697: PetscErrorCode RGGetComplement(RG rg,PetscBool *flg)
698: {
699: PetscFunctionBegin;
701: PetscAssertPointer(flg,2);
702: *flg = rg->complement;
703: PetscFunctionReturn(PETSC_SUCCESS);
704: }
706: /*@
707: RGSetScale - Sets the scaling factor to be used when checking that a
708: point is inside the region and when computing the contour.
710: Logically Collective
712: Input Parameters:
713: + rg - the region context
714: - sfactor - the scaling factor
716: Options Database Key:
717: . -rg_scale <real> - Sets the scaling factor
719: Level: intermediate
721: .seealso: [](sec:rg), `RGGetScale()`, `RGCheckInside()`, `RGComputeContour()`
722: @*/
723: PetscErrorCode RGSetScale(RG rg,PetscReal sfactor)
724: {
725: PetscFunctionBegin;
728: if (sfactor == (PetscReal)PETSC_DEFAULT || sfactor == (PetscReal)PETSC_DECIDE) sfactor = 1.0;
729: PetscCheck(sfactor>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of scaling factor. Must be > 0");
730: rg->sfactor = sfactor;
731: PetscFunctionReturn(PETSC_SUCCESS);
732: }
734: /*@
735: RGGetScale - Gets the scaling factor.
737: Not Collective
739: Input Parameter:
740: . rg - the region context
742: Output Parameter:
743: . sfactor - the scaling factor
745: Level: intermediate
747: .seealso: [](sec:rg), `RGSetScale()`
748: @*/
749: PetscErrorCode RGGetScale(RG rg,PetscReal *sfactor)
750: {
751: PetscFunctionBegin;
753: PetscAssertPointer(sfactor,2);
754: *sfactor = rg->sfactor;
755: PetscFunctionReturn(PETSC_SUCCESS);
756: }
758: /*@
759: RGPushScale - Sets an additional scaling factor, that will multiply the
760: user-defined scaling factor.
762: Logically Collective
764: Input Parameters:
765: + rg - the region context
766: - sfactor - the scaling factor
768: Notes:
769: The current implementation does not allow pushing several scaling factors.
771: This is intended for internal use, for instance in polynomial eigensolvers
772: that use parameter scaling.
774: Level: developer
776: .seealso: [](sec:rg), `RGPopScale()`, `RGSetScale()`
777: @*/
778: PetscErrorCode RGPushScale(RG rg,PetscReal sfactor)
779: {
780: PetscFunctionBegin;
783: PetscCheck(sfactor>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of scaling factor. Must be > 0");
784: PetscCheck(!rg->osfactor,PetscObjectComm((PetscObject)rg),PETSC_ERR_SUP,"Current implementation does not allow pushing several scaling factors");
785: rg->osfactor = rg->sfactor;
786: rg->sfactor *= sfactor;
787: PetscFunctionReturn(PETSC_SUCCESS);
788: }
790: /*@
791: RGPopScale - Pops the scaling factor set with `RGPushScale()`.
793: Logically Collective
795: Input Parameter:
796: . rg - the region context
798: Level: developer
800: .seealso: [](sec:rg), `RGPushScale()`
801: @*/
802: PetscErrorCode RGPopScale(RG rg)
803: {
804: PetscFunctionBegin;
806: PetscCheck(rg->osfactor,PetscObjectComm((PetscObject)rg),PETSC_ERR_ORDER,"Must call RGPushScale first");
807: rg->sfactor = rg->osfactor;
808: rg->osfactor = 0.0;
809: PetscFunctionReturn(PETSC_SUCCESS);
810: }
812: /*@
813: RGDestroy - Destroys `RG` context that was created with `RGCreate()`.
815: Collective
817: Input Parameter:
818: . rg - the region context
820: Level: beginner
822: .seealso: [](sec:rg), `RGCreate()`
823: @*/
824: PetscErrorCode RGDestroy(RG *rg)
825: {
826: PetscFunctionBegin;
827: if (!*rg) PetscFunctionReturn(PETSC_SUCCESS);
829: if (--((PetscObject)*rg)->refct > 0) { *rg = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
830: PetscTryTypeMethod(*rg,destroy);
831: PetscCall(PetscHeaderDestroy(rg));
832: PetscFunctionReturn(PETSC_SUCCESS);
833: }
835: /*@C
836: RGRegister - Adds a region to the `RG` package.
838: Not Collective
840: Input Parameters:
841: + name - name of a new user-defined RG
842: - function - routine to create context
844: Note:
845: `RGRegister()` may be called multiple times to add several user-defined regions.
847: Level: advanced
849: .seealso: [](sec:rg), `RGRegisterAll()`
850: @*/
851: PetscErrorCode RGRegister(const char *name,PetscErrorCode (*function)(RG))
852: {
853: PetscFunctionBegin;
854: PetscCall(RGInitializePackage());
855: PetscCall(PetscFunctionListAdd(&RGList,name,function));
856: PetscFunctionReturn(PETSC_SUCCESS);
857: }