Actual source code: test1.c
slepc-3.7.3 2016-09-29
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
8: SLEPc is free software: you can redistribute it and/or modify it under the
9: terms of version 3 of the GNU Lesser General Public License as published by
10: the Free Software Foundation.
12: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
13: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
15: more details.
17: You should have received a copy of the GNU Lesser General Public License
18: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
19: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: */
22: static char help[] = "Test RG interface functions.\n\n";
24: #include <slepcrg.h>
28: int main(int argc,char **argv)
29: {
31: RG rg;
32: PetscInt i,inside;
33: PetscReal re,im;
34: PetscScalar ar,ai,cr[10],ci[10],vr[7],vi[7];
36: SlepcInitialize(&argc,&argv,(char*)0,help);
37: RGCreate(PETSC_COMM_WORLD,&rg);
39: /* ellipse */
40: RGSetType(rg,RGELLIPSE);
41: RGEllipseSetParameters(rg,1.1,2,0.1);
42: RGSetFromOptions(rg);
43: RGView(rg,NULL);
44: re = 0.1; im = 0.3;
45: #if defined(PETSC_USE_COMPLEX)
46: ar = re+im*PETSC_i;
47: #else
48: ar = re; ai = im;
49: #endif
50: RGCheckInside(rg,1,&ar,&ai,&inside);
51: PetscPrintf(PETSC_COMM_WORLD,"Point (%g,%g) is %s the region\n",(double)re,(double)im,(inside>=0)?"inside":"outside");
53: PetscPrintf(PETSC_COMM_WORLD,"Contour points: ");
54: RGComputeContour(rg,10,cr,ci);
55: for (i=0;i<10;i++) {
56: #if defined(PETSC_USE_COMPLEX)
57: re = PetscRealPart(cr[i]);
58: im = PetscImaginaryPart(cr[i]);
59: #else
60: re = cr[i];
61: im = ci[i];
62: #endif
63: PetscPrintf(PETSC_COMM_WORLD,"(%.3g,%.3g) ",(double)re,(double)im);
64: }
65: PetscPrintf(PETSC_COMM_WORLD,"\n");
67: /* interval */
68: RGSetType(rg,RGINTERVAL);
69: RGIntervalSetEndpoints(rg,-1,1,-0.1,0.1);
70: RGSetFromOptions(rg);
71: RGView(rg,NULL);
72: re = 0.2; im = 0;
73: #if defined(PETSC_USE_COMPLEX)
74: ar = re+im*PETSC_i;
75: #else
76: ar = re; ai = im;
77: #endif
78: RGCheckInside(rg,1,&ar,&ai,&inside);
79: PetscPrintf(PETSC_COMM_WORLD,"Point (%g,%g) is %s the region\n",(double)re,(double)im,(inside>=0)?"inside":"outside");
81: /* polygon */
82: #if defined(PETSC_USE_COMPLEX)
83: vr[0] = 0.0+2.0*PETSC_i;
84: vr[1] = 1.0+4.0*PETSC_i;
85: vr[2] = 2.0+5.0*PETSC_i;
86: vr[3] = 4.0+3.0*PETSC_i;
87: vr[4] = 5.0+4.0*PETSC_i;
88: vr[5] = 6.0+1.0*PETSC_i;
89: vr[6] = 2.0+0.0*PETSC_i;
90: #else
91: vr[0] = 0.0; vi[0] = 2.0;
92: vr[1] = 1.0; vi[1] = 4.0;
93: vr[2] = 2.0; vi[2] = 5.0;
94: vr[3] = 4.0; vi[3] = 3.0;
95: vr[4] = 5.0; vi[4] = 4.0;
96: vr[5] = 6.0; vi[5] = 1.0;
97: vr[6] = 2.0; vi[6] = 0.0;
98: #endif
99: RGSetType(rg,RGPOLYGON);
100: RGPolygonSetVertices(rg,7,vr,vi);
101: RGSetFromOptions(rg);
102: RGView(rg,NULL);
103: re = 5; im = 0.9;
104: #if defined(PETSC_USE_COMPLEX)
105: ar = re+im*PETSC_i;
106: #else
107: ar = re; ai = im;
108: #endif
109: RGCheckInside(rg,1,&ar,&ai,&inside);
110: PetscPrintf(PETSC_COMM_WORLD,"Point (%g,%g) is %s the region\n",(double)re,(double)im,(inside>=0)?"inside":"outside");
112: RGDestroy(&rg);
113: SlepcFinalize();
114: return ierr;
115: }