LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ slasq1()

subroutine slasq1 ( integer  N,
real, dimension( * )  D,
real, dimension( * )  E,
real, dimension( * )  WORK,
integer  INFO 
)

SLASQ1 computes the singular values of a real square bidiagonal matrix. Used by sbdsqr.

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

Purpose:
 SLASQ1 computes the singular values of a real N-by-N bidiagonal
 matrix with diagonal D and off-diagonal E. The singular values
 are computed to high relative accuracy, in the absence of
 denormalization, underflow and overflow. The algorithm was first
 presented in

 "Accurate singular values and differential qd algorithms" by K. V.
 Fernando and B. N. Parlett, Numer. Math., Vol-67, No. 2, pp. 191-230,
 1994,

 and the present implementation is described in "An implementation of
 the dqds Algorithm (Positive Case)", LAPACK Working Note.
Parameters
[in]N
          N is INTEGER
        The number of rows and columns in the matrix. N >= 0.
[in,out]D
          D is REAL array, dimension (N)
        On entry, D contains the diagonal elements of the
        bidiagonal matrix whose SVD is desired. On normal exit,
        D contains the singular values in decreasing order.
[in,out]E
          E is REAL array, dimension (N)
        On entry, elements E(1:N-1) contain the off-diagonal elements
        of the bidiagonal matrix whose SVD is desired.
        On exit, E is overwritten.
[out]WORK
          WORK is REAL array, dimension (4*N)
[out]INFO
          INFO is INTEGER
        = 0: successful exit
        < 0: if INFO = -i, the i-th argument had an illegal value
        > 0: the algorithm failed
             = 1, a split was marked by a positive value in E
             = 2, current block of Z not diagonalized after 100*N
                  iterations (in inner while loop)  On exit D and E
                  represent a matrix with the same singular values
                  which the calling subroutine could use to finish the
                  computation, or even feed back into SLASQ1
             = 3, termination criterion of outer while loop not met
                  (program created more than N unreduced blocks)
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 107 of file slasq1.f.

108 *
109 * -- LAPACK computational routine --
110 * -- LAPACK is a software package provided by Univ. of Tennessee, --
111 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
112 *
113 * .. Scalar Arguments ..
114  INTEGER INFO, N
115 * ..
116 * .. Array Arguments ..
117  REAL D( * ), E( * ), WORK( * )
118 * ..
119 *
120 * =====================================================================
121 *
122 * .. Parameters ..
123  REAL ZERO
124  parameter( zero = 0.0e0 )
125 * ..
126 * .. Local Scalars ..
127  INTEGER I, IINFO
128  REAL EPS, SCALE, SAFMIN, SIGMN, SIGMX
129 * ..
130 * .. External Subroutines ..
131  EXTERNAL scopy, slas2, slascl, slasq2, slasrt, xerbla
132 * ..
133 * .. External Functions ..
134  REAL SLAMCH
135  EXTERNAL slamch
136 * ..
137 * .. Intrinsic Functions ..
138  INTRINSIC abs, max, sqrt
139 * ..
140 * .. Executable Statements ..
141 *
142  info = 0
143  IF( n.LT.0 ) THEN
144  info = -1
145  CALL xerbla( 'SLASQ1', -info )
146  RETURN
147  ELSE IF( n.EQ.0 ) THEN
148  RETURN
149  ELSE IF( n.EQ.1 ) THEN
150  d( 1 ) = abs( d( 1 ) )
151  RETURN
152  ELSE IF( n.EQ.2 ) THEN
153  CALL slas2( d( 1 ), e( 1 ), d( 2 ), sigmn, sigmx )
154  d( 1 ) = sigmx
155  d( 2 ) = sigmn
156  RETURN
157  END IF
158 *
159 * Estimate the largest singular value.
160 *
161  sigmx = zero
162  DO 10 i = 1, n - 1
163  d( i ) = abs( d( i ) )
164  sigmx = max( sigmx, abs( e( i ) ) )
165  10 CONTINUE
166  d( n ) = abs( d( n ) )
167 *
168 * Early return if SIGMX is zero (matrix is already diagonal).
169 *
170  IF( sigmx.EQ.zero ) THEN
171  CALL slasrt( 'D', n, d, iinfo )
172  RETURN
173  END IF
174 *
175  DO 20 i = 1, n
176  sigmx = max( sigmx, d( i ) )
177  20 CONTINUE
178 *
179 * Copy D and E into WORK (in the Z format) and scale (squaring the
180 * input data makes scaling by a power of the radix pointless).
181 *
182  eps = slamch( 'Precision' )
183  safmin = slamch( 'Safe minimum' )
184  scale = sqrt( eps / safmin )
185  CALL scopy( n, d, 1, work( 1 ), 2 )
186  CALL scopy( n-1, e, 1, work( 2 ), 2 )
187  CALL slascl( 'G', 0, 0, sigmx, scale, 2*n-1, 1, work, 2*n-1,
188  $ iinfo )
189 *
190 * Compute the q's and e's.
191 *
192  DO 30 i = 1, 2*n - 1
193  work( i ) = work( i )**2
194  30 CONTINUE
195  work( 2*n ) = zero
196 *
197  CALL slasq2( n, work, info )
198 *
199  IF( info.EQ.0 ) THEN
200  DO 40 i = 1, n
201  d( i ) = sqrt( work( i ) )
202  40 CONTINUE
203  CALL slascl( 'G', 0, 0, scale, sigmx, n, 1, d, n, iinfo )
204  ELSE IF( info.EQ.2 ) THEN
205 *
206 * Maximum number of iterations exceeded. Move data from WORK
207 * into D and E so the calling subroutine can try to finish
208 *
209  DO i = 1, n
210  d( i ) = sqrt( work( 2*i-1 ) )
211  e( i ) = sqrt( work( 2*i ) )
212  END DO
213  CALL slascl( 'G', 0, 0, scale, sigmx, n, 1, d, n, iinfo )
214  CALL slascl( 'G', 0, 0, scale, sigmx, n, 1, e, n, iinfo )
215  END IF
216 *
217  RETURN
218 *
219 * End of SLASQ1
220 *
subroutine slascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: slascl.f:143
subroutine slas2(F, G, H, SSMIN, SSMAX)
SLAS2 computes singular values of a 2-by-2 triangular matrix.
Definition: slas2.f:107
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine slasq2(N, Z, INFO)
SLASQ2 computes all the eigenvalues of the symmetric positive definite tridiagonal matrix associated ...
Definition: slasq2.f:112
subroutine slasrt(ID, N, D, INFO)
SLASRT sorts numbers in increasing or decreasing order.
Definition: slasrt.f:88
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:82
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:68
Here is the call graph for this function:
Here is the caller graph for this function: