LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine chptrd ( character  UPLO,
integer  N,
complex, dimension( * )  AP,
real, dimension( * )  D,
real, dimension( * )  E,
complex, dimension( * )  TAU,
integer  INFO 
)

CHPTRD

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

Purpose:
 CHPTRD reduces a complex Hermitian matrix A stored in packed form to
 real symmetric tridiagonal form T by a unitary similarity
 transformation: Q**H * A * Q = T.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]AP
          AP is COMPLEX array, dimension (N*(N+1)/2)
          On entry, the upper or lower triangle of the Hermitian matrix
          A, packed columnwise in a linear array.  The j-th column of A
          is stored in the array AP as follows:
          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
          if UPLO = 'L', AP(i + (j-1)*(2*n-j)/2) = A(i,j) for j<=i<=n.
          On exit, if UPLO = 'U', the diagonal and first superdiagonal
          of A are overwritten by the corresponding elements of the
          tridiagonal matrix T, and the elements above the first
          superdiagonal, with the array TAU, represent the unitary
          matrix Q as a product of elementary reflectors; if UPLO
          = 'L', the diagonal and first subdiagonal of A are over-
          written by the corresponding elements of the tridiagonal
          matrix T, and the elements below the first subdiagonal, with
          the array TAU, represent the unitary matrix Q as a product
          of elementary reflectors. See Further Details.
[out]D
          D is REAL array, dimension (N)
          The diagonal elements of the tridiagonal matrix T:
          D(i) = A(i,i).
[out]E
          E is REAL array, dimension (N-1)
          The off-diagonal elements of the tridiagonal matrix T:
          E(i) = A(i,i+1) if UPLO = 'U', E(i) = A(i+1,i) if UPLO = 'L'.
[out]TAU
          TAU is COMPLEX array, dimension (N-1)
          The scalar factors of the elementary reflectors (see Further
          Details).
[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
November 2011
Further Details:
  If UPLO = 'U', the matrix Q is represented as a product of elementary
  reflectors

     Q = H(n-1) . . . H(2) H(1).

  Each H(i) has the form

     H(i) = I - tau * v * v**H

  where tau is a complex scalar, and v is a complex vector with
  v(i+1:n) = 0 and v(i) = 1; v(1:i-1) is stored on exit in AP,
  overwriting A(1:i-1,i+1), and tau is stored in TAU(i).

  If UPLO = 'L', the matrix Q is represented as a product of elementary
  reflectors

     Q = H(1) H(2) . . . H(n-1).

  Each H(i) has the form

     H(i) = I - tau * v * v**H

  where tau is a complex scalar, and v is a complex vector with
  v(1:i) = 0 and v(i+1) = 1; v(i+2:n) is stored on exit in AP,
  overwriting A(i+2:n,i), and tau is stored in TAU(i).

Definition at line 153 of file chptrd.f.

153 *
154 * -- LAPACK computational routine (version 3.4.0) --
155 * -- LAPACK is a software package provided by Univ. of Tennessee, --
156 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
157 * November 2011
158 *
159 * .. Scalar Arguments ..
160  CHARACTER uplo
161  INTEGER info, n
162 * ..
163 * .. Array Arguments ..
164  REAL d( * ), e( * )
165  COMPLEX ap( * ), tau( * )
166 * ..
167 *
168 * =====================================================================
169 *
170 * .. Parameters ..
171  COMPLEX one, zero, half
172  parameter ( one = ( 1.0e+0, 0.0e+0 ),
173  $ zero = ( 0.0e+0, 0.0e+0 ),
174  $ half = ( 0.5e+0, 0.0e+0 ) )
175 * ..
176 * .. Local Scalars ..
177  LOGICAL upper
178  INTEGER i, i1, i1i1, ii
179  COMPLEX alpha, taui
180 * ..
181 * .. External Subroutines ..
182  EXTERNAL caxpy, chpmv, chpr2, clarfg, xerbla
183 * ..
184 * .. External Functions ..
185  LOGICAL lsame
186  COMPLEX cdotc
187  EXTERNAL lsame, cdotc
188 * ..
189 * .. Intrinsic Functions ..
190  INTRINSIC real
191 * ..
192 * .. Executable Statements ..
193 *
194 * Test the input parameters
195 *
196  info = 0
197  upper = lsame( uplo, 'U' )
198  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
199  info = -1
200  ELSE IF( n.LT.0 ) THEN
201  info = -2
202  END IF
203  IF( info.NE.0 ) THEN
204  CALL xerbla( 'CHPTRD', -info )
205  RETURN
206  END IF
207 *
208 * Quick return if possible
209 *
210  IF( n.LE.0 )
211  $ RETURN
212 *
213  IF( upper ) THEN
214 *
215 * Reduce the upper triangle of A.
216 * I1 is the index in AP of A(1,I+1).
217 *
218  i1 = n*( n-1 ) / 2 + 1
219  ap( i1+n-1 ) = REAL( AP( I1+N-1 ) )
220  DO 10 i = n - 1, 1, -1
221 *
222 * Generate elementary reflector H(i) = I - tau * v * v**H
223 * to annihilate A(1:i-1,i+1)
224 *
225  alpha = ap( i1+i-1 )
226  CALL clarfg( i, alpha, ap( i1 ), 1, taui )
227  e( i ) = alpha
228 *
229  IF( taui.NE.zero ) THEN
230 *
231 * Apply H(i) from both sides to A(1:i,1:i)
232 *
233  ap( i1+i-1 ) = one
234 *
235 * Compute y := tau * A * v storing y in TAU(1:i)
236 *
237  CALL chpmv( uplo, i, taui, ap, ap( i1 ), 1, zero, tau,
238  $ 1 )
239 *
240 * Compute w := y - 1/2 * tau * (y**H *v) * v
241 *
242  alpha = -half*taui*cdotc( i, tau, 1, ap( i1 ), 1 )
243  CALL caxpy( i, alpha, ap( i1 ), 1, tau, 1 )
244 *
245 * Apply the transformation as a rank-2 update:
246 * A := A - v * w**H - w * v**H
247 *
248  CALL chpr2( uplo, i, -one, ap( i1 ), 1, tau, 1, ap )
249 *
250  END IF
251  ap( i1+i-1 ) = e( i )
252  d( i+1 ) = ap( i1+i )
253  tau( i ) = taui
254  i1 = i1 - i
255  10 CONTINUE
256  d( 1 ) = ap( 1 )
257  ELSE
258 *
259 * Reduce the lower triangle of A. II is the index in AP of
260 * A(i,i) and I1I1 is the index of A(i+1,i+1).
261 *
262  ii = 1
263  ap( 1 ) = REAL( AP( 1 ) )
264  DO 20 i = 1, n - 1
265  i1i1 = ii + n - i + 1
266 *
267 * Generate elementary reflector H(i) = I - tau * v * v**H
268 * to annihilate A(i+2:n,i)
269 *
270  alpha = ap( ii+1 )
271  CALL clarfg( n-i, alpha, ap( ii+2 ), 1, taui )
272  e( i ) = alpha
273 *
274  IF( taui.NE.zero ) THEN
275 *
276 * Apply H(i) from both sides to A(i+1:n,i+1:n)
277 *
278  ap( ii+1 ) = one
279 *
280 * Compute y := tau * A * v storing y in TAU(i:n-1)
281 *
282  CALL chpmv( uplo, n-i, taui, ap( i1i1 ), ap( ii+1 ), 1,
283  $ zero, tau( i ), 1 )
284 *
285 * Compute w := y - 1/2 * tau * (y**H *v) * v
286 *
287  alpha = -half*taui*cdotc( n-i, tau( i ), 1, ap( ii+1 ),
288  $ 1 )
289  CALL caxpy( n-i, alpha, ap( ii+1 ), 1, tau( i ), 1 )
290 *
291 * Apply the transformation as a rank-2 update:
292 * A := A - v * w**H - w * v**H
293 *
294  CALL chpr2( uplo, n-i, -one, ap( ii+1 ), 1, tau( i ), 1,
295  $ ap( i1i1 ) )
296 *
297  END IF
298  ap( ii+1 ) = e( i )
299  d( i ) = ap( ii )
300  tau( i ) = taui
301  ii = i1i1
302  20 CONTINUE
303  d( n ) = ap( ii )
304  END IF
305 *
306  RETURN
307 *
308 * End of CHPTRD
309 *
subroutine chpmv(UPLO, N, ALPHA, AP, X, INCX, BETA, Y, INCY)
CHPMV
Definition: chpmv.f:151
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
complex function cdotc(N, CX, INCX, CY, INCY)
CDOTC
Definition: cdotc.f:54
subroutine chpr2(UPLO, N, ALPHA, X, INCX, Y, INCY, AP)
CHPR2
Definition: chpr2.f:147
subroutine caxpy(N, CA, CX, INCX, CY, INCY)
CAXPY
Definition: caxpy.f:53
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine clarfg(N, ALPHA, X, INCX, TAU)
CLARFG generates an elementary reflector (Householder matrix).
Definition: clarfg.f:108

Here is the call graph for this function:

Here is the caller graph for this function: