slaein.f.html
来自「famous linear algebra library (LAPACK) p」· HTML 代码 · 共 556 行 · 第 1/3 页
HTML
556 行
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html>
<head>
<title>slaein.f</title>
<meta name="generator" content="emacs 21.3.1; htmlfontify 0.20">
<style type="text/css"><!--
body { background: rgb(255, 255, 255); color: rgb(0, 0, 0); font-style: normal; font-weight: 500; font-stretch: normal; font-family: adobe-courier; font-size: 11pt; text-decoration: none; }
span.default { background: rgb(255, 255, 255); color: rgb(0, 0, 0); font-style: normal; font-weight: 500; font-stretch: normal; font-family: adobe-courier; font-size: 11pt; text-decoration: none; }
span.default a { background: rgb(255, 255, 255); color: rgb(0, 0, 0); font-style: normal; font-weight: 500; font-stretch: normal; font-family: adobe-courier; font-size: 11pt; text-decoration: underline; }
span.string { color: rgb(188, 143, 143); background: rgb(255, 255, 255); font-style: normal; font-weight: 500; font-stretch: normal; font-family: adobe-courier; font-size: 11pt; text-decoration: none; }
span.string a { color: rgb(188, 143, 143); background: rgb(255, 255, 255); font-style: normal; font-weight: 500; font-stretch: normal; font-family: adobe-courier; font-size: 11pt; text-decoration: underline; }
span.comment { color: rgb(178, 34, 34); background: rgb(255, 255, 255); font-style: normal; font-weight: 500; font-stretch: normal; font-family: adobe-courier; font-size: 11pt; text-decoration: none; }
span.comment a { color: rgb(178, 34, 34); background: rgb(255, 255, 255); font-style: normal; font-weight: 500; font-stretch: normal; font-family: adobe-courier; font-size: 11pt; text-decoration: underline; }
--></style>
</head>
<body>
<pre>
SUBROUTINE <a name="SLAEIN.1"></a><a href="slaein.f.html#SLAEIN.1">SLAEIN</a>( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, B,
$ LDB, WORK, EPS3, SMLNUM, BIGNUM, INFO )
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> -- LAPACK auxiliary routine (version 3.1) --
</span><span class="comment">*</span><span class="comment"> Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
</span><span class="comment">*</span><span class="comment"> November 2006
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> .. Scalar Arguments ..
</span> LOGICAL NOINIT, RIGHTV
INTEGER INFO, LDB, LDH, N
REAL BIGNUM, EPS3, SMLNUM, WI, WR
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment"> .. Array Arguments ..
</span> REAL B( LDB, * ), H( LDH, * ), VI( * ), VR( * ),
$ WORK( * )
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Purpose
</span><span class="comment">*</span><span class="comment"> =======
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> <a name="SLAEIN.21"></a><a href="slaein.f.html#SLAEIN.1">SLAEIN</a> uses inverse iteration to find a right or left eigenvector
</span><span class="comment">*</span><span class="comment"> corresponding to the eigenvalue (WR,WI) of a real upper Hessenberg
</span><span class="comment">*</span><span class="comment"> matrix H.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Arguments
</span><span class="comment">*</span><span class="comment"> =========
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> RIGHTV (input) LOGICAL
</span><span class="comment">*</span><span class="comment"> = .TRUE. : compute right eigenvector;
</span><span class="comment">*</span><span class="comment"> = .FALSE.: compute left eigenvector.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> NOINIT (input) LOGICAL
</span><span class="comment">*</span><span class="comment"> = .TRUE. : no initial vector supplied in (VR,VI).
</span><span class="comment">*</span><span class="comment"> = .FALSE.: initial vector supplied in (VR,VI).
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> N (input) INTEGER
</span><span class="comment">*</span><span class="comment"> The order of the matrix H. N >= 0.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> H (input) REAL array, dimension (LDH,N)
</span><span class="comment">*</span><span class="comment"> The upper Hessenberg matrix H.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> LDH (input) INTEGER
</span><span class="comment">*</span><span class="comment"> The leading dimension of the array H. LDH >= max(1,N).
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> WR (input) REAL
</span><span class="comment">*</span><span class="comment"> WI (input) REAL
</span><span class="comment">*</span><span class="comment"> The real and imaginary parts of the eigenvalue of H whose
</span><span class="comment">*</span><span class="comment"> corresponding right or left eigenvector is to be computed.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> VR (input/output) REAL array, dimension (N)
</span><span class="comment">*</span><span class="comment"> VI (input/output) REAL array, dimension (N)
</span><span class="comment">*</span><span class="comment"> On entry, if NOINIT = .FALSE. and WI = 0.0, VR must contain
</span><span class="comment">*</span><span class="comment"> a real starting vector for inverse iteration using the real
</span><span class="comment">*</span><span class="comment"> eigenvalue WR; if NOINIT = .FALSE. and WI.ne.0.0, VR and VI
</span><span class="comment">*</span><span class="comment"> must contain the real and imaginary parts of a complex
</span><span class="comment">*</span><span class="comment"> starting vector for inverse iteration using the complex
</span><span class="comment">*</span><span class="comment"> eigenvalue (WR,WI); otherwise VR and VI need not be set.
</span><span class="comment">*</span><span class="comment"> On exit, if WI = 0.0 (real eigenvalue), VR contains the
</span><span class="comment">*</span><span class="comment"> computed real eigenvector; if WI.ne.0.0 (complex eigenvalue),
</span><span class="comment">*</span><span class="comment"> VR and VI contain the real and imaginary parts of the
</span><span class="comment">*</span><span class="comment"> computed complex eigenvector. The eigenvector is normalized
</span><span class="comment">*</span><span class="comment"> so that the component of largest magnitude has magnitude 1;
</span><span class="comment">*</span><span class="comment"> here the magnitude of a complex number (x,y) is taken to be
</span><span class="comment">*</span><span class="comment"> |x| + |y|.
</span><span class="comment">*</span><span class="comment"> VI is not referenced if WI = 0.0.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> B (workspace) REAL array, dimension (LDB,N)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> LDB (input) INTEGER
</span><span class="comment">*</span><span class="comment"> The leading dimension of the array B. LDB >= N+1.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> WORK (workspace) REAL array, dimension (N)
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> EPS3 (input) REAL
</span><span class="comment">*</span><span class="comment"> A small machine-dependent value which is used to perturb
</span><span class="comment">*</span><span class="comment"> close eigenvalues, and to replace zero pivots.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> SMLNUM (input) REAL
</span><span class="comment">*</span><span class="comment"> A machine-dependent value close to the underflow threshold.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> BIGNUM (input) REAL
</span><span class="comment">*</span><span class="comment"> A machine-dependent value close to the overflow threshold.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> INFO (output) INTEGER
</span><span class="comment">*</span><span class="comment"> = 0: successful exit
</span><span class="comment">*</span><span class="comment"> = 1: inverse iteration did not converge; VR is set to the
</span><span class="comment">*</span><span class="comment"> last iterate, and so is VI if WI.ne.0.0.
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> =====================================================================
</span><span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> .. Parameters ..
</span> REAL ZERO, ONE, TENTH
PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0, TENTH = 1.0E-1 )
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment"> .. Local Scalars ..
</span> CHARACTER NORMIN, TRANS
INTEGER I, I1, I2, I3, IERR, ITS, J
REAL ABSBII, ABSBJJ, EI, EJ, GROWTO, NORM, NRMSML,
$ REC, ROOTN, SCALE, TEMP, VCRIT, VMAX, VNORM, W,
$ W1, X, XI, XR, Y
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment"> .. External Functions ..
</span> INTEGER ISAMAX
REAL SASUM, <a name="SLAPY2.104"></a><a href="slapy2.f.html#SLAPY2.1">SLAPY2</a>, SNRM2
EXTERNAL ISAMAX, SASUM, <a name="SLAPY2.105"></a><a href="slapy2.f.html#SLAPY2.1">SLAPY2</a>, SNRM2
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment"> .. External Subroutines ..
</span> EXTERNAL <a name="SLADIV.108"></a><a href="sladiv.f.html#SLADIV.1">SLADIV</a>, <a name="SLATRS.108"></a><a href="slatrs.f.html#SLATRS.1">SLATRS</a>, SSCAL
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment"> .. Intrinsic Functions ..
</span> INTRINSIC ABS, MAX, REAL, SQRT
<span class="comment">*</span><span class="comment"> ..
</span><span class="comment">*</span><span class="comment"> .. Executable Statements ..
</span><span class="comment">*</span><span class="comment">
</span> INFO = 0
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> GROWTO is the threshold used in the acceptance test for an
</span><span class="comment">*</span><span class="comment"> eigenvector.
</span><span class="comment">*</span><span class="comment">
</span> ROOTN = SQRT( REAL( N ) )
GROWTO = TENTH / ROOTN
NRMSML = MAX( ONE, EPS3*ROOTN )*SMLNUM
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Form B = H - (WR,WI)*I (except that the subdiagonal elements and
</span><span class="comment">*</span><span class="comment"> the imaginary parts of the diagonal elements are not stored).
</span><span class="comment">*</span><span class="comment">
</span> DO 20 J = 1, N
DO 10 I = 1, J - 1
B( I, J ) = H( I, J )
10 CONTINUE
B( J, J ) = H( J, J ) - WR
20 CONTINUE
<span class="comment">*</span><span class="comment">
</span> IF( WI.EQ.ZERO ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Real eigenvalue.
</span><span class="comment">*</span><span class="comment">
</span> IF( NOINIT ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Set initial vector.
</span><span class="comment">*</span><span class="comment">
</span> DO 30 I = 1, N
VR( I ) = EPS3
30 CONTINUE
ELSE
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Scale supplied initial vector.
</span><span class="comment">*</span><span class="comment">
</span> VNORM = SNRM2( N, VR, 1 )
CALL SSCAL( N, ( EPS3*ROOTN ) / MAX( VNORM, NRMSML ), VR,
$ 1 )
END IF
<span class="comment">*</span><span class="comment">
</span> IF( RIGHTV ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> LU decomposition with partial pivoting of B, replacing zero
</span><span class="comment">*</span><span class="comment"> pivots by EPS3.
</span><span class="comment">*</span><span class="comment">
</span> DO 60 I = 1, N - 1
EI = H( I+1, I )
IF( ABS( B( I, I ) ).LT.ABS( EI ) ) THEN
<span class="comment">*</span><span class="comment">
</span><span class="comment">*</span><span class="comment"> Interchange rows and eliminate.
</span><span class="comment">*</span><span class="comment">
</span> X = B( I, I ) / EI
B( I, I ) = EI
DO 40 J = I + 1, N
⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?