Actual source code: test1.c

slepc-3.8.3 2018-04-03
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2017, 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: int main(int argc,char **argv)
 16: {
 18:   RG             rg;
 19:   PetscInt       i,inside,nv;
 20:   PetscBool      triv;
 21:   PetscReal      re,im,a,b,c,d;
 22:   PetscScalar    ar,ai,cr[10],ci[10],vr[7],vi[7],*pr,*pi;

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

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

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

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

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

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

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

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

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

135:   PetscPrintf(PETSC_COMM_WORLD,"Contour points: ");
136:   RGComputeContour(rg,10,cr,ci);
137:   for (i=0;i<10;i++) {
138: #if defined(PETSC_USE_COMPLEX)
139:     re = PetscRealPart(cr[i]);
140:     im = PetscImaginaryPart(cr[i]);
141: #else
142:     re = cr[i];
143:     im = ci[i];
144: #endif
145:     PetscPrintf(PETSC_COMM_WORLD,"(%.3g,%.3g) ",(double)re,(double)im);
146:   }
147:   PetscPrintf(PETSC_COMM_WORLD,"\n");

149:   /* check vertices */
150:   RGPolygonGetVertices(rg,&nv,&pr,&pi);
151:   if (nv!=7) SETERRQ1(PETSC_COMM_WORLD,1,"Wrong number of vertices: %D",nv);
152:   for (i=0;i<nv;i++) {
153:     if (pr[i]!=vr[i]
154: #if !defined(PETSC_USE_COMPLEX)
155:         || pi[i]!=vi[i]
156: #endif
157:        ) SETERRQ1(PETSC_COMM_WORLD,1,"Vertex number %D does not match",i);
158:   }

160:   RGDestroy(&rg);
161:   SlepcFinalize();
162:   return ierr;
163: }