LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ zlaein()

subroutine zlaein ( logical rightv,
logical noinit,
integer n,
complex*16, dimension( ldh, * ) h,
integer ldh,
complex*16 w,
complex*16, dimension( * ) v,
complex*16, dimension( ldb, * ) b,
integer ldb,
double precision, dimension( * ) rwork,
double precision eps3,
double precision smlnum,
integer info )

ZLAEIN computes a specified right or left eigenvector of an upper Hessenberg matrix by inverse iteration.

Download ZLAEIN + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!> ZLAEIN uses inverse iteration to find a right or left eigenvector
!> corresponding to the eigenvalue W of a complex upper Hessenberg
!> matrix H.
!> 
Parameters
[in]RIGHTV
!>          RIGHTV is LOGICAL
!>          = .TRUE. : compute right eigenvector;
!>          = .FALSE.: compute left eigenvector.
!> 
[in]NOINIT
!>          NOINIT is LOGICAL
!>          = .TRUE. : no initial vector supplied in V
!>          = .FALSE.: initial vector supplied in V.
!> 
[in]N
!>          N is INTEGER
!>          The order of the matrix H.  N >= 0.
!> 
[in]H
!>          H is COMPLEX*16 array, dimension (LDH,N)
!>          The upper Hessenberg matrix H.
!> 
[in]LDH
!>          LDH is INTEGER
!>          The leading dimension of the array H.  LDH >= max(1,N).
!> 
[in]W
!>          W is COMPLEX*16
!>          The eigenvalue of H whose corresponding right or left
!>          eigenvector is to be computed.
!> 
[in,out]V
!>          V is COMPLEX*16 array, dimension (N)
!>          On entry, if NOINIT = .FALSE., V must contain a starting
!>          vector for inverse iteration; otherwise V need not be set.
!>          On exit, V contains the computed eigenvector, normalized so
!>          that the component of largest magnitude has magnitude 1; here
!>          the magnitude of a complex number (x,y) is taken to be
!>          |x| + |y|.
!> 
[out]B
!>          B is COMPLEX*16 array, dimension (LDB,N)
!> 
[in]LDB
!>          LDB is INTEGER
!>          The leading dimension of the array B.  LDB >= max(1,N).
!> 
[out]RWORK
!>          RWORK is DOUBLE PRECISION array, dimension (N)
!> 
[in]EPS3
!>          EPS3 is DOUBLE PRECISION
!>          A small machine-dependent value which is used to perturb
!>          close eigenvalues, and to replace zero pivots.
!> 
[in]SMLNUM
!>          SMLNUM is DOUBLE PRECISION
!>          A machine-dependent value close to the underflow threshold.
!> 
[out]INFO
!>          INFO is INTEGER
!>          = 0:  successful exit
!>          = 1:  inverse iteration did not converge; V is set to the
!>                last iterate.
!> 
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 145 of file zlaein.f.

148*
149* -- LAPACK auxiliary routine --
150* -- LAPACK is a software package provided by Univ. of Tennessee, --
151* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
152*
153* .. Scalar Arguments ..
154 LOGICAL NOINIT, RIGHTV
155 INTEGER INFO, LDB, LDH, N
156 DOUBLE PRECISION EPS3, SMLNUM
157 COMPLEX*16 W
158* ..
159* .. Array Arguments ..
160 DOUBLE PRECISION RWORK( * )
161 COMPLEX*16 B( LDB, * ), H( LDH, * ), V( * )
162* ..
163*
164* =====================================================================
165*
166* .. Parameters ..
167 DOUBLE PRECISION ONE, TENTH
168 parameter( one = 1.0d+0, tenth = 1.0d-1 )
169 COMPLEX*16 ZERO
170 parameter( zero = ( 0.0d+0, 0.0d+0 ) )
171* ..
172* .. Local Scalars ..
173 CHARACTER NORMIN, TRANS
174 INTEGER I, IERR, ITS, J
175 DOUBLE PRECISION GROWTO, NRMSML, ROOTN, RTEMP, SCALE, VNORM
176 COMPLEX*16 CDUM, EI, EJ, TEMP, X
177* ..
178* .. External Functions ..
179 INTEGER IZAMAX
180 DOUBLE PRECISION DZASUM, DZNRM2
181 COMPLEX*16 ZLADIV
182 EXTERNAL izamax, dzasum, dznrm2, zladiv
183* ..
184* .. External Subroutines ..
185 EXTERNAL zdscal, zlatrs
186* ..
187* .. Intrinsic Functions ..
188 INTRINSIC abs, dble, dimag, max, sqrt
189* ..
190* .. Statement Functions ..
191 DOUBLE PRECISION CABS1
192* ..
193* .. Statement Function definitions ..
194 cabs1( cdum ) = abs( dble( cdum ) ) + abs( dimag( cdum ) )
195* ..
196* .. Executable Statements ..
197*
198 info = 0
199*
200* GROWTO is the threshold used in the acceptance test for an
201* eigenvector.
202*
203 rootn = sqrt( dble( n ) )
204 growto = tenth / rootn
205 nrmsml = max( one, eps3*rootn )*smlnum
206*
207* Form B = H - W*I (except that the subdiagonal elements are not
208* stored).
209*
210 DO 20 j = 1, n
211 DO 10 i = 1, j - 1
212 b( i, j ) = h( i, j )
213 10 CONTINUE
214 b( j, j ) = h( j, j ) - w
215 20 CONTINUE
216*
217 IF( noinit ) THEN
218*
219* Initialize V.
220*
221 DO 30 i = 1, n
222 v( i ) = eps3
223 30 CONTINUE
224 ELSE
225*
226* Scale supplied initial vector.
227*
228 vnorm = dznrm2( n, v, 1 )
229 CALL zdscal( n, ( eps3*rootn ) / max( vnorm, nrmsml ), v,
230 $ 1 )
231 END IF
232*
233 IF( rightv ) THEN
234*
235* LU decomposition with partial pivoting of B, replacing zero
236* pivots by EPS3.
237*
238 DO 60 i = 1, n - 1
239 ei = h( i+1, i )
240 IF( cabs1( b( i, i ) ).LT.cabs1( ei ) ) THEN
241*
242* Interchange rows and eliminate.
243*
244 x = zladiv( b( i, i ), ei )
245 b( i, i ) = ei
246 DO 40 j = i + 1, n
247 temp = b( i+1, j )
248 b( i+1, j ) = b( i, j ) - x*temp
249 b( i, j ) = temp
250 40 CONTINUE
251 ELSE
252*
253* Eliminate without interchange.
254*
255 IF( b( i, i ).EQ.zero )
256 $ b( i, i ) = eps3
257 x = zladiv( ei, b( i, i ) )
258 IF( x.NE.zero ) THEN
259 DO 50 j = i + 1, n
260 b( i+1, j ) = b( i+1, j ) - x*b( i, j )
261 50 CONTINUE
262 END IF
263 END IF
264 60 CONTINUE
265 IF( b( n, n ).EQ.zero )
266 $ b( n, n ) = eps3
267*
268 trans = 'N'
269*
270 ELSE
271*
272* UL decomposition with partial pivoting of B, replacing zero
273* pivots by EPS3.
274*
275 DO 90 j = n, 2, -1
276 ej = h( j, j-1 )
277 IF( cabs1( b( j, j ) ).LT.cabs1( ej ) ) THEN
278*
279* Interchange columns and eliminate.
280*
281 x = zladiv( b( j, j ), ej )
282 b( j, j ) = ej
283 DO 70 i = 1, j - 1
284 temp = b( i, j-1 )
285 b( i, j-1 ) = b( i, j ) - x*temp
286 b( i, j ) = temp
287 70 CONTINUE
288 ELSE
289*
290* Eliminate without interchange.
291*
292 IF( b( j, j ).EQ.zero )
293 $ b( j, j ) = eps3
294 x = zladiv( ej, b( j, j ) )
295 IF( x.NE.zero ) THEN
296 DO 80 i = 1, j - 1
297 b( i, j-1 ) = b( i, j-1 ) - x*b( i, j )
298 80 CONTINUE
299 END IF
300 END IF
301 90 CONTINUE
302 IF( b( 1, 1 ).EQ.zero )
303 $ b( 1, 1 ) = eps3
304*
305 trans = 'C'
306*
307 END IF
308*
309 normin = 'N'
310 DO 110 its = 1, n
311*
312* Solve U*x = scale*v for a right eigenvector
313* or U**H *x = scale*v for a left eigenvector,
314* overwriting x on v.
315*
316 CALL zlatrs( 'Upper', trans, 'Nonunit', normin, n, b, ldb,
317 $ v,
318 $ scale, rwork, ierr )
319 normin = 'Y'
320*
321* Test for sufficient growth in the norm of v.
322*
323 vnorm = dzasum( n, v, 1 )
324 IF( vnorm.GE.growto*scale )
325 $ GO TO 120
326*
327* Choose new orthogonal starting vector and try again.
328*
329 rtemp = eps3 / ( rootn+one )
330 v( 1 ) = eps3
331 DO 100 i = 2, n
332 v( i ) = rtemp
333 100 CONTINUE
334 v( n-its+1 ) = v( n-its+1 ) - eps3*rootn
335 110 CONTINUE
336*
337* Failure to find eigenvector in N iterations.
338*
339 info = 1
340*
341 120 CONTINUE
342*
343* Normalize eigenvector.
344*
345 i = izamax( n, v, 1 )
346 CALL zdscal( n, one / cabs1( v( i ) ), v, 1 )
347*
348 RETURN
349*
350* End of ZLAEIN
351*
double precision function dzasum(n, zx, incx)
DZASUM
Definition dzasum.f:72
integer function izamax(n, zx, incx)
IZAMAX
Definition izamax.f:71
complex *16 function zladiv(x, y)
ZLADIV performs complex division in real arithmetic, avoiding unnecessary overflow.
Definition zladiv.f:62
subroutine zlatrs(uplo, trans, diag, normin, n, a, lda, x, scale, cnorm, info)
ZLATRS solves a triangular system of equations with the scale factor set to prevent overflow.
Definition zlatrs.f:238
real(wp) function dznrm2(n, x, incx)
DZNRM2
Definition dznrm2.f90:90
subroutine zdscal(n, da, zx, incx)
ZDSCAL
Definition zdscal.f:78
Here is the call graph for this function:
Here is the caller graph for this function: