Actual source code: test1f.F90
1: !
2: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: ! SLEPc - Scalable Library for Eigenvalue Problem Computations
4: ! Copyright (c) 2002-, 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,one
27: PetscInt inside(1)
28: PetscMPIInt rank
29: PetscErrorCode ierr
30: PetscReal re,im
31: PetscScalar ar,ai,cr(10),ci(10)
32: PetscScalar vr(7),vi(7)
33: PetscScalar center
34: PetscReal radius,vscale,a,b,c,d
36: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
37: ! Beginning of program
38: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
40: one = 1
41: PetscCallA(SlepcInitialize(PETSC_NULL_CHARACTER,ierr))
42: PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr))
43: PetscCallA(RGCreate(PETSC_COMM_WORLD,rg,ierr))
45: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
46: ! Ellipse
47: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
49: PetscCallA(RGSetType(rg,RGELLIPSE,ierr))
50: center = 1.1
51: radius = 2
52: vscale = 0.1
53: PetscCallA(RGEllipseSetParameters(rg,center,radius,vscale,ierr))
54: PetscCallA(RGSetFromOptions(rg,ierr))
55: PetscCallA(RGView(rg,PETSC_NULL_VIEWER,ierr))
56: re = 0.1
57: im = 0.3
58: #if defined(PETSC_USE_COMPLEX)
59: ar = re+im*PETSC_i
60: ai = 0.0
61: #else
62: ar = re
63: ai = im
64: #endif
65: PetscCallA(RGCheckInside(rg,one,[ar],[ai],inside,ierr))
66: if (rank .eq. 0) then
67: if (inside(1) >= 0) then
68: write(*,110) re, im, 'inside'
69: else
70: write(*,110) re, im, 'outside'
71: endif
72: endif
73: 110 format ('Point (',F4.1,',',F4.1,') is ',A7,' the region')
75: PetscCallA(RGComputeBoundingBox(rg,a,b,c,d,ierr))
76: if (rank .eq. 0) then
77: write(*,115) a, b, c, d
78: endif
79: 115 format ('Bounding box: [',F4.1,',',F4.1,']x[',F4.1,',',F4.1,']')
81: if (rank .eq. 0) then
82: write (*,*) 'Contour points:'
83: endif
84: n = 10
85: PetscCallA(RGComputeContour(rg,n,cr,ci,ierr))
86: do i=1,n
87: #if defined(PETSC_USE_COMPLEX)
88: re = PetscRealPart(cr(i))
89: im = PetscImaginaryPart(cr(i))
90: #else
91: re = cr(i)
92: im = ci(i)
93: #endif
94: if (rank .eq. 0) then
95: write(*,120) re, im
96: endif
97: enddo
98: 120 format ('(',F7.4,',',F7.4,')')
100: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
101: ! Interval
102: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
104: PetscCallA(RGSetType(rg,RGINTERVAL,ierr))
105: a = -1
106: b = 1
107: c = -0.1
108: d = 0.1
109: PetscCallA(RGIntervalSetEndpoints(rg,a,b,c,d,ierr))
110: PetscCallA(RGSetFromOptions(rg,ierr))
111: PetscCallA(RGView(rg,PETSC_NULL_VIEWER,ierr))
112: re = 0.2
113: im = 0
114: #if defined(PETSC_USE_COMPLEX)
115: ar = re+im*PETSC_i
116: ai = 0.0
117: #else
118: ar = re
119: ai = im
120: #endif
121: PetscCallA(RGCheckInside(rg,one,[ar],[ai],inside,ierr))
122: if (rank .eq. 0) then
123: if (inside(1) >= 0) then
124: write(*,110) re, im, 'inside'
125: else
126: write(*,110) re, im, 'outside'
127: endif
128: endif
130: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
131: ! Polygon
132: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
134: #if defined(PETSC_USE_COMPLEX)
135: vr(1) = (0.0,2.0)
136: vr(2) = (1.0,4.0)
137: vr(3) = (2.0,5.0)
138: vr(4) = (4.0,3.0)
139: vr(5) = (5.0,4.0)
140: vr(6) = (6.0,1.0)
141: vr(7) = (2.0,0.0)
142: #else
143: vr(1) = 0.0
144: vi(1) = 1.0
145: vr(2) = 0.0
146: vi(2) = -1.0
147: vr(3) = 0.6
148: vi(3) = -0.8
149: vr(4) = 1.0
150: vi(4) = -1.0
151: vr(5) = 2.0
152: vi(5) = 0.0
153: vr(6) = 1.0
154: vi(6) = 1.0
155: vr(7) = 0.6
156: vi(7) = 0.8
157: #endif
158: PetscCallA(RGSetType(rg,RGPOLYGON,ierr))
159: n = 7
160: PetscCallA(RGPolygonSetVertices(rg,n,vr,vi,ierr))
161: PetscCallA(RGSetFromOptions(rg,ierr))
162: PetscCallA(RGView(rg,PETSC_NULL_VIEWER,ierr))
163: re = 5
164: im = 0.9
165: #if defined(PETSC_USE_COMPLEX)
166: ar = re+im*PETSC_i
167: ai = 0.0
168: #else
169: ar = re
170: ai = im
171: #endif
172: PetscCallA(RGCheckInside(rg,one,[ar],[ai],inside,ierr))
173: if (rank .eq. 0) then
174: if (inside(1) >= 0) then
175: write(*,110) re, im, 'inside'
176: else
177: write(*,110) re, im, 'outside'
178: endif
179: endif
181: ! *** Clean up
182: PetscCallA(RGDestroy(rg,ierr))
183: PetscCallA(SlepcFinalize(ierr))
184: end
186: !/*TEST
187: !
188: ! test:
189: ! suffix: 1
190: ! requires: !complex
191: !
192: ! test:
193: ! suffix: 1_complex
194: ! requires: complex
195: !
196: !TEST*/