Actual source code: rgring.c

slepc-3.22.1 2024-10-28
Report Typos and Errors
  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:    Ring region, similar to the ellipse but with a start and end angle,
 12:    together with the width
 13: */

 15: #include <slepc/private/rgimpl.h>
 16: #include <petscdraw.h>

 18: typedef struct {
 19:   PetscScalar center;     /* center of the ellipse */
 20:   PetscReal   radius;     /* radius of the ellipse */
 21:   PetscReal   vscale;     /* vertical scale of the ellipse */
 22:   PetscReal   start_ang;  /* start angle */
 23:   PetscReal   end_ang;    /* end angle */
 24:   PetscReal   width;      /* ring width */
 25: } RG_RING;

 27: static PetscErrorCode RGRingSetParameters_Ring(RG rg,PetscScalar center,PetscReal radius,PetscReal vscale,PetscReal start_ang,PetscReal end_ang,PetscReal width)
 28: {
 29:   RG_RING *ctx = (RG_RING*)rg->data;

 31:   PetscFunctionBegin;
 32:   ctx->center = center;
 33:   if (radius == (PetscReal)PETSC_DETERMINE) {
 34:     ctx->radius = 1.0;
 35:   } else if (radius != (PetscReal)PETSC_CURRENT) {
 36:     PetscCheck(radius>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The radius argument must be > 0.0");
 37:     ctx->radius = radius;
 38:   }
 39:   if (vscale == (PetscReal)PETSC_DETERMINE) {
 40:     ctx->vscale = 1.0;
 41:   } else if (vscale != (PetscReal)PETSC_CURRENT) {
 42:     PetscCheck(vscale>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The vscale argument must be > 0.0");
 43:     ctx->vscale = vscale;
 44:   }
 45:   if (start_ang == (PetscReal)PETSC_DETERMINE) {
 46:     ctx->start_ang = 0.0;
 47:   } else if (start_ang != (PetscReal)PETSC_CURRENT) {
 48:     PetscCheck(start_ang>=0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The right-hand side angle argument must be >= 0.0");
 49:     PetscCheck(start_ang<=1.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The right-hand side angle argument must be <= 1.0");
 50:     ctx->start_ang = start_ang;
 51:   }
 52:   if (end_ang == (PetscReal)PETSC_DETERMINE) {
 53:     ctx->end_ang = 1.0;
 54:   } else if (end_ang != (PetscReal)PETSC_CURRENT) {
 55:     PetscCheck(end_ang>=0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The left-hand side angle argument must be >= 0.0");
 56:     PetscCheck(end_ang<=1.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The left-hand side angle argument must be <= 1.0");
 57:     ctx->end_ang = end_ang;
 58:   }
 59: #if !defined(PETSC_USE_COMPLEX)
 60:   PetscCheck(ctx->start_ang+ctx->end_ang==1.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_WRONG,"In real scalars the region must be symmetric wrt real axis");
 61: #endif
 62:   if (width == (PetscReal)PETSC_DETERMINE) {
 63:     ctx->width = 0.1;
 64:   } else if (width != (PetscReal)PETSC_CURRENT) {
 65:     PetscCheck(width>0.0,PetscObjectComm((PetscObject)rg),PETSC_ERR_ARG_OUTOFRANGE,"The width argument must be > 0.0");
 66:     ctx->width = width;
 67:   }
 68:   PetscFunctionReturn(PETSC_SUCCESS);
 69: }

 71: /*@
 72:    RGRingSetParameters - Sets the parameters defining the ring region.

 74:    Logically Collective

 76:    Input Parameters:
 77: +  rg        - the region context
 78: .  center    - center of the ellipse
 79: .  radius    - radius of the ellipse
 80: .  vscale    - vertical scale of the ellipse
 81: .  start_ang - the right-hand side angle
 82: .  end_ang   - the left-hand side angle
 83: -  width     - width of the ring

 85:    Options Database Keys:
 86: +  -rg_ring_center     - Sets the center
 87: .  -rg_ring_radius     - Sets the radius
 88: .  -rg_ring_vscale     - Sets the vertical scale
 89: .  -rg_ring_startangle - Sets the right-hand side angle
 90: .  -rg_ring_endangle   - Sets the left-hand side angle
 91: -  -rg_ring_width      - Sets the width of the ring

 93:    Notes:
 94:    The values of center, radius and vscale have the same meaning as in the
 95:    ellipse region. The startangle and endangle define the span of the ring
 96:    (by default it is the whole ring), while the width is the separation
 97:    between the two concentric ellipses (above and below the radius by
 98:    width/2).

100:    The start and end angles are expressed as a fraction of the circumference.
101:    The allowed range is [0..1], with 0 corresponding to 0 radians, 0.25 to
102:    pi/2 radians, and so on. It is allowed to have startangle>endangle, in
103:    which case the ring region crosses over the zero angle.

105:    In the case of complex scalars, a complex center can be provided in the
106:    command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
107:    -rg_ring_center 1.0+2.0i

109:    When PETSc is built with real scalars, the center is restricted to a real value,
110:    and the start and end angles must be such that the region is symmetric with
111:    respect to the real axis.

113:    For all arguments except center, you can use PETSC_CURRENT to keep the current
114:    value, and PETSC_DETERMINE to set them to a default value (1 for radius, vscale,
115:    end_ang, 0 for start_and, and 0.1 for width).

117:    Level: advanced

119: .seealso: RGRingGetParameters()
120: @*/
121: PetscErrorCode RGRingSetParameters(RG rg,PetscScalar center,PetscReal radius,PetscReal vscale,PetscReal start_ang,PetscReal end_ang,PetscReal width)
122: {
123:   PetscFunctionBegin;
131:   PetscTryMethod(rg,"RGRingSetParameters_C",(RG,PetscScalar,PetscReal,PetscReal,PetscReal,PetscReal,PetscReal),(rg,center,radius,vscale,start_ang,end_ang,width));
132:   PetscFunctionReturn(PETSC_SUCCESS);
133: }

135: static PetscErrorCode RGRingGetParameters_Ring(RG rg,PetscScalar *center,PetscReal *radius,PetscReal *vscale,PetscReal *start_ang,PetscReal *end_ang,PetscReal *width)
136: {
137:   RG_RING *ctx = (RG_RING*)rg->data;

139:   PetscFunctionBegin;
140:   if (center)    *center    = ctx->center;
141:   if (radius)    *radius    = ctx->radius;
142:   if (vscale)    *vscale    = ctx->vscale;
143:   if (start_ang) *start_ang = ctx->start_ang;
144:   if (end_ang)   *end_ang   = ctx->end_ang;
145:   if (width)     *width     = ctx->width;
146:   PetscFunctionReturn(PETSC_SUCCESS);
147: }

149: /*@
150:    RGRingGetParameters - Gets the parameters that define the ring region.

152:    Not Collective

154:    Input Parameter:
155: .  rg     - the region context

157:    Output Parameters:
158: +  center    - center of the region
159: .  radius    - radius of the region
160: .  vscale    - vertical scale of the region
161: .  start_ang - the right-hand side angle
162: .  end_ang   - the left-hand side angle
163: -  width     - width of the ring

165:    Level: advanced

167: .seealso: RGRingSetParameters()
168: @*/
169: PetscErrorCode RGRingGetParameters(RG rg,PetscScalar *center,PetscReal *radius,PetscReal *vscale,PetscReal *start_ang,PetscReal *end_ang,PetscReal *width)
170: {
171:   PetscFunctionBegin;
173:   PetscUseMethod(rg,"RGRingGetParameters_C",(RG,PetscScalar*,PetscReal*,PetscReal*,PetscReal*,PetscReal*,PetscReal*),(rg,center,radius,vscale,start_ang,end_ang,width));
174:   PetscFunctionReturn(PETSC_SUCCESS);
175: }

177: static PetscErrorCode RGView_Ring(RG rg,PetscViewer viewer)
178: {
179:   RG_RING        *ctx = (RG_RING*)rg->data;
180:   int            winw,winh;
181:   PetscBool      isdraw,isascii;
182:   PetscDraw      draw;
183:   PetscDrawAxis  axis;
184:   PetscReal      cx,cy,radius,width,ab,cd,lx,ly,w,end_ang,x1,y1,x2,y2,r,theta,scale=1.2;
185:   char           str[50];

187:   PetscFunctionBegin;
188:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
189:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
190:   if (isascii) {
191:     PetscCall(SlepcSNPrintfScalar(str,sizeof(str),ctx->center,PETSC_FALSE));
192:     PetscCall(PetscViewerASCIIPrintf(viewer,"  center: %s, radius: %g, vscale: %g, start angle: %g, end angle: %g, ring width: %g\n",str,RGShowReal(ctx->radius),RGShowReal(ctx->vscale),(double)ctx->start_ang,(double)ctx->end_ang,(double)ctx->width));
193:   } else if (isdraw) {
194:     PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
195:     PetscCall(PetscDrawCheckResizedWindow(draw));
196:     PetscCall(PetscDrawGetWindowSize(draw,&winw,&winh));
197:     winw = PetscMax(winw,1); winh = PetscMax(winh,1);
198:     PetscCall(PetscDrawClear(draw));
199:     PetscCall(PetscDrawSetTitle(draw,"Ring region"));
200:     PetscCall(PetscDrawAxisCreate(draw,&axis));
201:     cx = PetscRealPart(ctx->center)*rg->sfactor;
202:     cy = PetscImaginaryPart(ctx->center)*rg->sfactor;
203:     radius = ctx->radius*rg->sfactor;
204:     width  = ctx->width*rg->sfactor;
205:     lx = 2*(radius+width);
206:     ly = 2*(radius+width)*ctx->vscale;
207:     ab = cx;
208:     cd = cy;
209:     w  = scale*PetscMax(lx/winw,ly/winh)/2;
210:     PetscCall(PetscDrawAxisSetLimits(axis,ab-w*winw,ab+w*winw,cd-w*winh,cd+w*winh));
211:     PetscCall(PetscDrawAxisDraw(axis));
212:     PetscCall(PetscDrawAxisDestroy(&axis));
213:     /* draw outer ellipse */
214:     PetscCall(PetscDrawEllipse(draw,cx,cy,2*(radius+width),2*(radius+width)*ctx->vscale,PETSC_DRAW_ORANGE));
215:     /* remove inner part */
216:     PetscCall(PetscDrawEllipse(draw,cx,cy,2*(radius-width),2*(radius-width)*ctx->vscale,PETSC_DRAW_WHITE));
217:     if (ctx->start_ang!=ctx->end_ang) {
218:       /* remove section from end_ang to start_ang */
219:       end_ang = (ctx->start_ang<ctx->end_ang)? ctx->end_ang-1: ctx->end_ang;
220:       theta = end_ang;
221:       r = scale*(radius+width);
222:       if (ctx->vscale>1) r *= ctx->vscale;
223:       x1 = PetscMin(PetscMax(ab+r*PetscCosReal(2.0*PETSC_PI*theta),ab-w*winw),ab+w*winw);
224:       y1 = PetscMin(PetscMax(cd+r*PetscSinReal(2.0*PETSC_PI*theta),cd-w*winh),cd+w*winh);
225:       do {
226:         theta = PetscMin(PetscFloorReal(8*theta+1)/8,ctx->start_ang);
227:         x2 = PetscMin(PetscMax(ab+r*PetscCosReal(2.0*PETSC_PI*theta),ab-w*winw),ab+w*winw);
228:         y2 = PetscMin(PetscMax(cd+r*PetscSinReal(2.0*PETSC_PI*theta),cd-w*winh),cd+w*winh);
229:         PetscCall(PetscDrawTriangle(draw,cx,cy,x1,y1,x2,y2,PETSC_DRAW_WHITE,PETSC_DRAW_WHITE,PETSC_DRAW_WHITE));
230:         x1 = x2; y1 = y2;
231:       } while (theta<ctx->start_ang);
232:     }
233:     PetscCall(PetscDrawFlush(draw));
234:     PetscCall(PetscDrawSave(draw));
235:     PetscCall(PetscDrawPause(draw));
236:   }
237:   PetscFunctionReturn(PETSC_SUCCESS);
238: }

240: static PetscErrorCode RGIsTrivial_Ring(RG rg,PetscBool *trivial)
241: {
242:   RG_RING *ctx = (RG_RING*)rg->data;

244:   PetscFunctionBegin;
245:   if (rg->complement) *trivial = PetscNot(ctx->radius);
246:   else *trivial = PetscNot(ctx->radius<PETSC_MAX_REAL);
247:   PetscFunctionReturn(PETSC_SUCCESS);
248: }

250: static PetscErrorCode RGComputeContour_Ring(RG rg,PetscInt n,PetscScalar *cr,PetscScalar *ci)
251: {
252:   RG_RING   *ctx = (RG_RING*)rg->data;
253:   PetscReal theta,start_ang;
254:   PetscInt  i,n2=n/2;

256:   PetscFunctionBegin;
257:   start_ang = (ctx->start_ang>ctx->end_ang)? ctx->start_ang-1: ctx->start_ang;
258:   for (i=0;i<n;i++) {
259:     if (i < n2) {
260:       theta = ((ctx->end_ang-start_ang)*i/n2 + start_ang)*2.0*PETSC_PI;
261: #if defined(PETSC_USE_COMPLEX)
262:       cr[i] = ctx->center + (ctx->radius+ctx->width/2.0)*PetscCMPLX(PetscCosReal(theta),ctx->vscale*PetscSinReal(theta));
263: #else
264:       if (cr) cr[i] = ctx->center + (ctx->radius+ctx->width/2.0)*PetscCosReal(theta);
265:       if (ci) ci[i] = (ctx->radius+ctx->width/2.0)*ctx->vscale*PetscSinReal(theta);
266: #endif
267:     } else {
268:       theta = ((ctx->end_ang-start_ang)*(n-i)/n2 + start_ang)*2.0*PETSC_PI;
269: #if defined(PETSC_USE_COMPLEX)
270:       cr[i] = ctx->center + (ctx->radius-ctx->width/2.0)*PetscCMPLX(PetscCosReal(theta),ctx->vscale*PetscSinReal(theta));
271: #else
272:       if (cr) cr[i] = ctx->center + (ctx->radius-ctx->width/2.0)*PetscCosReal(theta);
273:       if (ci) ci[i] = (ctx->radius-ctx->width/2.0)*ctx->vscale*PetscSinReal(theta);
274: #endif
275:     }
276:   }
277:   PetscFunctionReturn(PETSC_SUCCESS);
278: }

280: static PetscErrorCode RGComputeBoundingBox_Ring(RG rg,PetscReal *a,PetscReal *b,PetscReal *c,PetscReal *d)
281: {
282:   RG_RING *ctx = (RG_RING*)rg->data;

284:   PetscFunctionBegin;
285:   /* current implementation does not return a tight bounding box */
286:   if (a) *a = PetscRealPart(ctx->center) - (ctx->radius+ctx->width/2.0);
287:   if (b) *b = PetscRealPart(ctx->center) + (ctx->radius+ctx->width/2.0);
288:   if (c) *c = PetscImaginaryPart(ctx->center) - (ctx->radius+ctx->width/2.0)*ctx->vscale;
289:   if (d) *d = PetscImaginaryPart(ctx->center) + (ctx->radius+ctx->width/2.0)*ctx->vscale;
290:   PetscFunctionReturn(PETSC_SUCCESS);
291: }

293: static PetscErrorCode RGComputeQuadrature_Ring(RG rg,RGQuadRule quad,PetscInt n,PetscScalar *z,PetscScalar *zn,PetscScalar *w)
294: {
295:   RG_RING     *ctx = (RG_RING*)rg->data;
296:   PetscReal   max_w=0.0;
297:   PetscScalar tmp,tmp2;
298:   PetscInt    i,j;

300:   PetscFunctionBegin;
301:   if (quad == RG_QUADRULE_CHEBYSHEV) {
302: #if defined(PETSC_USE_COMPLEX)
303:     PetscReal theta;
304:     for (i=0;i<n;i++) {
305:       theta = PETSC_PI*(i+0.5)/n;
306:       zn[i] = PetscCosReal(theta);
307:       w[i]  = PetscCosReal((n-1)*theta)/n;
308:       theta = (ctx->start_ang*2.0+(ctx->end_ang-ctx->start_ang)*(PetscRealPart(zn[i])+1.0))*PETSC_PI;
309:       z[i] = rg->sfactor*(ctx->center + ctx->radius*PetscCMPLX(PetscCosReal(theta),ctx->vscale*PetscSinReal(theta)));
310:     }
311: #else
312:     SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Integration points on a vertical line require complex arithmetic");
313: #endif
314:   } else {  /* RG_QUADRULE_TRAPEZOIDAL */
315:     for (i=0;i<n;i++) {
316:       zn[i] = (z[i]-rg->sfactor*ctx->center)/(rg->sfactor*ctx->radius);
317:       tmp = 1.0; tmp2 = 1.0;
318:       for (j=0;j<n;j++) {
319:         tmp *= z[j];
320:         if (i != j) tmp2 *= z[j]-z[i];
321:       }
322:       w[i] = tmp/tmp2;
323:       max_w = PetscMax(PetscAbsScalar(w[i]),max_w);
324:     }
325:     for (i=0;i<n;i++) w[i] /= (PetscScalar)max_w;
326:   }
327:   PetscFunctionReturn(PETSC_SUCCESS);
328: }

330: static PetscErrorCode RGCheckInside_Ring(RG rg,PetscReal px,PetscReal py,PetscInt *inside)
331: {
332:   RG_RING   *ctx = (RG_RING*)rg->data;
333:   PetscReal dx,dy,r;

335:   PetscFunctionBegin;
336:   /* outer ellipse */
337: #if defined(PETSC_USE_COMPLEX)
338:   dx = (px-PetscRealPart(ctx->center))/(ctx->radius+ctx->width/2.0);
339:   dy = (py-PetscImaginaryPart(ctx->center))/(ctx->radius+ctx->width/2.0);
340: #else
341:   dx = (px-ctx->center)/(ctx->radius+ctx->width/2.0);
342:   dy = py/(ctx->radius+ctx->width/2.0);
343: #endif
344:   r = 1.0-dx*dx-(dy*dy)/(ctx->vscale*ctx->vscale);
345:   *inside = PetscSign(r);
346:   /* inner ellipse */
347: #if defined(PETSC_USE_COMPLEX)
348:   dx = (px-PetscRealPart(ctx->center))/(ctx->radius-ctx->width/2.0);
349:   dy = (py-PetscImaginaryPart(ctx->center))/(ctx->radius-ctx->width/2.0);
350: #else
351:   dx = (px-ctx->center)/(ctx->radius-ctx->width/2.0);
352:   dy = py/(ctx->radius-ctx->width/2.0);
353: #endif
354:   r = -1.0+dx*dx+(dy*dy)/(ctx->vscale*ctx->vscale);
355:   *inside *= PetscSign(r);
356:   if (*inside == 1) {  /* check angles */
357: #if defined(PETSC_USE_COMPLEX)
358:     dx = (px-PetscRealPart(ctx->center));
359:     dy = (py-PetscImaginaryPart(ctx->center));
360: #else
361:     dx = px-ctx->center;
362:     dy = py;
363: #endif
364:     if (dx == 0) {
365:       if (dy == 0) r = -1;
366:       else if (dy > 0) r = 0.25;
367:       else r = 0.75;
368:     } else if (dx > 0) {
369:       r = PetscAtanReal((dy/ctx->vscale)/dx);
370:       if (dy >= 0) r /= 2*PETSC_PI;
371:       else r = r/(2*PETSC_PI)+1;
372:     } else r = PetscAtanReal((dy/ctx->vscale)/dx)/(2*PETSC_PI)+0.5;
373:     if (ctx->start_ang>ctx->end_ang) {
374:       if (r>ctx->end_ang && r<ctx->start_ang) *inside = -1;
375:     } else {
376:       if (r<ctx->start_ang || r>ctx->end_ang) *inside = -1;
377:     }
378:   }
379:   PetscFunctionReturn(PETSC_SUCCESS);
380: }

382: static PetscErrorCode RGIsAxisymmetric_Ring(RG rg,PetscBool vertical,PetscBool *symm)
383: {
384:   RG_RING *ctx = (RG_RING*)rg->data;

386:   PetscFunctionBegin;
387:   if (vertical) *symm = (PetscRealPart(ctx->center) == 0.0 && PetscAbs(ctx->start_ang+ctx->end_ang-PetscRealConstant(1.0)) == 0.5)? PETSC_TRUE: PETSC_FALSE;
388:   else *symm = (PetscImaginaryPart(ctx->center) == 0.0 && ctx->start_ang+ctx->end_ang == 1.0)? PETSC_TRUE: PETSC_FALSE;
389:   PetscFunctionReturn(PETSC_SUCCESS);
390: }

392: static PetscErrorCode RGSetFromOptions_Ring(RG rg,PetscOptionItems *PetscOptionsObject)
393: {
394:   PetscScalar    s;
395:   PetscReal      r1,r2,r3,r4,r5;
396:   PetscBool      flg1,flg2,flg3,flg4,flg5,flg6;

398:   PetscFunctionBegin;
399:   PetscOptionsHeadBegin(PetscOptionsObject,"RG Ring Options");

401:     PetscCall(RGRingGetParameters(rg,&s,&r1,&r2,&r3,&r4,&r5));
402:     PetscCall(PetscOptionsScalar("-rg_ring_center","Center of ellipse","RGRingSetParameters",s,&s,&flg1));
403:     PetscCall(PetscOptionsReal("-rg_ring_radius","Radius of ellipse","RGRingSetParameters",r1,&r1,&flg2));
404:     PetscCall(PetscOptionsReal("-rg_ring_vscale","Vertical scale of ellipse","RGRingSetParameters",r2,&r2,&flg3));
405:     PetscCall(PetscOptionsReal("-rg_ring_startangle","Right-hand side angle","RGRingSetParameters",r3,&r3,&flg4));
406:     PetscCall(PetscOptionsReal("-rg_ring_endangle","Left-hand side angle","RGRingSetParameters",r4,&r4,&flg5));
407:     PetscCall(PetscOptionsReal("-rg_ring_width","Width of ring","RGRingSetParameters",r5,&r5,&flg6));
408:     if (flg1 || flg2 || flg3 || flg4 || flg5 || flg6) PetscCall(RGRingSetParameters(rg,s,r1,r2,r3,r4,r5));

410:   PetscOptionsHeadEnd();
411:   PetscFunctionReturn(PETSC_SUCCESS);
412: }

414: static PetscErrorCode RGDestroy_Ring(RG rg)
415: {
416:   PetscFunctionBegin;
417:   PetscCall(PetscFree(rg->data));
418:   PetscCall(PetscObjectComposeFunction((PetscObject)rg,"RGRingSetParameters_C",NULL));
419:   PetscCall(PetscObjectComposeFunction((PetscObject)rg,"RGRingGetParameters_C",NULL));
420:   PetscFunctionReturn(PETSC_SUCCESS);
421: }

423: SLEPC_EXTERN PetscErrorCode RGCreate_Ring(RG rg)
424: {
425:   RG_RING        *ring;

427:   PetscFunctionBegin;
428:   PetscCall(PetscNew(&ring));
429:   ring->center    = 0.0;
430:   ring->radius    = PETSC_MAX_REAL;
431:   ring->vscale    = 1.0;
432:   ring->start_ang = 0.0;
433:   ring->end_ang   = 1.0;
434:   ring->width     = 0.1;
435:   rg->data = (void*)ring;

437:   rg->ops->istrivial         = RGIsTrivial_Ring;
438:   rg->ops->computecontour    = RGComputeContour_Ring;
439:   rg->ops->computebbox       = RGComputeBoundingBox_Ring;
440:   rg->ops->computequadrature = RGComputeQuadrature_Ring;
441:   rg->ops->checkinside       = RGCheckInside_Ring;
442:   rg->ops->isaxisymmetric    = RGIsAxisymmetric_Ring;
443:   rg->ops->setfromoptions    = RGSetFromOptions_Ring;
444:   rg->ops->view              = RGView_Ring;
445:   rg->ops->destroy           = RGDestroy_Ring;
446:   PetscCall(PetscObjectComposeFunction((PetscObject)rg,"RGRingSetParameters_C",RGRingSetParameters_Ring));
447:   PetscCall(PetscObjectComposeFunction((PetscObject)rg,"RGRingGetParameters_C",RGRingGetParameters_Ring));
448:   PetscFunctionReturn(PETSC_SUCCESS);
449: }