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: #include <slepc/finclude/slepcrg.h>
17: program test1f
18: use slepcrg
19: implicit none
21: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: ! Declarations
23: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
25: RG :: rg
26: PetscInt :: i, n
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: PetscCallA(SlepcInitialize(PETSC_NULL_CHARACTER, ierr))
41: PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rank, ierr))
42: PetscCallA(RGCreate(PETSC_COMM_WORLD, rg, ierr))
44: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
45: ! Ellipse
46: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
48: PetscCallA(RGSetType(rg, RGELLIPSE, ierr))
49: center = 1.1
50: radius = 2
51: vscale = 0.1
52: PetscCallA(RGEllipseSetParameters(rg, center, radius, vscale, ierr))
53: PetscCallA(RGSetFromOptions(rg, ierr))
54: PetscCallA(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: PetscCallA(RGCheckInside(rg, 1_PETSC_INT_KIND, [ar], [ai], inside, ierr))
65: if (rank == 0) then
66: if (inside(1) >= 0) then
67: write (*, '(a,f4.1,a,f4.1,a)') 'Point (', re, ',', im, ') is inside the region'
68: else
69: write (*, '(a,f4.1,a,f4.1,a)') 'Point (', re, ',', im, ') is outside the region'
70: end if
71: end if
73: PetscCallA(RGComputeBoundingBox(rg, a, b, c, d, ierr))
74: if (rank == 0) then
75: write (*, '(a,f4.1,a,f4.1,a,f4.1,a,f4.1,a)') 'Bounding box: [', a, ',', b, ']x[', c, ',', d, ']'
76: end if
78: if (rank == 0) then
79: write (*, *) 'Contour points:'
80: end if
81: n = 10
82: PetscCallA(RGComputeContour(rg, n, cr, ci, ierr))
83: do i = 1, n
84: #if defined(PETSC_USE_COMPLEX)
85: re = PetscRealPart(cr(i))
86: im = PetscImaginaryPart(cr(i))
87: #else
88: re = cr(i)
89: im = ci(i)
90: #endif
91: if (rank == 0) then
92: write (*, '(a,f7.4,a,f7.4,a)') '(', re, ',', im, ')'
93: end if
94: end do
96: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
97: ! Interval
98: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
100: PetscCallA(RGSetType(rg, RGINTERVAL, ierr))
101: a = -1
102: b = 1
103: c = -0.1
104: d = 0.1
105: PetscCallA(RGIntervalSetEndpoints(rg, a, b, c, d, ierr))
106: PetscCallA(RGSetFromOptions(rg, ierr))
107: PetscCallA(RGView(rg, PETSC_NULL_VIEWER, ierr))
108: re = 0.2
109: im = 0
110: #if defined(PETSC_USE_COMPLEX)
111: ar = re + im*PETSC_i
112: ai = 0.0
113: #else
114: ar = re
115: ai = im
116: #endif
117: PetscCallA(RGCheckInside(rg, 1_PETSC_INT_KIND, [ar], [ai], inside, ierr))
118: if (rank == 0) then
119: if (inside(1) >= 0) then
120: write (*, '(a,f4.1,a,f4.1,a)') 'Point (', re, ',', im, ') is inside the region'
121: else
122: write (*, '(a,f4.1,a,f4.1,a)') 'Point (', re, ',', im, ') is outside the region'
123: end if
124: end if
126: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
127: ! Polygon
128: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
130: #if defined(PETSC_USE_COMPLEX)
131: vr(1) = (0.0, 2.0)
132: vr(2) = (1.0, 4.0)
133: vr(3) = (2.0, 5.0)
134: vr(4) = (4.0, 3.0)
135: vr(5) = (5.0, 4.0)
136: vr(6) = (6.0, 1.0)
137: vr(7) = (2.0, 0.0)
138: #else
139: vr(1) = 0.0
140: vi(1) = 1.0
141: vr(2) = 0.0
142: vi(2) = -1.0
143: vr(3) = 0.6
144: vi(3) = -0.8
145: vr(4) = 1.0
146: vi(4) = -1.0
147: vr(5) = 2.0
148: vi(5) = 0.0
149: vr(6) = 1.0
150: vi(6) = 1.0
151: vr(7) = 0.6
152: vi(7) = 0.8
153: #endif
154: PetscCallA(RGSetType(rg, RGPOLYGON, ierr))
155: n = 7
156: PetscCallA(RGPolygonSetVertices(rg, n, vr, vi, ierr))
157: PetscCallA(RGSetFromOptions(rg, ierr))
158: PetscCallA(RGView(rg, PETSC_NULL_VIEWER, ierr))
159: re = 5
160: im = 0.9
161: #if defined(PETSC_USE_COMPLEX)
162: ar = re + im*PETSC_i
163: ai = 0.0
164: #else
165: ar = re
166: ai = im
167: #endif
168: PetscCallA(RGCheckInside(rg, 1_PETSC_INT_KIND, [ar], [ai], inside, ierr))
169: if (rank == 0) then
170: if (inside(1) >= 0) then
171: write (*, '(a,f4.1,a,f4.1,a)') 'Point (', re, ',', im, ') is inside the region'
172: else
173: write (*, '(a,f4.1,a,f4.1,a)') 'Point (', re, ',', im, ') is outside the region'
174: end if
175: end if
177: ! *** Clean up
178: PetscCallA(RGDestroy(rg, ierr))
179: PetscCallA(SlepcFinalize(ierr))
180: end program test1f
182: !/*TEST
183: !
184: ! test:
185: ! suffix: 1
186: ! requires: !complex
187: !
188: ! test:
189: ! suffix: 1_complex
190: ! requires: complex
191: !
192: !TEST*/