 LAPACK  3.10.0 LAPACK: Linear Algebra PACKage

## ◆ zlagsy()

 subroutine zlagsy ( integer N, integer K, double precision, dimension( * ) D, complex*16, dimension( lda, * ) A, integer LDA, integer, dimension( 4 ) ISEED, complex*16, dimension( * ) WORK, integer INFO )

ZLAGSY

Purpose:
``` ZLAGSY generates a complex symmetric matrix A, by pre- and post-
multiplying a real diagonal matrix D with a random unitary matrix:
A = U*D*U**T. The semi-bandwidth may then be reduced to k by
Parameters
 [in] N ``` N is INTEGER The order of the matrix A. N >= 0.``` [in] K ``` K is INTEGER The number of nonzero subdiagonals within the band of A. 0 <= K <= N-1.``` [in] D ``` D is DOUBLE PRECISION array, dimension (N) The diagonal elements of the diagonal matrix D.``` [out] A ``` A is COMPLEX*16 array, dimension (LDA,N) The generated n by n symmetric matrix A (the full matrix is stored).``` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. LDA >= N.``` [in,out] ISEED ``` ISEED is INTEGER array, dimension (4) On entry, the seed of the random number generator; the array elements must be between 0 and 4095, and ISEED(4) must be odd. On exit, the seed is updated.``` [out] WORK ` WORK is COMPLEX*16 array, dimension (2*N)` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value```

Definition at line 101 of file zlagsy.f.

102 *
103 * -- LAPACK auxiliary routine --
104 * -- LAPACK is a software package provided by Univ. of Tennessee, --
105 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
106 *
107 * .. Scalar Arguments ..
108  INTEGER INFO, K, LDA, N
109 * ..
110 * .. Array Arguments ..
111  INTEGER ISEED( 4 )
112  DOUBLE PRECISION D( * )
113  COMPLEX*16 A( LDA, * ), WORK( * )
114 * ..
115 *
116 * =====================================================================
117 *
118 * .. Parameters ..
119  COMPLEX*16 ZERO, ONE, HALF
120  parameter( zero = ( 0.0d+0, 0.0d+0 ),
121  \$ one = ( 1.0d+0, 0.0d+0 ),
122  \$ half = ( 0.5d+0, 0.0d+0 ) )
123 * ..
124 * .. Local Scalars ..
125  INTEGER I, II, J, JJ
126  DOUBLE PRECISION WN
127  COMPLEX*16 ALPHA, TAU, WA, WB
128 * ..
129 * .. External Subroutines ..
130  EXTERNAL xerbla, zaxpy, zgemv, zgerc, zlacgv, zlarnv,
131  \$ zscal, zsymv
132 * ..
133 * .. External Functions ..
134  DOUBLE PRECISION DZNRM2
135  COMPLEX*16 ZDOTC
136  EXTERNAL dznrm2, zdotc
137 * ..
138 * .. Intrinsic Functions ..
139  INTRINSIC abs, dble, max
140 * ..
141 * .. Executable Statements ..
142 *
143 * Test the input arguments
144 *
145  info = 0
146  IF( n.LT.0 ) THEN
147  info = -1
148  ELSE IF( k.LT.0 .OR. k.GT.n-1 ) THEN
149  info = -2
150  ELSE IF( lda.LT.max( 1, n ) ) THEN
151  info = -5
152  END IF
153  IF( info.LT.0 ) THEN
154  CALL xerbla( 'ZLAGSY', -info )
155  RETURN
156  END IF
157 *
158 * initialize lower triangle of A to diagonal matrix
159 *
160  DO 20 j = 1, n
161  DO 10 i = j + 1, n
162  a( i, j ) = zero
163  10 CONTINUE
164  20 CONTINUE
165  DO 30 i = 1, n
166  a( i, i ) = d( i )
167  30 CONTINUE
168 *
169 * Generate lower triangle of symmetric matrix
170 *
171  DO 60 i = n - 1, 1, -1
172 *
173 * generate random reflection
174 *
175  CALL zlarnv( 3, iseed, n-i+1, work )
176  wn = dznrm2( n-i+1, work, 1 )
177  wa = ( wn / abs( work( 1 ) ) )*work( 1 )
178  IF( wn.EQ.zero ) THEN
179  tau = zero
180  ELSE
181  wb = work( 1 ) + wa
182  CALL zscal( n-i, one / wb, work( 2 ), 1 )
183  work( 1 ) = one
184  tau = dble( wb / wa )
185  END IF
186 *
187 * apply random reflection to A(i:n,i:n) from the left
188 * and the right
189 *
190 * compute y := tau * A * conjg(u)
191 *
192  CALL zlacgv( n-i+1, work, 1 )
193  CALL zsymv( 'Lower', n-i+1, tau, a( i, i ), lda, work, 1, zero,
194  \$ work( n+1 ), 1 )
195  CALL zlacgv( n-i+1, work, 1 )
196 *
197 * compute v := y - 1/2 * tau * ( u, y ) * u
198 *
199  alpha = -half*tau*zdotc( n-i+1, work, 1, work( n+1 ), 1 )
200  CALL zaxpy( n-i+1, alpha, work, 1, work( n+1 ), 1 )
201 *
202 * apply the transformation as a rank-2 update to A(i:n,i:n)
203 *
204 * CALL ZSYR2( 'Lower', N-I+1, -ONE, WORK, 1, WORK( N+1 ), 1,
205 * \$ A( I, I ), LDA )
206 *
207  DO 50 jj = i, n
208  DO 40 ii = jj, n
209  a( ii, jj ) = a( ii, jj ) -
210  \$ work( ii-i+1 )*work( n+jj-i+1 ) -
211  \$ work( n+ii-i+1 )*work( jj-i+1 )
212  40 CONTINUE
213  50 CONTINUE
214  60 CONTINUE
215 *
216 * Reduce number of subdiagonals to K
217 *
218  DO 100 i = 1, n - 1 - k
219 *
220 * generate reflection to annihilate A(k+i+1:n,i)
221 *
222  wn = dznrm2( n-k-i+1, a( k+i, i ), 1 )
223  wa = ( wn / abs( a( k+i, i ) ) )*a( k+i, i )
224  IF( wn.EQ.zero ) THEN
225  tau = zero
226  ELSE
227  wb = a( k+i, i ) + wa
228  CALL zscal( n-k-i, one / wb, a( k+i+1, i ), 1 )
229  a( k+i, i ) = one
230  tau = dble( wb / wa )
231  END IF
232 *
233 * apply reflection to A(k+i:n,i+1:k+i-1) from the left
234 *
235  CALL zgemv( 'Conjugate transpose', n-k-i+1, k-1, one,
236  \$ a( k+i, i+1 ), lda, a( k+i, i ), 1, zero, work, 1 )
237  CALL zgerc( n-k-i+1, k-1, -tau, a( k+i, i ), 1, work, 1,
238  \$ a( k+i, i+1 ), lda )
239 *
240 * apply reflection to A(k+i:n,k+i:n) from the left and the right
241 *
242 * compute y := tau * A * conjg(u)
243 *
244  CALL zlacgv( n-k-i+1, a( k+i, i ), 1 )
245  CALL zsymv( 'Lower', n-k-i+1, tau, a( k+i, k+i ), lda,
246  \$ a( k+i, i ), 1, zero, work, 1 )
247  CALL zlacgv( n-k-i+1, a( k+i, i ), 1 )
248 *
249 * compute v := y - 1/2 * tau * ( u, y ) * u
250 *
251  alpha = -half*tau*zdotc( n-k-i+1, a( k+i, i ), 1, work, 1 )
252  CALL zaxpy( n-k-i+1, alpha, a( k+i, i ), 1, work, 1 )
253 *
254 * apply symmetric rank-2 update to A(k+i:n,k+i:n)
255 *
256 * CALL ZSYR2( 'Lower', N-K-I+1, -ONE, A( K+I, I ), 1, WORK, 1,
257 * \$ A( K+I, K+I ), LDA )
258 *
259  DO 80 jj = k + i, n
260  DO 70 ii = jj, n
261  a( ii, jj ) = a( ii, jj ) - a( ii, i )*work( jj-k-i+1 ) -
262  \$ work( ii-k-i+1 )*a( jj, i )
263  70 CONTINUE
264  80 CONTINUE
265 *
266  a( k+i, i ) = -wa
267  DO 90 j = k + i + 1, n
268  a( j, i ) = zero
269  90 CONTINUE
270  100 CONTINUE
271 *
272 * Store full symmetric matrix
273 *
274  DO 120 j = 1, n
275  DO 110 i = j + 1, n
276  a( j, i ) = a( i, j )
277  110 CONTINUE
278  120 CONTINUE
279  RETURN
280 *
281 * End of ZLAGSY
282 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
complex *16 function zdotc(N, ZX, INCX, ZY, INCY)
ZDOTC
Definition: zdotc.f:83
subroutine zaxpy(N, ZA, ZX, INCX, ZY, INCY)
ZAXPY
Definition: zaxpy.f:88
subroutine zscal(N, ZA, ZX, INCX)
ZSCAL
Definition: zscal.f:78
subroutine zgerc(M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
ZGERC
Definition: zgerc.f:130
subroutine zgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
ZGEMV
Definition: zgemv.f:158
subroutine zlarnv(IDIST, ISEED, N, X)
ZLARNV returns a vector of random numbers from a uniform or normal distribution.
Definition: zlarnv.f:99
subroutine zlacgv(N, X, INCX)
ZLACGV conjugates a complex vector.
Definition: zlacgv.f:74
subroutine zsymv(UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
ZSYMV computes a matrix-vector product for a complex symmetric matrix.
Definition: zsymv.f:157
real(wp) function dznrm2(n, x, incx)
DZNRM2
Definition: dznrm2.f90:90
Here is the call graph for this function:
Here is the caller graph for this function: