Actual source code: rgellipse.c
slepc-3.5.0 2014-07-29
1: /*
2: Region enclosed in an ellipse (aligned with the coordinate axes).
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2014, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include <slepc-private/rgimpl.h> /*I "slepcrg.h" I*/
26: typedef struct {
27: PetscScalar center; /* center of the ellipse */
28: PetscReal radius; /* radius of the ellipse */
29: PetscReal vscale; /* vertical scale of the ellipse */
30: } RG_ELLIPSE;
34: static PetscErrorCode RGEllipseSetParameters_Ellipse(RG rg,PetscScalar center,PetscReal radius,PetscReal vscale)
35: {
36: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
39: ctx->center = center;
40: if (radius == PETSC_DEFAULT) {
41: ctx->radius = 1.0;
42: } else {
43: if (radius<=0.0) SETERRQ(PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The radius argument must be > 0.0");
44: ctx->radius = radius;
45: }
46: if (vscale<=0.0) SETERRQ(PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The vscale argument must be > 0.0");
47: ctx->vscale = vscale;
48: return(0);
49: }
53: /*@
54: RGEllipseSetParameters - Sets the parameters defining the ellipse region.
56: Logically Collective on RG
58: Input Parameters:
59: + rg - the region context
60: . center - center of the ellipse
61: . radius - radius of the ellipse
62: - vscale - vertical scale of the ellipse
64: Options Database Keys:
65: + -rg_ellipse_center - Sets the center
66: . -rg_ellipse_radius - Sets the radius
67: - -rg_ellipse_vscale - Sets the vertical scale
69: Note:
70: When PETSc is built with real scalars, the center is restricted to a real value.
72: Level: advanced
74: .seealso: RGEllipseGetParameters()
75: @*/
76: PetscErrorCode RGEllipseSetParameters(RG rg,PetscScalar center,PetscReal radius,PetscReal vscale)
77: {
85: PetscTryMethod(rg,"RGEllipseSetParameters_C",(RG,PetscScalar,PetscReal,PetscReal),(rg,center,radius,vscale));
86: return(0);
87: }
91: static PetscErrorCode RGEllipseGetParameters_Ellipse(RG rg,PetscScalar *center,PetscReal *radius,PetscReal *vscale)
92: {
93: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
96: if (center) *center = ctx->center;
97: if (radius) *radius = ctx->radius;
98: if (vscale) *vscale = ctx->vscale;
99: return(0);
100: }
104: /*@
105: RGEllipseGetParameters - Gets the parameters that define the ellipse region.
107: Not Collective
109: Input Parameter:
110: . rg - the region context
112: Output Parameters:
113: + center - center of the region
114: . radius - radius of the region
115: - vscale - vertical scale of the region
117: Level: advanced
119: .seealso: RGEllipseSetParameters()
120: @*/
121: PetscErrorCode RGEllipseGetParameters(RG rg,PetscScalar *center,PetscReal *radius,PetscReal *vscale)
122: {
127: PetscTryMethod(rg,"RGEllipseGetParameters_C",(RG,PetscScalar*,PetscReal*,PetscReal*),(rg,center,radius,vscale));
128: return(0);
129: }
133: PetscErrorCode RGView_Ellipse(RG rg,PetscViewer viewer)
134: {
136: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
137: PetscBool isascii;
138: char str[50];
141: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
142: if (isascii) {
143: SlepcSNPrintfScalar(str,50,ctx->center,PETSC_FALSE);
144: PetscViewerASCIIPrintf(viewer,"center: %s, radius: %g, vscale: %g\n",str,RGShowReal(ctx->radius),RGShowReal(ctx->vscale));
145: }
146: return(0);
147: }
151: PetscErrorCode RGIsTrivial_Ellipse(RG rg,PetscBool *trivial)
152: {
153: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
156: if (rg->complement) *trivial = (ctx->radius==0.0)? PETSC_TRUE: PETSC_FALSE;
157: else *trivial = (ctx->radius>=PETSC_MAX_REAL)? PETSC_TRUE: PETSC_FALSE;
158: return(0);
159: }
163: PetscErrorCode RGComputeContour_Ellipse(RG rg,PetscInt n,PetscScalar *cr,PetscScalar *ci)
164: {
165: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
166: PetscReal theta;
167: PetscInt i;
170: for (i=0;i<n;i++) {
171: theta = 2.0*PETSC_PI*(i+0.5)/n;
172: #if defined(PETSC_USE_COMPLEX)
173: cr[i] = ctx->center + ctx->radius*(PetscCosReal(theta)+ctx->vscale*PetscSinReal(theta)*PETSC_i);
174: #else
175: cr[i] = ctx->center + ctx->radius*PetscCosReal(theta);
176: ci[i] = ctx->radius*ctx->vscale*PetscSinReal(theta);
177: #endif
178: }
179: return(0);
180: }
184: PetscErrorCode RGCheckInside_Ellipse(RG rg,PetscInt n,PetscScalar *ar,PetscScalar *ai,PetscInt *inside)
185: {
186: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
187: PetscInt i;
188: PetscReal dx,dy,r;
189: #if defined(PETSC_USE_COMPLEX)
190: PetscScalar d;
191: #endif
194: for (i=0;i<n;i++) {
195: #if defined(PETSC_USE_COMPLEX)
196: d = (ar[i]-ctx->center)/ctx->radius;
197: dx = PetscRealPart(d);
198: dy = PetscImaginaryPart(d);
199: #else
200: dx = (ar[i]-ctx->center)/ctx->radius;
201: dy = ai[i]/ctx->radius;
202: #endif
203: r = 1.0-dx*dx-(dy*dy)/(ctx->vscale*ctx->vscale);
204: inside[i] = PetscSign(r);
205: }
206: return(0);
207: }
211: PetscErrorCode RGSetFromOptions_Ellipse(RG rg)
212: {
214: PetscScalar s;
215: PetscReal r1,r2;
216: PetscBool flg1,flg2,flg3;
219: PetscOptionsHead("RG Ellipse Options");
221: RGEllipseGetParameters(rg,&s,&r1,&r2);
222: PetscOptionsScalar("-rg_ellipse_center","Center of ellipse","RGEllipseSetParameters",s,&s,&flg1);
223: PetscOptionsReal("-rg_ellipse_radius","Radius of ellipse","RGEllipseSetParameters",r1,&r1,&flg2);
224: PetscOptionsReal("-rg_ellipse_vscale","Vertical scale of ellipse","RGEllipseSetParameters",r2,&r2,&flg3);
225: if (flg1 || flg2 || flg3) {
226: RGEllipseSetParameters(rg,s,r1,r2);
227: }
229: PetscOptionsTail();
230: return(0);
231: }
235: PetscErrorCode RGDestroy_Ellipse(RG rg)
236: {
240: PetscFree(rg->data);
241: PetscObjectComposeFunction((PetscObject)rg,"RGEllipseSetParameters_C",NULL);
242: PetscObjectComposeFunction((PetscObject)rg,"RGEllipseGetParameters_C",NULL);
243: return(0);
244: }
248: PETSC_EXTERN PetscErrorCode RGCreate_Ellipse(RG rg)
249: {
250: RG_ELLIPSE *ellipse;
254: PetscNewLog(rg,&ellipse);
255: ellipse->center = 0.0;
256: ellipse->radius = 1.0;
257: ellipse->vscale = 1.0;
258: rg->data = (void*)ellipse;
260: rg->ops->istrivial = RGIsTrivial_Ellipse;
261: rg->ops->computecontour = RGComputeContour_Ellipse;
262: rg->ops->checkinside = RGCheckInside_Ellipse;
263: rg->ops->setfromoptions = RGSetFromOptions_Ellipse;
264: rg->ops->view = RGView_Ellipse;
265: rg->ops->destroy = RGDestroy_Ellipse;
266: PetscObjectComposeFunction((PetscObject)rg,"RGEllipseSetParameters_C",RGEllipseSetParameters_Ellipse);
267: PetscObjectComposeFunction((PetscObject)rg,"RGEllipseGetParameters_C",RGEllipseGetParameters_Ellipse);
268: return(0);
269: }