Actual source code: test2.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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: {
20: PetscInt inside;
21: PetscScalar ar,ai;
24: #if defined(PETSC_USE_COMPLEX)
25: ar = re+im*PETSC_i;
26: #else
27: ar = re; ai = im;
28: #endif
29: RGCheckInside(rg,1,&ar,&ai,&inside);
30: PetscPrintf(PETSC_COMM_WORLD,"Point (%g,%g) is %s the region\n",(double)re,(double)im,(inside>=0)?"inside":"outside");
31: return(0);
32: }
34: int main(int argc,char **argv)
35: {
37: RG rg;
38: RGType rtype;
39: PetscInt i;
40: PetscBool triv;
41: PetscReal re,im,radius,vscale,start_ang,end_ang,width;
42: PetscScalar center,cr[NPOINTS],ci[NPOINTS];
44: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
45: RGCreate(PETSC_COMM_WORLD,&rg);
47: RGSetType(rg,RGRING);
48: RGIsTrivial(rg,&triv);
49: if (!triv) SETERRQ(PETSC_COMM_WORLD,1,"Region should be trivial before setting parameters");
50: RGRingSetParameters(rg,2,PETSC_DEFAULT,0.5,0.25,0.75,0.1);
51: RGSetFromOptions(rg);
52: RGIsTrivial(rg,&triv);
53: if (triv) SETERRQ(PETSC_COMM_WORLD,1,"Region should be non-trivial after setting parameters");
54: RGView(rg,NULL);
55: RGViewFromOptions(rg,NULL,"-rg_view");
57: RGGetType(rg,&rtype);
58: RGRingGetParameters(rg,¢er,&radius,&vscale,&start_ang,&end_ang,&width);
59: 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);
61: CheckPoint(rg,3.0,0.3);
62: CheckPoint(rg,1.1747,0.28253);
64: PetscPrintf(PETSC_COMM_WORLD,"\nContour points: ");
65: RGComputeContour(rg,NPOINTS,cr,ci);
66: for (i=0;i<NPOINTS;i++) {
67: #if defined(PETSC_USE_COMPLEX)
68: re = PetscRealPart(cr[i]);
69: im = PetscImaginaryPart(cr[i]);
70: #else
71: re = cr[i];
72: im = ci[i];
73: #endif
74: PetscPrintf(PETSC_COMM_WORLD,"(%.3g,%.3g) ",(double)re,(double)im);
75: }
76: PetscPrintf(PETSC_COMM_WORLD,"\n");
78: RGDestroy(&rg);
79: SlepcFinalize();
80: return ierr;
81: }
83: /*TEST
85: test:
86: suffix: 1
87: args: -rg_ring_width 0.015
89: test:
90: suffix: 2
91: args: -rg_ring_width 0.015 -rg_scale 1.5
93: test:
94: suffix: 3
95: args: -rg_view draw:tikz:test2_3_ring.tikz
96: filter: cat - test2_3_ring.tikz
97: requires: !single
99: TEST*/