195 SUBROUTINE dsgesv( N, NRHS, A, LDA, IPIV, B, LDB, X, LDX, WORK,
196 $ swork, iter, info )
204 INTEGER INFO, ITER, LDA, LDB, LDX, N, NRHS
209 DOUBLE PRECISION A( lda, * ), B( ldb, * ), WORK( n, * ),
217 parameter( doitref = .true. )
220 parameter( itermax = 30 )
222 DOUBLE PRECISION BWDMAX
223 parameter( bwdmax = 1.0e+00 )
225 DOUBLE PRECISION NEGONE, ONE
226 parameter( negone = -1.0d+0, one = 1.0d+0 )
229 INTEGER I, IITER, PTSA, PTSX
230 DOUBLE PRECISION ANRM, CTE, EPS, RNRM, XNRM
238 DOUBLE PRECISION DLAMCH, DLANGE
239 EXTERNAL idamax, dlamch, dlange
242 INTRINSIC abs, dble, max, sqrt
253 ELSE IF( nrhs.LT.0 )
THEN
255 ELSE IF( lda.LT.max( 1, n ) )
THEN
257 ELSE IF( ldb.LT.max( 1, n ) )
THEN
259 ELSE IF( ldx.LT.max( 1, n ) )
THEN
263 CALL
xerbla(
'DSGESV', -info )
275 IF( .NOT.doitref )
THEN
282 anrm = dlange(
'I', n, n, a, lda, work )
283 eps = dlamch(
'Epsilon' )
284 cte = anrm*eps*sqrt( dble( n ) )*bwdmax
294 CALL
dlag2s( n, nrhs, b, ldb, swork( ptsx ), n, info )
304 CALL
dlag2s( n, n, a, lda, swork( ptsa ), n, info )
313 CALL
sgetrf( n, n, swork( ptsa ), n, ipiv, info )
322 CALL
sgetrs(
'No transpose', n, nrhs, swork( ptsa ), n, ipiv,
323 $ swork( ptsx ), n, info )
327 CALL
slag2d( n, nrhs, swork( ptsx ), n, x, ldx, info )
331 CALL
dlacpy(
'All', n, nrhs, b, ldb, work, n )
333 CALL
dgemm(
'No Transpose',
'No Transpose', n, nrhs, n, negone, a,
334 $ lda, x, ldx, one, work, n )
340 xnrm = abs( x( idamax( n, x( 1, i ), 1 ), i ) )
341 rnrm = abs( work( idamax( n, work( 1, i ), 1 ), i ) )
342 IF( rnrm.GT.xnrm*cte )
354 DO 30 iiter = 1, itermax
359 CALL
dlag2s( n, nrhs, work, n, swork( ptsx ), n, info )
368 CALL
sgetrs(
'No transpose', n, nrhs, swork( ptsa ), n, ipiv,
369 $ swork( ptsx ), n, info )
374 CALL
slag2d( n, nrhs, swork( ptsx ), n, work, n, info )
377 CALL
daxpy( n, one, work( 1, i ), 1, x( 1, i ), 1 )
382 CALL
dlacpy(
'All', n, nrhs, b, ldb, work, n )
384 CALL
dgemm(
'No Transpose',
'No Transpose', n, nrhs, n, negone,
385 $ a, lda, x, ldx, one, work, n )
391 xnrm = abs( x( idamax( n, x( 1, i ), 1 ), i ) )
392 rnrm = abs( work( idamax( n, work( 1, i ), 1 ), i ) )
393 IF( rnrm.GT.xnrm*cte )
420 CALL
dgetrf( n, n, a, lda, ipiv, info )
425 CALL
dlacpy(
'All', n, nrhs, b, ldb, x, ldx )
426 CALL
dgetrs(
'No transpose', n, nrhs, a, lda, ipiv, x, ldx,
subroutine sgetrs(TRANS, N, NRHS, A, LDA, IPIV, B, LDB, INFO)
SGETRS
subroutine dlag2s(M, N, A, LDA, SA, LDSA, INFO)
DLAG2S converts a double precision matrix to a single precision matrix.
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
subroutine sgetrf(M, N, A, LDA, IPIV, INFO)
SGETRF
subroutine dgetrs(TRANS, N, NRHS, A, LDA, IPIV, B, LDB, INFO)
DGETRS
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine slag2d(M, N, SA, LDSA, A, LDA, INFO)
SLAG2D converts a single precision matrix to a double precision matrix.
subroutine daxpy(N, DA, DX, INCX, DY, INCY)
DAXPY
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
subroutine dsgesv(N, NRHS, A, LDA, IPIV, B, LDB, X, LDX, WORK, SWORK, ITER, INFO)
DSGESV computes the solution to system of linear equations A * X = B for GE matrices (mixed precisio...
subroutine dgetrf(M, N, A, LDA, IPIV, INFO)
DGETRF