1 SUBROUTINE ZLATSP( UPLO, N, X, ISEED )
2 *
3 * -- LAPACK auxiliary test routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8 CHARACTER UPLO
9 INTEGER N
10 * ..
11 * .. Array Arguments ..
12 INTEGER ISEED( * )
13 COMPLEX*16 X( * )
14 * ..
15 *
16 * Purpose
17 * =======
18 *
19 * ZLATSP generates a special test matrix for the complex symmetric
20 * (indefinite) factorization for packed matrices. The pivot blocks of
21 * the generated matrix will be in the following order:
22 * 2x2 pivot block, non diagonalizable
23 * 1x1 pivot block
24 * 2x2 pivot block, diagonalizable
25 * (cycle repeats)
26 * A row interchange is required for each non-diagonalizable 2x2 block.
27 *
28 * Arguments
29 * =========
30 *
31 * UPLO (input) CHARACTER
32 * Specifies whether the generated matrix is to be upper or
33 * lower triangular.
34 * = 'U': Upper triangular
35 * = 'L': Lower triangular
36 *
37 * N (input) INTEGER
38 * The dimension of the matrix to be generated.
39 *
40 * X (output) COMPLEX*16 array, dimension (N*(N+1)/2)
41 * The generated matrix in packed storage format. The matrix
42 * consists of 3x3 and 2x2 diagonal blocks which result in the
43 * pivot sequence given above. The matrix outside these
44 * diagonal blocks is zero.
45 *
46 * ISEED (input/output) INTEGER array, dimension (4)
47 * On entry, the seed for the random number generator. The last
48 * of the four integers must be odd. (modified on exit)
49 *
50 * =====================================================================
51 *
52 * .. Parameters ..
53 COMPLEX*16 EYE
54 PARAMETER ( EYE = ( 0.0D0, 1.0D0 ) )
55 * ..
56 * .. Local Scalars ..
57 INTEGER J, JJ, N5
58 DOUBLE PRECISION ALPHA, ALPHA3, BETA
59 COMPLEX*16 A, B, C, R
60 * ..
61 * .. External Functions ..
62 COMPLEX*16 ZLARND
63 EXTERNAL ZLARND
64 * ..
65 * .. Intrinsic Functions ..
66 INTRINSIC ABS, SQRT
67 * ..
68 * .. Executable Statements ..
69 *
70 * Initialize constants
71 *
72 ALPHA = ( 1.D0+SQRT( 17.D0 ) ) / 8.D0
73 BETA = ALPHA - 1.D0 / 1000.D0
74 ALPHA3 = ALPHA*ALPHA*ALPHA
75 *
76 * Fill the matrix with zeros.
77 *
78 DO 10 J = 1, N*( N+1 ) / 2
79 X( J ) = 0.0D0
80 10 CONTINUE
81 *
82 * UPLO = 'U': Upper triangular storage
83 *
84 IF( UPLO.EQ.'U' ) THEN
85 N5 = N / 5
86 N5 = N - 5*N5 + 1
87 *
88 JJ = N*( N+1 ) / 2
89 DO 20 J = N, N5, -5
90 A = ALPHA3*ZLARND( 5, ISEED )
91 B = ZLARND( 5, ISEED ) / ALPHA
92 C = A - 2.D0*B*EYE
93 R = C / BETA
94 X( JJ ) = A
95 X( JJ-2 ) = B
96 JJ = JJ - J
97 X( JJ ) = ZLARND( 2, ISEED )
98 X( JJ-1 ) = R
99 JJ = JJ - ( J-1 )
100 X( JJ ) = C
101 JJ = JJ - ( J-2 )
102 X( JJ ) = ZLARND( 2, ISEED )
103 JJ = JJ - ( J-3 )
104 X( JJ ) = ZLARND( 2, ISEED )
105 IF( ABS( X( JJ+( J-3 ) ) ).GT.ABS( X( JJ ) ) ) THEN
106 X( JJ+( J-4 ) ) = 2.0D0*X( JJ+( J-3 ) )
107 ELSE
108 X( JJ+( J-4 ) ) = 2.0D0*X( JJ )
109 END IF
110 JJ = JJ - ( J-4 )
111 20 CONTINUE
112 *
113 * Clean-up for N not a multiple of 5.
114 *
115 J = N5 - 1
116 IF( J.GT.2 ) THEN
117 A = ALPHA3*ZLARND( 5, ISEED )
118 B = ZLARND( 5, ISEED ) / ALPHA
119 C = A - 2.D0*B*EYE
120 R = C / BETA
121 X( JJ ) = A
122 X( JJ-2 ) = B
123 JJ = JJ - J
124 X( JJ ) = ZLARND( 2, ISEED )
125 X( JJ-1 ) = R
126 JJ = JJ - ( J-1 )
127 X( JJ ) = C
128 JJ = JJ - ( J-2 )
129 J = J - 3
130 END IF
131 IF( J.GT.1 ) THEN
132 X( JJ ) = ZLARND( 2, ISEED )
133 X( JJ-J ) = ZLARND( 2, ISEED )
134 IF( ABS( X( JJ ) ).GT.ABS( X( JJ-J ) ) ) THEN
135 X( JJ-1 ) = 2.0D0*X( JJ )
136 ELSE
137 X( JJ-1 ) = 2.0D0*X( JJ-J )
138 END IF
139 JJ = JJ - J - ( J-1 )
140 J = J - 2
141 ELSE IF( J.EQ.1 ) THEN
142 X( JJ ) = ZLARND( 2, ISEED )
143 J = J - 1
144 END IF
145 *
146 * UPLO = 'L': Lower triangular storage
147 *
148 ELSE
149 N5 = N / 5
150 N5 = N5*5
151 *
152 JJ = 1
153 DO 30 J = 1, N5, 5
154 A = ALPHA3*ZLARND( 5, ISEED )
155 B = ZLARND( 5, ISEED ) / ALPHA
156 C = A - 2.D0*B*EYE
157 R = C / BETA
158 X( JJ ) = A
159 X( JJ+2 ) = B
160 JJ = JJ + ( N-J+1 )
161 X( JJ ) = ZLARND( 2, ISEED )
162 X( JJ+1 ) = R
163 JJ = JJ + ( N-J )
164 X( JJ ) = C
165 JJ = JJ + ( N-J-1 )
166 X( JJ ) = ZLARND( 2, ISEED )
167 JJ = JJ + ( N-J-2 )
168 X( JJ ) = ZLARND( 2, ISEED )
169 IF( ABS( X( JJ-( N-J-2 ) ) ).GT.ABS( X( JJ ) ) ) THEN
170 X( JJ-( N-J-2 )+1 ) = 2.0D0*X( JJ-( N-J-2 ) )
171 ELSE
172 X( JJ-( N-J-2 )+1 ) = 2.0D0*X( JJ )
173 END IF
174 JJ = JJ + ( N-J-3 )
175 30 CONTINUE
176 *
177 * Clean-up for N not a multiple of 5.
178 *
179 J = N5 + 1
180 IF( J.LT.N-1 ) THEN
181 A = ALPHA3*ZLARND( 5, ISEED )
182 B = ZLARND( 5, ISEED ) / ALPHA
183 C = A - 2.D0*B*EYE
184 R = C / BETA
185 X( JJ ) = A
186 X( JJ+2 ) = B
187 JJ = JJ + ( N-J+1 )
188 X( JJ ) = ZLARND( 2, ISEED )
189 X( JJ+1 ) = R
190 JJ = JJ + ( N-J )
191 X( JJ ) = C
192 JJ = JJ + ( N-J-1 )
193 J = J + 3
194 END IF
195 IF( J.LT.N ) THEN
196 X( JJ ) = ZLARND( 2, ISEED )
197 X( JJ+( N-J+1 ) ) = ZLARND( 2, ISEED )
198 IF( ABS( X( JJ ) ).GT.ABS( X( JJ+( N-J+1 ) ) ) ) THEN
199 X( JJ+1 ) = 2.0D0*X( JJ )
200 ELSE
201 X( JJ+1 ) = 2.0D0*X( JJ+( N-J+1 ) )
202 END IF
203 JJ = JJ + ( N-J+1 ) + ( N-J )
204 J = J + 2
205 ELSE IF( J.EQ.N ) THEN
206 X( JJ ) = ZLARND( 2, ISEED )
207 JJ = JJ + ( N-J+1 )
208 J = J + 1
209 END IF
210 END IF
211 *
212 RETURN
213 *
214 * End of ZLATSP
215 *
216 END
2 *
3 * -- LAPACK auxiliary test routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8 CHARACTER UPLO
9 INTEGER N
10 * ..
11 * .. Array Arguments ..
12 INTEGER ISEED( * )
13 COMPLEX*16 X( * )
14 * ..
15 *
16 * Purpose
17 * =======
18 *
19 * ZLATSP generates a special test matrix for the complex symmetric
20 * (indefinite) factorization for packed matrices. The pivot blocks of
21 * the generated matrix will be in the following order:
22 * 2x2 pivot block, non diagonalizable
23 * 1x1 pivot block
24 * 2x2 pivot block, diagonalizable
25 * (cycle repeats)
26 * A row interchange is required for each non-diagonalizable 2x2 block.
27 *
28 * Arguments
29 * =========
30 *
31 * UPLO (input) CHARACTER
32 * Specifies whether the generated matrix is to be upper or
33 * lower triangular.
34 * = 'U': Upper triangular
35 * = 'L': Lower triangular
36 *
37 * N (input) INTEGER
38 * The dimension of the matrix to be generated.
39 *
40 * X (output) COMPLEX*16 array, dimension (N*(N+1)/2)
41 * The generated matrix in packed storage format. The matrix
42 * consists of 3x3 and 2x2 diagonal blocks which result in the
43 * pivot sequence given above. The matrix outside these
44 * diagonal blocks is zero.
45 *
46 * ISEED (input/output) INTEGER array, dimension (4)
47 * On entry, the seed for the random number generator. The last
48 * of the four integers must be odd. (modified on exit)
49 *
50 * =====================================================================
51 *
52 * .. Parameters ..
53 COMPLEX*16 EYE
54 PARAMETER ( EYE = ( 0.0D0, 1.0D0 ) )
55 * ..
56 * .. Local Scalars ..
57 INTEGER J, JJ, N5
58 DOUBLE PRECISION ALPHA, ALPHA3, BETA
59 COMPLEX*16 A, B, C, R
60 * ..
61 * .. External Functions ..
62 COMPLEX*16 ZLARND
63 EXTERNAL ZLARND
64 * ..
65 * .. Intrinsic Functions ..
66 INTRINSIC ABS, SQRT
67 * ..
68 * .. Executable Statements ..
69 *
70 * Initialize constants
71 *
72 ALPHA = ( 1.D0+SQRT( 17.D0 ) ) / 8.D0
73 BETA = ALPHA - 1.D0 / 1000.D0
74 ALPHA3 = ALPHA*ALPHA*ALPHA
75 *
76 * Fill the matrix with zeros.
77 *
78 DO 10 J = 1, N*( N+1 ) / 2
79 X( J ) = 0.0D0
80 10 CONTINUE
81 *
82 * UPLO = 'U': Upper triangular storage
83 *
84 IF( UPLO.EQ.'U' ) THEN
85 N5 = N / 5
86 N5 = N - 5*N5 + 1
87 *
88 JJ = N*( N+1 ) / 2
89 DO 20 J = N, N5, -5
90 A = ALPHA3*ZLARND( 5, ISEED )
91 B = ZLARND( 5, ISEED ) / ALPHA
92 C = A - 2.D0*B*EYE
93 R = C / BETA
94 X( JJ ) = A
95 X( JJ-2 ) = B
96 JJ = JJ - J
97 X( JJ ) = ZLARND( 2, ISEED )
98 X( JJ-1 ) = R
99 JJ = JJ - ( J-1 )
100 X( JJ ) = C
101 JJ = JJ - ( J-2 )
102 X( JJ ) = ZLARND( 2, ISEED )
103 JJ = JJ - ( J-3 )
104 X( JJ ) = ZLARND( 2, ISEED )
105 IF( ABS( X( JJ+( J-3 ) ) ).GT.ABS( X( JJ ) ) ) THEN
106 X( JJ+( J-4 ) ) = 2.0D0*X( JJ+( J-3 ) )
107 ELSE
108 X( JJ+( J-4 ) ) = 2.0D0*X( JJ )
109 END IF
110 JJ = JJ - ( J-4 )
111 20 CONTINUE
112 *
113 * Clean-up for N not a multiple of 5.
114 *
115 J = N5 - 1
116 IF( J.GT.2 ) THEN
117 A = ALPHA3*ZLARND( 5, ISEED )
118 B = ZLARND( 5, ISEED ) / ALPHA
119 C = A - 2.D0*B*EYE
120 R = C / BETA
121 X( JJ ) = A
122 X( JJ-2 ) = B
123 JJ = JJ - J
124 X( JJ ) = ZLARND( 2, ISEED )
125 X( JJ-1 ) = R
126 JJ = JJ - ( J-1 )
127 X( JJ ) = C
128 JJ = JJ - ( J-2 )
129 J = J - 3
130 END IF
131 IF( J.GT.1 ) THEN
132 X( JJ ) = ZLARND( 2, ISEED )
133 X( JJ-J ) = ZLARND( 2, ISEED )
134 IF( ABS( X( JJ ) ).GT.ABS( X( JJ-J ) ) ) THEN
135 X( JJ-1 ) = 2.0D0*X( JJ )
136 ELSE
137 X( JJ-1 ) = 2.0D0*X( JJ-J )
138 END IF
139 JJ = JJ - J - ( J-1 )
140 J = J - 2
141 ELSE IF( J.EQ.1 ) THEN
142 X( JJ ) = ZLARND( 2, ISEED )
143 J = J - 1
144 END IF
145 *
146 * UPLO = 'L': Lower triangular storage
147 *
148 ELSE
149 N5 = N / 5
150 N5 = N5*5
151 *
152 JJ = 1
153 DO 30 J = 1, N5, 5
154 A = ALPHA3*ZLARND( 5, ISEED )
155 B = ZLARND( 5, ISEED ) / ALPHA
156 C = A - 2.D0*B*EYE
157 R = C / BETA
158 X( JJ ) = A
159 X( JJ+2 ) = B
160 JJ = JJ + ( N-J+1 )
161 X( JJ ) = ZLARND( 2, ISEED )
162 X( JJ+1 ) = R
163 JJ = JJ + ( N-J )
164 X( JJ ) = C
165 JJ = JJ + ( N-J-1 )
166 X( JJ ) = ZLARND( 2, ISEED )
167 JJ = JJ + ( N-J-2 )
168 X( JJ ) = ZLARND( 2, ISEED )
169 IF( ABS( X( JJ-( N-J-2 ) ) ).GT.ABS( X( JJ ) ) ) THEN
170 X( JJ-( N-J-2 )+1 ) = 2.0D0*X( JJ-( N-J-2 ) )
171 ELSE
172 X( JJ-( N-J-2 )+1 ) = 2.0D0*X( JJ )
173 END IF
174 JJ = JJ + ( N-J-3 )
175 30 CONTINUE
176 *
177 * Clean-up for N not a multiple of 5.
178 *
179 J = N5 + 1
180 IF( J.LT.N-1 ) THEN
181 A = ALPHA3*ZLARND( 5, ISEED )
182 B = ZLARND( 5, ISEED ) / ALPHA
183 C = A - 2.D0*B*EYE
184 R = C / BETA
185 X( JJ ) = A
186 X( JJ+2 ) = B
187 JJ = JJ + ( N-J+1 )
188 X( JJ ) = ZLARND( 2, ISEED )
189 X( JJ+1 ) = R
190 JJ = JJ + ( N-J )
191 X( JJ ) = C
192 JJ = JJ + ( N-J-1 )
193 J = J + 3
194 END IF
195 IF( J.LT.N ) THEN
196 X( JJ ) = ZLARND( 2, ISEED )
197 X( JJ+( N-J+1 ) ) = ZLARND( 2, ISEED )
198 IF( ABS( X( JJ ) ).GT.ABS( X( JJ+( N-J+1 ) ) ) ) THEN
199 X( JJ+1 ) = 2.0D0*X( JJ )
200 ELSE
201 X( JJ+1 ) = 2.0D0*X( JJ+( N-J+1 ) )
202 END IF
203 JJ = JJ + ( N-J+1 ) + ( N-J )
204 J = J + 2
205 ELSE IF( J.EQ.N ) THEN
206 X( JJ ) = ZLARND( 2, ISEED )
207 JJ = JJ + ( N-J+1 )
208 J = J + 1
209 END IF
210 END IF
211 *
212 RETURN
213 *
214 * End of ZLATSP
215 *
216 END