LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
cckgsv.f
Go to the documentation of this file.
1*> \brief \b CCKGSV
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 CCKGSV( NM, MVAL, PVAL, NVAL, NMATS, ISEED, THRESH,
12* NMAX, A, AF, B, BF, U, V, Q, ALPHA, BETA, R,
13* IWORK, WORK, RWORK, NIN, NOUT, INFO )
14*
15* .. Scalar Arguments ..
16* INTEGER INFO, NIN, NM, NMATS, NMAX, NOUT
17* REAL THRESH
18* ..
19* .. Array Arguments ..
20* INTEGER ISEED( 4 ), IWORK( * ), MVAL( * ), NVAL( * ),
21* $ PVAL( * )
22* REAL ALPHA( * ), BETA( * ), RWORK( * )
23* COMPLEX A( * ), AF( * ), B( * ), BF( * ), Q( * ),
24* $ R( * ), U( * ), V( * ), WORK( * )
25* ..
26*
27*
28*> \par Purpose:
29* =============
30*>
31*> \verbatim
32*>
33*> CCKGSV tests CGGSVD:
34*> the GSVD for M-by-N matrix A and P-by-N matrix B.
35*> \endverbatim
36*
37* Arguments:
38* ==========
39*
40*> \param[in] NM
41*> \verbatim
42*> NM is INTEGER
43*> The number of values of M contained in the vector MVAL.
44*> \endverbatim
45*>
46*> \param[in] MVAL
47*> \verbatim
48*> MVAL is INTEGER array, dimension (NM)
49*> The values of the matrix row dimension M.
50*> \endverbatim
51*>
52*> \param[in] PVAL
53*> \verbatim
54*> PVAL is INTEGER array, dimension (NP)
55*> The values of the matrix row dimension P.
56*> \endverbatim
57*>
58*> \param[in] NVAL
59*> \verbatim
60*> NVAL is INTEGER array, dimension (NN)
61*> The values of the matrix column dimension N.
62*> \endverbatim
63*>
64*> \param[in] NMATS
65*> \verbatim
66*> NMATS is INTEGER
67*> The number of matrix types to be tested for each combination
68*> of matrix dimensions. If NMATS >= NTYPES (the maximum
69*> number of matrix types), then all the different types are
70*> generated for testing. If NMATS < NTYPES, another input line
71*> is read to get the numbers of the matrix types to be used.
72*> \endverbatim
73*>
74*> \param[in,out] ISEED
75*> \verbatim
76*> ISEED is INTEGER array, dimension (4)
77*> On entry, the seed of the random number generator. The array
78*> elements should be between 0 and 4095, otherwise they will be
79*> reduced mod 4096, and ISEED(4) must be odd.
80*> On exit, the next seed in the random number sequence after
81*> all the test matrices have been generated.
82*> \endverbatim
83*>
84*> \param[in] THRESH
85*> \verbatim
86*> THRESH is REAL
87*> The threshold value for the test ratios. A result is
88*> included in the output file if RESULT >= THRESH. To have
89*> every test ratio printed, use THRESH = 0.
90*> \endverbatim
91*>
92*> \param[in] NMAX
93*> \verbatim
94*> NMAX is INTEGER
95*> The maximum value permitted for M or N, used in dimensioning
96*> the work arrays.
97*> \endverbatim
98*>
99*> \param[out] A
100*> \verbatim
101*> A is COMPLEX array, dimension (NMAX*NMAX)
102*> \endverbatim
103*>
104*> \param[out] AF
105*> \verbatim
106*> AF is COMPLEX array, dimension (NMAX*NMAX)
107*> \endverbatim
108*>
109*> \param[out] B
110*> \verbatim
111*> B is COMPLEX array, dimension (NMAX*NMAX)
112*> \endverbatim
113*>
114*> \param[out] BF
115*> \verbatim
116*> BF is COMPLEX array, dimension (NMAX*NMAX)
117*> \endverbatim
118*>
119*> \param[out] U
120*> \verbatim
121*> U is COMPLEX array, dimension (NMAX*NMAX)
122*> \endverbatim
123*>
124*> \param[out] V
125*> \verbatim
126*> V is COMPLEX array, dimension (NMAX*NMAX)
127*> \endverbatim
128*>
129*> \param[out] Q
130*> \verbatim
131*> Q is COMPLEX array, dimension (NMAX*NMAX)
132*> \endverbatim
133*>
134*> \param[out] ALPHA
135*> \verbatim
136*> ALPHA is REAL array, dimension (NMAX)
137*> \endverbatim
138*>
139*> \param[out] BETA
140*> \verbatim
141*> BETA is REAL array, dimension (NMAX)
142*> \endverbatim
143*>
144*> \param[out] R
145*> \verbatim
146*> R is COMPLEX array, dimension (NMAX*NMAX)
147*> \endverbatim
148*>
149*> \param[out] IWORK
150*> \verbatim
151*> IWORK is INTEGER array, dimension (NMAX)
152*> \endverbatim
153*>
154*> \param[out] WORK
155*> \verbatim
156*> WORK is COMPLEX array, dimension (NMAX*NMAX)
157*> \endverbatim
158*>
159*> \param[out] RWORK
160*> \verbatim
161*> RWORK is REAL array, dimension (NMAX)
162*> \endverbatim
163*>
164*> \param[in] NIN
165*> \verbatim
166*> NIN is INTEGER
167*> The unit number for input.
168*> \endverbatim
169*>
170*> \param[in] NOUT
171*> \verbatim
172*> NOUT is INTEGER
173*> The unit number for output.
174*> \endverbatim
175*>
176*> \param[out] INFO
177*> \verbatim
178*> INFO is INTEGER
179*> = 0 : successful exit
180*> > 0 : If CLATMS returns an error code, the absolute value
181*> of it is returned.
182*> \endverbatim
183*
184* Authors:
185* ========
186*
187*> \author Univ. of Tennessee
188*> \author Univ. of California Berkeley
189*> \author Univ. of Colorado Denver
190*> \author NAG Ltd.
191*
192*> \ingroup complex_eig
193*
194* =====================================================================
195 SUBROUTINE cckgsv( NM, MVAL, PVAL, NVAL, NMATS, ISEED, THRESH,
196 $ NMAX, A, AF, B, BF, U, V, Q, ALPHA, BETA, R,
197 $ IWORK, WORK, RWORK, NIN, NOUT, INFO )
198*
199* -- LAPACK test routine --
200* -- LAPACK is a software package provided by Univ. of Tennessee, --
201* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
202*
203* .. Scalar Arguments ..
204 INTEGER INFO, NIN, NM, NMATS, NMAX, NOUT
205 REAL THRESH
206* ..
207* .. Array Arguments ..
208 INTEGER ISEED( 4 ), IWORK( * ), MVAL( * ), NVAL( * ),
209 $ PVAL( * )
210 REAL ALPHA( * ), BETA( * ), RWORK( * )
211 COMPLEX A( * ), AF( * ), B( * ), BF( * ), Q( * ),
212 $ r( * ), u( * ), v( * ), work( * )
213* ..
214*
215* =====================================================================
216*
217* .. Parameters ..
218 INTEGER NTESTS
219 PARAMETER ( NTESTS = 12 )
220 INTEGER NTYPES
221 parameter( ntypes = 8 )
222* ..
223* .. Local Scalars ..
224 LOGICAL FIRSTT
225 CHARACTER DISTA, DISTB, TYPE
226 CHARACTER*3 PATH
227 INTEGER I, IINFO, IM, IMAT, KLA, KLB, KUA, KUB, LDA,
228 $ ldb, ldq, ldr, ldu, ldv, lwork, m, modea,
229 $ modeb, n, nfail, nrun, nt, p, k, l
230 REAL ANORM, BNORM, CNDNMA, CNDNMB
231* ..
232* .. Local Arrays ..
233 LOGICAL DOTYPE( NTYPES )
234 REAL RESULT( NTESTS )
235* ..
236* .. External Subroutines ..
238* ..
239* .. Intrinsic Functions ..
240 INTRINSIC abs
241* ..
242* .. Executable Statements ..
243*
244* Initialize constants and the random number seed.
245*
246 path( 1: 3 ) = 'GSV'
247 info = 0
248 nrun = 0
249 nfail = 0
250 firstt = .true.
251 CALL alareq( path, nmats, dotype, ntypes, nin, nout )
252 lda = nmax
253 ldb = nmax
254 ldu = nmax
255 ldv = nmax
256 ldq = nmax
257 ldr = nmax
258 lwork = nmax*nmax
259*
260* Specific cases
261*
262* Test: https://github.com/Reference-LAPACK/lapack/issues/411#issue-608776973
263*
264 m = 6
265 p = 6
266 n = 6
267 a(1:m*n) = cmplx(1.e0, 0.e0)
268 b(1:m*n) = cmplx(0.e0, 0.e0)
269 b(1+0*m) = cmplx(9.e19, 0.e0)
270 b(2+1*m) = cmplx(9.e18, 0.e0)
271 b(3+2*m) = cmplx(9.e17, 0.e0)
272 b(4+3*m) = cmplx(9.e16, 0.e0)
273 b(5+4*m) = cmplx(9.e15, 0.e0)
274 b(6+5*m) = cmplx(9.e14, 0.e0)
275 CALL cggsvd3('N','N','N', m, p, n, k, l, a, m, b, m,
276 $ alpha, beta, u, 1, v, 1, q, 1,
277 $ work, m*n, rwork, iwork, info)
278*
279* Print information there is a NAN in BETA
280 DO 40 i = 1, l
281 IF( beta(i).NE.beta(i) ) THEN
282 info = -i
283 EXIT
284 END IF
285 40 CONTINUE
286 IF( info.LT.0 ) THEN
287 IF( nfail.EQ.0 .AND. firstt ) THEN
288 firstt = .false.
289 CALL alahdg( nout, path )
290 END IF
291 WRITE( nout, fmt = 9997 ) -info
292 nfail = nfail + 1
293 END IF
294 nrun = nrun + 1
295 info = 0
296*
297* Do for each value of M in MVAL.
298*
299 DO 30 im = 1, nm
300 m = mval( im )
301 p = pval( im )
302 n = nval( im )
303*
304 DO 20 imat = 1, ntypes
305*
306* Do the tests only if DOTYPE( IMAT ) is true.
307*
308 IF( .NOT.dotype( imat ) )
309 $ GO TO 20
310*
311* Set up parameters with SLATB9 and generate test
312* matrices A and B with CLATMS.
313*
314 CALL slatb9( path, imat, m, p, n, TYPE, kla, kua, klb, kub,
315 $ anorm, bnorm, modea, modeb, cndnma, cndnmb,
316 $ dista, distb )
317*
318* Generate M by N matrix A
319*
320 CALL clatms( m, n, dista, iseed, TYPE, rwork, modea, cndnma,
321 $ anorm, kla, kua, 'No packing', a, lda, work,
322 $ iinfo )
323 IF( iinfo.NE.0 ) THEN
324 WRITE( nout, fmt = 9999 )iinfo
325 info = abs( iinfo )
326 GO TO 20
327 END IF
328*
329* Generate P by N matrix B
330*
331 CALL clatms( p, n, distb, iseed, TYPE, rwork, modeb, cndnmb,
332 $ bnorm, klb, kub, 'No packing', b, ldb, work,
333 $ iinfo )
334 IF( iinfo.NE.0 ) THEN
335 WRITE( nout, fmt = 9999 )iinfo
336 info = abs( iinfo )
337 GO TO 20
338 END IF
339*
340 nt = 6
341*
342 CALL cgsvts3( m, p, n, a, af, lda, b, bf, ldb, u, ldu, v,
343 $ ldv, q, ldq, alpha, beta, r, ldr, iwork, work,
344 $ lwork, rwork, result )
345*
346* Print information about the tests that did not
347* pass the threshold.
348*
349 DO 10 i = 1, nt
350 IF( result( i ).GE.thresh ) THEN
351 IF( nfail.EQ.0 .AND. firstt ) THEN
352 firstt = .false.
353 CALL alahdg( nout, path )
354 END IF
355 WRITE( nout, fmt = 9998 )m, p, n, imat, i,
356 $ result( i )
357 nfail = nfail + 1
358 END IF
359 10 CONTINUE
360 nrun = nrun + nt
361*
362 20 CONTINUE
363 30 CONTINUE
364*
365* Print a summary of the results.
366*
367 CALL alasum( path, nout, nfail, nrun, 0 )
368*
369 9999 FORMAT( ' CLATMS in CCKGSV INFO = ', i5 )
370 9998 FORMAT( ' M=', i4, ' P=', i4, ', N=', i4, ', type ', i2,
371 $ ', test ', i2, ', ratio=', g13.6 )
372 9997 FORMAT( ' FOUND NaN in BETA(', i4,')' )
373 RETURN
374*
375* End of CCKGSV
376*
377 END
subroutine alareq(path, nmats, dotype, ntypes, nin, nout)
ALAREQ
Definition alareq.f:90
subroutine alasum(type, nout, nfail, nrun, nerrs)
ALASUM
Definition alasum.f:73
subroutine alahdg(iounit, path)
ALAHDG
Definition alahdg.f:62
subroutine cckgsv(nm, mval, pval, nval, nmats, iseed, thresh, nmax, a, af, b, bf, u, v, q, alpha, beta, r, iwork, work, rwork, nin, nout, info)
CCKGSV
Definition cckgsv.f:198
subroutine cgsvts3(m, p, n, a, af, lda, b, bf, ldb, u, ldu, v, ldv, q, ldq, alpha, beta, r, ldr, iwork, work, lwork, rwork, result)
CGSVTS3
Definition cgsvts3.f:209
subroutine clatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
CLATMS
Definition clatms.f:332
subroutine cggsvd3(jobu, jobv, jobq, m, n, p, k, l, a, lda, b, ldb, alpha, beta, u, ldu, v, ldv, q, ldq, work, lwork, rwork, iwork, info)
CGGSVD3 computes the singular value decomposition (SVD) for OTHER matrices
Definition cggsvd3.f:354
subroutine slatb9(path, imat, m, p, n, type, kla, kua, klb, kub, anorm, bnorm, modea, modeb, cndnma, cndnmb, dista, distb)
SLATB9
Definition slatb9.f:170