Actual source code: test1f.F

slepc-3.11.2 2019-07-30
Report Typos and Errors
  1: !
  2: !  - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3: !  SLEPc - Scalable Library for Eigenvalue Problem Computations
  4: !  Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain
  5: !
  6: !  This file is part of SLEPc.
  7: !  SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8: !  - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: !
 10: !  Program usage: mpiexec -n <np> ./test1f [-help]
 11: !
 12: !  Description: Simple example that tests RG interface functions.
 13: !
 14: ! ----------------------------------------------------------------------
 15: !
 16:       program main
 17: #include <slepc/finclude/slepcrg.h>
 18:       use slepcrg
 19:       implicit none

 21: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 22: !     Declarations
 23: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 25:       RG             rg
 26:       PetscInt       i,n,inside,one
 27:       PetscMPIInt    rank
 28:       PetscErrorCode ierr
 29:       PetscReal      re,im
 30:       PetscScalar    ar,ai,cr(10),ci(10)
 31:       PetscScalar    vr(7),vi(7)
 32:       PetscScalar    center
 33:       PetscReal      radius,vscale,a,b,c,d

 35: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 36: !     Beginning of program
 37: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 39:       one = 1
 40:       call SlepcInitialize(PETSC_NULL_CHARACTER,ierr)
 41:       call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
 42:       call RGCreate(PETSC_COMM_WORLD,rg,ierr)

 44: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 45: !     Ellipse
 46: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

 48:       call RGSetType(rg,RGELLIPSE,ierr)
 49:       center = 1.1
 50:       radius = 2
 51:       vscale = 0.1
 52:       call RGEllipseSetParameters(rg,center,radius,vscale,ierr)
 53:       call RGSetFromOptions(rg,ierr)
 54:       call RGView(rg,PETSC_NULL_VIEWER,ierr)
 55:       re = 0.1
 56:       im = 0.3
 57: #if defined(PETSC_USE_COMPLEX)
 58:       ar = re+im*PETSC_i
 59:       ai = 0.0
 60: #else
 61:       ar = re
 62:       ai = im
 63: #endif
 64:       call RGCheckInside(rg,one,ar,ai,inside,ierr)
 65:       if (rank .eq. 0) then
 66:         if (inside >= 0) then
 67:           write(*,110) re, im, 'inside'
 68:         else
 69:           write(*,110) re, im, 'outside'
 70:         endif
 71:       endif
 72:  110  format ('Point (',F4.1,',',F4.1,') is ',A7,' the region')

 74:       call RGComputeBoundingBox(rg,a,b,c,d,ierr)
 75:       if (rank .eq. 0) then
 76:         write(*,115) a, b, c, d
 77:       endif
 78:  115  format ('Bounding box: [',F4.1,',',F4.1,']x[',F4.1,',',F4.1,']')

 80:       if (rank .eq. 0) then
 81:         write (*,*) 'Contour points:'
 82:       endif
 83:       n = 10
 84:       call RGComputeContour(rg,n,cr,ci,ierr)
 85:       do i=1,n
 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:         if (rank .eq. 0) then
 94:           write(*,120) re, im
 95:         endif
 96:       enddo
 97:  120  format ('(',F7.4,',',F7.4,')')

 99: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
100: !     Interval
101: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

103:       call RGSetType(rg,RGINTERVAL,ierr)
104:       a = -1
105:       b =  1
106:       c = -0.1
107:       d =  0.1
108:       call RGIntervalSetEndpoints(rg,a,b,c,d,ierr)
109:       call RGSetFromOptions(rg,ierr)
110:       call RGView(rg,PETSC_NULL_VIEWER,ierr)
111:       re = 0.2
112:       im = 0
113: #if defined(PETSC_USE_COMPLEX)
114:       ar = re+im*PETSC_i
115:       ai = 0.0
116: #else
117:       ar = re
118:       ai = im
119: #endif
120:       call RGCheckInside(rg,one,ar,ai,inside,ierr)
121:       if (rank .eq. 0) then
122:         if (inside >= 0) then
123:           write(*,110) re, im, 'inside'
124:         else
125:           write(*,110) re, im, 'outside'
126:         endif
127:       endif

129: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
130: !     Polygon
131: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

133: #if defined(PETSC_USE_COMPLEX)
134:       vr(1) = 0.0+2.0*PETSC_i
135:       vr(2) = 1.0+4.0*PETSC_i
136:       vr(3) = 2.0+5.0*PETSC_i
137:       vr(4) = 4.0+3.0*PETSC_i
138:       vr(5) = 5.0+4.0*PETSC_i
139:       vr(6) = 6.0+1.0*PETSC_i
140:       vr(7) = 2.0+0.0*PETSC_i
141: #else
142:       vr(1) = 0.0
143:       vi(1) = 1.0
144:       vr(2) = 0.0
145:       vi(2) = -1.0
146:       vr(3) = 0.6
147:       vi(3) = -0.8
148:       vr(4) = 1.0
149:       vi(4) = -1.0
150:       vr(5) = 2.0
151:       vi(5) = 0.0
152:       vr(6) = 1.0
153:       vi(6) = 1.0
154:       vr(7) = 0.6
155:       vi(7) = 0.8
156: #endif
157:       call RGSetType(rg,RGPOLYGON,ierr)
158:       n = 7
159:       call RGPolygonSetVertices(rg,n,vr,vi,ierr)
160:       call RGSetFromOptions(rg,ierr)
161:       call RGView(rg,PETSC_NULL_VIEWER,ierr)
162:       re = 5
163:       im = 0.9
164: #if defined(PETSC_USE_COMPLEX)
165:       ar = re+im*PETSC_i
166:       ai = 0.0
167: #else
168:       ar = re
169:       ai = im
170: #endif
171:       call RGCheckInside(rg,one,ar,ai,inside,ierr)
172:       if (rank .eq. 0) then
173:         if (inside >= 0) then
174:           write(*,110) re, im, 'inside'
175:         else
176:           write(*,110) re, im, 'outside'
177:         endif
178:       endif

180: !     *** Clean up
181:       call RGDestroy(rg,ierr)
182:       call SlepcFinalize(ierr)
183:       end

185: !/*TEST
186: !
187: !   test:
188: !      suffix: 1
189: !      requires: !complex
190: !
191: !   test:
192: !      suffix: 1_complex
193: !      requires: complex
194: !
195: !TEST*/