Routine SSTEQR2() Documented By Universal Report

Routine: SSTEQR2()

# lines:	498
# code:	498
# comment:	0
# blank:	0
# Variables:	32
# Callers:	1
# Callings:	0
# Words:	188
# Keywords:	100

..
     .. Local Scalars ..
     ..
     .. External Functions ..
     ..
     .. External Subroutines ..
     ..
     .. Intrinsic Functions ..
     ..
     .. Executable Statements ..
     Test the input parameters.
     Quick return if possible
     Determine the unit roundoff and over/underflow thresholds.
     Compute the eigenvalues and eigenvectors of the tridiagonal
     matrix.
     Determine where the matrix splits and choose QL or QR iteration
     for each block, according to whether top or bottom diagonal
     element is smaller.
     Scale submatrix in rows and columns L to LEND
     Choose between QL and QR iteration
        QL Iteration
        Look for small subdiagonal element.
        If remaining matrix is 2-by-2, use SLAE2 or SLAEV2
        to compute its eigensystem.
        Form shift.
        Inner loop
           If eigenvectors are desired, then save rotations.
        If eigenvectors are desired, then apply saved rotations.
        Eigenvalue found.

Display dynamic version

FontSize: - +

001        SUBROUTINE SSTEQR2( COMPZ , N , D , E , Z , LDZ , NR , WORK , INFO )
002  
003  *     -- LAPACK routine(version 2.0) --
004  *     Univ. of Tennessee , Univ. of California Berkeley , NAG Ltd. ,
005  *     Courant Institute , Argonne National Lab , and Rice University
006  *     September 30 , 1994
007  
008  *     .. Scalar Arguments ..
009        CHARACTER COMPZ
010        INTEGER INFO , LDZ , N , NR
011  *     ..
012  *     .. Array Arguments ..
013        REAL D( * ) , E( * ) , WORK( * ) , Z( LDZ , * )
014  *     ..
015  
016  *     Purpose
017  *     === ====
018  
019  *     SSTEQR2 is a modified version of LAPACK routine SSTEQR.
020  *     SSTEQR2 computes all eigenvalues and , optionally , eigenvectors of a
021  *     symmetric tridiagonal matrix using the implicit QL or QR method.
022  *     running SSTEQR2 to perform updates on a distributed matrix Q.
023  *     Proper usage of SSTEQR2 can be gleaned from examination of ScaLAPACK's
024  *     PSSYEV.
025  
026  *     Arguments
027  *     === ======
028  
029  *     COMPZ(input) CHARACTER*1
030  *     = 'N' : Compute eigenvalues only.
031  *     = 'I' : Compute eigenvalues and eigenvectors of the
032  *     tridiagonal matrix. Z must be initialized to the
033  *     identity matrix by PDLASET or DLASET prior to entering
034  *     this subroutine.
035  
036  *     N(input) INTEGER
037  *     The order of the matrix. N >= 0.
038  
039  *     D(input / output) REAL array , dimension(N)
040  *     On entry , the diagonal elements of the tridiagonal matrix.
041  *     On exit , if INFO = 0 , the eigenvalues in ascending order.
042  
043  *     E(input / output) REAL array , dimension(N - 1)
044  *     On entry , the(n - 1) subdiagonal elements of the tridiagonal
045  *     matrix.
046  *     On exit , E has been destroyed.
047  
048  *     Z(local input / local output) REAL array , global
049  *     dimension(N , N) , local dimension(LDZ , NR).
050  *     On entry , if COMPZ = 'V' , then Z contains the orthogonal
051  *     matrix used in the reduction to tridiagonal form.
052  *     On exit , if INFO = 0 , then if COMPZ = 'V' , Z contains the
053  *     orthonormal eigenvectors of the original symmetric matrix ,
054  *     and if COMPZ = 'I' , Z contains the orthonormal eigenvectors
055  *     of the symmetric tridiagonal matrix.
056  *     If COMPZ = 'N' , then Z is not referenced.
057  
058  *     LDZ(input) INTEGER
059  *     The leading dimension of the array Z. LDZ >= 1 , and if
060  *     eigenvectors are desired , then LDZ >= max(1 , N).
061  
062  *     NR(input) INTEGER
063  *     NR = MAX(1 , NUMROC( N , NB , MYPROW , 0 , NPROCS ) ).
064  *     If COMPZ = 'N' , then NR is not referenced.
065  
066  *     WORK(workspace) REAL array , dimension(max(1 , 2*N - 2))
067  *     If COMPZ = 'N' , then WORK is not referenced.
068  
069  *     INFO(output) INTEGER
070  *     = 0 : successful exit
071  *     < 0 : if INFO = - i , the i - th argument had an illegal value
072  *     > 0 : the algorithm has failed to find all the eigenvalues in
073  *     a total of 30*N iterations ; if INFO = i , then i
074  *     elements of E have not converged to zero ; on exit , D
075  *     and E contain the elements of a symmetric tridiagonal
076  *     matrix which is orthogonally similar to the original
077  *     matrix.
078  
079  *     === ==================================================================
080  
081  *     .. Parameters ..
082        REAL ZERO , ONE , TWO , THREE
083        PARAMETER( ZERO = 0.0E0 , ONE = 1.0E0 , TWO = 2.0E0 ,
084       $THREE = 3.0E0 )
085        INTEGER MAXIT
086        PARAMETER( MAXIT = 30 )
087        ELSE
088  
089  *         QR Iteration
090  
091  *         Look for small superdiagonal element.
092  
093     90 CONTINUE
094        IF( L.NE.LEND ) THEN
094               
095            LENDP1 = LEND + 1
096            DO 100 M = L , LENDP1 , - 1
096                   
097                TST = ABS( E( M - 1 ) )**2
098                IF( TST.LE.( EPS2*ABS( D( M ) ) )*ABS( D( M - 1 ) ) +
099       $        SAFMIN )GO TO 110
100    100     CONTINUE
101        END IF
102  
103        M = LEND
104  
105    110 CONTINUE
106        IF( M.GT.LEND )
106               
107       $    E( M - 1 ) = ZERO
108            P = D( L )
109            IF( M.EQ.L )
109                   
110       $        GO TO 130
111  
112  *             If remaining matrix is 2 - by - 2 , use SLAE2 or SLAEV2
113  *             to compute its eigensystem.
114  
115                IF( M.EQ.L - 1 ) THEN
115                       
116                    IF( ICOMPZ.GT.0 ) THEN
116                           
117                        CALL SLAEV2( D( L - 1 ) , E( L - 1 ) , D( L ) , RT1 , RT2 , C , S )
118                        WORK( M ) = C
119                        WORK( N - 1 + M ) = S
120                        CALL SLASR( 'R' , 'V' , 'F' , NR , 2 , WORK( M ) ,
121       $                WORK( N - 1 + M ) , Z( 1 , L - 1 ) , LDZ )
122                    ELSE
122                           
123                        CALL SLAE2( D( L - 1 ) , E( L - 1 ) , D( L ) , RT1 , RT2 )
124                    END IF
125                    D( L - 1 ) = RT1
126                    D( L ) = RT2
127                    E( L - 1 ) = ZERO
128                    L = L - 2
129                    IF( L.GE.LEND )
129                           
130       $                GO TO 90
131                        GO TO 140
132                    END IF
133  
134                    IF( JTOT.EQ.NMAXIT )
134                           
135       $                GO TO 140
136                        JTOT = JTOT + 1
137  
138  *                     Form shift.
139  
140                        G =( D( L - 1 ) - P ) / ( TWO*E( L - 1 ) )
141                        R = SLAPY2( G , ONE )
142                        G = D( M ) - P + ( E( L - 1 ) / ( G + SIGN( R , G ) ) )
143  
144                        S = ONE
145                        C = ONE
146                        P = ZERO
147  
148  *                     Inner loop
149  
150                        LM1 = L - 1
151                        DO 120 I = M , LM1
151                               
152                            F = S*E( I )
153                            B = C*E( I )
154                            CALL SLARTG( G , F , C , S , R )
155                            IF( I.NE.M )
155                                   
156       $                        E( I - 1 ) = R
157                                G = D( I ) - P
158                                R =( D( I + 1 ) - G )*S + TWO*C*B
159                                P = S*R
160                                D( I ) = G + P
161                                G = C*R - B
162  
163  *                             If eigenvectors are desired , then save rotations.
164  
165                                IF( ICOMPZ.GT.0 ) THEN
165                                       
166                                    WORK( I ) = C
167                                    WORK( N - 1 + I ) = S
168                                END IF
169  
170    120                 CONTINUE
171  
172  *                     If eigenvectors are desired , then apply saved rotations.
173  
173                           
174                        IF( ICOMPZ.GT.0 ) THEN
174                               
175                            MM = L - M + 1
176                            CALL SLASR( 'R' , 'V' , 'F' , NR , MM , WORK( M ) , WORK( N - 1 + M ) ,
177       $                    Z( 1 , M ) , LDZ )
178                        END IF
179  
180                        D( L ) = D( L ) - P
181                        E( LM1 ) = G
182                        GO TO 90
183  
184  *                     Eigenvalue found.
185  
186    130 CONTINUE
187        D( L ) = P
188  
189        L = L - 1
190        IF( L.GE.LEND )
190               
191       $    GO TO 90
192            GO TO 140
193  
194        END IF
195  
196  *     Undo scaling if necessary
197  
198    140 CONTINUE
199        IF( ISCALE.EQ.1 ) THEN
199               
200            CALL SLASCL( 'G' , 0 , 0 , SSFMAX , ANORM , LENDSV - LSV + 1 , 1 ,
201       $    D( LSV ) , N , INFO )
202            CALL SLASCL( 'G' , 0 , 0 , SSFMAX , ANORM , LENDSV - LSV , 1 , E( LSV ) ,
203       $    N , INFO )
204        ELSE IF( ISCALE.EQ.2 ) THEN
204               
205            CALL SLASCL( 'G' , 0 , 0 , SSFMIN , ANORM , LENDSV - LSV + 1 , 1 ,
206       $    D( LSV ) , N , INFO )
207            CALL SLASCL( 'G' , 0 , 0 , SSFMIN , ANORM , LENDSV - LSV , 1 , E( LSV ) ,
208       $    N , INFO )
209        END IF
210  
211  *     Check for no convergence to an eigenvalue after a total
212  *     of N*MAXIT iterations.
213  
214        IF( JTOT.LT.NMAXIT )
214               
215       $    GO TO 10
216            DO 150 I = 1 , N - 1
216                   
217                IF( E( I ).NE.ZERO )
217                       
218       $            INFO = INFO + 1
219    150     CONTINUE
219               
220            GO TO 190
221  
222  *         Order eigenvalues and eigenvectors.
223  
224    160 CONTINUE
225        IF( ICOMPZ.EQ.0 ) THEN
226  
227  *         Use Quick Sort
228  
228               
229            CALL SLASRT( 'I' , N , D , INFO )
230  
231        ELSE
232  
233  *         Use Selection Sort to minimize swaps of eigenvectors
234  
234               
235            DO 180 II = 2 , N
235                   
236                I = II - 1
237                K = I
238                P = D( I )
239                DO 170 J = II , N
239                       
240                    IF( D( J ).LT.P ) THEN
240                           
241                        K = J
242                        P = D( J )
243                    END IF
244    170         CONTINUE
244                   
245                IF( K.NE.I ) THEN
245                       
246                    D( K ) = D( I )
247                    D( I ) = P
248                    CALL SSWAP( NR , Z( 1 , I ) , 1 , Z( 1 , K ) , 1 )
249                END IF
250    180     CONTINUE
251        END IF
252  
253    190 CONTINUE
254        RETURN
255  
256  *     End of SSTEQR2
257  
258        END7753

Variables in Routine SSTEQR2()

Summary Report
Data Type	Quantity	Size(byte)
CHARACTER	1	1
INTEGER	15	60
REAL	16	64
TOTAL	32	125

List of Variables

CHARACTER

COMPZ

INTEGER

I II INFO J JTOT
K L LDZ LENDP1 LM1
M MAXIT MM N NR

REAL

B C D( * ) E( * ) F
G ONE P R S
THREE TST TWO WORK( * ) Z( LDZ, * )
ZERO

Variables Dependence Graph
Put the mouse over a right hand side variable to display the corresponding line of the dependence
-
- -
B <--- CB = C*E( I ), EB = C*E( I ), IB = C*E( I )

C <--- ONEC = ONE

D <--- LD( L ) = D( L ) - P, PD( I ) = G + P{2D( L ) = D( L ) - P, 3D( L ) = P, 4D( I ) = P}, DD( L ) = D( L ) - P{2D( K ) = D( I )}, GD( I ) = G + P, ID( K ) = D( I )

E <--- ZEROE( L-1 ) = ZERO, GE( LM1 ) = G

F <--- SF = S*E( I ), EF = S*E( I ), IF = S*E( I )

G <--- BG = C*R - B, LG = ( D( L-1 )-P ) / ( TWO*E( L-1 ) ){2}, M, CG = C*R - B, PG = ( D( L-1 )-P ) / ( TWO*E( L-1 ) ){2, 3G = D( I ) - P}, R{2G = C*R - B}, TWOG = ( D( L-1 )-P ) / ( TWO*E( L-1 ) ), DG = ( D( L-1 )-P ) / ( TWO*E( L-1 ) ){2, 3G = D( I ) - P}, EG = ( D( L-1 )-P ) / ( TWO*E( L-1 ) ){2}, G, IG = D( I ) - P

I <--- LM1DO 120 I = M, LM1, MDO 120 I = M, LM1, NDO 150 I = 1, N - 1, III = II - 1

II <--- NDO 180 II = 2, N

J <--- NDO 170 J = II, N, IIDO 170 J = II, N

JTOT <--- JTOTJTOT = JTOT + 1

K <--- JK = J, IK = I

L <--- LL = L - 2{2L = L - 1}

LM1 <--- LLM1 = L - 1

M <--- LDO 100 M = L, LENDP1, -1, LENDP1DO 100 M = L, LENDP1, -1

MM <--- LMM = L - M + 1, MMM = L - M + 1

P <--- JP = D( J ), LP = D( L ), RP = S*R, SP = S*R, ZEROP = ZERO, DP = D( L ){2P = D( I ), 3P = D( J )}, IP = D( I )

R <--- BR = ( D( I+1 )-G )*S + TWO*C*B, CR = ( D( I+1 )-G )*S + TWO*C*B, ONER = SLAPY2( G, ONE ), SR = ( D( I+1 )-G )*S + TWO*C*B, TWOR = ( D( I+1 )-G )*S + TWO*C*B, DR = ( D( I+1 )-G )*S + TWO*C*B, GR = SLAPY2( G, ONE ){2R = ( D( I+1 )-G )*S + TWO*C*B}, IR = ( D( I+1 )-G )*S + TWO*C*B

S <--- ONES = ONE

TST <--- MTST = ABS( E( M-1 ) )**2, ETST = ABS( E( M-1 ) )**2

WORK <--- CWORK( M ) = C{2WORK( I ) = C}, SWORK( N-1+M ) = S{2WORK( N-1+I ) = S}

Analysis elements of the routine SSTEQR2()
Put the mouse over each element to display detailed matching information

Assigned variables
		B , C , COMPZ , F , G , I , II , INFO , J , JTOT , K , L , LDZ , LENDP1 , LM1 , M , MAXIT , MM , N , NR , ONE , P , R , S , THREE , TST , TWO , ZERO

Active variables
		B , C , COMPZ , D , E , F , G , I , II , INFO , J , JTOT , K , L , LDZ , LENDP1 , LM1 , M , MAXIT , MM , N , NR , ONE , P , R , s , THREE , TST , TWO , WORK , Z , zero

Accessed arrays [ array name : associated index ]
	D	: * , I , I , I , I , I , I+1 , J , J , K , L , L , L , L , L , L , L-1 , L-1 , L-1 , L-1 , LSV , LSV , M , M , M-1
	E	: * , I , I , I , I-1 , L-1 , L-1 , L-1 , L-1 , L-1 , LM1 , LSV , LSV , M-1 , M-1
	WORK	: * , I , M , M , M , N-1+I , N-1+M , N-1+M , N-1+M
	Z	: 1, I , 1, K , 1, L-1 , 1, M , LDZ, *

Conditional statements [ statement : associated predicate ]
	do	: `(` 100 M = L , LENDP1 , - 1 `)` , `(` 120 I = M , LM1 `)` , `(` 150 I = 1 , N - 1 `)` , `(` 180 II = 2 , N `)` , `(` 170 J = II , N `)`
	for	: `(` small superdiagonal element. `)` , `(` no convergence to an eigenvalue after a total `)`
	if	: `(` INFO = 0 , the eigenvalues in ascending order. `)` , `(` COMPZ = 'V' , `)` , `(` INFO = 0 , `)` , `(` COMPZ = 'V' , Z contains the `)` , `(` COMPZ = 'I' , Z contains the orthonormal eigenvectors `)` , `(` COMPZ = 'N' , `)` , `(` COMPZ = 'N' , `)` , `(` COMPZ = 'N' , `)` , `(` INFO = - i , the i - th argument had an illegal value `)` , `(` INFO = i , `)` , `(` L.NE.LEND `)` , `(` (TST.LE.( EPS2ABS( D( M ) ) )ABS( D( M - 1 ) ) + `)` , `(` M.GT.LEND `)` , `(` M.EQ.L `)` , `(` remaining matrix is 2 - by - 2 , use SLAE2 or SLAEV2 `)` , `(` M.EQ.L - 1 `)` , `(` ICOMPZ.GT.0 `)` , `(` L.GE.LEND `)` , `(` JTOT.EQ.NMAXIT `)` , `(` I.NE.M `)` , `(` eigenvectors are desired , `)` , `(` ICOMPZ.GT.0 `)` , `(` eigenvectors are desired , `)` , `(` ICOMPZ.GT.0 `)` , `(` L.GE.LEND `)` , `(` necessary `)` , `(` ISCALE.EQ.1 `)` , `(` ISCALE.EQ.2 `)` , `(` JTOT.LT.NMAXIT `)` , `(` (E( I ).NE.ZERO ) `)` , `(` ICOMPZ.EQ.0 `)` , `(` (D( J ).LT.P ) `)` , `(` K.NE.I `)`