Actual source code: rgellipse.c

slepc-3.8.0 2017-10-20
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2017, 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>      /*I "slepcrg.h" I*/

 16: typedef struct {
 17:   PetscScalar center;     /* center of the ellipse */
 18:   PetscReal   radius;     /* radius of the ellipse */
 19:   PetscReal   vscale;     /* vertical scale of the ellipse */
 20: } RG_ELLIPSE;

 22: static PetscErrorCode RGEllipseSetParameters_Ellipse(RG rg,PetscScalar center,PetscReal radius,PetscReal vscale)
 23: {
 24:   RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;

 27:   ctx->center = center;
 28:   if (radius == PETSC_DEFAULT) {
 29:     ctx->radius = 1.0;
 30:   } else {
 31:     if (radius<=0.0) SETERRQ(PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The radius argument must be > 0.0");
 32:     ctx->radius = radius;
 33:   }
 34:   if (vscale<=0.0) SETERRQ(PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The vscale argument must be > 0.0");
 35:   ctx->vscale = vscale;
 36:   return(0);
 37: }

 39: /*@
 40:    RGEllipseSetParameters - Sets the parameters defining the ellipse region.

 42:    Logically Collective on RG

 44:    Input Parameters:
 45: +  rg     - the region context
 46: .  center - center of the ellipse
 47: .  radius - radius of the ellipse
 48: -  vscale - vertical scale of the ellipse

 50:    Options Database Keys:
 51: +  -rg_ellipse_center - Sets the center
 52: .  -rg_ellipse_radius - Sets the radius
 53: -  -rg_ellipse_vscale - Sets the vertical scale

 55:    Notes:
 56:    In the case of complex scalars, a complex center can be provided in the
 57:    command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
 58:    -rg_ellipse_center 1.0+2.0i

 60:    When PETSc is built with real scalars, the center is restricted to a real value.

 62:    Level: advanced

 64: .seealso: RGEllipseGetParameters()
 65: @*/
 66: PetscErrorCode RGEllipseSetParameters(RG rg,PetscScalar center,PetscReal radius,PetscReal vscale)
 67: {

 75:   PetscTryMethod(rg,"RGEllipseSetParameters_C",(RG,PetscScalar,PetscReal,PetscReal),(rg,center,radius,vscale));
 76:   return(0);
 77: }

 79: static PetscErrorCode RGEllipseGetParameters_Ellipse(RG rg,PetscScalar *center,PetscReal *radius,PetscReal *vscale)
 80: {
 81:   RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;

 84:   if (center) *center = ctx->center;
 85:   if (radius) *radius = ctx->radius;
 86:   if (vscale) *vscale = ctx->vscale;
 87:   return(0);
 88: }

 90: /*@
 91:    RGEllipseGetParameters - Gets the parameters that define the ellipse region.

 93:    Not Collective

 95:    Input Parameter:
 96: .  rg     - the region context

 98:    Output Parameters:
 99: +  center - center of the region
100: .  radius - radius of the region
101: -  vscale - vertical scale of the region

103:    Level: advanced

105: .seealso: RGEllipseSetParameters()
106: @*/
107: PetscErrorCode RGEllipseGetParameters(RG rg,PetscScalar *center,PetscReal *radius,PetscReal *vscale)
108: {

113:   PetscUseMethod(rg,"RGEllipseGetParameters_C",(RG,PetscScalar*,PetscReal*,PetscReal*),(rg,center,radius,vscale));
114:   return(0);
115: }

117: PetscErrorCode RGView_Ellipse(RG rg,PetscViewer viewer)
118: {
120:   RG_ELLIPSE     *ctx = (RG_ELLIPSE*)rg->data;
121:   PetscBool      isascii;
122:   char           str[50];

125:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
126:   if (isascii) {
127:     SlepcSNPrintfScalar(str,50,ctx->center,PETSC_FALSE);
128:     PetscViewerASCIIPrintf(viewer,"  center: %s, radius: %g, vscale: %g\n",str,RGShowReal(ctx->radius),RGShowReal(ctx->vscale));
129:   }
130:   return(0);
131: }

133: PetscErrorCode RGIsTrivial_Ellipse(RG rg,PetscBool *trivial)
134: {
135:   RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;

138:   if (rg->complement) *trivial = PetscNot(ctx->radius);
139:   else *trivial = PetscNot(ctx->radius<PETSC_MAX_REAL);
140:   return(0);
141: }

143: PetscErrorCode RGComputeContour_Ellipse(RG rg,PetscInt n,PetscScalar *cr,PetscScalar *ci)
144: {
145:   RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
146:   PetscReal  theta;
147:   PetscInt   i;

150:   for (i=0;i<n;i++) {
151:     theta = 2.0*PETSC_PI*(i+0.5)/n;
152: #if defined(PETSC_USE_COMPLEX)
153:     cr[i] = ctx->center + ctx->radius*(PetscCosReal(theta)+ctx->vscale*PetscSinReal(theta)*PETSC_i);
154: #else
155:     cr[i] = ctx->center + ctx->radius*PetscCosReal(theta);
156:     ci[i] = ctx->radius*ctx->vscale*PetscSinReal(theta);
157: #endif
158:   }
159:   return(0);
160: }

162: PetscErrorCode RGComputeBoundingBox_Ellipse(RG rg,PetscReal *a,PetscReal *b,PetscReal *c,PetscReal *d)
163: {
164:   RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;

167:   if (a) *a = PetscRealPart(ctx->center) - ctx->radius;
168:   if (b) *b = PetscRealPart(ctx->center) + ctx->radius;
169:   if (c) *c = PetscImaginaryPart(ctx->center) - ctx->radius*ctx->vscale;
170:   if (d) *d = PetscImaginaryPart(ctx->center) + ctx->radius*ctx->vscale;
171:   return(0);
172: }

174: PetscErrorCode RGCheckInside_Ellipse(RG rg,PetscReal px,PetscReal py,PetscInt *inside)
175: {
176:   RG_ELLIPSE *ctx = (RG_ELLIPSE*)rg->data;
177:   PetscReal  dx,dy,r;

180: #if defined(PETSC_USE_COMPLEX)
181:   dx = (px-PetscRealPart(ctx->center))/ctx->radius;
182:   dy = (py-PetscImaginaryPart(ctx->center))/ctx->radius;
183: #else
184:   dx = (px-ctx->center)/ctx->radius;
185:   dy = py/ctx->radius;
186: #endif
187:   r = 1.0-dx*dx-(dy*dy)/(ctx->vscale*ctx->vscale);
188:   *inside = PetscSign(r);
189:   return(0);
190: }

192: PetscErrorCode RGSetFromOptions_Ellipse(PetscOptionItems *PetscOptionsObject,RG rg)
193: {
195:   PetscScalar    s;
196:   PetscReal      r1,r2;
197:   PetscBool      flg1,flg2,flg3;

200:   PetscOptionsHead(PetscOptionsObject,"RG Ellipse Options");

202:     RGEllipseGetParameters(rg,&s,&r1,&r2);
203:     PetscOptionsScalar("-rg_ellipse_center","Center of ellipse","RGEllipseSetParameters",s,&s,&flg1);
204:     PetscOptionsReal("-rg_ellipse_radius","Radius of ellipse","RGEllipseSetParameters",r1,&r1,&flg2);
205:     PetscOptionsReal("-rg_ellipse_vscale","Vertical scale of ellipse","RGEllipseSetParameters",r2,&r2,&flg3);
206:     if (flg1 || flg2 || flg3) { RGEllipseSetParameters(rg,s,r1,r2); }

208:   PetscOptionsTail();
209:   return(0);
210: }

212: PetscErrorCode RGDestroy_Ellipse(RG rg)
213: {

217:   PetscFree(rg->data);
218:   PetscObjectComposeFunction((PetscObject)rg,"RGEllipseSetParameters_C",NULL);
219:   PetscObjectComposeFunction((PetscObject)rg,"RGEllipseGetParameters_C",NULL);
220:   return(0);
221: }

223: PETSC_EXTERN PetscErrorCode RGCreate_Ellipse(RG rg)
224: {
225:   RG_ELLIPSE     *ellipse;

229:   PetscNewLog(rg,&ellipse);
230:   ellipse->center = 0.0;
231:   ellipse->radius = PETSC_MAX_REAL;
232:   ellipse->vscale = 1.0;
233:   rg->data = (void*)ellipse;

235:   rg->ops->istrivial      = RGIsTrivial_Ellipse;
236:   rg->ops->computecontour = RGComputeContour_Ellipse;
237:   rg->ops->computebbox    = RGComputeBoundingBox_Ellipse;
238:   rg->ops->checkinside    = RGCheckInside_Ellipse;
239:   rg->ops->setfromoptions = RGSetFromOptions_Ellipse;
240:   rg->ops->view           = RGView_Ellipse;
241:   rg->ops->destroy        = RGDestroy_Ellipse;
242:   PetscObjectComposeFunction((PetscObject)rg,"RGEllipseSetParameters_C",RGEllipseSetParameters_Ellipse);
243:   PetscObjectComposeFunction((PetscObject)rg,"RGEllipseGetParameters_C",RGEllipseGetParameters_Ellipse);
244:   return(0);
245: }