97
98
99
100
101
102
103 INTEGER LDA, LWORK, M, N
104
105
106 DOUBLE PRECISION RWORK( * ), S( * )
107 COMPLEX*16 A( LDA, * ), WORK( LWORK )
108
109
110
111
112
113 DOUBLE PRECISION ZERO, ONE
114 parameter( zero = 0.0d0, one = 1.0d0 )
115
116
117 INTEGER I, INFO, ISCL, J, MN
118 DOUBLE PRECISION ANRM, BIGNUM, NRMSVL, SMLNUM
119
120
121 DOUBLE PRECISION DUMMY( 1 )
122
123
124 DOUBLE PRECISION DASUM, DLAMCH, DNRM2, ZLANGE
126
127
130
131
132 INTRINSIC dble, dcmplx, max, min
133
134
135
137
138
139
140 IF( lwork.LT.m*n+2*min( m, n )+max( m, n ) ) THEN
141 CALL xerbla(
'ZQRT12', 7 )
142 RETURN
143 END IF
144
145
146
147 mn = min( m, n )
148 IF( mn.LE.zero )
149 $ RETURN
150
151 nrmsvl =
dnrm2( mn, s, 1 )
152
153
154
155 CALL zlaset(
'Full', m, n, dcmplx( zero ), dcmplx( zero ), work,
156 $ m )
157 DO j = 1, n
158 DO i = 1, min( j, m )
159 work( ( j-1 )*m+i ) = a( i, j )
160 END DO
161 END DO
162
163
164
166 bignum = one / smlnum
167
168
169
170 anrm =
zlange(
'M', m, n, work, m, dummy )
171 iscl = 0
172 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
173
174
175
176 CALL zlascl(
'G', 0, 0, anrm, smlnum, m, n, work, m, info )
177 iscl = 1
178 ELSE IF( anrm.GT.bignum ) THEN
179
180
181
182 CALL zlascl(
'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
189
190 CALL zgebd2( 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 dbdsqr(
'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 dlascl(
'G', 0, 0, bignum, anrm, mn, 1, rwork( 1 ),
200 $ mn, info )
201 END IF
202 IF( anrm.LT.smlnum ) THEN
203 CALL dlascl(
'G', 0, 0, smlnum, anrm, mn, 1, rwork( 1 ),
204 $ mn, info )
205 END IF
206 END IF
207
208 ELSE
209
210 DO i = 1, mn
211 rwork( i ) = zero
212 END DO
213 END IF
214
215
216
217 CALL daxpy( mn, -one, s, 1, rwork( 1 ), 1 )
219 $ (
dlamch(
'Epsilon' )*dble( max( m, n ) ) )
220
221 IF( nrmsvl.NE.zero )
223
224 RETURN
225
226
227
subroutine xerbla(srname, info)
double precision function dasum(n, dx, incx)
DASUM
subroutine daxpy(n, da, dx, incx, dy, incy)
DAXPY
subroutine dbdsqr(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, work, info)
DBDSQR
subroutine zgebd2(m, n, a, lda, d, e, tauq, taup, work, info)
ZGEBD2 reduces a general matrix to bidiagonal form using an unblocked algorithm.
double precision function dlamch(cmach)
DLAMCH
double precision function zlange(norm, m, n, a, lda, work)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
subroutine zlascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
ZLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
subroutine dlascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
subroutine zlaset(uplo, m, n, alpha, beta, a, lda)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
real(wp) function dnrm2(n, x, incx)
DNRM2
double precision function zqrt12(m, n, a, lda, s, work, lwork, rwork)
ZQRT12