Actual source code: test1.c

slepc-3.16.0 2021-09-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2021, 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 RG interface functions.\n\n";

 13: #include <slepcrg.h>

 15: #define NPOINTS 10
 16: #define NVERTEX 7

 18: int main(int argc,char **argv)
 19: {
 21:   RG             rg;
 22:   PetscInt       i,inside,nv;
 23:   PetscBool      triv;
 24:   PetscReal      re,im,a,b,c,d;
 25:   PetscScalar    ar,ai,cr[NPOINTS],ci[NPOINTS],vr[NVERTEX],vi[NVERTEX],*pr,*pi;

 27:   SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
 28:   RGCreate(PETSC_COMM_WORLD,&rg);

 30:   /* ellipse */
 31:   RGSetType(rg,RGELLIPSE);
 32:   RGIsTrivial(rg,&triv);
 33:   if (!triv) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Region should be trivial before setting parameters");
 34:   RGEllipseSetParameters(rg,1.1,2,0.1);
 35:   RGSetFromOptions(rg);
 36:   RGIsTrivial(rg,&triv);
 37:   if (triv) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Region should be non-trivial after setting parameters");
 38:   RGView(rg,NULL);
 39:   RGViewFromOptions(rg,NULL,"-rg_ellipse_view");
 40:   re = 0.1; im = 0.3;
 41: #if defined(PETSC_USE_COMPLEX)
 42:   ar = PetscCMPLX(re,im);
 43: #else
 44:   ar = re; ai = im;
 45: #endif
 46:   RGCheckInside(rg,1,&ar,&ai,&inside);
 47:   PetscPrintf(PETSC_COMM_WORLD,"Point (%g,%g) is %s the region\n",(double)re,(double)im,(inside>=0)?"inside":"outside");

 49:   RGComputeBoundingBox(rg,&a,&b,&c,&d);
 50:   PetscPrintf(PETSC_COMM_WORLD,"The bounding box is [%g,%g]x[%g,%g]\n",(double)a,(double)b,(double)c,(double)d);

 52:   PetscPrintf(PETSC_COMM_WORLD,"Contour points: ");
 53:   RGComputeContour(rg,NPOINTS,cr,ci);
 54:   for (i=0;i<NPOINTS;i++) {
 55: #if defined(PETSC_USE_COMPLEX)
 56:     re = PetscRealPart(cr[i]);
 57:     im = PetscImaginaryPart(cr[i]);
 58: #else
 59:     re = cr[i];
 60:     im = ci[i];
 61: #endif
 62:     PetscPrintf(PETSC_COMM_WORLD,"(%.3g,%.3g) ",(double)re,(double)im);
 63:   }
 64:   PetscPrintf(PETSC_COMM_WORLD,"\n");

 66:   /* interval */
 67:   RGSetType(rg,RGINTERVAL);
 68:   RGIsTrivial(rg,&triv);
 69:   if (!triv) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Region should be trivial before setting parameters");
 70:   RGIntervalSetEndpoints(rg,-1,1,-0.1,0.1);
 71:   RGSetFromOptions(rg);
 72:   RGIsTrivial(rg,&triv);
 73:   if (triv) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Region should be non-trivial after setting parameters");
 74:   RGView(rg,NULL);
 75:   RGViewFromOptions(rg,NULL,"-rg_interval_view");
 76:   re = 0.2; im = 0;
 77: #if defined(PETSC_USE_COMPLEX)
 78:   ar = PetscCMPLX(re,im);
 79: #else
 80:   ar = re; ai = im;
 81: #endif
 82:   RGCheckInside(rg,1,&ar,&ai,&inside);
 83:   PetscPrintf(PETSC_COMM_WORLD,"Point (%g,%g) is %s the region\n",(double)re,(double)im,(inside>=0)?"inside":"outside");

 85:   RGComputeBoundingBox(rg,&a,&b,&c,&d);
 86:   PetscPrintf(PETSC_COMM_WORLD,"The bounding box is [%g,%g]x[%g,%g]\n",(double)a,(double)b,(double)c,(double)d);

 88:   PetscPrintf(PETSC_COMM_WORLD,"Contour points: ");
 89:   RGComputeContour(rg,NPOINTS,cr,ci);
 90:   for (i=0;i<NPOINTS;i++) {
 91: #if defined(PETSC_USE_COMPLEX)
 92:     re = PetscRealPart(cr[i]);
 93:     im = PetscImaginaryPart(cr[i]);
 94: #else
 95:     re = cr[i];
 96:     im = ci[i];
 97: #endif
 98:     PetscPrintf(PETSC_COMM_WORLD,"(%.3g,%.3g) ",(double)re,(double)im);
 99:   }
100:   PetscPrintf(PETSC_COMM_WORLD,"\n");

102:   /* polygon */
103: #if defined(PETSC_USE_COMPLEX)
104:   vr[0] = PetscCMPLX(0.0,2.0);
105:   vr[1] = PetscCMPLX(1.0,4.0);
106:   vr[2] = PetscCMPLX(2.0,5.0);
107:   vr[3] = PetscCMPLX(4.0,3.0);
108:   vr[4] = PetscCMPLX(5.0,4.0);
109:   vr[5] = PetscCMPLX(6.0,1.0);
110:   vr[6] = PetscCMPLX(2.0,0.0);
111: #else
112:   vr[0] = 0.0; vi[0] = 1.0;
113:   vr[1] = 0.0; vi[1] = -1.0;
114:   vr[2] = 0.6; vi[2] = -0.8;
115:   vr[3] = 1.0; vi[3] = -1.0;
116:   vr[4] = 2.0; vi[4] = 0.0;
117:   vr[5] = 1.0; vi[5] = 1.0;
118:   vr[6] = 0.6; vi[6] = 0.8;
119: #endif
120:   RGSetType(rg,RGPOLYGON);
121:   RGIsTrivial(rg,&triv);
122:   if (!triv) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Region should be trivial before setting parameters");
123:   RGPolygonSetVertices(rg,NVERTEX,vr,vi);
124:   RGSetFromOptions(rg);
125:   RGIsTrivial(rg,&triv);
126:   if (triv) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Region should be non-trivial after setting parameters");
127:   RGView(rg,NULL);
128:   RGViewFromOptions(rg,NULL,"-rg_polygon_view");
129:   re = 5; im = 0.9;
130: #if defined(PETSC_USE_COMPLEX)
131:   ar = PetscCMPLX(re,im);
132: #else
133:   ar = re; ai = im;
134: #endif
135:   RGCheckInside(rg,1,&ar,&ai,&inside);
136:   PetscPrintf(PETSC_COMM_WORLD,"Point (%g,%g) is %s the region\n",(double)re,(double)im,(inside>=0)?"inside":"outside");

138:   RGComputeBoundingBox(rg,&a,&b,&c,&d);
139:   PetscPrintf(PETSC_COMM_WORLD,"The bounding box is [%g,%g]x[%g,%g]\n",(double)a,(double)b,(double)c,(double)d);

141:   PetscPrintf(PETSC_COMM_WORLD,"Contour points: ");
142:   RGComputeContour(rg,NPOINTS,cr,ci);
143:   for (i=0;i<NPOINTS;i++) {
144: #if defined(PETSC_USE_COMPLEX)
145:     re = PetscRealPart(cr[i]);
146:     im = PetscImaginaryPart(cr[i]);
147: #else
148:     re = cr[i];
149:     im = ci[i];
150: #endif
151:     PetscPrintf(PETSC_COMM_WORLD,"(%.3g,%.3g) ",(double)re,(double)im);
152:   }
153:   PetscPrintf(PETSC_COMM_WORLD,"\n");

155:   /* check vertices */
156:   RGPolygonGetVertices(rg,&nv,&pr,&pi);
157:   if (nv!=NVERTEX) SETERRQ1(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Wrong number of vertices: %D",nv);
158:   for (i=0;i<nv;i++) {
159: #if !defined(PETSC_USE_COMPLEX)
160:     if (pr[i]!=vr[i] || pi[i]!=vi[i])
161: #else
162:     if (pr[i]!=vr[i])
163: #endif
164:        SETERRQ1(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Vertex number %D does not match",i);
165:   }

167:   PetscFree(pr);
168: #if !defined(PETSC_USE_COMPLEX)
169:   PetscFree(pi);
170: #endif
171:   RGDestroy(&rg);
172:   SlepcFinalize();
173:   return ierr;
174: }

176: /*TEST

178:    test:
179:       suffix: 1
180:       requires: !complex

182:    test:
183:       suffix: 1_complex
184:       requires: complex

186:    test:
187:       suffix: 2
188:       requires: !complex
189:       args: -rg_ellipse_view draw:tikz:ellipse.tikz -rg_interval_view draw:tikz:interval.tikz -rg_polygon_view draw:tikz:polygon.tikz
190:       filter: cat - ellipse.tikz interval.tikz polygon.tikz
191:       requires: !single

193:    test:
194:       suffix: 2_complex
195:       requires: complex !single
196:       args: -rg_ellipse_view draw:tikz:ellipse.tikz -rg_interval_view draw:tikz:interval.tikz -rg_polygon_view draw:tikz:polygon.tikz
197:       filter: cat - ellipse.tikz interval.tikz polygon.tikz

199: TEST*/