Actual source code: rgellipse.c
slepc-3.21.0 2024-03-30
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: Region enclosed in an ellipse (aligned with the coordinate axes)
12: */
14: #include <slepc/private/rgimpl.h>
15: #include <petscdraw.h>
17: typedef struct {
18: PetscScalar center; /* center of the ellipse */
19: PetscReal radius; /* radius of the ellipse */
20: PetscReal vscale; /* vertical scale of the ellipse */
21: } RG_ELLIPSE;
23: static PetscErrorCode RGEllipseSetParameters_Ellipse(RG rg,PetscScalar center,PetscReal radius,PetscReal vscale)
24: {
25: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
27: PetscFunctionBegin;
28: ctx->center = center;
29: if (radius == (PetscReal)PETSC_DEFAULT) {
30: ctx->radius = 1.0;
31: } else {
32: PetscCheck(radius>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The radius argument must be > 0.0");
33: ctx->radius = radius;
34: }
35: PetscCheck(vscale>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The vscale argument must be > 0.0");
36: ctx->vscale = vscale;
37: PetscFunctionReturn(PETSC_SUCCESS);
38: }
40: /*@
41: RGEllipseSetParameters - Sets the parameters defining the ellipse region.
43: Logically Collective
45: Input Parameters:
46: + rg - the region context
47: . center - center of the ellipse
48: . radius - radius of the ellipse
49: - vscale - vertical scale of the ellipse
51: Options Database Keys:
52: + -rg_ellipse_center - Sets the center
53: . -rg_ellipse_radius - Sets the radius
54: - -rg_ellipse_vscale - Sets the vertical scale
56: Notes:
57: In the case of complex scalars, a complex center can be provided in the
58: command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
59: -rg_ellipse_center 1.0+2.0i
61: When PETSc is built with real scalars, the center is restricted to a real value.
63: Level: advanced
65: .seealso: RGEllipseGetParameters()
66: @*/
67: PetscErrorCode RGEllipseSetParameters(RG rg,PetscScalar center,PetscReal radius,PetscReal vscale)
68: {
69: PetscFunctionBegin;
74: PetscTryMethod(rg,"RGEllipseSetParameters_C",(RG,PetscScalar,PetscReal,PetscReal),(rg,center,radius,vscale));
75: PetscFunctionReturn(PETSC_SUCCESS);
76: }
78: static PetscErrorCode RGEllipseGetParameters_Ellipse(RG rg,PetscScalar *center,PetscReal *radius,PetscReal *vscale)
79: {
80: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
82: PetscFunctionBegin;
83: if (center) *center = ctx->center;
84: if (radius) *radius = ctx->radius;
85: if (vscale) *vscale = ctx->vscale;
86: PetscFunctionReturn(PETSC_SUCCESS);
87: }
89: /*@C
90: RGEllipseGetParameters - Gets the parameters that define the ellipse region.
92: Not Collective
94: Input Parameter:
95: . rg - the region context
97: Output Parameters:
98: + center - center of the region
99: . radius - radius of the region
100: - vscale - vertical scale of the region
102: Level: advanced
104: .seealso: RGEllipseSetParameters()
105: @*/
106: PetscErrorCode RGEllipseGetParameters(RG rg,PetscScalar *center,PetscReal *radius,PetscReal *vscale)
107: {
108: PetscFunctionBegin;
110: PetscUseMethod(rg,"RGEllipseGetParameters_C",(RG,PetscScalar*,PetscReal*,PetscReal*),(rg,center,radius,vscale));
111: PetscFunctionReturn(PETSC_SUCCESS);
112: }
114: static PetscErrorCode RGView_Ellipse(RG rg,PetscViewer viewer)
115: {
116: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
117: PetscBool isdraw,isascii;
118: int winw,winh;
119: PetscDraw draw;
120: PetscDrawAxis axis;
121: PetscReal cx,cy,r,ab,cd,lx,ly,w,scale=1.2;
122: char str[50];
124: PetscFunctionBegin;
125: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
126: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
127: if (isascii) {
128: PetscCall(SlepcSNPrintfScalar(str,sizeof(str),ctx->center,PETSC_FALSE));
129: PetscCall(PetscViewerASCIIPrintf(viewer," center: %s, radius: %g, vscale: %g\n",str,RGShowReal(ctx->radius),RGShowReal(ctx->vscale)));
130: } else if (isdraw) {
131: PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
132: PetscCall(PetscDrawCheckResizedWindow(draw));
133: PetscCall(PetscDrawGetWindowSize(draw,&winw,&winh));
134: winw = PetscMax(winw,1); winh = PetscMax(winh,1);
135: PetscCall(PetscDrawClear(draw));
136: PetscCall(PetscDrawSetTitle(draw,"Ellipse region"));
137: PetscCall(PetscDrawAxisCreate(draw,&axis));
138: cx = PetscRealPart(ctx->center)*rg->sfactor;
139: cy = PetscImaginaryPart(ctx->center)*rg->sfactor;
140: r = ctx->radius*rg->sfactor;
141: lx = 2*r;
142: ly = 2*r*ctx->vscale;
143: ab = cx;
144: cd = cy;
145: w = scale*PetscMax(lx/winw,ly/winh)/2;
146: PetscCall(PetscDrawAxisSetLimits(axis,ab-w*winw,ab+w*winw,cd-w*winh,cd+w*winh));
147: PetscCall(PetscDrawAxisDraw(axis));
148: PetscCall(PetscDrawAxisDestroy(&axis));
149: PetscCall(PetscDrawEllipse(draw,cx,cy,2*r,2*r*ctx->vscale,PETSC_DRAW_RED));
150: PetscCall(PetscDrawFlush(draw));
151: PetscCall(PetscDrawSave(draw));
152: PetscCall(PetscDrawPause(draw));
153: }
154: PetscFunctionReturn(PETSC_SUCCESS);
155: }
157: static PetscErrorCode RGIsTrivial_Ellipse(RG rg,PetscBool *trivial)
158: {
159: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
161: PetscFunctionBegin;
162: if (rg->complement) *trivial = PetscNot(ctx->radius);
163: else *trivial = PetscNot(ctx->radius<PETSC_MAX_REAL);
164: PetscFunctionReturn(PETSC_SUCCESS);
165: }
167: static PetscErrorCode RGComputeContour_Ellipse(RG rg,PetscInt n,PetscScalar *cr,PetscScalar *ci)
168: {
169: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
170: PetscReal theta;
171: PetscInt i;
173: PetscFunctionBegin;
174: for (i=0;i<n;i++) {
175: theta = 2.0*PETSC_PI*(i+0.5)/n;
176: #if defined(PETSC_USE_COMPLEX)
177: cr[i] = ctx->center + ctx->radius*PetscCMPLX(PetscCosReal(theta),ctx->vscale*PetscSinReal(theta));
178: #else
179: if (cr) cr[i] = ctx->center + ctx->radius*PetscCosReal(theta);
180: if (ci) ci[i] = ctx->radius*ctx->vscale*PetscSinReal(theta);
181: #endif
182: }
183: PetscFunctionReturn(PETSC_SUCCESS);
184: }
186: static PetscErrorCode RGComputeBoundingBox_Ellipse(RG rg,PetscReal *a,PetscReal *b,PetscReal *c,PetscReal *d)
187: {
188: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
190: PetscFunctionBegin;
191: if (a) *a = PetscRealPart(ctx->center) - ctx->radius;
192: if (b) *b = PetscRealPart(ctx->center) + ctx->radius;
193: if (c) *c = PetscImaginaryPart(ctx->center) - ctx->radius*ctx->vscale;
194: if (d) *d = PetscImaginaryPart(ctx->center) + ctx->radius*ctx->vscale;
195: PetscFunctionReturn(PETSC_SUCCESS);
196: }
198: static PetscErrorCode RGComputeQuadrature_Ellipse(RG rg,RGQuadRule quad,PetscInt n,PetscScalar *z,PetscScalar *zn,PetscScalar *w)
199: {
200: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
201: PetscReal theta;
202: PetscInt i;
204: PetscFunctionBegin;
205: for (i=0;i<n;i++) {
206: #if defined(PETSC_USE_COMPLEX)
207: theta = 2.0*PETSC_PI*(i+0.5)/n;
208: zn[i] = PetscCMPLX(PetscCosReal(theta),ctx->vscale*PetscSinReal(theta));
209: w[i] = rg->sfactor*ctx->radius*(PetscCMPLX(ctx->vscale*PetscCosReal(theta),PetscSinReal(theta)))/n;
210: #else
211: theta = PETSC_PI*(i+0.5)/n;
212: zn[i] = PetscCosReal(theta);
213: w[i] = PetscCosReal((n-1)*theta)/n;
214: z[i] = rg->sfactor*(ctx->center + ctx->radius*zn[i]);
215: #endif
216: }
217: PetscFunctionReturn(PETSC_SUCCESS);
218: }
220: static PetscErrorCode RGCheckInside_Ellipse(RG rg,PetscReal px,PetscReal py,PetscInt *inside)
221: {
222: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
223: PetscReal dx,dy,r;
225: PetscFunctionBegin;
226: #if defined(PETSC_USE_COMPLEX)
227: dx = (px-PetscRealPart(ctx->center))/ctx->radius;
228: dy = (py-PetscImaginaryPart(ctx->center))/ctx->radius;
229: #else
230: dx = (px-ctx->center)/ctx->radius;
231: dy = py/ctx->radius;
232: #endif
233: r = 1.0-dx*dx-(dy*dy)/(ctx->vscale*ctx->vscale);
234: *inside = PetscSign(r);
235: PetscFunctionReturn(PETSC_SUCCESS);
236: }
238: static PetscErrorCode RGIsAxisymmetric_Ellipse(RG rg,PetscBool vertical,PetscBool *symm)
239: {
240: RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
242: PetscFunctionBegin;
243: if (vertical) *symm = (PetscRealPart(ctx->center) == 0.0)? PETSC_TRUE: PETSC_FALSE;
244: else *symm = (PetscImaginaryPart(ctx->center) == 0.0)? PETSC_TRUE: PETSC_FALSE;
245: PetscFunctionReturn(PETSC_SUCCESS);
246: }
248: static PetscErrorCode RGSetFromOptions_Ellipse(RG rg,PetscOptionItems *PetscOptionsObject)
249: {
250: PetscScalar s;
251: PetscReal r1,r2;
252: PetscBool flg1,flg2,flg3;
254: PetscFunctionBegin;
255: PetscOptionsHeadBegin(PetscOptionsObject,"RG Ellipse Options");
257: PetscCall(RGEllipseGetParameters(rg,&s,&r1,&r2));
258: PetscCall(PetscOptionsScalar("-rg_ellipse_center","Center of ellipse","RGEllipseSetParameters",s,&s,&flg1));
259: PetscCall(PetscOptionsReal("-rg_ellipse_radius","Radius of ellipse","RGEllipseSetParameters",r1,&r1,&flg2));
260: PetscCall(PetscOptionsReal("-rg_ellipse_vscale","Vertical scale of ellipse","RGEllipseSetParameters",r2,&r2,&flg3));
261: if (flg1 || flg2 || flg3) PetscCall(RGEllipseSetParameters(rg,s,r1,r2));
263: PetscOptionsHeadEnd();
264: PetscFunctionReturn(PETSC_SUCCESS);
265: }
267: static PetscErrorCode RGDestroy_Ellipse(RG rg)
268: {
269: PetscFunctionBegin;
270: PetscCall(PetscFree(rg->data));
271: PetscCall(PetscObjectComposeFunction((PetscObject)rg,"RGEllipseSetParameters_C",NULL));
272: PetscCall(PetscObjectComposeFunction((PetscObject)rg,"RGEllipseGetParameters_C",NULL));
273: PetscFunctionReturn(PETSC_SUCCESS);
274: }
276: SLEPC_EXTERN PetscErrorCode RGCreate_Ellipse(RG rg)
277: {
278: RG_ELLIPSE *ellipse;
280: PetscFunctionBegin;
281: PetscCall(PetscNew(&ellipse));
282: ellipse->center = 0.0;
283: ellipse->radius = PETSC_MAX_REAL;
284: ellipse->vscale = 1.0;
285: rg->data = (void*)ellipse;
287: rg->ops->istrivial = RGIsTrivial_Ellipse;
288: rg->ops->computecontour = RGComputeContour_Ellipse;
289: rg->ops->computebbox = RGComputeBoundingBox_Ellipse;
290: rg->ops->computequadrature = RGComputeQuadrature_Ellipse;
291: rg->ops->checkinside = RGCheckInside_Ellipse;
292: rg->ops->isaxisymmetric = RGIsAxisymmetric_Ellipse;
293: rg->ops->setfromoptions = RGSetFromOptions_Ellipse;
294: rg->ops->view = RGView_Ellipse;
295: rg->ops->destroy = RGDestroy_Ellipse;
296: PetscCall(PetscObjectComposeFunction((PetscObject)rg,"RGEllipseSetParameters_C",RGEllipseSetParameters_Ellipse));
297: PetscCall(PetscObjectComposeFunction((PetscObject)rg,"RGEllipseGetParameters_C",RGEllipseGetParameters_Ellipse));
298: PetscFunctionReturn(PETSC_SUCCESS);
299: }