LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine dlaeda ( integer  N,
integer  TLVLS,
integer  CURLVL,
integer  CURPBM,
integer, dimension( * )  PRMPTR,
integer, dimension( * )  PERM,
integer, dimension( * )  GIVPTR,
integer, dimension( 2, * )  GIVCOL,
double precision, dimension( 2, * )  GIVNUM,
double precision, dimension( * )  Q,
integer, dimension( * )  QPTR,
double precision, dimension( * )  Z,
double precision, dimension( * )  ZTEMP,
integer  INFO 
)

DLAEDA used by sstedc. Computes the Z vector determining the rank-one modification of the diagonal matrix. Used when the original matrix is dense.

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

Purpose: 
 DLAEDA computes the Z vector corresponding to the merge step in the
 CURLVLth step of the merge process with TLVLS steps for the CURPBMth
 problem.
Parameters
[in]N
          N is INTEGER
         The dimension of the symmetric tridiagonal matrix.  N >= 0.
[in]TLVLS
          TLVLS is INTEGER
         The total number of merging levels in the overall divide and
         conquer tree.
[in]CURLVL
          CURLVL is INTEGER
         The current level in the overall merge routine,
         0 <= curlvl <= tlvls.
[in]CURPBM
          CURPBM is INTEGER
         The current problem in the current level in the overall
         merge routine (counting from upper left to lower right).
[in]PRMPTR
          PRMPTR is INTEGER array, dimension (N lg N)
         Contains a list of pointers which indicate where in PERM a
         level's permutation is stored.  PRMPTR(i+1) - PRMPTR(i)
         indicates the size of the permutation and incidentally the
         size of the full, non-deflated problem.
[in]PERM
          PERM is INTEGER array, dimension (N lg N)
         Contains the permutations (from deflation and sorting) to be
         applied to each eigenblock.
[in]GIVPTR
          GIVPTR is INTEGER array, dimension (N lg N)
         Contains a list of pointers which indicate where in GIVCOL a
         level's Givens rotations are stored.  GIVPTR(i+1) - GIVPTR(i)
         indicates the number of Givens rotations.
[in]GIVCOL
          GIVCOL is INTEGER array, dimension (2, N lg N)
         Each pair of numbers indicates a pair of columns to take place
         in a Givens rotation.
[in]GIVNUM
          GIVNUM is DOUBLE PRECISION array, dimension (2, N lg N)
         Each number indicates the S value to be used in the
         corresponding Givens rotation.
[in]Q
          Q is DOUBLE PRECISION array, dimension (N**2)
         Contains the square eigenblocks from previous levels, the
         starting positions for blocks are given by QPTR.
[in]QPTR
          QPTR is INTEGER array, dimension (N+2)
         Contains a list of pointers which indicate where in Q an
         eigenblock is stored.  SQRT( QPTR(i+1) - QPTR(i) ) indicates
         the size of the block.
[out]Z
          Z is DOUBLE PRECISION array, dimension (N)
         On output this vector contains the updating vector (the last
         row of the first sub-eigenvector matrix and the first row of
         the second sub-eigenvector matrix).
[out]ZTEMP
          ZTEMP is DOUBLE PRECISION array, dimension (N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit.
          < 0:  if INFO = -i, the i-th argument had an illegal value.
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 168 of file dlaeda.f.

168 *
169 * -- LAPACK computational routine (version 3.4.2) --
170 * -- LAPACK is a software package provided by Univ. of Tennessee, --
171 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
172 * September 2012
173 *
174 * .. Scalar Arguments ..
175  INTEGER curlvl, curpbm, info, n, tlvls
176 * ..
177 * .. Array Arguments ..
178  INTEGER givcol( 2, * ), givptr( * ), perm( * ),
179  $ prmptr( * ), qptr( * )
180  DOUBLE PRECISION givnum( 2, * ), q( * ), z( * ), ztemp( * )
181 * ..
182 *
183 * =====================================================================
184 *
185 * .. Parameters ..
186  DOUBLE PRECISION zero, half, one
187  parameter ( zero = 0.0d0, half = 0.5d0, one = 1.0d0 )
188 * ..
189 * .. Local Scalars ..
190  INTEGER bsiz1, bsiz2, curr, i, k, mid, psiz1, psiz2,
191  $ ptr, zptr1
192 * ..
193 * .. External Subroutines ..
194  EXTERNAL dcopy, dgemv, drot, xerbla
195 * ..
196 * .. Intrinsic Functions ..
197  INTRINSIC dble, int, sqrt
198 * ..
199 * .. Executable Statements ..
200 *
201 * Test the input parameters.
202 *
203  info = 0
204 *
205  IF( n.LT.0 ) THEN
206  info = -1
207  END IF
208  IF( info.NE.0 ) THEN
209  CALL xerbla( 'DLAEDA', -info )
210  RETURN
211  END IF
212 *
213 * Quick return if possible
214 *
215  IF( n.EQ.0 )
216  $ RETURN
217 *
218 * Determine location of first number in second half.
219 *
220  mid = n / 2 + 1
221 *
222 * Gather last/first rows of appropriate eigenblocks into center of Z
223 *
224  ptr = 1
225 *
226 * Determine location of lowest level subproblem in the full storage
227 * scheme
228 *
229  curr = ptr + curpbm*2**curlvl + 2**( curlvl-1 ) - 1
230 *
231 * Determine size of these matrices. We add HALF to the value of
232 * the SQRT in case the machine underestimates one of these square
233 * roots.
234 *
235  bsiz1 = int( half+sqrt( dble( qptr( curr+1 )-qptr( curr ) ) ) )
236  bsiz2 = int( half+sqrt( dble( qptr( curr+2 )-qptr( curr+1 ) ) ) )
237  DO 10 k = 1, mid - bsiz1 - 1
238  z( k ) = zero
239  10 CONTINUE
240  CALL dcopy( bsiz1, q( qptr( curr )+bsiz1-1 ), bsiz1,
241  $ z( mid-bsiz1 ), 1 )
242  CALL dcopy( bsiz2, q( qptr( curr+1 ) ), bsiz2, z( mid ), 1 )
243  DO 20 k = mid + bsiz2, n
244  z( k ) = zero
245  20 CONTINUE
246 *
247 * Loop through remaining levels 1 -> CURLVL applying the Givens
248 * rotations and permutation and then multiplying the center matrices
249 * against the current Z.
250 *
251  ptr = 2**tlvls + 1
252  DO 70 k = 1, curlvl - 1
253  curr = ptr + curpbm*2**( curlvl-k ) + 2**( curlvl-k-1 ) - 1
254  psiz1 = prmptr( curr+1 ) - prmptr( curr )
255  psiz2 = prmptr( curr+2 ) - prmptr( curr+1 )
256  zptr1 = mid - psiz1
257 *
258 * Apply Givens at CURR and CURR+1
259 *
260  DO 30 i = givptr( curr ), givptr( curr+1 ) - 1
261  CALL drot( 1, z( zptr1+givcol( 1, i )-1 ), 1,
262  $ z( zptr1+givcol( 2, i )-1 ), 1, givnum( 1, i ),
263  $ givnum( 2, i ) )
264  30 CONTINUE
265  DO 40 i = givptr( curr+1 ), givptr( curr+2 ) - 1
266  CALL drot( 1, z( mid-1+givcol( 1, i ) ), 1,
267  $ z( mid-1+givcol( 2, i ) ), 1, givnum( 1, i ),
268  $ givnum( 2, i ) )
269  40 CONTINUE
270  psiz1 = prmptr( curr+1 ) - prmptr( curr )
271  psiz2 = prmptr( curr+2 ) - prmptr( curr+1 )
272  DO 50 i = 0, psiz1 - 1
273  ztemp( i+1 ) = z( zptr1+perm( prmptr( curr )+i )-1 )
274  50 CONTINUE
275  DO 60 i = 0, psiz2 - 1
276  ztemp( psiz1+i+1 ) = z( mid+perm( prmptr( curr+1 )+i )-1 )
277  60 CONTINUE
278 *
279 * Multiply Blocks at CURR and CURR+1
280 *
281 * Determine size of these matrices. We add HALF to the value of
282 * the SQRT in case the machine underestimates one of these
283 * square roots.
284 *
285  bsiz1 = int( half+sqrt( dble( qptr( curr+1 )-qptr( curr ) ) ) )
286  bsiz2 = int( half+sqrt( dble( qptr( curr+2 )-qptr( curr+
287  $ 1 ) ) ) )
288  IF( bsiz1.GT.0 ) THEN
289  CALL dgemv( 'T', bsiz1, bsiz1, one, q( qptr( curr ) ),
290  $ bsiz1, ztemp( 1 ), 1, zero, z( zptr1 ), 1 )
291  END IF
292  CALL dcopy( psiz1-bsiz1, ztemp( bsiz1+1 ), 1, z( zptr1+bsiz1 ),
293  $ 1 )
294  IF( bsiz2.GT.0 ) THEN
295  CALL dgemv( 'T', bsiz2, bsiz2, one, q( qptr( curr+1 ) ),
296  $ bsiz2, ztemp( psiz1+1 ), 1, zero, z( mid ), 1 )
297  END IF
298  CALL dcopy( psiz2-bsiz2, ztemp( psiz1+bsiz2+1 ), 1,
299  $ z( mid+bsiz2 ), 1 )
300 *
301  ptr = ptr + 2**( tlvls-k )
302  70 CONTINUE
303 *
304  RETURN
305 *
306 * End of DLAEDA
307 *
dcopy
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:53
dgemv
subroutine dgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
DGEMV
Definition: dgemv.f:158
drot
subroutine drot(N, DX, INCX, DY, INCY, C, S)
DROT
Definition: drot.f:53
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62

Here is the call graph for this function:

Here is the caller graph for this function: