Actual source code: test2.c
slepc-3.21.1 2024-04-26
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: */
11: static char help[] = "Test the ring region.\n\n";
13: #include <slepcrg.h>
15: #define NPOINTS 11
17: PetscErrorCode CheckPoint(RG rg,PetscReal re,PetscReal im)
18: {
19: PetscInt inside;
20: PetscScalar ar,ai;
22: PetscFunctionBeginUser;
23: #if defined(PETSC_USE_COMPLEX)
24: ar = PetscCMPLX(re,im);
25: #else
26: ar = re; ai = im;
27: #endif
28: PetscCall(RGCheckInside(rg,1,&ar,&ai,&inside));
29: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Point (%g,%g) is %s the region\n",(double)re,(double)im,(inside>=0)?"inside":"outside"));
30: PetscFunctionReturn(PETSC_SUCCESS);
31: }
33: int main(int argc,char **argv)
34: {
35: RG rg;
36: RGType rtype;
37: PetscInt i;
38: PetscBool triv;
39: PetscReal re,im,radius,vscale,start_ang,end_ang,width,a,b,c,d;
40: PetscScalar center,cr[NPOINTS],ci[NPOINTS];
42: PetscFunctionBeginUser;
43: PetscCall(SlepcInitialize(&argc,&argv,(char*)0,help));
44: PetscCall(RGCreate(PETSC_COMM_WORLD,&rg));
46: PetscCall(RGSetType(rg,RGRING));
47: PetscCall(RGIsTrivial(rg,&triv));
48: PetscCheck(triv,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Region should be trivial before setting parameters");
49: PetscCall(RGRingSetParameters(rg,2,PETSC_DEFAULT,0.5,0.25,0.75,0.1));
50: PetscCall(RGSetFromOptions(rg));
51: PetscCall(RGIsTrivial(rg,&triv));
52: PetscCheck(!triv,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Region should be non-trivial after setting parameters");
53: PetscCall(RGView(rg,NULL));
54: PetscCall(RGViewFromOptions(rg,NULL,"-rg_view"));
56: PetscCall(RGGetType(rg,&rtype));
57: PetscCall(RGRingGetParameters(rg,¢er,&radius,&vscale,&start_ang,&end_ang,&width));
58: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"%s region: \n center=%g, radius=%g, vscale=%g\n start angle=%g, end angle=%g, width=%g\n\n",rtype,(double)PetscRealPart(center),(double)radius,(double)vscale,(double)start_ang,(double)end_ang,(double)width));
60: PetscCall(CheckPoint(rg,3.0,0.3));
61: PetscCall(CheckPoint(rg,1.1747,0.28253));
63: PetscCall(RGComputeBoundingBox(rg,&a,&b,&c,&d));
64: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"The bounding box is [%g,%g]x[%g,%g]\n",(double)a,(double)b,(double)c,(double)d));
66: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Contour points: "));
67: PetscCall(RGComputeContour(rg,NPOINTS,cr,ci));
68: for (i=0;i<NPOINTS;i++) {
69: #if defined(PETSC_USE_COMPLEX)
70: re = PetscRealPart(cr[i]);
71: im = PetscImaginaryPart(cr[i]);
72: #else
73: re = cr[i];
74: im = ci[i];
75: #endif
76: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"(%.3g,%.3g) ",(double)re,(double)im));
77: }
78: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\n"));
80: PetscCall(RGDestroy(&rg));
81: PetscCall(SlepcFinalize());
82: return 0;
83: }
85: /*TEST
87: test:
88: suffix: 1
89: args: -rg_ring_width 0.015
91: test:
92: suffix: 2
93: args: -rg_ring_width 0.015 -rg_scale 1.5
95: test:
96: suffix: 3
97: args: -rg_view draw:tikz:test2_3_ring.tikz
98: filter: cat - test2_3_ring.tikz
99: requires: !single
101: test:
102: suffix: 4
103: args: -rg_ring_width 0.015 -rg_ring_center 3 -rg_ring_radius 0.3 -rg_ring_vscale 1
104: requires: !single
106: test:
107: suffix: 5
108: args: -rg_ring_width 0.1 -rg_ring_center 0.35 -rg_ring_radius 0.86 -rg_ring_vscale 1 -rg_ring_startangle 0.75 -rg_ring_endangle 0.25
110: TEST*/