LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine dlaed0 ( integer  ICOMPQ,
integer  QSIZ,
integer  N,
double precision, dimension( * )  D,
double precision, dimension( * )  E,
double precision, dimension( ldq, * )  Q,
integer  LDQ,
double precision, dimension( ldqs, * )  QSTORE,
integer  LDQS,
double precision, dimension( * )  WORK,
integer, dimension( * )  IWORK,
integer  INFO 
)

DLAED0 used by sstedc. Computes all eigenvalues and corresponding eigenvectors of an unreduced symmetric tridiagonal matrix using the divide and conquer method.

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

Purpose: 
 DLAED0 computes all eigenvalues and corresponding eigenvectors of a
 symmetric tridiagonal matrix using the divide and conquer method.
Parameters
[in]ICOMPQ
          ICOMPQ is INTEGER
          = 0:  Compute eigenvalues only.
          = 1:  Compute eigenvectors of original dense symmetric matrix
                also.  On entry, Q contains the orthogonal matrix used
                to reduce the original matrix to tridiagonal form.
          = 2:  Compute eigenvalues and eigenvectors of tridiagonal
                matrix.
[in]QSIZ
          QSIZ is INTEGER
         The dimension of the orthogonal matrix used to reduce
         the full matrix to tridiagonal form.  QSIZ >= N if ICOMPQ = 1.
[in]N
          N is INTEGER
         The dimension of the symmetric tridiagonal matrix.  N >= 0.
[in,out]D
          D is DOUBLE PRECISION array, dimension (N)
         On entry, the main diagonal of the tridiagonal matrix.
         On exit, its eigenvalues.
[in]E
          E is DOUBLE PRECISION array, dimension (N-1)
         The off-diagonal elements of the tridiagonal matrix.
         On exit, E has been destroyed.
[in,out]Q
          Q is DOUBLE PRECISION array, dimension (LDQ, N)
         On entry, Q must contain an N-by-N orthogonal matrix.
         If ICOMPQ = 0    Q is not referenced.
         If ICOMPQ = 1    On entry, Q is a subset of the columns of the
                          orthogonal matrix used to reduce the full
                          matrix to tridiagonal form corresponding to
                          the subset of the full matrix which is being
                          decomposed at this time.
         If ICOMPQ = 2    On entry, Q will be the identity matrix.
                          On exit, Q contains the eigenvectors of the
                          tridiagonal matrix.
[in]LDQ
          LDQ is INTEGER
         The leading dimension of the array Q.  If eigenvectors are
         desired, then  LDQ >= max(1,N).  In any case,  LDQ >= 1.
[out]QSTORE
          QSTORE is DOUBLE PRECISION array, dimension (LDQS, N)
         Referenced only when ICOMPQ = 1.  Used to store parts of
         the eigenvector matrix when the updating matrix multiplies
         take place.
[in]LDQS
          LDQS is INTEGER
         The leading dimension of the array QSTORE.  If ICOMPQ = 1,
         then  LDQS >= max(1,N).  In any case,  LDQS >= 1.
[out]WORK
          WORK is DOUBLE PRECISION array,
         If ICOMPQ = 0 or 1, the dimension of WORK must be at least
                     1 + 3*N + 2*N*lg N + 3*N**2
                     ( lg( N ) = smallest integer k
                                 such that 2^k >= N )
         If ICOMPQ = 2, the dimension of WORK must be at least
                     4*N + N**2.
[out]IWORK
          IWORK is INTEGER array,
         If ICOMPQ = 0 or 1, the dimension of IWORK must be at least
                        6 + 6*N + 5*N*lg N.
                        ( lg( N ) = smallest integer k
                                    such that 2^k >= N )
         If ICOMPQ = 2, the dimension of IWORK must be at least
                        3 + 5*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 to compute an eigenvalue while
                working on the submatrix lying in rows and columns
                INFO/(N+1) through mod(INFO,N+1).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012
Contributors:
Jeff Rutter, Computer Science Division, University of California at Berkeley, USA

Definition at line 174 of file dlaed0.f.

174 *
175 * -- LAPACK computational routine (version 3.4.2) --
176 * -- LAPACK is a software package provided by Univ. of Tennessee, --
177 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
178 * September 2012
179 *
180 * .. Scalar Arguments ..
181  INTEGER icompq, info, ldq, ldqs, n, qsiz
182 * ..
183 * .. Array Arguments ..
184  INTEGER iwork( * )
185  DOUBLE PRECISION d( * ), e( * ), q( ldq, * ), qstore( ldqs, * ),
186  $ work( * )
187 * ..
188 *
189 * =====================================================================
190 *
191 * .. Parameters ..
192  DOUBLE PRECISION zero, one, two
193  parameter ( zero = 0.d0, one = 1.d0, two = 2.d0 )
194 * ..
195 * .. Local Scalars ..
196  INTEGER curlvl, curprb, curr, i, igivcl, igivnm,
197  $ igivpt, indxq, iperm, iprmpt, iq, iqptr, iwrem,
198  $ j, k, lgn, matsiz, msd2, smlsiz, smm1, spm1,
199  $ spm2, submat, subpbs, tlvls
200  DOUBLE PRECISION temp
201 * ..
202 * .. External Subroutines ..
203  EXTERNAL dcopy, dgemm, dlacpy, dlaed1, dlaed7, dsteqr,
204  $ xerbla
205 * ..
206 * .. External Functions ..
207  INTEGER ilaenv
208  EXTERNAL ilaenv
209 * ..
210 * .. Intrinsic Functions ..
211  INTRINSIC abs, dble, int, log, max
212 * ..
213 * .. Executable Statements ..
214 *
215 * Test the input parameters.
216 *
217  info = 0
218 *
219  IF( icompq.LT.0 .OR. icompq.GT.2 ) THEN
220  info = -1
221  ELSE IF( ( icompq.EQ.1 ) .AND. ( qsiz.LT.max( 0, n ) ) ) THEN
222  info = -2
223  ELSE IF( n.LT.0 ) THEN
224  info = -3
225  ELSE IF( ldq.LT.max( 1, n ) ) THEN
226  info = -7
227  ELSE IF( ldqs.LT.max( 1, n ) ) THEN
228  info = -9
229  END IF
230  IF( info.NE.0 ) THEN
231  CALL xerbla( 'DLAED0', -info )
232  RETURN
233  END IF
234 *
235 * Quick return if possible
236 *
237  IF( n.EQ.0 )
238  $ RETURN
239 *
240  smlsiz = ilaenv( 9, 'DLAED0', ' ', 0, 0, 0, 0 )
241 *
242 * Determine the size and placement of the submatrices, and save in
243 * the leading elements of IWORK.
244 *
245  iwork( 1 ) = n
246  subpbs = 1
247  tlvls = 0
248  10 CONTINUE
249  IF( iwork( subpbs ).GT.smlsiz ) THEN
250  DO 20 j = subpbs, 1, -1
251  iwork( 2*j ) = ( iwork( j )+1 ) / 2
252  iwork( 2*j-1 ) = iwork( j ) / 2
253  20 CONTINUE
254  tlvls = tlvls + 1
255  subpbs = 2*subpbs
256  GO TO 10
257  END IF
258  DO 30 j = 2, subpbs
259  iwork( j ) = iwork( j ) + iwork( j-1 )
260  30 CONTINUE
261 *
262 * Divide the matrix into SUBPBS submatrices of size at most SMLSIZ+1
263 * using rank-1 modifications (cuts).
264 *
265  spm1 = subpbs - 1
266  DO 40 i = 1, spm1
267  submat = iwork( i ) + 1
268  smm1 = submat - 1
269  d( smm1 ) = d( smm1 ) - abs( e( smm1 ) )
270  d( submat ) = d( submat ) - abs( e( smm1 ) )
271  40 CONTINUE
272 *
273  indxq = 4*n + 3
274  IF( icompq.NE.2 ) THEN
275 *
276 * Set up workspaces for eigenvalues only/accumulate new vectors
277 * routine
278 *
279  temp = log( dble( n ) ) / log( two )
280  lgn = int( temp )
281  IF( 2**lgn.LT.n )
282  $ lgn = lgn + 1
283  IF( 2**lgn.LT.n )
284  $ lgn = lgn + 1
285  iprmpt = indxq + n + 1
286  iperm = iprmpt + n*lgn
287  iqptr = iperm + n*lgn
288  igivpt = iqptr + n + 2
289  igivcl = igivpt + n*lgn
290 *
291  igivnm = 1
292  iq = igivnm + 2*n*lgn
293  iwrem = iq + n**2 + 1
294 *
295 * Initialize pointers
296 *
297  DO 50 i = 0, subpbs
298  iwork( iprmpt+i ) = 1
299  iwork( igivpt+i ) = 1
300  50 CONTINUE
301  iwork( iqptr ) = 1
302  END IF
303 *
304 * Solve each submatrix eigenproblem at the bottom of the divide and
305 * conquer tree.
306 *
307  curr = 0
308  DO 70 i = 0, spm1
309  IF( i.EQ.0 ) THEN
310  submat = 1
311  matsiz = iwork( 1 )
312  ELSE
313  submat = iwork( i ) + 1
314  matsiz = iwork( i+1 ) - iwork( i )
315  END IF
316  IF( icompq.EQ.2 ) THEN
317  CALL dsteqr( 'I', matsiz, d( submat ), e( submat ),
318  $ q( submat, submat ), ldq, work, info )
319  IF( info.NE.0 )
320  $ GO TO 130
321  ELSE
322  CALL dsteqr( 'I', matsiz, d( submat ), e( submat ),
323  $ work( iq-1+iwork( iqptr+curr ) ), matsiz, work,
324  $ info )
325  IF( info.NE.0 )
326  $ GO TO 130
327  IF( icompq.EQ.1 ) THEN
328  CALL dgemm( 'N', 'N', qsiz, matsiz, matsiz, one,
329  $ q( 1, submat ), ldq, work( iq-1+iwork( iqptr+
330  $ curr ) ), matsiz, zero, qstore( 1, submat ),
331  $ ldqs )
332  END IF
333  iwork( iqptr+curr+1 ) = iwork( iqptr+curr ) + matsiz**2
334  curr = curr + 1
335  END IF
336  k = 1
337  DO 60 j = submat, iwork( i+1 )
338  iwork( indxq+j ) = k
339  k = k + 1
340  60 CONTINUE
341  70 CONTINUE
342 *
343 * Successively merge eigensystems of adjacent submatrices
344 * into eigensystem for the corresponding larger matrix.
345 *
346 * while ( SUBPBS > 1 )
347 *
348  curlvl = 1
349  80 CONTINUE
350  IF( subpbs.GT.1 ) THEN
351  spm2 = subpbs - 2
352  DO 90 i = 0, spm2, 2
353  IF( i.EQ.0 ) THEN
354  submat = 1
355  matsiz = iwork( 2 )
356  msd2 = iwork( 1 )
357  curprb = 0
358  ELSE
359  submat = iwork( i ) + 1
360  matsiz = iwork( i+2 ) - iwork( i )
361  msd2 = matsiz / 2
362  curprb = curprb + 1
363  END IF
364 *
365 * Merge lower order eigensystems (of size MSD2 and MATSIZ - MSD2)
366 * into an eigensystem of size MATSIZ.
367 * DLAED1 is used only for the full eigensystem of a tridiagonal
368 * matrix.
369 * DLAED7 handles the cases in which eigenvalues only or eigenvalues
370 * and eigenvectors of a full symmetric matrix (which was reduced to
371 * tridiagonal form) are desired.
372 *
373  IF( icompq.EQ.2 ) THEN
374  CALL dlaed1( matsiz, d( submat ), q( submat, submat ),
375  $ ldq, iwork( indxq+submat ),
376  $ e( submat+msd2-1 ), msd2, work,
377  $ iwork( subpbs+1 ), info )
378  ELSE
379  CALL dlaed7( icompq, matsiz, qsiz, tlvls, curlvl, curprb,
380  $ d( submat ), qstore( 1, submat ), ldqs,
381  $ iwork( indxq+submat ), e( submat+msd2-1 ),
382  $ msd2, work( iq ), iwork( iqptr ),
383  $ iwork( iprmpt ), iwork( iperm ),
384  $ iwork( igivpt ), iwork( igivcl ),
385  $ work( igivnm ), work( iwrem ),
386  $ iwork( subpbs+1 ), info )
387  END IF
388  IF( info.NE.0 )
389  $ GO TO 130
390  iwork( i / 2+1 ) = iwork( i+2 )
391  90 CONTINUE
392  subpbs = subpbs / 2
393  curlvl = curlvl + 1
394  GO TO 80
395  END IF
396 *
397 * end while
398 *
399 * Re-merge the eigenvalues/vectors which were deflated at the final
400 * merge step.
401 *
402  IF( icompq.EQ.1 ) THEN
403  DO 100 i = 1, n
404  j = iwork( indxq+i )
405  work( i ) = d( j )
406  CALL dcopy( qsiz, qstore( 1, j ), 1, q( 1, i ), 1 )
407  100 CONTINUE
408  CALL dcopy( n, work, 1, d, 1 )
409  ELSE IF( icompq.EQ.2 ) THEN
410  DO 110 i = 1, n
411  j = iwork( indxq+i )
412  work( i ) = d( j )
413  CALL dcopy( n, q( 1, j ), 1, work( n*i+1 ), 1 )
414  110 CONTINUE
415  CALL dcopy( n, work, 1, d, 1 )
416  CALL dlacpy( 'A', n, n, work( n+1 ), n, q, ldq )
417  ELSE
418  DO 120 i = 1, n
419  j = iwork( indxq+i )
420  work( i ) = d( j )
421  120 CONTINUE
422  CALL dcopy( n, work, 1, d, 1 )
423  END IF
424  GO TO 140
425 *
426  130 CONTINUE
427  info = submat*( n+1 ) + submat + matsiz - 1
428 *
429  140 CONTINUE
430  RETURN
431 *
432 * End of DLAED0
433 *
dcopy
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:53
dlaed1
subroutine dlaed1(N, D, Q, LDQ, INDXQ, RHO, CUTPNT, WORK, IWORK, INFO)
DLAED1 used by sstedc. Computes the updated eigensystem of a diagonal matrix after modification by a ...
Definition: dlaed1.f:165
dlacpy
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:105
dgemm
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
Definition: dgemm.f:189
dsteqr
subroutine dsteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
DSTEQR
Definition: dsteqr.f:133
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
ilaenv
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
dlaed7
subroutine dlaed7(ICOMPQ, N, QSIZ, TLVLS, CURLVL, CURPBM, D, Q, LDQ, INDXQ, RHO, CUTPNT, QSTORE, QPTR, PRMPTR, PERM, GIVPTR, GIVCOL, GIVNUM, WORK, IWORK, INFO)
DLAED7 used by sstedc. Computes the updated eigensystem of a diagonal matrix after modification by a ...
Definition: dlaed7.f:262

Here is the call graph for this function:

Here is the caller graph for this function: