LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
sckcsd.f
Go to the documentation of this file.
1 *> \brief \b SCKCSD
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 SCKCSD( NM, MVAL, PVAL, QVAL, NMATS, ISEED, THRESH,
12 * MMAX, X, XF, U1, U2, V1T, V2T, THETA, IWORK,
13 * WORK, RWORK, NIN, NOUT, INFO )
14 *
15 * .. Scalar Arguments ..
16 * INTEGER INFO, NIN, NM, NMATS, MMAX, NOUT
17 * REAL THRESH
18 * ..
19 * .. Array Arguments ..
20 * INTEGER ISEED( 4 ), IWORK( * ), MVAL( * ), PVAL( * ),
21 * $ QVAL( * )
22 * REAL RWORK( * ), THETA( * )
23 * REAL U1( * ), U2( * ), V1T( * ), V2T( * ),
24 * $ WORK( * ), X( * ), XF( * )
25 * ..
26 *
27 *
28 *> \par Purpose:
29 * =============
30 *>
31 *> \verbatim
32 *>
33 *> SCKCSD tests SORCSD:
34 *> the CSD for an M-by-M orthogonal matrix X partitioned as
35 *> [ X11 X12; X21 X22 ]. X11 is P-by-Q.
36 *> \endverbatim
37 *
38 * Arguments:
39 * ==========
40 *
41 *> \param[in] NM
42 *> \verbatim
43 *> NM is INTEGER
44 *> The number of values of M contained in the vector MVAL.
45 *> \endverbatim
46 *>
47 *> \param[in] MVAL
48 *> \verbatim
49 *> MVAL is INTEGER array, dimension (NM)
50 *> The values of the matrix row dimension M.
51 *> \endverbatim
52 *>
53 *> \param[in] PVAL
54 *> \verbatim
55 *> PVAL is INTEGER array, dimension (NM)
56 *> The values of the matrix row dimension P.
57 *> \endverbatim
58 *>
59 *> \param[in] QVAL
60 *> \verbatim
61 *> QVAL is INTEGER array, dimension (NM)
62 *> The values of the matrix column dimension Q.
63 *> \endverbatim
64 *>
65 *> \param[in] NMATS
66 *> \verbatim
67 *> NMATS is INTEGER
68 *> The number of matrix types to be tested for each combination
69 *> of matrix dimensions. If NMATS >= NTYPES (the maximum
70 *> number of matrix types), then all the different types are
71 *> generated for testing. If NMATS < NTYPES, another input line
72 *> is read to get the numbers of the matrix types to be used.
73 *> \endverbatim
74 *>
75 *> \param[in,out] ISEED
76 *> \verbatim
77 *> ISEED is INTEGER array, dimension (4)
78 *> On entry, the seed of the random number generator. The array
79 *> elements should be between 0 and 4095, otherwise they will be
80 *> reduced mod 4096, and ISEED(4) must be odd.
81 *> On exit, the next seed in the random number sequence after
82 *> all the test matrices have been generated.
83 *> \endverbatim
84 *>
85 *> \param[in] THRESH
86 *> \verbatim
87 *> THRESH is REAL
88 *> The threshold value for the test ratios. A result is
89 *> included in the output file if RESULT >= THRESH. To have
90 *> every test ratio printed, use THRESH = 0.
91 *> \endverbatim
92 *>
93 *> \param[in] MMAX
94 *> \verbatim
95 *> MMAX is INTEGER
96 *> The maximum value permitted for M, used in dimensioning the
97 *> work arrays.
98 *> \endverbatim
99 *>
100 *> \param[out] X
101 *> \verbatim
102 *> X is REAL array, dimension (MMAX*MMAX)
103 *> \endverbatim
104 *>
105 *> \param[out] XF
106 *> \verbatim
107 *> XF is REAL array, dimension (MMAX*MMAX)
108 *> \endverbatim
109 *>
110 *> \param[out] U1
111 *> \verbatim
112 *> U1 is REAL array, dimension (MMAX*MMAX)
113 *> \endverbatim
114 *>
115 *> \param[out] U2
116 *> \verbatim
117 *> U2 is REAL array, dimension (MMAX*MMAX)
118 *> \endverbatim
119 *>
120 *> \param[out] V1T
121 *> \verbatim
122 *> V1T is REAL array, dimension (MMAX*MMAX)
123 *> \endverbatim
124 *>
125 *> \param[out] V2T
126 *> \verbatim
127 *> V2T is REAL array, dimension (MMAX*MMAX)
128 *> \endverbatim
129 *>
130 *> \param[out] THETA
131 *> \verbatim
132 *> THETA is REAL array, dimension (MMAX)
133 *> \endverbatim
134 *>
135 *> \param[out] IWORK
136 *> \verbatim
137 *> IWORK is INTEGER array, dimension (MMAX)
138 *> \endverbatim
139 *>
140 *> \param[out] WORK
141 *> \verbatim
142 *> WORK is REAL array
143 *> \endverbatim
144 *>
145 *> \param[out] RWORK
146 *> \verbatim
147 *> RWORK is REAL array
148 *> \endverbatim
149 *>
150 *> \param[in] NIN
151 *> \verbatim
152 *> NIN is INTEGER
153 *> The unit number for input.
154 *> \endverbatim
155 *>
156 *> \param[in] NOUT
157 *> \verbatim
158 *> NOUT is INTEGER
159 *> The unit number for output.
160 *> \endverbatim
161 *>
162 *> \param[out] INFO
163 *> \verbatim
164 *> INFO is INTEGER
165 *> = 0 : successful exit
166 *> > 0 : If SLAROR returns an error code, the absolute value
167 *> of it is returned.
168 *> \endverbatim
169 *
170 * Authors:
171 * ========
172 *
173 *> \author Univ. of Tennessee
174 *> \author Univ. of California Berkeley
175 *> \author Univ. of Colorado Denver
176 *> \author NAG Ltd.
177 *
178 *> \ingroup single_eig
179 *
180 * =====================================================================
181  SUBROUTINE sckcsd( NM, MVAL, PVAL, QVAL, NMATS, ISEED, THRESH,
182  $ MMAX, X, XF, U1, U2, V1T, V2T, THETA, IWORK,
183  $ WORK, RWORK, NIN, NOUT, INFO )
184 *
185 * -- LAPACK test routine --
186 * -- LAPACK is a software package provided by Univ. of Tennessee, --
187 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
188 *
189 * .. Scalar Arguments ..
190  INTEGER INFO, NIN, NM, NMATS, MMAX, NOUT
191  REAL THRESH
192 * ..
193 * .. Array Arguments ..
194  INTEGER ISEED( 4 ), IWORK( * ), MVAL( * ), PVAL( * ),
195  $ QVAL( * )
196  REAL RWORK( * ), THETA( * )
197  REAL U1( * ), U2( * ), V1T( * ), V2T( * ),
198  $ work( * ), x( * ), xf( * )
199 * ..
200 *
201 * =====================================================================
202 *
203 * .. Parameters ..
204  INTEGER NTESTS
205  PARAMETER ( NTESTS = 15 )
206  INTEGER NTYPES
207  parameter( ntypes = 4 )
208  REAL GAPDIGIT, ONE, ORTH, TEN, ZERO
209  parameter( gapdigit = 10.0e0, one = 1.0e0,
210  $ orth = 1.0e-4, ten = 10.0e0, zero = 0.0e0 )
211  REAL PIOVER2
212  PARAMETER ( PIOVER2 = 1.57079632679489661923132169163975144210e0 )
213 * ..
214 * .. Local Scalars ..
215  LOGICAL FIRSTT
216  CHARACTER*3 PATH
217  INTEGER I, IINFO, IM, IMAT, J, LDU1, LDU2, LDV1T,
218  $ ldv2t, ldx, lwork, m, nfail, nrun, nt, p, q, r
219 * ..
220 * .. Local Arrays ..
221  LOGICAL DOTYPE( NTYPES )
222  REAL RESULT( NTESTS )
223 * ..
224 * .. External Subroutines ..
225  EXTERNAL alahdg, alareq, alasum, scsdts, slacsg, slaror,
226  $ slaset, srot
227 * ..
228 * .. Intrinsic Functions ..
229  INTRINSIC abs, min
230 * ..
231 * .. External Functions ..
232  REAL SLARAN, SLARND
233  EXTERNAL SLARAN, SLARND
234 * ..
235 * .. Executable Statements ..
236 *
237 * Initialize constants and the random number seed.
238 *
239  path( 1: 3 ) = 'CSD'
240  info = 0
241  nrun = 0
242  nfail = 0
243  firstt = .true.
244  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
245  ldx = mmax
246  ldu1 = mmax
247  ldu2 = mmax
248  ldv1t = mmax
249  ldv2t = mmax
250  lwork = mmax*mmax
251 *
252 * Do for each value of M in MVAL.
253 *
254  DO 30 im = 1, nm
255  m = mval( im )
256  p = pval( im )
257  q = qval( im )
258 *
259  DO 20 imat = 1, ntypes
260 *
261 * Do the tests only if DOTYPE( IMAT ) is true.
262 *
263  IF( .NOT.dotype( imat ) )
264  $ GO TO 20
265 *
266 * Generate X
267 *
268  IF( imat.EQ.1 ) THEN
269  CALL slaror( 'L', 'I', m, m, x, ldx, iseed, work, iinfo )
270  IF( m .NE. 0 .AND. iinfo .NE. 0 ) THEN
271  WRITE( nout, fmt = 9999 ) m, iinfo
272  info = abs( iinfo )
273  GO TO 20
274  END IF
275  ELSE IF( imat.EQ.2 ) THEN
276  r = min( p, m-p, q, m-q )
277  DO i = 1, r
278  theta(i) = piover2 * slarnd( 1, iseed )
279  END DO
280  CALL slacsg( m, p, q, theta, iseed, x, ldx, work )
281  DO i = 1, m
282  DO j = 1, m
283  x(i+(j-1)*ldx) = x(i+(j-1)*ldx) +
284  $ orth*slarnd(2,iseed)
285  END DO
286  END DO
287  ELSE IF( imat.EQ.3 ) THEN
288  r = min( p, m-p, q, m-q )
289  DO i = 1, r+1
290  theta(i) = ten**(-slarnd(1,iseed)*gapdigit)
291  END DO
292  DO i = 2, r+1
293  theta(i) = theta(i-1) + theta(i)
294  END DO
295  DO i = 1, r
296  theta(i) = piover2 * theta(i) / theta(r+1)
297  END DO
298  CALL slacsg( m, p, q, theta, iseed, x, ldx, work )
299  ELSE
300  CALL slaset( 'F', m, m, zero, one, x, ldx )
301  DO i = 1, m
302  j = int( slaran( iseed ) * m ) + 1
303  IF( j .NE. i ) THEN
304  CALL srot( m, x(1+(i-1)*ldx), 1, x(1+(j-1)*ldx), 1,
305  $ zero, one )
306  END IF
307  END DO
308  END IF
309 *
310  nt = 15
311 *
312  CALL scsdts( m, p, q, x, xf, ldx, u1, ldu1, u2, ldu2, v1t,
313  $ ldv1t, v2t, ldv2t, theta, iwork, work, lwork,
314  $ rwork, result )
315 *
316 * Print information about the tests that did not
317 * pass the threshold.
318 *
319  DO 10 i = 1, nt
320  IF( result( i ).GE.thresh ) THEN
321  IF( nfail.EQ.0 .AND. firstt ) THEN
322  firstt = .false.
323  CALL alahdg( nout, path )
324  END IF
325  WRITE( nout, fmt = 9998 )m, p, q, imat, i,
326  $ result( i )
327  nfail = nfail + 1
328  END IF
329  10 CONTINUE
330  nrun = nrun + nt
331  20 CONTINUE
332  30 CONTINUE
333 *
334 * Print a summary of the results.
335 *
336  CALL alasum( path, nout, nfail, nrun, 0 )
337 *
338  9999 FORMAT( ' SLAROR in SCKCSD: M = ', i5, ', INFO = ', i15 )
339  9998 FORMAT( ' M=', i4, ' P=', i4, ', Q=', i4, ', type ', i2,
340  $ ', test ', i2, ', ratio=', g13.6 )
341  RETURN
342 *
343 * End of SCKCSD
344 *
345  END
346 *
347 *
348 *
349  SUBROUTINE slacsg( M, P, Q, THETA, ISEED, X, LDX, WORK )
350  IMPLICIT NONE
351 *
352  INTEGER LDX, M, P, Q
353  INTEGER ISEED( 4 )
354  REAL THETA( * )
355  REAL WORK( * ), X( LDX, * )
356 *
357  REAL ONE, ZERO
358  PARAMETER ( ONE = 1.0e0, zero = 0.0e0 )
359 *
360  INTEGER I, INFO, R
361 *
362  r = min( p, m-p, q, m-q )
363 *
364  CALL slaset( 'Full', m, m, zero, zero, x, ldx )
365 *
366  DO i = 1, min(p,q)-r
367  x(i,i) = one
368  END DO
369  DO i = 1, r
370  x(min(p,q)-r+i,min(p,q)-r+i) = cos(theta(i))
371  END DO
372  DO i = 1, min(p,m-q)-r
373  x(p-i+1,m-i+1) = -one
374  END DO
375  DO i = 1, r
376  x(p-(min(p,m-q)-r)+1-i,m-(min(p,m-q)-r)+1-i) =
377  $ -sin(theta(r-i+1))
378  END DO
379  DO i = 1, min(m-p,q)-r
380  x(m-i+1,q-i+1) = one
381  END DO
382  DO i = 1, r
383  x(m-(min(m-p,q)-r)+1-i,q-(min(m-p,q)-r)+1-i) =
384  $ sin(theta(r-i+1))
385  END DO
386  DO i = 1, min(m-p,m-q)-r
387  x(p+i,q+i) = one
388  END DO
389  DO i = 1, r
390  x(p+(min(m-p,m-q)-r)+i,q+(min(m-p,m-q)-r)+i) =
391  $ cos(theta(i))
392  END DO
393  CALL slaror( 'Left', 'No init', p, m, x, ldx, iseed, work, info )
394  CALL slaror( 'Left', 'No init', m-p, m, x(p+1,1), ldx,
395  $ iseed, work, info )
396  CALL slaror( 'Right', 'No init', m, q, x, ldx, iseed,
397  $ work, info )
398  CALL slaror( 'Right', 'No init', m, m-q,
399  $ x(1,q+1), ldx, iseed, work, info )
400 *
401  END
402 
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: slaset.f:110
subroutine alareq(PATH, NMATS, DOTYPE, NTYPES, NIN, NOUT)
ALAREQ
Definition: alareq.f:90
subroutine alahdg(IOUNIT, PATH)
ALAHDG
Definition: alahdg.f:62
subroutine alasum(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASUM
Definition: alasum.f:73
subroutine slaror(SIDE, INIT, M, N, A, LDA, ISEED, X, INFO)
SLAROR
Definition: slaror.f:146
subroutine srot(N, SX, INCX, SY, INCY, C, S)
SROT
Definition: srot.f:92
subroutine sckcsd(NM, MVAL, PVAL, QVAL, NMATS, ISEED, THRESH, MMAX, X, XF, U1, U2, V1T, V2T, THETA, IWORK, WORK, RWORK, NIN, NOUT, INFO)
SCKCSD
Definition: sckcsd.f:184
subroutine scsdts(M, P, Q, X, XF, LDX, U1, LDU1, U2, LDU2, V1T, LDV1T, V2T, LDV2T, THETA, IWORK, WORK, LWORK, RWORK, RESULT)
SCSDTS
Definition: scsdts.f:229
subroutine slacsg(M, P, Q, THETA, ISEED, X, LDX, WORK)
Definition: sckcsd.f:350