Book a Demo!
CoCalc Logo Icon
StoreFeaturesDocsShareSupportNewsAboutPoliciesSign UpSign In
ElmerCSC
GitHub Repository: ElmerCSC/elmerfem
Path: blob/devel/mathlibs/src/arpack/cneigh.f
5218 views
1
c\BeginDoc
2
c
3
c\Name: cneigh
4
c
5
c\Description:
6
c Compute the eigenvalues of the current upper Hessenberg matrix
7
c and the corresponding Ritz estimates given the current residual norm.
8
c
9
c\Usage:
10
c call cneigh
11
c ( RNORM, N, H, LDH, RITZ, BOUNDS, Q, LDQ, WORKL, RWORK, IERR )
12
c
13
c\Arguments
14
c RNORM Real scalar. (INPUT)
15
c Residual norm corresponding to the current upper Hessenberg
16
c matrix H.
17
c
18
c N Integer. (INPUT)
19
c Size of the matrix H.
20
c
21
c H Complex N by N array. (INPUT)
22
c H contains the current upper Hessenberg matrix.
23
c
24
c LDH Integer. (INPUT)
25
c Leading dimension of H exactly as declared in the calling
26
c program.
27
c
28
c RITZ Complex array of length N. (OUTPUT)
29
c On output, RITZ(1:N) contains the eigenvalues of H.
30
c
31
c BOUNDS Complex array of length N. (OUTPUT)
32
c On output, BOUNDS contains the Ritz estimates associated with
33
c the eigenvalues held in RITZ. This is equal to RNORM
34
c times the last components of the eigenvectors corresponding
35
c to the eigenvalues in RITZ.
36
c
37
c Q Complex N by N array. (WORKSPACE)
38
c Workspace needed to store the eigenvectors of H.
39
c
40
c LDQ Integer. (INPUT)
41
c Leading dimension of Q exactly as declared in the calling
42
c program.
43
c
44
c WORKL Complex work array of length N**2 + 3*N. (WORKSPACE)
45
c Private (replicated) array on each PE or array allocated on
46
c the front end. This is needed to keep the full Schur form
47
c of H and also in the calculation of the eigenvectors of H.
48
c
49
c RWORK Real work array of length N (WORKSPACE)
50
c Private (replicated) array on each PE or array allocated on
51
c the front end.
52
c
53
c IERR Integer. (OUTPUT)
54
c Error exit flag from clahqr or ctrevc.
55
c
56
c\EndDoc
57
c
58
c-----------------------------------------------------------------------
59
c
60
c\BeginLib
61
c
62
c\Local variables:
63
c xxxxxx Complex
64
c
65
c\Routines called:
66
c ivout ARPACK utility routine that prints integers.
67
c second ARPACK utility routine for timing.
68
c cmout ARPACK utility routine that prints matrices
69
c cvout ARPACK utility routine that prints vectors.
70
c svout ARPACK utility routine that prints vectors.
71
c clacpy LAPACK matrix copy routine.
72
c clahqr LAPACK routine to compute the Schur form of an
73
c upper Hessenberg matrix.
74
c claset LAPACK matrix initialization routine.
75
c ctrevc LAPACK routine to compute the eigenvectors of a matrix
76
c in upper triangular form
77
c ccopy Level 1 BLAS that copies one vector to another.
78
c csscal Level 1 BLAS that scales a complex vector by a real number.
79
c scnrm2 Level 1 BLAS that computes the norm of a vector.
80
c
81
c
82
c\Author
83
c Danny Sorensen Phuong Vu
84
c Richard Lehoucq CRPC / Rice University
85
c Dept. of Computational & Houston, Texas
86
c Applied Mathematics
87
c Rice University
88
c Houston, Texas
89
c
90
c\SCCS Information: @(#)
91
c FILE: neigh.F SID: 2.2 DATE OF SID: 4/20/96 RELEASE: 2
92
c
93
c\Remarks
94
c None
95
c
96
c\EndLib
97
c
98
c-----------------------------------------------------------------------
99
c
100
subroutine cneigh (rnorm, n, h, ldh, ritz, bounds,
101
& q, ldq, workl, rwork, ierr)
102
c
103
c %----------------------------------------------------%
104
c | Include files for debugging and timing information |
105
c %----------------------------------------------------%
106
c
107
include 'debug.h'
108
include 'stat.h'
109
c
110
c %------------------%
111
c | Scalar Arguments |
112
c %------------------%
113
c
114
integer ierr, n, ldh, ldq
115
Real
116
& rnorm
117
c
118
c %-----------------%
119
c | Array Arguments |
120
c %-----------------%
121
c
122
Complex
123
& bounds(n), h(ldh,n), q(ldq,n), ritz(n),
124
& workl(n*(n+3))
125
Real
126
& rwork(n)
127
c
128
c %------------%
129
c | Parameters |
130
c %------------%
131
c
132
Complex
133
& one, zero
134
Real
135
& rone
136
parameter (one = (1.0E+0, 0.0E+0), zero = (0.0E+0, 0.0E+0),
137
& rone = 1.0E+0)
138
c
139
c %------------------------%
140
c | Local Scalars & Arrays |
141
c %------------------------%
142
c
143
logical select(1)
144
integer j, msglvl
145
Complex
146
& vl(1)
147
Real
148
& temp
149
c
150
c %----------------------%
151
c | External Subroutines |
152
c %----------------------%
153
c
154
external clacpy, clahqr, ctrevc, ccopy,
155
& csscal, cmout, cvout, second
156
c
157
c %--------------------%
158
c | External Functions |
159
c %--------------------%
160
c
161
Real
162
& scnrm2
163
external scnrm2
164
c
165
c %-----------------------%
166
c | Executable Statements |
167
c %-----------------------%
168
c
169
c %-------------------------------%
170
c | Initialize timing statistics |
171
c | & message level for debugging |
172
c %-------------------------------%
173
c
174
call second (t0)
175
msglvl = mceigh
176
c
177
if (msglvl .gt. 2) then
178
call cmout (logfil, n, n, h, ldh, ndigit,
179
& '_neigh: Entering upper Hessenberg matrix H ')
180
end if
181
c
182
c %----------------------------------------------------------%
183
c | 1. Compute the eigenvalues, the last components of the |
184
c | corresponding Schur vectors and the full Schur form T |
185
c | of the current upper Hessenberg matrix H. |
186
c | clahqr returns the full Schur form of H |
187
c | in WORKL(1:N**2), and the Schur vectors in q. |
188
c %----------------------------------------------------------%
189
c
190
call clacpy ('All', n, n, h, ldh, workl, n)
191
call claset ('All', n, n, zero, one, q, ldq)
192
call clahqr (.true., .true., n, 1, n, workl, ldh, ritz,
193
& 1, n, q, ldq, ierr)
194
if (ierr .ne. 0) go to 9000
195
c
196
call ccopy (n, q(n-1,1), ldq, bounds, 1)
197
if (msglvl .gt. 1) then
198
call cvout (logfil, n, bounds, ndigit,
199
& '_neigh: last row of the Schur matrix for H')
200
end if
201
c
202
c %----------------------------------------------------------%
203
c | 2. Compute the eigenvectors of the full Schur form T and |
204
c | apply the Schur vectors to get the corresponding |
205
c | eigenvectors. |
206
c %----------------------------------------------------------%
207
c
208
call ctrevc ('Right', 'Back', select, n, workl, n, vl, n, q,
209
& ldq, n, n, workl(n*n+1), rwork, ierr)
210
c
211
if (ierr .ne. 0) go to 9000
212
c
213
c %------------------------------------------------%
214
c | Scale the returning eigenvectors so that their |
215
c | Euclidean norms are all one. LAPACK subroutine |
216
c | ctrevc returns each eigenvector normalized so |
217
c | that the element of largest magnitude has |
218
c | magnitude 1; here the magnitude of a complex |
219
c | number (x,y) is taken to be |x| + |y|. |
220
c %------------------------------------------------%
221
c
222
do 10 j=1, n
223
temp = scnrm2( n, q(1,j), 1 )
224
call csscal ( n, rone / temp, q(1,j), 1 )
225
10 continue
226
c
227
if (msglvl .gt. 1) then
228
call ccopy(n, q(n,1), ldq, workl, 1)
229
call cvout (logfil, n, workl, ndigit,
230
& '_neigh: Last row of the eigenvector matrix for H')
231
end if
232
c
233
c %----------------------------%
234
c | Compute the Ritz estimates |
235
c %----------------------------%
236
c
237
call ccopy(n, q(n,1), n, bounds, 1)
238
call csscal(n, rnorm, bounds, 1)
239
c
240
if (msglvl .gt. 2) then
241
call cvout (logfil, n, ritz, ndigit,
242
& '_neigh: The eigenvalues of H')
243
call cvout (logfil, n, bounds, ndigit,
244
& '_neigh: Ritz estimates for the eigenvalues of H')
245
end if
246
c
247
call second(t1)
248
tceigh = tceigh + (t1 - t0)
249
c
250
9000 continue
251
return
252
c
253
c %---------------%
254
c | End of cneigh |
255
c %---------------%
256
c
257
end
258
259