1 SUBROUTINE SCKCSD( NM, MVAL, PVAL, QVAL, NMATS, ISEED, THRESH,
2 $ MMAX, X, XF, U1, U2, V1T, V2T, THETA, IWORK,
3 $ WORK, RWORK, NIN, NOUT, INFO )
4 IMPLICIT NONE
5 *
6 * Originally SCKGSV
7 * -- LAPACK test routine (version 3.1) --
8 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
9 * November 2006
10 *
11 * Adapted to SCKCSD
12 * July 2010
13 *
14 * .. Scalar Arguments ..
15 INTEGER INFO, NIN, NM, NMATS, MMAX, NOUT
16 REAL THRESH
17 * ..
18 * .. Array Arguments ..
19 INTEGER ISEED( 4 ), IWORK( * ), MVAL( * ), PVAL( * ),
20 $ QVAL( * )
21 REAL RWORK( * ), THETA( * )
22 REAL U1( * ), U2( * ), V1T( * ), V2T( * ),
23 $ WORK( * ), X( * ), XF( * )
24 * ..
25 *
26 * Purpose
27 * =======
28 *
29 * SCKCSD tests SORCSD:
30 * the CSD for an M-by-M orthogonal matrix X partitioned as
31 * [ X11 X12; X21 X22 ]. X11 is P-by-Q.
32 *
33 * Arguments
34 * =========
35 *
36 * NM (input) INTEGER
37 * The number of values of M contained in the vector MVAL.
38 *
39 * MVAL (input) INTEGER array, dimension (NM)
40 * The values of the matrix row dimension M.
41 *
42 * PVAL (input) INTEGER array, dimension (NM)
43 * The values of the matrix row dimension P.
44 *
45 * QVAL (input) INTEGER array, dimension (NM)
46 * The values of the matrix column dimension Q.
47 *
48 * NMATS (input) INTEGER
49 * The number of matrix types to be tested for each combination
50 * of matrix dimensions. If NMATS >= NTYPES (the maximum
51 * number of matrix types), then all the different types are
52 * generated for testing. If NMATS < NTYPES, another input line
53 * is read to get the numbers of the matrix types to be used.
54 *
55 * ISEED (input/output) INTEGER array, dimension (4)
56 * On entry, the seed of the random number generator. The array
57 * elements should be between 0 and 4095, otherwise they will be
58 * reduced mod 4096, and ISEED(4) must be odd.
59 * On exit, the next seed in the random number sequence after
60 * all the test matrices have been generated.
61 *
62 * THRESH (input) REAL
63 * The threshold value for the test ratios. A result is
64 * included in the output file if RESULT >= THRESH. To have
65 * every test ratio printed, use THRESH = 0.
66 *
67 * MMAX (input) INTEGER
68 * The maximum value permitted for M, used in dimensioning the
69 * work arrays.
70 *
71 * X (workspace) REAL array, dimension (MMAX*MMAX)
72 *
73 * XF (workspace) REAL array, dimension (MMAX*MMAX)
74 *
75 * U1 (workspace) REAL array, dimension (MMAX*MMAX)
76 *
77 * U2 (workspace) REAL array, dimension (MMAX*MMAX)
78 *
79 * V1T (workspace) REAL array, dimension (MMAX*MMAX)
80 *
81 * V2T (workspace) REAL array, dimension (MMAX*MMAX)
82 *
83 * THETA (workspace) REAL array, dimension (MMAX)
84 *
85 * IWORK (workspace) INTEGER array, dimension (MMAX)
86 *
87 * WORK (workspace) REAL array
88 *
89 * RWORK (workspace) REAL array
90 *
91 * NIN (input) INTEGER
92 * The unit number for input.
93 *
94 * NOUT (input) INTEGER
95 * The unit number for output.
96 *
97 * INFO (output) INTEGER
98 * = 0 : successful exit
99 * > 0 : If SLAROR returns an error code, the absolute value
100 * of it is returned.
101 *
102 * =====================================================================
103 *
104 * .. Parameters ..
105 INTEGER NTESTS
106 PARAMETER ( NTESTS = 9 )
107 INTEGER NTYPES
108 PARAMETER ( NTYPES = 3 )
109 REAL GAPDIGIT, ORTH, PIOVER2, TEN
110 PARAMETER ( GAPDIGIT = 10.0E0, ORTH = 1.0E-4,
111 $ PIOVER2 = 1.57079632679489662E0,
112 $ TEN = 10.0D0 )
113 REAL ZERO, ONE
114 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0 )
115 * ..
116 * .. Local Scalars ..
117 LOGICAL FIRSTT
118 CHARACTER*3 PATH
119 INTEGER I, IINFO, IM, IMAT, J, LDU1, LDU2, LDV1T,
120 $ LDV2T, LDX, LWORK, M, NFAIL, NRUN, NT, P, Q, R
121 * ..
122 * .. Local Arrays ..
123 LOGICAL DOTYPE( NTYPES )
124 REAL RESULT( NTESTS )
125 * ..
126 * .. External Subroutines ..
127 EXTERNAL ALAHDG, ALAREQ, ALASUM, SCSDTS, SLACSG, SLAROR,
128 $ SLASET
129 * ..
130 * .. Intrinsic Functions ..
131 INTRINSIC ABS, COS, MIN, SIN
132 * ..
133 * .. External Functions ..
134 REAL SLANGE, SLARND
135 EXTERNAL SLANGE, SLARND
136 * ..
137 * .. Executable Statements ..
138 *
139 * Initialize constants and the random number seed.
140 *
141 PATH( 1: 3 ) = 'CSD'
142 INFO = 0
143 NRUN = 0
144 NFAIL = 0
145 FIRSTT = .TRUE.
146 CALL ALAREQ( PATH, NMATS, DOTYPE, NTYPES, NIN, NOUT )
147 LDX = MMAX
148 LDU1 = MMAX
149 LDU2 = MMAX
150 LDV1T = MMAX
151 LDV2T = MMAX
152 LWORK = MMAX*MMAX
153 *
154 * Do for each value of M in MVAL.
155 *
156 DO 30 IM = 1, NM
157 M = MVAL( IM )
158 P = PVAL( IM )
159 Q = QVAL( IM )
160 *
161 DO 20 IMAT = 1, NTYPES
162 *
163 * Do the tests only if DOTYPE( IMAT ) is true.
164 *
165 IF( .NOT.DOTYPE( IMAT ) )
166 $ GO TO 20
167 *
168 * Generate X
169 *
170 IF( IMAT.EQ.1 ) THEN
171 CALL SLAROR( 'L', 'I', M, M, X, LDX, ISEED, WORK, IINFO )
172 IF( M .NE. 0 .AND. IINFO .NE. 0 ) THEN
173 WRITE( NOUT, FMT = 9999 ) M, IINFO
174 INFO = ABS( IINFO )
175 GO TO 20
176 END IF
177 ELSE IF( IMAT.EQ.2 ) THEN
178 R = MIN( P, M-P, Q, M-Q )
179 DO I = 1, R
180 THETA(I) = PIOVER2 * SLARND( 1, ISEED )
181 END DO
182 CALL SLACSG( M, P, Q, THETA, ISEED, X, LDX, WORK )
183 DO I = 1, M
184 DO J = 1, M
185 X(I+(J-1)*LDX) = X(I+(J-1)*LDX) +
186 $ ORTH*SLARND(2,ISEED)
187 END DO
188 END DO
189 ELSE
190 R = MIN( P, M-P, Q, M-Q )
191 DO I = 1, R+1
192 THETA(I) = TEN**(-SLARND(1,ISEED)*GAPDIGIT)
193 END DO
194 DO I = 2, R+1
195 THETA(I) = THETA(I-1) + THETA(I)
196 END DO
197 DO I = 1, R
198 THETA(I) = PIOVER2 * THETA(I) / THETA(R+1)
199 END DO
200 CALL SLACSG( M, P, Q, THETA, ISEED, X, LDX, WORK )
201 END IF
202 *
203 NT = 9
204 *
205 CALL SCSDTS( M, P, Q, X, XF, LDX, U1, LDU1, U2, LDU2, V1T,
206 $ LDV1T, V2T, LDV2T, THETA, IWORK, WORK, LWORK,
207 $ RWORK, RESULT )
208 *
209 * Print information about the tests that did not
210 * pass the threshold.
211 *
212 DO 10 I = 1, NT
213 IF( RESULT( I ).GE.THRESH ) THEN
214 IF( NFAIL.EQ.0 .AND. FIRSTT ) THEN
215 FIRSTT = .FALSE.
216 CALL ALAHDG( NOUT, PATH )
217 END IF
218 WRITE( NOUT, FMT = 9998 )M, P, Q, IMAT, I,
219 $ RESULT( I )
220 NFAIL = NFAIL + 1
221 END IF
222 10 CONTINUE
223 NRUN = NRUN + NT
224 20 CONTINUE
225 30 CONTINUE
226 *
227 * Print a summary of the results.
228 *
229 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, 0 )
230 *
231 9999 FORMAT( ' SLAROR in SCKCSD: M = ', I5, ', INFO = ', I15 )
232 9998 FORMAT( ' M=', I4, ' P=', I4, ', Q=', I4, ', type ', I2,
233 $ ', test ', I2, ', ratio=', G13.6 )
234 RETURN
235 *
236 * End of SCKCSD
237 *
238 END
239 *
240 *
241 *
242 SUBROUTINE SLACSG( M, P, Q, THETA, ISEED, X, LDX, WORK )
243 IMPLICIT NONE
244 *
245 INTEGER LDX, M, P, Q
246 INTEGER ISEED( 4 )
247 REAL THETA( * )
248 REAL WORK( * ), X( LDX, * )
249 *
250 REAL ONE, ZERO
251 PARAMETER ( ONE = 1.0E0, ZERO = 0.0E0 )
252 *
253 INTEGER I, INFO, R
254 *
255 R = MIN( P, M-P, Q, M-Q )
256 *
257 CALL SLASET( 'Full', M, M, ZERO, ZERO, X, LDX )
258 *
259 DO I = 1, MIN(P,Q)-R
260 X(I,I) = ONE
261 END DO
262 DO I = 1, R
263 X(MIN(P,Q)-R+I,MIN(P,Q)-R+I) = COS(THETA(I))
264 END DO
265 DO I = 1, MIN(P,M-Q)-R
266 X(P-I+1,M-I+1) = -ONE
267 END DO
268 DO I = 1, R
269 X(P-(MIN(P,M-Q)-R)+1-I,M-(MIN(P,M-Q)-R)+1-I) =
270 $ -SIN(THETA(R-I+1))
271 END DO
272 DO I = 1, MIN(M-P,Q)-R
273 X(M-I+1,Q-I+1) = ONE
274 END DO
275 DO I = 1, R
276 X(M-(MIN(M-P,Q)-R)+1-I,Q-(MIN(M-P,Q)-R)+1-I) =
277 $ SIN(THETA(R-I+1))
278 END DO
279 DO I = 1, MIN(M-P,M-Q)-R
280 X(P+I,Q+I) = ONE
281 END DO
282 DO I = 1, R
283 X(P+(MIN(M-P,M-Q)-R)+I,Q+(MIN(M-P,M-Q)-R)+I) =
284 $ COS(THETA(I))
285 END DO
286 CALL SLAROR( 'Left', 'No init', P, M, X, LDX, ISEED, WORK, INFO )
287 CALL SLAROR( 'Left', 'No init', M-P, M, X(P+1,1), LDX,
288 $ ISEED, WORK, INFO )
289 CALL SLAROR( 'Right', 'No init', M, Q, X, LDX, ISEED,
290 $ WORK, INFO )
291 CALL SLAROR( 'Right', 'No init', M, M-Q,
292 $ X(1,Q+1), LDX, ISEED, WORK, INFO )
293 *
294 END
295
2 $ MMAX, X, XF, U1, U2, V1T, V2T, THETA, IWORK,
3 $ WORK, RWORK, NIN, NOUT, INFO )
4 IMPLICIT NONE
5 *
6 * Originally SCKGSV
7 * -- LAPACK test routine (version 3.1) --
8 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
9 * November 2006
10 *
11 * Adapted to SCKCSD
12 * July 2010
13 *
14 * .. Scalar Arguments ..
15 INTEGER INFO, NIN, NM, NMATS, MMAX, NOUT
16 REAL THRESH
17 * ..
18 * .. Array Arguments ..
19 INTEGER ISEED( 4 ), IWORK( * ), MVAL( * ), PVAL( * ),
20 $ QVAL( * )
21 REAL RWORK( * ), THETA( * )
22 REAL U1( * ), U2( * ), V1T( * ), V2T( * ),
23 $ WORK( * ), X( * ), XF( * )
24 * ..
25 *
26 * Purpose
27 * =======
28 *
29 * SCKCSD tests SORCSD:
30 * the CSD for an M-by-M orthogonal matrix X partitioned as
31 * [ X11 X12; X21 X22 ]. X11 is P-by-Q.
32 *
33 * Arguments
34 * =========
35 *
36 * NM (input) INTEGER
37 * The number of values of M contained in the vector MVAL.
38 *
39 * MVAL (input) INTEGER array, dimension (NM)
40 * The values of the matrix row dimension M.
41 *
42 * PVAL (input) INTEGER array, dimension (NM)
43 * The values of the matrix row dimension P.
44 *
45 * QVAL (input) INTEGER array, dimension (NM)
46 * The values of the matrix column dimension Q.
47 *
48 * NMATS (input) INTEGER
49 * The number of matrix types to be tested for each combination
50 * of matrix dimensions. If NMATS >= NTYPES (the maximum
51 * number of matrix types), then all the different types are
52 * generated for testing. If NMATS < NTYPES, another input line
53 * is read to get the numbers of the matrix types to be used.
54 *
55 * ISEED (input/output) INTEGER array, dimension (4)
56 * On entry, the seed of the random number generator. The array
57 * elements should be between 0 and 4095, otherwise they will be
58 * reduced mod 4096, and ISEED(4) must be odd.
59 * On exit, the next seed in the random number sequence after
60 * all the test matrices have been generated.
61 *
62 * THRESH (input) REAL
63 * The threshold value for the test ratios. A result is
64 * included in the output file if RESULT >= THRESH. To have
65 * every test ratio printed, use THRESH = 0.
66 *
67 * MMAX (input) INTEGER
68 * The maximum value permitted for M, used in dimensioning the
69 * work arrays.
70 *
71 * X (workspace) REAL array, dimension (MMAX*MMAX)
72 *
73 * XF (workspace) REAL array, dimension (MMAX*MMAX)
74 *
75 * U1 (workspace) REAL array, dimension (MMAX*MMAX)
76 *
77 * U2 (workspace) REAL array, dimension (MMAX*MMAX)
78 *
79 * V1T (workspace) REAL array, dimension (MMAX*MMAX)
80 *
81 * V2T (workspace) REAL array, dimension (MMAX*MMAX)
82 *
83 * THETA (workspace) REAL array, dimension (MMAX)
84 *
85 * IWORK (workspace) INTEGER array, dimension (MMAX)
86 *
87 * WORK (workspace) REAL array
88 *
89 * RWORK (workspace) REAL array
90 *
91 * NIN (input) INTEGER
92 * The unit number for input.
93 *
94 * NOUT (input) INTEGER
95 * The unit number for output.
96 *
97 * INFO (output) INTEGER
98 * = 0 : successful exit
99 * > 0 : If SLAROR returns an error code, the absolute value
100 * of it is returned.
101 *
102 * =====================================================================
103 *
104 * .. Parameters ..
105 INTEGER NTESTS
106 PARAMETER ( NTESTS = 9 )
107 INTEGER NTYPES
108 PARAMETER ( NTYPES = 3 )
109 REAL GAPDIGIT, ORTH, PIOVER2, TEN
110 PARAMETER ( GAPDIGIT = 10.0E0, ORTH = 1.0E-4,
111 $ PIOVER2 = 1.57079632679489662E0,
112 $ TEN = 10.0D0 )
113 REAL ZERO, ONE
114 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0 )
115 * ..
116 * .. Local Scalars ..
117 LOGICAL FIRSTT
118 CHARACTER*3 PATH
119 INTEGER I, IINFO, IM, IMAT, J, LDU1, LDU2, LDV1T,
120 $ LDV2T, LDX, LWORK, M, NFAIL, NRUN, NT, P, Q, R
121 * ..
122 * .. Local Arrays ..
123 LOGICAL DOTYPE( NTYPES )
124 REAL RESULT( NTESTS )
125 * ..
126 * .. External Subroutines ..
127 EXTERNAL ALAHDG, ALAREQ, ALASUM, SCSDTS, SLACSG, SLAROR,
128 $ SLASET
129 * ..
130 * .. Intrinsic Functions ..
131 INTRINSIC ABS, COS, MIN, SIN
132 * ..
133 * .. External Functions ..
134 REAL SLANGE, SLARND
135 EXTERNAL SLANGE, SLARND
136 * ..
137 * .. Executable Statements ..
138 *
139 * Initialize constants and the random number seed.
140 *
141 PATH( 1: 3 ) = 'CSD'
142 INFO = 0
143 NRUN = 0
144 NFAIL = 0
145 FIRSTT = .TRUE.
146 CALL ALAREQ( PATH, NMATS, DOTYPE, NTYPES, NIN, NOUT )
147 LDX = MMAX
148 LDU1 = MMAX
149 LDU2 = MMAX
150 LDV1T = MMAX
151 LDV2T = MMAX
152 LWORK = MMAX*MMAX
153 *
154 * Do for each value of M in MVAL.
155 *
156 DO 30 IM = 1, NM
157 M = MVAL( IM )
158 P = PVAL( IM )
159 Q = QVAL( IM )
160 *
161 DO 20 IMAT = 1, NTYPES
162 *
163 * Do the tests only if DOTYPE( IMAT ) is true.
164 *
165 IF( .NOT.DOTYPE( IMAT ) )
166 $ GO TO 20
167 *
168 * Generate X
169 *
170 IF( IMAT.EQ.1 ) THEN
171 CALL SLAROR( 'L', 'I', M, M, X, LDX, ISEED, WORK, IINFO )
172 IF( M .NE. 0 .AND. IINFO .NE. 0 ) THEN
173 WRITE( NOUT, FMT = 9999 ) M, IINFO
174 INFO = ABS( IINFO )
175 GO TO 20
176 END IF
177 ELSE IF( IMAT.EQ.2 ) THEN
178 R = MIN( P, M-P, Q, M-Q )
179 DO I = 1, R
180 THETA(I) = PIOVER2 * SLARND( 1, ISEED )
181 END DO
182 CALL SLACSG( M, P, Q, THETA, ISEED, X, LDX, WORK )
183 DO I = 1, M
184 DO J = 1, M
185 X(I+(J-1)*LDX) = X(I+(J-1)*LDX) +
186 $ ORTH*SLARND(2,ISEED)
187 END DO
188 END DO
189 ELSE
190 R = MIN( P, M-P, Q, M-Q )
191 DO I = 1, R+1
192 THETA(I) = TEN**(-SLARND(1,ISEED)*GAPDIGIT)
193 END DO
194 DO I = 2, R+1
195 THETA(I) = THETA(I-1) + THETA(I)
196 END DO
197 DO I = 1, R
198 THETA(I) = PIOVER2 * THETA(I) / THETA(R+1)
199 END DO
200 CALL SLACSG( M, P, Q, THETA, ISEED, X, LDX, WORK )
201 END IF
202 *
203 NT = 9
204 *
205 CALL SCSDTS( M, P, Q, X, XF, LDX, U1, LDU1, U2, LDU2, V1T,
206 $ LDV1T, V2T, LDV2T, THETA, IWORK, WORK, LWORK,
207 $ RWORK, RESULT )
208 *
209 * Print information about the tests that did not
210 * pass the threshold.
211 *
212 DO 10 I = 1, NT
213 IF( RESULT( I ).GE.THRESH ) THEN
214 IF( NFAIL.EQ.0 .AND. FIRSTT ) THEN
215 FIRSTT = .FALSE.
216 CALL ALAHDG( NOUT, PATH )
217 END IF
218 WRITE( NOUT, FMT = 9998 )M, P, Q, IMAT, I,
219 $ RESULT( I )
220 NFAIL = NFAIL + 1
221 END IF
222 10 CONTINUE
223 NRUN = NRUN + NT
224 20 CONTINUE
225 30 CONTINUE
226 *
227 * Print a summary of the results.
228 *
229 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, 0 )
230 *
231 9999 FORMAT( ' SLAROR in SCKCSD: M = ', I5, ', INFO = ', I15 )
232 9998 FORMAT( ' M=', I4, ' P=', I4, ', Q=', I4, ', type ', I2,
233 $ ', test ', I2, ', ratio=', G13.6 )
234 RETURN
235 *
236 * End of SCKCSD
237 *
238 END
239 *
240 *
241 *
242 SUBROUTINE SLACSG( M, P, Q, THETA, ISEED, X, LDX, WORK )
243 IMPLICIT NONE
244 *
245 INTEGER LDX, M, P, Q
246 INTEGER ISEED( 4 )
247 REAL THETA( * )
248 REAL WORK( * ), X( LDX, * )
249 *
250 REAL ONE, ZERO
251 PARAMETER ( ONE = 1.0E0, ZERO = 0.0E0 )
252 *
253 INTEGER I, INFO, R
254 *
255 R = MIN( P, M-P, Q, M-Q )
256 *
257 CALL SLASET( 'Full', M, M, ZERO, ZERO, X, LDX )
258 *
259 DO I = 1, MIN(P,Q)-R
260 X(I,I) = ONE
261 END DO
262 DO I = 1, R
263 X(MIN(P,Q)-R+I,MIN(P,Q)-R+I) = COS(THETA(I))
264 END DO
265 DO I = 1, MIN(P,M-Q)-R
266 X(P-I+1,M-I+1) = -ONE
267 END DO
268 DO I = 1, R
269 X(P-(MIN(P,M-Q)-R)+1-I,M-(MIN(P,M-Q)-R)+1-I) =
270 $ -SIN(THETA(R-I+1))
271 END DO
272 DO I = 1, MIN(M-P,Q)-R
273 X(M-I+1,Q-I+1) = ONE
274 END DO
275 DO I = 1, R
276 X(M-(MIN(M-P,Q)-R)+1-I,Q-(MIN(M-P,Q)-R)+1-I) =
277 $ SIN(THETA(R-I+1))
278 END DO
279 DO I = 1, MIN(M-P,M-Q)-R
280 X(P+I,Q+I) = ONE
281 END DO
282 DO I = 1, R
283 X(P+(MIN(M-P,M-Q)-R)+I,Q+(MIN(M-P,M-Q)-R)+I) =
284 $ COS(THETA(I))
285 END DO
286 CALL SLAROR( 'Left', 'No init', P, M, X, LDX, ISEED, WORK, INFO )
287 CALL SLAROR( 'Left', 'No init', M-P, M, X(P+1,1), LDX,
288 $ ISEED, WORK, INFO )
289 CALL SLAROR( 'Right', 'No init', M, Q, X, LDX, ISEED,
290 $ WORK, INFO )
291 CALL SLAROR( 'Right', 'No init', M, M-Q,
292 $ X(1,Q+1), LDX, ISEED, WORK, INFO )
293 *
294 END
295