1*> \brief \b CQRT12
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*  Definition:
9*  ===========
10*
11*       REAL             FUNCTION CQRT12( M, N, A, LDA, S, WORK, LWORK,
12*                        RWORK )
13*
14*       .. Scalar Arguments ..
15*       INTEGER            LDA, LWORK, M, N
16*       ..
17*       .. Array Arguments ..
18*       REAL               RWORK( * ), S( * )
19*       COMPLEX            A( LDA, * ), WORK( LWORK )
20*       ..
21*
22*
23*> \par Purpose:
24*  =============
25*>
26*> \verbatim
27*>
28*> CQRT12 computes the singular values `svlues' of the upper trapezoid
29*> of A(1:M,1:N) and returns the ratio
30*>
31*>      || s - svlues||/(||svlues||*eps*max(M,N))
32*> \endverbatim
33*
34*  Arguments:
35*  ==========
36*
37*> \param[in] M
38*> \verbatim
39*>          M is INTEGER
40*>          The number of rows of the matrix A.
41*> \endverbatim
42*>
43*> \param[in] N
44*> \verbatim
45*>          N is INTEGER
46*>          The number of columns of the matrix A.
47*> \endverbatim
48*>
49*> \param[in] A
50*> \verbatim
51*>          A is COMPLEX array, dimension (LDA,N)
52*>          The M-by-N matrix A. Only the upper trapezoid is referenced.
53*> \endverbatim
54*>
55*> \param[in] LDA
56*> \verbatim
57*>          LDA is INTEGER
58*>          The leading dimension of the array A.
59*> \endverbatim
60*>
61*> \param[in] S
62*> \verbatim
63*>          S is REAL array, dimension (min(M,N))
64*>          The singular values of the matrix A.
65*> \endverbatim
66*>
67*> \param[out] WORK
68*> \verbatim
69*>          WORK is COMPLEX array, dimension (LWORK)
70*> \endverbatim
71*>
72*> \param[in] LWORK
73*> \verbatim
74*>          LWORK is INTEGER
75*>          The length of the array WORK. LWORK >= M*N + 2*min(M,N) +
76*>          max(M,N).
77*> \endverbatim
78*>
79*> \param[out] RWORK
80*> \verbatim
81*>          RWORK is REAL array, dimension (4*min(M,N))
82*> \endverbatim
83*
84*  Authors:
85*  ========
86*
87*> \author Univ. of Tennessee
88*> \author Univ. of California Berkeley
89*> \author Univ. of Colorado Denver
90*> \author NAG Ltd.
91*
92*> \ingroup complex_lin
93*
94*  =====================================================================
95      REAL             FUNCTION CQRT12( M, N, A, LDA, S, WORK, LWORK,
96     $                 RWORK )
97*
98*  -- LAPACK test routine --
99*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
100*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
101*
102*     .. Scalar Arguments ..
103      INTEGER            LDA, LWORK, M, N
104*     ..
105*     .. Array Arguments ..
106      REAL               RWORK( * ), S( * )
107      COMPLEX            A( LDA, * ), WORK( LWORK )
108*     ..
109*
110*  =====================================================================
111*
112*     .. Parameters ..
113      REAL               ZERO, ONE
114      PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0 )
115*     ..
116*     .. Local Scalars ..
117      INTEGER            I, INFO, ISCL, J, MN
118      REAL               ANRM, BIGNUM, NRMSVL, SMLNUM
119*     ..
120*     .. Local Arrays ..
121      REAL               DUMMY( 1 )
122*     ..
123*     .. External Functions ..
124      REAL               CLANGE, SASUM, SLAMCH, SNRM2
125      EXTERNAL           CLANGE, SASUM, SLAMCH, SNRM2
126*     ..
127*     .. External Subroutines ..
128      EXTERNAL           CGEBD2, CLASCL, CLASET, SAXPY, SBDSQR, SLABAD,
129     $                   SLASCL, XERBLA
130*     ..
131*     .. Intrinsic Functions ..
132      INTRINSIC          CMPLX, MAX, MIN, REAL
133*     ..
134*     .. Executable Statements ..
135*
136      CQRT12 = ZERO
137*
138*     Test that enough workspace is supplied
139*
140      IF( LWORK.LT.M*N+2*MIN( M, N )+MAX( M, N ) ) THEN
141         CALL XERBLA( 'CQRT12', 7 )
142         RETURN
143      END IF
144*
145*     Quick return if possible
146*
147      MN = MIN( M, N )
148      IF( MN.LE.ZERO )
149     $   RETURN
150*
151      NRMSVL = SNRM2( MN, S, 1 )
152*
153*     Copy upper triangle of A into work
154*
155      CALL CLASET( 'Full', M, N, CMPLX( ZERO ), CMPLX( ZERO ), WORK, M )
156      DO 20 J = 1, N
157         DO 10 I = 1, MIN( J, M )
158            WORK( ( J-1 )*M+I ) = A( I, J )
159   10    CONTINUE
160   20 CONTINUE
161*
162*     Get machine parameters
163*
164      SMLNUM = SLAMCH( 'S' ) / SLAMCH( 'P' )
165      BIGNUM = ONE / SMLNUM
166      CALL SLABAD( SMLNUM, BIGNUM )
167*
168*     Scale work if max entry outside range [SMLNUM,BIGNUM]
169*
170      ANRM = CLANGE( 'M', M, N, WORK, M, DUMMY )
171      ISCL = 0
172      IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN
173*
174*        Scale matrix norm up to SMLNUM
175*
176         CALL CLASCL( 'G', 0, 0, ANRM, SMLNUM, M, N, WORK, M, INFO )
177         ISCL = 1
178      ELSE IF( ANRM.GT.BIGNUM ) THEN
179*
180*        Scale matrix norm down to BIGNUM
181*
182         CALL CLASCL( 'G', 0, 0, ANRM, BIGNUM, M, N, WORK, M, INFO )
183         ISCL = 1
184      END IF
185*
186      IF( ANRM.NE.ZERO ) THEN
187*
188*        Compute SVD of work
189*
190         CALL CGEBD2( M, N, WORK, M, RWORK( 1 ), RWORK( MN+1 ),
191     $                WORK( M*N+1 ), WORK( M*N+MN+1 ),
192     $                WORK( M*N+2*MN+1 ), INFO )
193         CALL SBDSQR( 'Upper', MN, 0, 0, 0, RWORK( 1 ), RWORK( MN+1 ),
194     $                DUMMY, MN, DUMMY, 1, DUMMY, MN, RWORK( 2*MN+1 ),
195     $                INFO )
196*
197         IF( ISCL.EQ.1 ) THEN
198            IF( ANRM.GT.BIGNUM ) THEN
199               CALL SLASCL( 'G', 0, 0, BIGNUM, ANRM, MN, 1, RWORK( 1 ),
200     $                      MN, INFO )
201            END IF
202            IF( ANRM.LT.SMLNUM ) THEN
203               CALL SLASCL( 'G', 0, 0, SMLNUM, ANRM, MN, 1, RWORK( 1 ),
204     $                      MN, INFO )
205            END IF
206         END IF
207*
208      ELSE
209*
210         DO 30 I = 1, MN
211            RWORK( I ) = ZERO
212   30    CONTINUE
213      END IF
214*
215*     Compare s and singular values of work
216*
217      CALL SAXPY( MN, -ONE, S, 1, RWORK( 1 ), 1 )
218      CQRT12 = SASUM( MN, RWORK( 1 ), 1 ) /
219     $         ( SLAMCH( 'Epsilon' )*REAL( MAX( M, N ) ) )
220      IF( NRMSVL.NE.ZERO )
221     $   CQRT12 = CQRT12 / NRMSVL
222*
223      RETURN
224*
225*     End of CQRT12
226*
227      END
228