LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
clatsy.f
Go to the documentation of this file.
1 *> \brief \b CLATSY
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * SUBROUTINE CLATSY( UPLO, N, X, LDX, ISEED )
12 *
13 * .. Scalar Arguments ..
14 * CHARACTER UPLO
15 * INTEGER LDX, N
16 * ..
17 * .. Array Arguments ..
18 * INTEGER ISEED( * )
19 * COMPLEX X( LDX, * )
20 * ..
21 *
22 *
23 *> \par Purpose:
24 * =============
25 *>
26 *> \verbatim
27 *>
28 *> CLATSY generates a special test matrix for the complex symmetric
29 *> (indefinite) factorization. The pivot blocks of the generated matrix
30 *> will be in the following order:
31 *> 2x2 pivot block, non diagonalizable
32 *> 1x1 pivot block
33 *> 2x2 pivot block, diagonalizable
34 *> (cycle repeats)
35 *> A row interchange is required for each non-diagonalizable 2x2 block.
36 *> \endverbatim
37 *
38 * Arguments:
39 * ==========
40 *
41 *> \param[in] UPLO
42 *> \verbatim
43 *> UPLO is CHARACTER
44 *> Specifies whether the generated matrix is to be upper or
45 *> lower triangular.
46 *> = 'U': Upper triangular
47 *> = 'L': Lower triangular
48 *> \endverbatim
49 *>
50 *> \param[in] N
51 *> \verbatim
52 *> N is INTEGER
53 *> The dimension of the matrix to be generated.
54 *> \endverbatim
55 *>
56 *> \param[out] X
57 *> \verbatim
58 *> X is COMPLEX array, dimension (LDX,N)
59 *> The generated matrix, consisting of 3x3 and 2x2 diagonal
60 *> blocks which result in the pivot sequence given above.
61 *> The matrix outside of these diagonal blocks is zero.
62 *> \endverbatim
63 *>
64 *> \param[in] LDX
65 *> \verbatim
66 *> LDX is INTEGER
67 *> The leading dimension of the array X.
68 *> \endverbatim
69 *>
70 *> \param[in,out] ISEED
71 *> \verbatim
72 *> ISEED is INTEGER array, dimension (4)
73 *> On entry, the seed for the random number generator. The last
74 *> of the four integers must be odd. (modified on exit)
75 *> \endverbatim
76 *
77 * Authors:
78 * ========
79 *
80 *> \author Univ. of Tennessee
81 *> \author Univ. of California Berkeley
82 *> \author Univ. of Colorado Denver
83 *> \author NAG Ltd.
84 *
85 *> \ingroup complex_lin
86 *
87 * =====================================================================
88  SUBROUTINE clatsy( UPLO, N, X, LDX, ISEED )
89 *
90 * -- LAPACK test routine --
91 * -- LAPACK is a software package provided by Univ. of Tennessee, --
92 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
93 *
94 * .. Scalar Arguments ..
95  CHARACTER UPLO
96  INTEGER LDX, N
97 * ..
98 * .. Array Arguments ..
99  INTEGER ISEED( * )
100  COMPLEX X( LDX, * )
101 * ..
102 *
103 * =====================================================================
104 *
105 * .. Parameters ..
106  COMPLEX EYE
107  parameter( eye = ( 0.0, 1.0 ) )
108 * ..
109 * .. Local Scalars ..
110  INTEGER I, J, N5
111  REAL ALPHA, ALPHA3, BETA
112  COMPLEX A, B, C, R
113 * ..
114 * .. External Functions ..
115  COMPLEX CLARND
116  EXTERNAL clarnd
117 * ..
118 * .. Intrinsic Functions ..
119  INTRINSIC abs, sqrt
120 * ..
121 * .. Executable Statements ..
122 *
123 * Initialize constants
124 *
125  alpha = ( 1.+sqrt( 17. ) ) / 8.
126  beta = alpha - 1. / 1000.
127  alpha3 = alpha*alpha*alpha
128 *
129 * UPLO = 'U': Upper triangular storage
130 *
131  IF( uplo.EQ.'U' ) THEN
132 *
133 * Fill the upper triangle of the matrix with zeros.
134 *
135  DO 20 j = 1, n
136  DO 10 i = 1, j
137  x( i, j ) = 0.0
138  10 CONTINUE
139  20 CONTINUE
140  n5 = n / 5
141  n5 = n - 5*n5 + 1
142 *
143  DO 30 i = n, n5, -5
144  a = alpha3*clarnd( 5, iseed )
145  b = clarnd( 5, iseed ) / alpha
146  c = a - 2.*b*eye
147  r = c / beta
148  x( i, i ) = a
149  x( i-2, i ) = b
150  x( i-2, i-1 ) = r
151  x( i-2, i-2 ) = c
152  x( i-1, i-1 ) = clarnd( 2, iseed )
153  x( i-3, i-3 ) = clarnd( 2, iseed )
154  x( i-4, i-4 ) = clarnd( 2, iseed )
155  IF( abs( x( i-3, i-3 ) ).GT.abs( x( i-4, i-4 ) ) ) THEN
156  x( i-4, i-3 ) = 2.0*x( i-3, i-3 )
157  ELSE
158  x( i-4, i-3 ) = 2.0*x( i-4, i-4 )
159  END IF
160  30 CONTINUE
161 *
162 * Clean-up for N not a multiple of 5.
163 *
164  i = n5 - 1
165  IF( i.GT.2 ) THEN
166  a = alpha3*clarnd( 5, iseed )
167  b = clarnd( 5, iseed ) / alpha
168  c = a - 2.*b*eye
169  r = c / beta
170  x( i, i ) = a
171  x( i-2, i ) = b
172  x( i-2, i-1 ) = r
173  x( i-2, i-2 ) = c
174  x( i-1, i-1 ) = clarnd( 2, iseed )
175  i = i - 3
176  END IF
177  IF( i.GT.1 ) THEN
178  x( i, i ) = clarnd( 2, iseed )
179  x( i-1, i-1 ) = clarnd( 2, iseed )
180  IF( abs( x( i, i ) ).GT.abs( x( i-1, i-1 ) ) ) THEN
181  x( i-1, i ) = 2.0*x( i, i )
182  ELSE
183  x( i-1, i ) = 2.0*x( i-1, i-1 )
184  END IF
185  i = i - 2
186  ELSE IF( i.EQ.1 ) THEN
187  x( i, i ) = clarnd( 2, iseed )
188  i = i - 1
189  END IF
190 *
191 * UPLO = 'L': Lower triangular storage
192 *
193  ELSE
194 *
195 * Fill the lower triangle of the matrix with zeros.
196 *
197  DO 50 j = 1, n
198  DO 40 i = j, n
199  x( i, j ) = 0.0
200  40 CONTINUE
201  50 CONTINUE
202  n5 = n / 5
203  n5 = n5*5
204 *
205  DO 60 i = 1, n5, 5
206  a = alpha3*clarnd( 5, iseed )
207  b = clarnd( 5, iseed ) / alpha
208  c = a - 2.*b*eye
209  r = c / beta
210  x( i, i ) = a
211  x( i+2, i ) = b
212  x( i+2, i+1 ) = r
213  x( i+2, i+2 ) = c
214  x( i+1, i+1 ) = clarnd( 2, iseed )
215  x( i+3, i+3 ) = clarnd( 2, iseed )
216  x( i+4, i+4 ) = clarnd( 2, iseed )
217  IF( abs( x( i+3, i+3 ) ).GT.abs( x( i+4, i+4 ) ) ) THEN
218  x( i+4, i+3 ) = 2.0*x( i+3, i+3 )
219  ELSE
220  x( i+4, i+3 ) = 2.0*x( i+4, i+4 )
221  END IF
222  60 CONTINUE
223 *
224 * Clean-up for N not a multiple of 5.
225 *
226  i = n5 + 1
227  IF( i.LT.n-1 ) THEN
228  a = alpha3*clarnd( 5, iseed )
229  b = clarnd( 5, iseed ) / alpha
230  c = a - 2.*b*eye
231  r = c / beta
232  x( i, i ) = a
233  x( i+2, i ) = b
234  x( i+2, i+1 ) = r
235  x( i+2, i+2 ) = c
236  x( i+1, i+1 ) = clarnd( 2, iseed )
237  i = i + 3
238  END IF
239  IF( i.LT.n ) THEN
240  x( i, i ) = clarnd( 2, iseed )
241  x( i+1, i+1 ) = clarnd( 2, iseed )
242  IF( abs( x( i, i ) ).GT.abs( x( i+1, i+1 ) ) ) THEN
243  x( i+1, i ) = 2.0*x( i, i )
244  ELSE
245  x( i+1, i ) = 2.0*x( i+1, i+1 )
246  END IF
247  i = i + 2
248  ELSE IF( i.EQ.n ) THEN
249  x( i, i ) = clarnd( 2, iseed )
250  i = i + 1
251  END IF
252  END IF
253 *
254  RETURN
255 *
256 * End of CLATSY
257 *
258  END
subroutine clatsy(UPLO, N, X, LDX, ISEED)
CLATSY
Definition: clatsy.f:89