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: Note:
238: `type` should not be retained for later use as it will be an invalid pointer
239: if the `RGType` of `rg` is changed.
241: Level: beginner
243: .seealso: [](sec:rg), `RGSetType()`, `PetscObjectTypeCompare()`, `PetscObjectTypeCompareAny()`
244: @*/
245: PetscErrorCode RGGetType(RG rg,RGType *type)
246: {
247: PetscFunctionBegin;
249: PetscAssertPointer(type,2);
250: *type = ((PetscObject)rg)->type_name;
251: PetscFunctionReturn(PETSC_SUCCESS);
252: }
254: /*@
255: RGSetFromOptions - Sets `RG` options from the options database.
257: Collective
259: Input Parameter:
260: . rg - the region context
262: Note:
263: To see all options, run your program with the `-help` option.
265: Level: beginner
267: .seealso: [](sec:rg), `RGSetOptionsPrefix()`
268: @*/
269: PetscErrorCode RGSetFromOptions(RG rg)
270: {
271: char type[256];
272: PetscBool flg;
273: PetscReal sfactor;
275: PetscFunctionBegin;
277: PetscCall(RGRegisterAll());
278: PetscObjectOptionsBegin((PetscObject)rg);
279: PetscCall(PetscOptionsFList("-rg_type","Region type","RGSetType",RGList,(char*)(((PetscObject)rg)->type_name?((PetscObject)rg)->type_name:RGINTERVAL),type,sizeof(type),&flg));
280: if (flg) PetscCall(RGSetType(rg,type));
281: else if (!((PetscObject)rg)->type_name) PetscCall(RGSetType(rg,RGINTERVAL));
283: PetscCall(PetscOptionsBool("-rg_complement","Whether region is complemented or not","RGSetComplement",rg->complement,&rg->complement,NULL));
285: PetscCall(PetscOptionsReal("-rg_scale","Scaling factor","RGSetScale",1.0,&sfactor,&flg));
286: if (flg) PetscCall(RGSetScale(rg,sfactor));
288: PetscTryTypeMethod(rg,setfromoptions,PetscOptionsObject);
289: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)rg,PetscOptionsObject));
290: PetscOptionsEnd();
291: PetscFunctionReturn(PETSC_SUCCESS);
292: }
294: /*@
295: RGView - Prints the `RG` data structure.
297: Collective
299: Input Parameters:
300: + rg - the region context
301: - viewer - optional visualization context
303: Notes:
304: The available visualization contexts include
305: + `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
306: - `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard output where only the
307: first process opens the file; all other processes send their data to the
308: first one to print
310: The user can open an alternative visualization context with `PetscViewerASCIIOpen()`
311: to output to a specified file.
313: Use `RGViewFromOptions()` to allow the user to select many different `PetscViewerType`
314: and formats from the options database.
316: Level: beginner
318: .seealso: [](sec:rg), `RGCreate()`, `RGViewFromOptions()`
319: @*/
320: PetscErrorCode RGView(RG rg,PetscViewer viewer)
321: {
322: PetscBool isdraw,isascii;
324: PetscFunctionBegin;
326: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)rg),&viewer));
328: PetscCheckSameComm(rg,1,viewer,2);
329: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
330: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
331: if (isascii) {
332: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)rg,viewer));
333: PetscCall(PetscViewerASCIIPushTab(viewer));
334: PetscTryTypeMethod(rg,view,viewer);
335: PetscCall(PetscViewerASCIIPopTab(viewer));
336: if (rg->complement) PetscCall(PetscViewerASCIIPrintf(viewer," selected region is the complement of the specified one\n"));
337: if (rg->sfactor!=1.0) PetscCall(PetscViewerASCIIPrintf(viewer," scaling factor = %g\n",(double)rg->sfactor));
338: } else if (isdraw) PetscTryTypeMethod(rg,view,viewer);
339: PetscFunctionReturn(PETSC_SUCCESS);
340: }
342: /*@
343: RGViewFromOptions - View (print) an `RG` object based on values in the options database.
345: Collective
347: Input Parameters:
348: + rg - the region context
349: . obj - optional object that provides the options prefix used to query the options database
350: - name - command line option
352: Level: intermediate
354: .seealso: [](sec:rg), `RGView()`, `RGCreate()`, `PetscObjectViewFromOptions()`
355: @*/
356: PetscErrorCode RGViewFromOptions(RG rg,PetscObject obj,const char name[])
357: {
358: PetscFunctionBegin;
360: PetscCall(PetscObjectViewFromOptions((PetscObject)rg,obj,name));
361: PetscFunctionReturn(PETSC_SUCCESS);
362: }
364: /*@
365: RGIsTrivial - Whether it is the trivial region (whole complex plane).
367: Not Collective
369: Input Parameter:
370: . rg - the region context
372: Output Parameter:
373: . trivial - true if the region is equal to the whole complex plane, e.g.,
374: an interval region with all four endpoints unbounded or an
375: ellipse with infinite radius.
377: Level: beginner
379: .seealso: [](sec:rg), `RGCheckInside()`
380: @*/
381: PetscErrorCode RGIsTrivial(RG rg,PetscBool *trivial)
382: {
383: PetscFunctionBegin;
386: PetscAssertPointer(trivial,2);
387: *trivial = PETSC_FALSE;
388: PetscTryTypeMethod(rg,istrivial,trivial);
389: PetscFunctionReturn(PETSC_SUCCESS);
390: }
392: /*@
393: RGCheckInside - Determines if a set of given points are inside the region or not.
395: Not Collective
397: Input Parameters:
398: + rg - the region context
399: . n - number of points to check
400: . ar - array of real parts
401: - ai - array of imaginary parts
403: Output Parameter:
404: . inside - array of results (1=inside, 0=on the contour, -1=outside)
406: Notes:
407: The point `a` is expressed as a couple of `PetscScalar` variables `ar`, `ai`.
408: If built with complex scalars, the point is supposed to be stored in `ar`,
409: otherwise `ar`, `ai` contain the real and imaginary parts, respectively.
411: If a scaling factor was set, the points are scaled before checking.
413: Level: intermediate
415: .seealso: [](sec:rg), `RGSetScale()`, `RGSetComplement()`
416: @*/
417: PetscErrorCode RGCheckInside(RG rg,PetscInt n,PetscScalar ar[],PetscScalar ai[],PetscInt inside[])
418: {
419: PetscReal px,py;
420: PetscInt i;
422: PetscFunctionBegin;
425: PetscAssertPointer(ar,3);
426: #if !defined(PETSC_USE_COMPLEX)
427: PetscAssertPointer(ai,4);
428: #endif
429: PetscAssertPointer(inside,5);
431: for (i=0;i<n;i++) {
432: #if defined(PETSC_USE_COMPLEX)
433: px = PetscRealPart(ar[i]);
434: py = PetscImaginaryPart(ar[i]);
435: #else
436: px = ar[i];
437: py = ai[i];
438: #endif
439: if (PetscUnlikely(rg->sfactor != 1.0)) {
440: px /= rg->sfactor;
441: py /= rg->sfactor;
442: }
443: PetscUseTypeMethod(rg,checkinside,px,py,inside+i);
444: if (PetscUnlikely(rg->complement)) inside[i] = -inside[i];
445: }
446: PetscFunctionReturn(PETSC_SUCCESS);
447: }
449: /*@
450: RGIsAxisymmetric - Determines if the region is symmetric with respect
451: to the real or imaginary axis.
453: Not Collective
455: Input Parameters:
456: + rg - the region context
457: - vertical - `PETSC_TRUE` if symmetry must be checked against the vertical axis
459: Output Parameter:
460: . symm - `PETSC_TRUE` if the region is axisymmetric
462: Note:
463: If the vertical argument is true, symmetry is checked with respect to
464: the vertical axis, otherwise with respect to the horizontal axis.
466: Level: intermediate
468: .seealso: [](sec:rg), `RGCanUseConjugates()`
469: @*/
470: PetscErrorCode RGIsAxisymmetric(RG rg,PetscBool vertical,PetscBool *symm)
471: {
472: PetscFunctionBegin;
475: PetscAssertPointer(symm,3);
476: *symm = PETSC_FALSE;
477: PetscTryTypeMethod(rg,isaxisymmetric,vertical,symm);
478: PetscFunctionReturn(PETSC_SUCCESS);
479: }
481: /*@
482: RGCanUseConjugates - Used in contour integral methods to determine whether
483: half of integration points can be avoided (use their conjugates).
485: Not Collective
487: Input Parameters:
488: + rg - the region context
489: - realmats - `PETSC_TRUE` if the problem matrices are real
491: Output Parameter:
492: . useconj - whether it is possible to use conjugates
494: Notes:
495: If some integration points are the conjugates of other points, then the
496: associated computational cost can be saved. This depends on the problem
497: matrices being real and also the region being symmetric with respect to
498: the horizontal axis. The result is false if using real arithmetic or
499: in the case of a flat region (height equal to zero).
501: Level: developer
503: .seealso: [](sec:rg), `RGIsAxisymmetric()`
504: @*/
505: PetscErrorCode RGCanUseConjugates(RG rg,PetscBool realmats,PetscBool *useconj)
506: {
507: #if defined(PETSC_USE_COMPLEX)
508: PetscReal c,d;
509: PetscBool isaxisymm;
510: #endif
512: PetscFunctionBegin;
515: PetscAssertPointer(useconj,3);
516: *useconj = PETSC_FALSE;
517: #if defined(PETSC_USE_COMPLEX)
518: if (realmats) {
519: PetscCall(RGIsAxisymmetric(rg,PETSC_FALSE,&isaxisymm));
520: if (isaxisymm) {
521: PetscCall(RGComputeBoundingBox(rg,NULL,NULL,&c,&d));
522: if (c!=d) *useconj = PETSC_TRUE;
523: }
524: }
525: #endif
526: PetscFunctionReturn(PETSC_SUCCESS);
527: }
529: /*@
530: RGComputeContour - Computes the coordinates of several points lying on the
531: contour of the region.
533: Not Collective
535: Input Parameters:
536: + rg - the region context
537: - n - number of points to compute
539: Output Parameters:
540: + cr - location to store real parts
541: - ci - location to store imaginary parts
543: Note:
544: In real scalars, either `cr` or `ci` can be `NULL` (but not both). In complex
545: scalars, the coordinates are stored in `cr`, which cannot be `NULL` (`ci` is
546: not referenced).
548: Level: intermediate
550: .seealso: [](sec:rg), `RGComputeBoundingBox()`, `RGSetScale()`
551: @*/
552: PetscErrorCode RGComputeContour(RG rg,PetscInt n,PetscScalar cr[],PetscScalar ci[])
553: {
554: PetscInt i;
556: PetscFunctionBegin;
559: #if defined(PETSC_USE_COMPLEX)
560: PetscAssertPointer(cr,3);
561: #else
562: PetscCheck(cr || ci,PetscObjectComm((PetscObject)rg),PETSC_ERR_SUP,"cr and ci cannot be NULL simultaneously");
563: #endif
564: PetscCheck(!rg->complement,PetscObjectComm((PetscObject)rg),PETSC_ERR_SUP,"Cannot compute contour of region with complement flag set");
565: PetscUseTypeMethod(rg,computecontour,n,cr,ci);
566: for (i=0;i<n;i++) {
567: if (cr) cr[i] *= rg->sfactor;
568: if (ci) ci[i] *= rg->sfactor;
569: }
570: PetscFunctionReturn(PETSC_SUCCESS);
571: }
573: /*@
574: RGComputeBoundingBox - Determines the endpoints of a rectangle in the complex plane that
575: contains the region.
577: Not Collective
579: Input Parameter:
580: . rg - the region context
582: Output Parameters:
583: + a - left endpoint of the bounding box in the real axis
584: . b - right endpoint of the bounding box in the real axis
585: . c - bottom endpoint of the bounding box in the imaginary axis
586: - d - top endpoint of the bounding box in the imaginary axis
588: Note:
589: The bounding box is defined as $[a,b]\times[c,d]$. In regions that are not bounded (e.g., an
590: open interval) or with the complement flag set, it makes no sense to compute a bounding
591: box, so the return values are infinite.
593: Level: intermediate
595: .seealso: [](sec:rg), `RGComputeContour()`, `RGSetScale()`, `RGSetComplement()`
596: @*/
597: PetscErrorCode RGComputeBoundingBox(RG rg,PetscReal *a,PetscReal *b,PetscReal *c,PetscReal *d)
598: {
599: PetscFunctionBegin;
603: if (rg->complement) { /* cannot compute bounding box */
604: if (a) *a = -PETSC_MAX_REAL;
605: if (b) *b = PETSC_MAX_REAL;
606: if (c) *c = -PETSC_MAX_REAL;
607: if (d) *d = PETSC_MAX_REAL;
608: } else {
609: PetscUseTypeMethod(rg,computebbox,a,b,c,d);
610: if (a && *a!=-PETSC_MAX_REAL) *a *= rg->sfactor;
611: if (b && *b!= PETSC_MAX_REAL) *b *= rg->sfactor;
612: if (c && *c!=-PETSC_MAX_REAL) *c *= rg->sfactor;
613: if (d && *d!= PETSC_MAX_REAL) *d *= rg->sfactor;
614: }
615: PetscFunctionReturn(PETSC_SUCCESS);
616: }
618: /*@
619: RGComputeQuadrature - Computes the values of the parameters used in a
620: quadrature rule for a contour integral around the boundary of the region.
622: Not Collective
624: Input Parameters:
625: + rg - the region context
626: . quad - the type of quadrature
627: - n - number of quadrature points to compute
629: Output Parameters:
630: + z - quadrature points
631: . zn - normalized quadrature points
632: - w - quadrature weights
634: Notes:
635: In complex scalars, the values returned in `z` are often the same as those
636: computed by `RGComputeContour()`, but this is not the case in real scalars
637: where all output arguments are real.
639: The computed values change for different quadrature rules.
641: Level: advanced
643: .seealso: [](sec:rg), `RGComputeContour()`
644: @*/
645: PetscErrorCode RGComputeQuadrature(RG rg,RGQuadRule quad,PetscInt n,PetscScalar z[],PetscScalar zn[],PetscScalar w[])
646: {
647: PetscFunctionBegin;
650: PetscAssertPointer(z,4);
651: PetscAssertPointer(zn,5);
652: PetscAssertPointer(w,6);
654: PetscCall(RGComputeContour(rg,n,z,NULL));
655: PetscUseTypeMethod(rg,computequadrature,quad,n,z,zn,w);
656: PetscFunctionReturn(PETSC_SUCCESS);
657: }
659: /*@
660: RGSetComplement - Sets a flag to indicate that the region is the complement
661: of the specified one.
663: Logically Collective
665: Input Parameters:
666: + rg - the region context
667: - flg - the boolean flag
669: Options Database Key:
670: . -rg_complement (true|false) - activate/deactivate complementing of the region
672: Level: intermediate
674: .seealso: [](sec:rg), `RGGetComplement()`
675: @*/
676: PetscErrorCode RGSetComplement(RG rg,PetscBool flg)
677: {
678: PetscFunctionBegin;
681: rg->complement = flg;
682: PetscFunctionReturn(PETSC_SUCCESS);
683: }
685: /*@
686: RGGetComplement - Gets a flag that indicates whether the region
687: is complemented or not.
689: Not Collective
691: Input Parameter:
692: . rg - the region context
694: Output Parameter:
695: . flg - the flag
697: Level: intermediate
699: .seealso: [](sec:rg), `RGSetComplement()`
700: @*/
701: PetscErrorCode RGGetComplement(RG rg,PetscBool *flg)
702: {
703: PetscFunctionBegin;
705: PetscAssertPointer(flg,2);
706: *flg = rg->complement;
707: PetscFunctionReturn(PETSC_SUCCESS);
708: }
710: /*@
711: RGSetScale - Sets the scaling factor to be used when checking that a
712: point is inside the region and when computing the contour.
714: Logically Collective
716: Input Parameters:
717: + rg - the region context
718: - sfactor - the scaling factor
720: Options Database Key:
721: . -rg_scale sfactor - sets the scaling factor
723: Level: intermediate
725: .seealso: [](sec:rg), `RGGetScale()`, `RGCheckInside()`, `RGComputeContour()`
726: @*/
727: PetscErrorCode RGSetScale(RG rg,PetscReal sfactor)
728: {
729: PetscFunctionBegin;
732: if (sfactor == (PetscReal)PETSC_DEFAULT || sfactor == (PetscReal)PETSC_DECIDE) sfactor = 1.0;
733: PetscCheck(sfactor>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of scaling factor. Must be > 0");
734: rg->sfactor = sfactor;
735: PetscFunctionReturn(PETSC_SUCCESS);
736: }
738: /*@
739: RGGetScale - Gets the scaling factor.
741: Not Collective
743: Input Parameter:
744: . rg - the region context
746: Output Parameter:
747: . sfactor - the scaling factor
749: Level: intermediate
751: .seealso: [](sec:rg), `RGSetScale()`
752: @*/
753: PetscErrorCode RGGetScale(RG rg,PetscReal *sfactor)
754: {
755: PetscFunctionBegin;
757: PetscAssertPointer(sfactor,2);
758: *sfactor = rg->sfactor;
759: PetscFunctionReturn(PETSC_SUCCESS);
760: }
762: /*@
763: RGPushScale - Sets an additional scaling factor, that will multiply the
764: user-defined scaling factor.
766: Logically Collective
768: Input Parameters:
769: + rg - the region context
770: - sfactor - the scaling factor
772: Notes:
773: The current implementation does not allow pushing several scaling factors.
775: This is intended for internal use, for instance in polynomial eigensolvers
776: that use parameter scaling.
778: Level: developer
780: .seealso: [](sec:rg), `RGPopScale()`, `RGSetScale()`
781: @*/
782: PetscErrorCode RGPushScale(RG rg,PetscReal sfactor)
783: {
784: PetscFunctionBegin;
787: PetscCheck(sfactor>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of scaling factor. Must be > 0");
788: PetscCheck(!rg->osfactor,PetscObjectComm((PetscObject)rg),PETSC_ERR_SUP,"Current implementation does not allow pushing several scaling factors");
789: rg->osfactor = rg->sfactor;
790: rg->sfactor *= sfactor;
791: PetscFunctionReturn(PETSC_SUCCESS);
792: }
794: /*@
795: RGPopScale - Pops the scaling factor set with `RGPushScale()`.
797: Logically Collective
799: Input Parameter:
800: . rg - the region context
802: Level: developer
804: .seealso: [](sec:rg), `RGPushScale()`
805: @*/
806: PetscErrorCode RGPopScale(RG rg)
807: {
808: PetscFunctionBegin;
810: PetscCheck(rg->osfactor,PetscObjectComm((PetscObject)rg),PETSC_ERR_ORDER,"Must call RGPushScale first");
811: rg->sfactor = rg->osfactor;
812: rg->osfactor = 0.0;
813: PetscFunctionReturn(PETSC_SUCCESS);
814: }
816: /*@
817: RGDestroy - Destroys `RG` context that was created with `RGCreate()`.
819: Collective
821: Input Parameter:
822: . rg - the region context
824: Level: beginner
826: .seealso: [](sec:rg), `RGCreate()`
827: @*/
828: PetscErrorCode RGDestroy(RG *rg)
829: {
830: PetscFunctionBegin;
831: if (!*rg) PetscFunctionReturn(PETSC_SUCCESS);
833: if (--((PetscObject)*rg)->refct > 0) { *rg = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
834: PetscTryTypeMethod(*rg,destroy);
835: PetscCall(PetscHeaderDestroy(rg));
836: PetscFunctionReturn(PETSC_SUCCESS);
837: }
839: /*@C
840: RGRegister - Adds a region to the `RG` package.
842: Not Collective
844: Input Parameters:
845: + name - name of a new user-defined RG
846: - function - routine to create context
848: Note:
849: `RGRegister()` may be called multiple times to add several user-defined regions.
851: Level: advanced
853: .seealso: [](sec:rg), `RGRegisterAll()`
854: @*/
855: PetscErrorCode RGRegister(const char *name,PetscErrorCode (*function)(RG))
856: {
857: PetscFunctionBegin;
858: PetscCall(RGInitializePackage());
859: PetscCall(PetscFunctionListAdd(&RGList,name,function));
860: PetscFunctionReturn(PETSC_SUCCESS);
861: }