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

◆ slasd3()

subroutine slasd3 ( integer nl,
integer nr,
integer sqre,
integer k,
real, dimension( * ) d,
real, dimension( ldq, * ) q,
integer ldq,
real, dimension( * ) dsigma,
real, dimension( ldu, * ) u,
integer ldu,
real, dimension( ldu2, * ) u2,
integer ldu2,
real, dimension( ldvt, * ) vt,
integer ldvt,
real, dimension( ldvt2, * ) vt2,
integer ldvt2,
integer, dimension( * ) idxc,
integer, dimension( * ) ctot,
real, dimension( * ) z,
integer info )

SLASD3 finds all square roots of the roots of the secular equation, as defined by the values in D and Z, and then updates the singular vectors by matrix multiplication. Used by sbdsdc.

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

Purpose:
!>
!> SLASD3 finds all the square roots of the roots of the secular
!> equation, as defined by the values in D and Z.  It makes the
!> appropriate calls to SLASD4 and then updates the singular
!> vectors by matrix multiplication.
!>
!> SLASD3 is called from SLASD1.
!>
Parameters
[in]NL
!>          NL is INTEGER
!>         The row dimension of the upper block.  NL >= 1.
!>
[in]NR
!>          NR is INTEGER
!>         The row dimension of the lower block.  NR >= 1.
!>
[in]SQRE
!>          SQRE is INTEGER
!>         = 0: the lower block is an NR-by-NR square matrix.
!>         = 1: the lower block is an NR-by-(NR+1) rectangular matrix.
!>
!>         The bidiagonal matrix has N = NL + NR + 1 rows and
!>         M = N + SQRE >= N columns.
!>
[in]K
!>          K is INTEGER
!>         The size of the secular equation, 1 =< K = < N.
!>
[out]D
!>          D is REAL array, dimension(K)
!>         On exit the square roots of the roots of the secular equation,
!>         in ascending order.
!>
[out]Q
!>          Q is REAL array, dimension (LDQ,K)
!>
[in]LDQ
!>          LDQ is INTEGER
!>         The leading dimension of the array Q.  LDQ >= K.
!>
[in]DSIGMA
!>          DSIGMA is REAL array, dimension(K)
!>         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]U
!>          U is REAL array, dimension (LDU, N)
!>         The last N - K columns of this matrix contain the deflated
!>         left singular vectors.
!>
[in]LDU
!>          LDU is INTEGER
!>         The leading dimension of the array U.  LDU >= N.
!>
[in]U2
!>          U2 is REAL array, dimension (LDU2, N)
!>         The first K columns of this matrix contain the non-deflated
!>         left singular vectors for the split problem.
!>
[in]LDU2
!>          LDU2 is INTEGER
!>         The leading dimension of the array U2.  LDU2 >= N.
!>
[out]VT
!>          VT is REAL array, dimension (LDVT, M)
!>         The last M - K columns of VT**T contain the deflated
!>         right singular vectors.
!>
[in]LDVT
!>          LDVT is INTEGER
!>         The leading dimension of the array VT.  LDVT >= N.
!>
[in,out]VT2
!>          VT2 is REAL array, dimension (LDVT2, N)
!>         The first K columns of VT2**T contain the non-deflated
!>         right singular vectors for the split problem.
!>
[in]LDVT2
!>          LDVT2 is INTEGER
!>         The leading dimension of the array VT2.  LDVT2 >= N.
!>
[in]IDXC
!>          IDXC is INTEGER array, dimension (N)
!>         The permutation used to arrange the columns of U (and rows of
!>         VT) into three groups:  the first group contains non-zero
!>         entries only at and above (or before) NL +1; the second
!>         contains non-zero entries only at and below (or after) NL+2;
!>         and the third is dense. The first column of U and the row of
!>         VT are treated separately, however.
!>
!>         The rows of the singular vectors found by SLASD4
!>         must be likewise permuted before the matrix multiplies can
!>         take place.
!>
[in]CTOT
!>          CTOT is INTEGER array, dimension (4)
!>         A count of the total number of the various types of columns
!>         in U (or rows in VT), as described in IDXC. The fourth column
!>         type is any column which has been deflated.
!>
[in,out]Z
!>          Z is REAL array, dimension (K)
!>         The first K elements of this array contain the components
!>         of the deflation-adjusted updating row vector.
!>
[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 212 of file slasd3.f.

216*
217* -- LAPACK auxiliary routine --
218* -- LAPACK is a software package provided by Univ. of Tennessee, --
219* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
220*
221* .. Scalar Arguments ..
222 INTEGER INFO, K, LDQ, LDU, LDU2, LDVT, LDVT2, NL, NR,
223 $ SQRE
224* ..
225* .. Array Arguments ..
226 INTEGER CTOT( * ), IDXC( * )
227 REAL D( * ), DSIGMA( * ), Q( LDQ, * ), U( LDU, * ),
228 $ U2( LDU2, * ), VT( LDVT, * ), VT2( LDVT2, * ),
229 $ Z( * )
230* ..
231*
232* =====================================================================
233*
234* .. Parameters ..
235 REAL ONE, ZERO, NEGONE
236 parameter( one = 1.0e+0, zero = 0.0e+0,
237 $ negone = -1.0e+0 )
238* ..
239* .. Local Scalars ..
240 INTEGER CTEMP, I, J, JC, KTEMP, M, N, NLP1, NLP2, NRP1
241 REAL RHO, TEMP
242* ..
243* .. External Functions ..
244 REAL SNRM2
245 EXTERNAL snrm2
246* ..
247* .. External Subroutines ..
248 EXTERNAL scopy, sgemm, slacpy, slascl, slasd4,
249 $ xerbla
250* ..
251* .. Intrinsic Functions ..
252 INTRINSIC abs, sign, sqrt
253* ..
254* .. Executable Statements ..
255*
256* Test the input parameters.
257*
258 info = 0
259*
260 IF( nl.LT.1 ) THEN
261 info = -1
262 ELSE IF( nr.LT.1 ) THEN
263 info = -2
264 ELSE IF( ( sqre.NE.1 ) .AND. ( sqre.NE.0 ) ) THEN
265 info = -3
266 END IF
267*
268 n = nl + nr + 1
269 m = n + sqre
270 nlp1 = nl + 1
271 nlp2 = nl + 2
272*
273 IF( ( k.LT.1 ) .OR. ( k.GT.n ) ) THEN
274 info = -4
275 ELSE IF( ldq.LT.k ) THEN
276 info = -7
277 ELSE IF( ldu.LT.n ) THEN
278 info = -10
279 ELSE IF( ldu2.LT.n ) THEN
280 info = -12
281 ELSE IF( ldvt.LT.m ) THEN
282 info = -14
283 ELSE IF( ldvt2.LT.m ) THEN
284 info = -16
285 END IF
286 IF( info.NE.0 ) THEN
287 CALL xerbla( 'SLASD3', -info )
288 RETURN
289 END IF
290*
291* Quick return if possible
292*
293 IF( k.EQ.1 ) THEN
294 d( 1 ) = abs( z( 1 ) )
295 CALL scopy( m, vt2( 1, 1 ), ldvt2, vt( 1, 1 ), ldvt )
296 IF( z( 1 ).GT.zero ) THEN
297 CALL scopy( n, u2( 1, 1 ), 1, u( 1, 1 ), 1 )
298 ELSE
299 DO 10 i = 1, n
300 u( i, 1 ) = -u2( i, 1 )
301 10 CONTINUE
302 END IF
303 RETURN
304 END IF
305*
306* Keep a copy of Z.
307*
308 CALL scopy( k, z, 1, q, 1 )
309*
310* Normalize Z.
311*
312 rho = snrm2( k, z, 1 )
313 CALL slascl( 'G', 0, 0, rho, one, k, 1, z, k, info )
314 rho = rho*rho
315*
316* Find the new singular values.
317*
318 DO 30 j = 1, k
319 CALL slasd4( k, j, dsigma, z, u( 1, j ), rho, d( j ),
320 $ vt( 1, j ), info )
321*
322* If the zero finder fails, report the convergence failure.
323*
324 IF( info.NE.0 ) THEN
325 RETURN
326 END IF
327 30 CONTINUE
328*
329* Compute updated Z.
330*
331 DO 60 i = 1, k
332 z( i ) = u( i, k )*vt( i, k )
333 DO 40 j = 1, i - 1
334 z( i ) = z( i )*( u( i, j )*vt( i, j ) /
335 $ ( dsigma( i )-dsigma( j ) ) /
336 $ ( dsigma( i )+dsigma( j ) ) )
337 40 CONTINUE
338 DO 50 j = i, k - 1
339 z( i ) = z( i )*( u( i, j )*vt( i, j ) /
340 $ ( dsigma( i )-dsigma( j+1 ) ) /
341 $ ( dsigma( i )+dsigma( j+1 ) ) )
342 50 CONTINUE
343 z( i ) = sign( sqrt( abs( z( i ) ) ), q( i, 1 ) )
344 60 CONTINUE
345*
346* Compute left singular vectors of the modified diagonal matrix,
347* and store related information for the right singular vectors.
348*
349 DO 90 i = 1, k
350 vt( 1, i ) = z( 1 ) / u( 1, i ) / vt( 1, i )
351 u( 1, i ) = negone
352 DO 70 j = 2, k
353 vt( j, i ) = z( j ) / u( j, i ) / vt( j, i )
354 u( j, i ) = dsigma( j )*vt( j, i )
355 70 CONTINUE
356 temp = snrm2( k, u( 1, i ), 1 )
357 q( 1, i ) = u( 1, i ) / temp
358 DO 80 j = 2, k
359 jc = idxc( j )
360 q( j, i ) = u( jc, i ) / temp
361 80 CONTINUE
362 90 CONTINUE
363*
364* Update the left singular vector matrix.
365*
366 IF( k.EQ.2 ) THEN
367 CALL sgemm( 'N', 'N', n, k, k, one, u2, ldu2, q, ldq, zero,
368 $ u,
369 $ ldu )
370 GO TO 100
371 END IF
372 IF( ctot( 1 ).GT.0 ) THEN
373 CALL sgemm( 'N', 'N', nl, k, ctot( 1 ), one, u2( 1, 2 ),
374 $ ldu2,
375 $ q( 2, 1 ), ldq, zero, u( 1, 1 ), ldu )
376 IF( ctot( 3 ).GT.0 ) THEN
377 ktemp = 2 + ctot( 1 ) + ctot( 2 )
378 CALL sgemm( 'N', 'N', nl, k, ctot( 3 ), one, u2( 1,
379 $ ktemp ),
380 $ ldu2, q( ktemp, 1 ), ldq, one, u( 1, 1 ), ldu )
381 END IF
382 ELSE IF( ctot( 3 ).GT.0 ) THEN
383 ktemp = 2 + ctot( 1 ) + ctot( 2 )
384 CALL sgemm( 'N', 'N', nl, k, ctot( 3 ), one, u2( 1, ktemp ),
385 $ ldu2, q( ktemp, 1 ), ldq, zero, u( 1, 1 ), ldu )
386 ELSE
387 CALL slacpy( 'F', nl, k, u2, ldu2, u, ldu )
388 END IF
389 CALL scopy( k, q( 1, 1 ), ldq, u( nlp1, 1 ), ldu )
390 ktemp = 2 + ctot( 1 )
391 ctemp = ctot( 2 ) + ctot( 3 )
392 CALL sgemm( 'N', 'N', nr, k, ctemp, one, u2( nlp2, ktemp ),
393 $ ldu2,
394 $ q( ktemp, 1 ), ldq, zero, u( nlp2, 1 ), ldu )
395*
396* Generate the right singular vectors.
397*
398 100 CONTINUE
399 DO 120 i = 1, k
400 temp = snrm2( k, vt( 1, i ), 1 )
401 q( i, 1 ) = vt( 1, i ) / temp
402 DO 110 j = 2, k
403 jc = idxc( j )
404 q( i, j ) = vt( jc, i ) / temp
405 110 CONTINUE
406 120 CONTINUE
407*
408* Update the right singular vector matrix.
409*
410 IF( k.EQ.2 ) THEN
411 CALL sgemm( 'N', 'N', k, m, k, one, q, ldq, vt2, ldvt2,
412 $ zero,
413 $ vt, ldvt )
414 RETURN
415 END IF
416 ktemp = 1 + ctot( 1 )
417 CALL sgemm( 'N', 'N', k, nlp1, ktemp, one, q( 1, 1 ), ldq,
418 $ vt2( 1, 1 ), ldvt2, zero, vt( 1, 1 ), ldvt )
419 ktemp = 2 + ctot( 1 ) + ctot( 2 )
420 IF( ktemp.LE.ldvt2 )
421 $ CALL sgemm( 'N', 'N', k, nlp1, ctot( 3 ), one, q( 1,
422 $ ktemp ),
423 $ ldq, vt2( ktemp, 1 ), ldvt2, one, vt( 1, 1 ),
424 $ ldvt )
425*
426 ktemp = ctot( 1 ) + 1
427 nrp1 = nr + sqre
428 IF( ktemp.GT.1 ) THEN
429 DO 130 i = 1, k
430 q( i, ktemp ) = q( i, 1 )
431 130 CONTINUE
432 DO 140 i = nlp2, m
433 vt2( ktemp, i ) = vt2( 1, i )
434 140 CONTINUE
435 END IF
436 ctemp = 1 + ctot( 2 ) + ctot( 3 )
437 CALL sgemm( 'N', 'N', k, nrp1, ctemp, one, q( 1, ktemp ), ldq,
438 $ vt2( ktemp, nlp2 ), ldvt2, zero, vt( 1, nlp2 ), ldvt )
439*
440 RETURN
441*
442* End of SLASD3
443*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
scopy
subroutine scopy(n, sx, incx, sy, incy)
SCOPY
Definition scopy.f:82
sgemm
subroutine sgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
SGEMM
Definition sgemm.f:188
slacpy
subroutine slacpy(uplo, m, n, a, lda, b, ldb)
SLACPY copies all or part of one two-dimensional array to another.
Definition slacpy.f:101
slascl
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:142
slasd4
subroutine slasd4(n, i, d, z, delta, rho, sigma, work, info)
SLASD4 computes the square root of the i-th updated eigenvalue of a positive symmetric rank-one modif...
Definition slasd4.f:151
snrm2
real(wp) function snrm2(n, x, incx)
SNRM2
Definition snrm2.f90:89
Here is the call graph for this function:
Here is the caller graph for this function: