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

◆ dlasd8()

subroutine dlasd8 ( integer icompq,
integer k,
double precision, dimension( * ) d,
double precision, dimension( * ) z,
double precision, dimension( * ) vf,
double precision, dimension( * ) vl,
double precision, dimension( * ) difl,
double precision, dimension( lddifr, * ) difr,
integer lddifr,
double precision, dimension( * ) dsigma,
double precision, dimension( * ) work,
integer info )

DLASD8 finds the square roots of the roots of the secular equation, and stores, for each element in D, the distance to its two nearest poles. Used by sbdsdc.

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

Purpose:
!>
!> DLASD8 finds the square roots of the roots of the secular equation,
!> as defined by the values in DSIGMA and Z. It makes the appropriate
!> calls to DLASD4, and stores, for each  element in D, the distance
!> to its two nearest poles (elements in DSIGMA). It also updates
!> the arrays VF and VL, the first and last components of all the
!> right singular vectors of the original bidiagonal matrix.
!>
!> DLASD8 is called from DLASD6.
!>
Parameters
[in]ICOMPQ
!>          ICOMPQ is INTEGER
!>          Specifies whether singular vectors are to be computed in
!>          factored form in the calling routine:
!>          = 0: Compute singular values only.
!>          = 1: Compute singular vectors in factored form as well.
!>
[in]K
!>          K is INTEGER
!>          The number of terms in the rational function to be solved
!>          by DLASD4.  K >= 1.
!>
[out]D
!>          D is DOUBLE PRECISION array, dimension ( K )
!>          On output, D contains the updated singular values.
!>
[in,out]Z
!>          Z is DOUBLE PRECISION array, dimension ( K )
!>          On entry, the first K elements of this array contain the
!>          components of the deflation-adjusted updating row vector.
!>          On exit, Z is updated.
!>
[in,out]VF
!>          VF is DOUBLE PRECISION array, dimension ( K )
!>          On entry, VF contains  information passed through DBEDE8.
!>          On exit, VF contains the first K components of the first
!>          components of all right singular vectors of the bidiagonal
!>          matrix.
!>
[in,out]VL
!>          VL is DOUBLE PRECISION array, dimension ( K )
!>          On entry, VL contains  information passed through DBEDE8.
!>          On exit, VL contains the first K components of the last
!>          components of all right singular vectors of the bidiagonal
!>          matrix.
!>
[out]DIFL
!>          DIFL is DOUBLE PRECISION array, dimension ( K )
!>          On exit, DIFL(I) = D(I) - DSIGMA(I).
!>
[out]DIFR
!>          DIFR is DOUBLE PRECISION array,
!>                   dimension ( LDDIFR, 2 ) if ICOMPQ = 1 and
!>                   dimension ( K ) if ICOMPQ = 0.
!>          On exit, DIFR(I,1) = D(I) - DSIGMA(I+1), DIFR(K,1) is not
!>          defined and will not be referenced.
!>
!>          If ICOMPQ = 1, DIFR(1:K,2) is an array containing the
!>          normalizing factors for the right singular vector matrix.
!>
[in]LDDIFR
!>          LDDIFR is INTEGER
!>          The leading dimension of DIFR, must be at least K.
!>
[in]DSIGMA
!>          DSIGMA is DOUBLE PRECISION array, dimension ( K )
!>          On entry, the first K elements of this array contain the old
!>          roots of the deflated updating problem.  These are the poles
!>          of the secular equation.
!>
[out]WORK
!>          WORK is DOUBLE PRECISION array, dimension (3*K)
!>
[out]INFO
!>          INFO is INTEGER
!>          = 0:  successful exit.
!>          < 0:  if INFO = -i, the i-th argument had an illegal value.
!>          > 0:  if INFO = 1, a singular value did not converge
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Contributors:
Ming Gu and Huan Ren, Computer Science Division, University of California at Berkeley, USA

Definition at line 160 of file dlasd8.f.

162*
163* -- LAPACK auxiliary routine --
164* -- LAPACK is a software package provided by Univ. of Tennessee, --
165* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
166*
167* .. Scalar Arguments ..
168 INTEGER ICOMPQ, INFO, K, LDDIFR
169* ..
170* .. Array Arguments ..
171 DOUBLE PRECISION D( * ), DIFL( * ), DIFR( LDDIFR, * ),
172 $ DSIGMA( * ), VF( * ), VL( * ), WORK( * ),
173 $ Z( * )
174* ..
175*
176* =====================================================================
177*
178* .. Parameters ..
179 DOUBLE PRECISION ONE
180 parameter( one = 1.0d+0 )
181* ..
182* .. Local Scalars ..
183 INTEGER I, IWK1, IWK2, IWK2I, IWK3, IWK3I, J
184 DOUBLE PRECISION DIFLJ, DIFRJ, DJ, DSIGJ, DSIGJP, RHO, TEMP
185* ..
186* .. External Subroutines ..
187 EXTERNAL dcopy, dlascl, dlasd4, dlaset,
188 $ xerbla
189* ..
190* .. External Functions ..
191 DOUBLE PRECISION DDOT, DLAMC3, DNRM2
192 EXTERNAL ddot, dlamc3, dnrm2
193* ..
194* .. Intrinsic Functions ..
195 INTRINSIC abs, sign, sqrt
196* ..
197* .. Executable Statements ..
198*
199* Test the input parameters.
200*
201 info = 0
202*
203 IF( ( icompq.LT.0 ) .OR. ( icompq.GT.1 ) ) THEN
204 info = -1
205 ELSE IF( k.LT.1 ) THEN
206 info = -2
207 ELSE IF( lddifr.LT.k ) THEN
208 info = -9
209 END IF
210 IF( info.NE.0 ) THEN
211 CALL xerbla( 'DLASD8', -info )
212 RETURN
213 END IF
214*
215* Quick return if possible
216*
217 IF( k.EQ.1 ) THEN
218 d( 1 ) = abs( z( 1 ) )
219 difl( 1 ) = d( 1 )
220 IF( icompq.EQ.1 ) THEN
221 difl( 2 ) = one
222 difr( 1, 2 ) = one
223 END IF
224 RETURN
225 END IF
226*
227* Book keeping.
228*
229 iwk1 = 1
230 iwk2 = iwk1 + k
231 iwk3 = iwk2 + k
232 iwk2i = iwk2 - 1
233 iwk3i = iwk3 - 1
234*
235* Normalize Z.
236*
237 rho = dnrm2( k, z, 1 )
238 CALL dlascl( 'G', 0, 0, rho, one, k, 1, z, k, info )
239 rho = rho*rho
240*
241* Initialize WORK(IWK3).
242*
243 CALL dlaset( 'A', k, 1, one, one, work( iwk3 ), k )
244*
245* Compute the updated singular values, the arrays DIFL, DIFR,
246* and the updated Z.
247*
248 DO 40 j = 1, k
249 CALL dlasd4( k, j, dsigma, z, work( iwk1 ), rho, d( j ),
250 $ work( iwk2 ), info )
251*
252* If the root finder fails, report the convergence failure.
253*
254 IF( info.NE.0 ) THEN
255 RETURN
256 END IF
257 work( iwk3i+j ) = work( iwk3i+j )*work( j )*work( iwk2i+j )
258 difl( j ) = -work( j )
259 difr( j, 1 ) = -work( j+1 )
260 DO 20 i = 1, j - 1
261 work( iwk3i+i ) = work( iwk3i+i )*work( i )*
262 $ work( iwk2i+i ) / ( dsigma( i )-
263 $ dsigma( j ) ) / ( dsigma( i )+
264 $ dsigma( j ) )
265 20 CONTINUE
266 DO 30 i = j + 1, k
267 work( iwk3i+i ) = work( iwk3i+i )*work( i )*
268 $ work( iwk2i+i ) / ( dsigma( i )-
269 $ dsigma( j ) ) / ( dsigma( i )+
270 $ dsigma( j ) )
271 30 CONTINUE
272 40 CONTINUE
273*
274* Compute updated Z.
275*
276 DO 50 i = 1, k
277 z( i ) = sign( sqrt( abs( work( iwk3i+i ) ) ), z( i ) )
278 50 CONTINUE
279*
280* Update VF and VL.
281*
282 DO 80 j = 1, k
283 diflj = difl( j )
284 dj = d( j )
285 dsigj = -dsigma( j )
286 IF( j.LT.k ) THEN
287 difrj = -difr( j, 1 )
288 dsigjp = -dsigma( j+1 )
289 END IF
290 work( j ) = -z( j ) / diflj / ( dsigma( j )+dj )
291*
292* Use calls to the subroutine DLAMC3 to enforce the parentheses
293* (x+y)+z. The goal is to prevent optimizing compilers
294* from doing x+(y+z).
295*
296 DO 60 i = 1, j - 1
297 work( i ) = z( i ) / ( dlamc3( dsigma( i ),
298 $ dsigj )-diflj )
299 $ / ( dsigma( i )+dj )
300 60 CONTINUE
301 DO 70 i = j + 1, k
302 work( i ) = z( i ) / ( dlamc3( dsigma( i ),
303 $ dsigjp )+difrj )
304 $ / ( dsigma( i )+dj )
305 70 CONTINUE
306 temp = dnrm2( k, work, 1 )
307 work( iwk2i+j ) = ddot( k, work, 1, vf, 1 ) / temp
308 work( iwk3i+j ) = ddot( k, work, 1, vl, 1 ) / temp
309 IF( icompq.EQ.1 ) THEN
310 difr( j, 2 ) = temp
311 END IF
312 80 CONTINUE
313*
314 CALL dcopy( k, work( iwk2 ), 1, vf, 1 )
315 CALL dcopy( k, work( iwk3 ), 1, vl, 1 )
316*
317 RETURN
318*
319* End of DLASD8
320*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
dcopy
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
Definition dcopy.f:82
ddot
double precision function ddot(n, dx, incx, dy, incy)
DDOT
Definition ddot.f:82
dlamc3
double precision function dlamc3(a, b)
DLAMC3
Definition dlamch.f:172
dlascl
subroutine dlascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition dlascl.f:142
dlasd4
subroutine dlasd4(n, i, d, z, delta, rho, sigma, work, info)
DLASD4 computes the square root of the i-th updated eigenvalue of a positive symmetric rank-one modif...
Definition dlasd4.f:151
dlaset
subroutine dlaset(uplo, m, n, alpha, beta, a, lda)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition dlaset.f:108
dnrm2
real(wp) function dnrm2(n, x, incx)
DNRM2
Definition dnrm2.f90:89
Here is the call graph for this function:
Here is the caller graph for this function: