 LAPACK 3.11.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ dlaed2()

 subroutine dlaed2 ( integer K, integer N, integer N1, double precision, dimension( * ) D, double precision, dimension( ldq, * ) Q, integer LDQ, integer, dimension( * ) INDXQ, double precision RHO, double precision, dimension( * ) Z, double precision, dimension( * ) DLAMDA, double precision, dimension( * ) W, double precision, dimension( * ) Q2, integer, dimension( * ) INDX, integer, dimension( * ) INDXC, integer, dimension( * ) INDXP, integer, dimension( * ) COLTYP, integer INFO )

DLAED2 used by DSTEDC. Merges eigenvalues and deflates secular equation. Used when the original matrix is tridiagonal.

Purpose:
``` DLAED2 merges the two sets of eigenvalues together into a single
sorted set.  Then it tries to deflate the size of the problem.
There are two ways in which deflation can occur:  when two or more
eigenvalues are close together or if there is a tiny entry in the
Z vector.  For each such occurrence the order of the related secular
equation problem is reduced by one.```
Parameters
 [out] K ``` K is INTEGER The number of non-deflated eigenvalues, and the order of the related secular equation. 0 <= K <=N.``` [in] N ``` N is INTEGER The dimension of the symmetric tridiagonal matrix. N >= 0.``` [in] N1 ``` N1 is INTEGER The location of the last eigenvalue in the leading sub-matrix. min(1,N) <= N1 <= N/2.``` [in,out] D ``` D is DOUBLE PRECISION array, dimension (N) On entry, D contains the eigenvalues of the two submatrices to be combined. On exit, D contains the trailing (N-K) updated eigenvalues (those which were deflated) sorted into increasing order.``` [in,out] Q ``` Q is DOUBLE PRECISION array, dimension (LDQ, N) On entry, Q contains the eigenvectors of two submatrices in the two square blocks with corners at (1,1), (N1,N1) and (N1+1, N1+1), (N,N). On exit, Q contains the trailing (N-K) updated eigenvectors (those which were deflated) in its last N-K columns.``` [in] LDQ ``` LDQ is INTEGER The leading dimension of the array Q. LDQ >= max(1,N).``` [in,out] INDXQ ``` INDXQ is INTEGER array, dimension (N) The permutation which separately sorts the two sub-problems in D into ascending order. Note that elements in the second half of this permutation must first have N1 added to their values. Destroyed on exit.``` [in,out] RHO ``` RHO is DOUBLE PRECISION On entry, the off-diagonal element associated with the rank-1 cut which originally split the two submatrices which are now being recombined. On exit, RHO has been modified to the value required by DLAED3.``` [in] Z ``` Z is DOUBLE PRECISION array, dimension (N) On entry, Z contains the updating vector (the last row of the first sub-eigenvector matrix and the first row of the second sub-eigenvector matrix). On exit, the contents of Z have been destroyed by the updating process.``` [out] DLAMDA ``` DLAMDA is DOUBLE PRECISION array, dimension (N) A copy of the first K eigenvalues which will be used by DLAED3 to form the secular equation.``` [out] W ``` W is DOUBLE PRECISION array, dimension (N) The first k values of the final deflation-altered z-vector which will be passed to DLAED3.``` [out] Q2 ``` Q2 is DOUBLE PRECISION array, dimension (N1**2+(N-N1)**2) A copy of the first K eigenvectors which will be used by DLAED3 in a matrix multiply (DGEMM) to solve for the new eigenvectors.``` [out] INDX ``` INDX is INTEGER array, dimension (N) The permutation used to sort the contents of DLAMDA into ascending order.``` [out] INDXC ``` INDXC is INTEGER array, dimension (N) The permutation used to arrange the columns of the deflated Q matrix into three groups: the first group contains non-zero elements only at and above N1, the second contains non-zero elements only below N1, and the third is dense.``` [out] INDXP ``` INDXP is INTEGER array, dimension (N) The permutation used to place deflated values of D at the end of the array. INDXP(1:K) points to the nondeflated D-values and INDXP(K+1:N) points to the deflated eigenvalues.``` [out] COLTYP ``` COLTYP is INTEGER array, dimension (N) During execution, a label which will indicate which of the following types a column in the Q2 matrix is: 1 : non-zero in the upper half only; 2 : dense; 3 : non-zero in the lower half only; 4 : deflated. On exit, COLTYP(i) is the number of columns of type i, for i=1 to 4 only.``` [out] INFO ``` INFO is INTEGER = 0: successful exit. < 0: if INFO = -i, the i-th argument had an illegal value.```
Contributors:
Jeff Rutter, Computer Science Division, University of California at Berkeley, USA
Modified by Francoise Tisseur, University of Tennessee

Definition at line 210 of file dlaed2.f.

212*
213* -- LAPACK computational routine --
214* -- LAPACK is a software package provided by Univ. of Tennessee, --
215* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
216*
217* .. Scalar Arguments ..
218 INTEGER INFO, K, LDQ, N, N1
219 DOUBLE PRECISION RHO
220* ..
221* .. Array Arguments ..
222 INTEGER COLTYP( * ), INDX( * ), INDXC( * ), INDXP( * ),
223 \$ INDXQ( * )
224 DOUBLE PRECISION D( * ), DLAMDA( * ), Q( LDQ, * ), Q2( * ),
225 \$ W( * ), Z( * )
226* ..
227*
228* =====================================================================
229*
230* .. Parameters ..
231 DOUBLE PRECISION MONE, ZERO, ONE, TWO, EIGHT
232 parameter( mone = -1.0d0, zero = 0.0d0, one = 1.0d0,
233 \$ two = 2.0d0, eight = 8.0d0 )
234* ..
235* .. Local Arrays ..
236 INTEGER CTOT( 4 ), PSM( 4 )
237* ..
238* .. Local Scalars ..
239 INTEGER CT, I, IMAX, IQ1, IQ2, J, JMAX, JS, K2, N1P1,
240 \$ N2, NJ, PJ
241 DOUBLE PRECISION C, EPS, S, T, TAU, TOL
242* ..
243* .. External Functions ..
244 INTEGER IDAMAX
245 DOUBLE PRECISION DLAMCH, DLAPY2
246 EXTERNAL idamax, dlamch, dlapy2
247* ..
248* .. External Subroutines ..
249 EXTERNAL dcopy, dlacpy, dlamrg, drot, dscal, xerbla
250* ..
251* .. Intrinsic Functions ..
252 INTRINSIC abs, max, min, sqrt
253* ..
254* .. Executable Statements ..
255*
256* Test the input parameters.
257*
258 info = 0
259*
260 IF( n.LT.0 ) THEN
261 info = -2
262 ELSE IF( ldq.LT.max( 1, n ) ) THEN
263 info = -6
264 ELSE IF( min( 1, ( n / 2 ) ).GT.n1 .OR. ( n / 2 ).LT.n1 ) THEN
265 info = -3
266 END IF
267 IF( info.NE.0 ) THEN
268 CALL xerbla( 'DLAED2', -info )
269 RETURN
270 END IF
271*
272* Quick return if possible
273*
274 IF( n.EQ.0 )
275 \$ RETURN
276*
277 n2 = n - n1
278 n1p1 = n1 + 1
279*
280 IF( rho.LT.zero ) THEN
281 CALL dscal( n2, mone, z( n1p1 ), 1 )
282 END IF
283*
284* Normalize z so that norm(z) = 1. Since z is the concatenation of
285* two normalized vectors, norm2(z) = sqrt(2).
286*
287 t = one / sqrt( two )
288 CALL dscal( n, t, z, 1 )
289*
290* RHO = ABS( norm(z)**2 * RHO )
291*
292 rho = abs( two*rho )
293*
294* Sort the eigenvalues into increasing order
295*
296 DO 10 i = n1p1, n
297 indxq( i ) = indxq( i ) + n1
298 10 CONTINUE
299*
300* re-integrate the deflated parts from the last pass
301*
302 DO 20 i = 1, n
303 dlamda( i ) = d( indxq( i ) )
304 20 CONTINUE
305 CALL dlamrg( n1, n2, dlamda, 1, 1, indxc )
306 DO 30 i = 1, n
307 indx( i ) = indxq( indxc( i ) )
308 30 CONTINUE
309*
310* Calculate the allowable deflation tolerance
311*
312 imax = idamax( n, z, 1 )
313 jmax = idamax( n, d, 1 )
314 eps = dlamch( 'Epsilon' )
315 tol = eight*eps*max( abs( d( jmax ) ), abs( z( imax ) ) )
316*
317* If the rank-1 modifier is small enough, no more needs to be done
318* except to reorganize Q so that its columns correspond with the
319* elements in D.
320*
321 IF( rho*abs( z( imax ) ).LE.tol ) THEN
322 k = 0
323 iq2 = 1
324 DO 40 j = 1, n
325 i = indx( j )
326 CALL dcopy( n, q( 1, i ), 1, q2( iq2 ), 1 )
327 dlamda( j ) = d( i )
328 iq2 = iq2 + n
329 40 CONTINUE
330 CALL dlacpy( 'A', n, n, q2, n, q, ldq )
331 CALL dcopy( n, dlamda, 1, d, 1 )
332 GO TO 190
333 END IF
334*
335* If there are multiple eigenvalues then the problem deflates. Here
336* the number of equal eigenvalues are found. As each equal
337* eigenvalue is found, an elementary reflector is computed to rotate
338* the corresponding eigensubspace so that the corresponding
339* components of Z are zero in this new basis.
340*
341 DO 50 i = 1, n1
342 coltyp( i ) = 1
343 50 CONTINUE
344 DO 60 i = n1p1, n
345 coltyp( i ) = 3
346 60 CONTINUE
347*
348*
349 k = 0
350 k2 = n + 1
351 DO 70 j = 1, n
352 nj = indx( j )
353 IF( rho*abs( z( nj ) ).LE.tol ) THEN
354*
355* Deflate due to small z component.
356*
357 k2 = k2 - 1
358 coltyp( nj ) = 4
359 indxp( k2 ) = nj
360 IF( j.EQ.n )
361 \$ GO TO 100
362 ELSE
363 pj = nj
364 GO TO 80
365 END IF
366 70 CONTINUE
367 80 CONTINUE
368 j = j + 1
369 nj = indx( j )
370 IF( j.GT.n )
371 \$ GO TO 100
372 IF( rho*abs( z( nj ) ).LE.tol ) THEN
373*
374* Deflate due to small z component.
375*
376 k2 = k2 - 1
377 coltyp( nj ) = 4
378 indxp( k2 ) = nj
379 ELSE
380*
381* Check if eigenvalues are close enough to allow deflation.
382*
383 s = z( pj )
384 c = z( nj )
385*
386* Find sqrt(a**2+b**2) without overflow or
387* destructive underflow.
388*
389 tau = dlapy2( c, s )
390 t = d( nj ) - d( pj )
391 c = c / tau
392 s = -s / tau
393 IF( abs( t*c*s ).LE.tol ) THEN
394*
395* Deflation is possible.
396*
397 z( nj ) = tau
398 z( pj ) = zero
399 IF( coltyp( nj ).NE.coltyp( pj ) )
400 \$ coltyp( nj ) = 2
401 coltyp( pj ) = 4
402 CALL drot( n, q( 1, pj ), 1, q( 1, nj ), 1, c, s )
403 t = d( pj )*c**2 + d( nj )*s**2
404 d( nj ) = d( pj )*s**2 + d( nj )*c**2
405 d( pj ) = t
406 k2 = k2 - 1
407 i = 1
408 90 CONTINUE
409 IF( k2+i.LE.n ) THEN
410 IF( d( pj ).LT.d( indxp( k2+i ) ) ) THEN
411 indxp( k2+i-1 ) = indxp( k2+i )
412 indxp( k2+i ) = pj
413 i = i + 1
414 GO TO 90
415 ELSE
416 indxp( k2+i-1 ) = pj
417 END IF
418 ELSE
419 indxp( k2+i-1 ) = pj
420 END IF
421 pj = nj
422 ELSE
423 k = k + 1
424 dlamda( k ) = d( pj )
425 w( k ) = z( pj )
426 indxp( k ) = pj
427 pj = nj
428 END IF
429 END IF
430 GO TO 80
431 100 CONTINUE
432*
433* Record the last eigenvalue.
434*
435 k = k + 1
436 dlamda( k ) = d( pj )
437 w( k ) = z( pj )
438 indxp( k ) = pj
439*
440* Count up the total number of the various types of columns, then
441* form a permutation which positions the four column types into
442* four uniform groups (although one or more of these groups may be
443* empty).
444*
445 DO 110 j = 1, 4
446 ctot( j ) = 0
447 110 CONTINUE
448 DO 120 j = 1, n
449 ct = coltyp( j )
450 ctot( ct ) = ctot( ct ) + 1
451 120 CONTINUE
452*
453* PSM(*) = Position in SubMatrix (of types 1 through 4)
454*
455 psm( 1 ) = 1
456 psm( 2 ) = 1 + ctot( 1 )
457 psm( 3 ) = psm( 2 ) + ctot( 2 )
458 psm( 4 ) = psm( 3 ) + ctot( 3 )
459 k = n - ctot( 4 )
460*
461* Fill out the INDXC array so that the permutation which it induces
462* will place all type-1 columns first, all type-2 columns next,
463* then all type-3's, and finally all type-4's.
464*
465 DO 130 j = 1, n
466 js = indxp( j )
467 ct = coltyp( js )
468 indx( psm( ct ) ) = js
469 indxc( psm( ct ) ) = j
470 psm( ct ) = psm( ct ) + 1
471 130 CONTINUE
472*
473* Sort the eigenvalues and corresponding eigenvectors into DLAMDA
474* and Q2 respectively. The eigenvalues/vectors which were not
475* deflated go into the first K slots of DLAMDA and Q2 respectively,
476* while those which were deflated go into the last N - K slots.
477*
478 i = 1
479 iq1 = 1
480 iq2 = 1 + ( ctot( 1 )+ctot( 2 ) )*n1
481 DO 140 j = 1, ctot( 1 )
482 js = indx( i )
483 CALL dcopy( n1, q( 1, js ), 1, q2( iq1 ), 1 )
484 z( i ) = d( js )
485 i = i + 1
486 iq1 = iq1 + n1
487 140 CONTINUE
488*
489 DO 150 j = 1, ctot( 2 )
490 js = indx( i )
491 CALL dcopy( n1, q( 1, js ), 1, q2( iq1 ), 1 )
492 CALL dcopy( n2, q( n1+1, js ), 1, q2( iq2 ), 1 )
493 z( i ) = d( js )
494 i = i + 1
495 iq1 = iq1 + n1
496 iq2 = iq2 + n2
497 150 CONTINUE
498*
499 DO 160 j = 1, ctot( 3 )
500 js = indx( i )
501 CALL dcopy( n2, q( n1+1, js ), 1, q2( iq2 ), 1 )
502 z( i ) = d( js )
503 i = i + 1
504 iq2 = iq2 + n2
505 160 CONTINUE
506*
507 iq1 = iq2
508 DO 170 j = 1, ctot( 4 )
509 js = indx( i )
510 CALL dcopy( n, q( 1, js ), 1, q2( iq2 ), 1 )
511 iq2 = iq2 + n
512 z( i ) = d( js )
513 i = i + 1
514 170 CONTINUE
515*
516* The deflated eigenvalues and their corresponding vectors go back
517* into the last N - K slots of D and Q respectively.
518*
519 IF( k.LT.n ) THEN
520 CALL dlacpy( 'A', n, ctot( 4 ), q2( iq1 ), n,
521 \$ q( 1, k+1 ), ldq )
522 CALL dcopy( n-k, z( k+1 ), 1, d( k+1 ), 1 )
523 END IF
524*
525* Copy CTOT into COLTYP for referencing in DLAED3.
526*
527 DO 180 j = 1, 4
528 coltyp( j ) = ctot( j )
529 180 CONTINUE
530*
531 190 CONTINUE
532 RETURN
533*
534* End of DLAED2
535*
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:103
double precision function dlapy2(X, Y)
DLAPY2 returns sqrt(x2+y2).
Definition: dlapy2.f:63
integer function idamax(N, DX, INCX)
IDAMAX
Definition: idamax.f:71
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine dlamrg(N1, N2, A, DTRD1, DTRD2, INDEX)
DLAMRG creates a permutation list to merge the entries of two independently sorted sets into a single...
Definition: dlamrg.f:99
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:82
subroutine drot(N, DX, INCX, DY, INCY, C, S)
DROT
Definition: drot.f:92
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:79
Here is the call graph for this function:
Here is the caller graph for this function: