89
90
91
92
93
94
95 INTEGER LDA, LWORK, M, N
96
97
98 REAL A( LDA, * ), S( * ), WORK( LWORK )
99
100
101
102
103
104 REAL ZERO, ONE
105 parameter( zero = 0.0e0, one = 1.0e0 )
106
107
108 INTEGER I, INFO, ISCL, J, MN
109 REAL ANRM, BIGNUM, NRMSVL, SMLNUM
110
111
112 REAL SASUM, SLAMCH, SLANGE, SNRM2
114
115
117
118
119 INTRINSIC max, min, real
120
121
122 REAL DUMMY( 1 )
123
124
125
127
128
129
130 IF( lwork.LT.max( m*n+4*min( m, n )+max( m, n ),
131 $ m*n+2*min( m, n )+4*n) ) THEN
132 CALL xerbla(
'SQRT12', 7 )
133 RETURN
134 END IF
135
136
137
138 mn = min( m, n )
139 IF( mn.LE.zero )
140 $ RETURN
141
142 nrmsvl =
snrm2( mn, s, 1 )
143
144
145
146 CALL slaset(
'Full', m, n, zero, zero, work, m )
147 DO j = 1, n
148 DO i = 1, min( j, m )
149 work( ( j-1 )*m+i ) = a( i, j )
150 END DO
151 END DO
152
153
154
156 bignum = one / smlnum
157
158
159
160 anrm =
slange(
'M', m, n, work, m, dummy )
161 iscl = 0
162 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
163
164
165
166 CALL slascl(
'G', 0, 0, anrm, smlnum, m, n, work, m, info )
167 iscl = 1
168 ELSE IF( anrm.GT.bignum ) THEN
169
170
171
172 CALL slascl(
'G', 0, 0, anrm, bignum, m, n, work, m, info )
173 iscl = 1
174 END IF
175
176 IF( anrm.NE.zero ) THEN
177
178
179
180 CALL sgebd2( m, n, work, m, work( m*n+1 ), work( m*n+mn+1 ),
181 $ work( m*n+2*mn+1 ), work( m*n+3*mn+1 ),
182 $ work( m*n+4*mn+1 ), info )
183 CALL sbdsqr(
'Upper', mn, 0, 0, 0, work( m*n+1 ),
184 $ work( m*n+mn+1 ), dummy, mn, dummy, 1, dummy, mn,
185 $ work( m*n+2*mn+1 ), info )
186
187 IF( iscl.EQ.1 ) THEN
188 IF( anrm.GT.bignum ) THEN
189 CALL slascl(
'G', 0, 0, bignum, anrm, mn, 1,
190 $ work( m*n+1 ), mn, info )
191 END IF
192 IF( anrm.LT.smlnum ) THEN
193 CALL slascl(
'G', 0, 0, smlnum, anrm, mn, 1,
194 $ work( m*n+1 ), mn, info )
195 END IF
196 END IF
197
198 ELSE
199
200 DO i = 1, mn
201 work( m*n+i ) = zero
202 END DO
203 END IF
204
205
206
207 CALL saxpy( mn, -one, s, 1, work( m*n+1 ), 1 )
209 $ (
slamch(
'Epsilon' )*real( max( m, n ) ) )
210 IF( nrmsvl.NE.zero )
212
213 RETURN
214
215
216
subroutine xerbla(srname, info)
real function sasum(n, sx, incx)
SASUM
subroutine saxpy(n, sa, sx, incx, sy, incy)
SAXPY
subroutine sbdsqr(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, work, info)
SBDSQR
subroutine sgebd2(m, n, a, lda, d, e, tauq, taup, work, info)
SGEBD2 reduces a general matrix to bidiagonal form using an unblocked algorithm.
real function slamch(cmach)
SLAMCH
real function slange(norm, m, n, a, lda, work)
SLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
subroutine slascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
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.
real(wp) function snrm2(n, x, incx)
SNRM2
real function sqrt12(m, n, a, lda, s, work, lwork)
SQRT12