LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dgesvx.f
Go to the documentation of this file.
1*> \brief <b> DGESVX computes the solution to system of linear equations A * X = B for GE matrices</b>
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download DGESVX + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgesvx.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgesvx.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgesvx.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE DGESVX( FACT, TRANS, N, NRHS, A, LDA, AF, LDAF, IPIV,
22* EQUED, R, C, B, LDB, X, LDX, RCOND, FERR, BERR,
23* WORK, IWORK, INFO )
24*
25* .. Scalar Arguments ..
26* CHARACTER EQUED, FACT, TRANS
27* INTEGER INFO, LDA, LDAF, LDB, LDX, N, NRHS
28* DOUBLE PRECISION RCOND
29* ..
30* .. Array Arguments ..
31* INTEGER IPIV( * ), IWORK( * )
32* DOUBLE PRECISION A( LDA, * ), AF( LDAF, * ), B( LDB, * ),
33* $ BERR( * ), C( * ), FERR( * ), R( * ),
34* $ WORK( * ), X( LDX, * )
35* ..
36*
37*
38*> \par Purpose:
39* =============
40*>
41*> \verbatim
42*>
43*> DGESVX uses the LU factorization to compute the solution to a real
44*> system of linear equations
45*> A * X = B,
46*> where A is an N-by-N matrix and X and B are N-by-NRHS matrices.
47*>
48*> Error bounds on the solution and a condition estimate are also
49*> provided.
50*> \endverbatim
51*
52*> \par Description:
53* =================
54*>
55*> \verbatim
56*>
57*> The following steps are performed:
58*>
59*> 1. If FACT = 'E', real scaling factors are computed to equilibrate
60*> the system:
61*> TRANS = 'N': diag(R)*A*diag(C) *inv(diag(C))*X = diag(R)*B
62*> TRANS = 'T': (diag(R)*A*diag(C))**T *inv(diag(R))*X = diag(C)*B
63*> TRANS = 'C': (diag(R)*A*diag(C))**H *inv(diag(R))*X = diag(C)*B
64*> Whether or not the system will be equilibrated depends on the
65*> scaling of the matrix A, but if equilibration is used, A is
66*> overwritten by diag(R)*A*diag(C) and B by diag(R)*B (if TRANS='N')
67*> or diag(C)*B (if TRANS = 'T' or 'C').
68*>
69*> 2. If FACT = 'N' or 'E', the LU decomposition is used to factor the
70*> matrix A (after equilibration if FACT = 'E') as
71*> A = P * L * U,
72*> where P is a permutation matrix, L is a unit lower triangular
73*> matrix, and U is upper triangular.
74*>
75*> 3. If some U(i,i)=0, so that U is exactly singular, then the routine
76*> returns with INFO = i. Otherwise, the factored form of A is used
77*> to estimate the condition number of the matrix A. If the
78*> reciprocal of the condition number is less than machine precision,
79*> INFO = N+1 is returned as a warning, but the routine still goes on
80*> to solve for X and compute error bounds as described below.
81*>
82*> 4. The system of equations is solved for X using the factored form
83*> of A.
84*>
85*> 5. Iterative refinement is applied to improve the computed solution
86*> matrix and calculate error bounds and backward error estimates
87*> for it.
88*>
89*> 6. If equilibration was used, the matrix X is premultiplied by
90*> diag(C) (if TRANS = 'N') or diag(R) (if TRANS = 'T' or 'C') so
91*> that it solves the original system before equilibration.
92*> \endverbatim
93*
94* Arguments:
95* ==========
96*
97*> \param[in] FACT
98*> \verbatim
99*> FACT is CHARACTER*1
100*> Specifies whether or not the factored form of the matrix A is
101*> supplied on entry, and if not, whether the matrix A should be
102*> equilibrated before it is factored.
103*> = 'F': On entry, AF and IPIV contain the factored form of A.
104*> If EQUED is not 'N', the matrix A has been
105*> equilibrated with scaling factors given by R and C.
106*> A, AF, and IPIV are not modified.
107*> = 'N': The matrix A will be copied to AF and factored.
108*> = 'E': The matrix A will be equilibrated if necessary, then
109*> copied to AF and factored.
110*> \endverbatim
111*>
112*> \param[in] TRANS
113*> \verbatim
114*> TRANS is CHARACTER*1
115*> Specifies the form of the system of equations:
116*> = 'N': A * X = B (No transpose)
117*> = 'T': A**T * X = B (Transpose)
118*> = 'C': A**H * X = B (Transpose)
119*> \endverbatim
120*>
121*> \param[in] N
122*> \verbatim
123*> N is INTEGER
124*> The number of linear equations, i.e., the order of the
125*> matrix A. N >= 0.
126*> \endverbatim
127*>
128*> \param[in] NRHS
129*> \verbatim
130*> NRHS is INTEGER
131*> The number of right hand sides, i.e., the number of columns
132*> of the matrices B and X. NRHS >= 0.
133*> \endverbatim
134*>
135*> \param[in,out] A
136*> \verbatim
137*> A is DOUBLE PRECISION array, dimension (LDA,N)
138*> On entry, the N-by-N matrix A. If FACT = 'F' and EQUED is
139*> not 'N', then A must have been equilibrated by the scaling
140*> factors in R and/or C. A is not modified if FACT = 'F' or
141*> 'N', or if FACT = 'E' and EQUED = 'N' on exit.
142*>
143*> On exit, if EQUED .ne. 'N', A is scaled as follows:
144*> EQUED = 'R': A := diag(R) * A
145*> EQUED = 'C': A := A * diag(C)
146*> EQUED = 'B': A := diag(R) * A * diag(C).
147*> \endverbatim
148*>
149*> \param[in] LDA
150*> \verbatim
151*> LDA is INTEGER
152*> The leading dimension of the array A. LDA >= max(1,N).
153*> \endverbatim
154*>
155*> \param[in,out] AF
156*> \verbatim
157*> AF is DOUBLE PRECISION array, dimension (LDAF,N)
158*> If FACT = 'F', then AF is an input argument and on entry
159*> contains the factors L and U from the factorization
160*> A = P*L*U as computed by DGETRF. If EQUED .ne. 'N', then
161*> AF is the factored form of the equilibrated matrix A.
162*>
163*> If FACT = 'N', then AF is an output argument and on exit
164*> returns the factors L and U from the factorization A = P*L*U
165*> of the original matrix A.
166*>
167*> If FACT = 'E', then AF is an output argument and on exit
168*> returns the factors L and U from the factorization A = P*L*U
169*> of the equilibrated matrix A (see the description of A for
170*> the form of the equilibrated matrix).
171*> \endverbatim
172*>
173*> \param[in] LDAF
174*> \verbatim
175*> LDAF is INTEGER
176*> The leading dimension of the array AF. LDAF >= max(1,N).
177*> \endverbatim
178*>
179*> \param[in,out] IPIV
180*> \verbatim
181*> IPIV is INTEGER array, dimension (N)
182*> If FACT = 'F', then IPIV is an input argument and on entry
183*> contains the pivot indices from the factorization A = P*L*U
184*> as computed by DGETRF; row i of the matrix was interchanged
185*> with row IPIV(i).
186*>
187*> If FACT = 'N', then IPIV is an output argument and on exit
188*> contains the pivot indices from the factorization A = P*L*U
189*> of the original matrix A.
190*>
191*> If FACT = 'E', then IPIV is an output argument and on exit
192*> contains the pivot indices from the factorization A = P*L*U
193*> of the equilibrated matrix A.
194*> \endverbatim
195*>
196*> \param[in,out] EQUED
197*> \verbatim
198*> EQUED is CHARACTER*1
199*> Specifies the form of equilibration that was done.
200*> = 'N': No equilibration (always true if FACT = 'N').
201*> = 'R': Row equilibration, i.e., A has been premultiplied by
202*> diag(R).
203*> = 'C': Column equilibration, i.e., A has been postmultiplied
204*> by diag(C).
205*> = 'B': Both row and column equilibration, i.e., A has been
206*> replaced by diag(R) * A * diag(C).
207*> EQUED is an input argument if FACT = 'F'; otherwise, it is an
208*> output argument.
209*> \endverbatim
210*>
211*> \param[in,out] R
212*> \verbatim
213*> R is DOUBLE PRECISION array, dimension (N)
214*> The row scale factors for A. If EQUED = 'R' or 'B', A is
215*> multiplied on the left by diag(R); if EQUED = 'N' or 'C', R
216*> is not accessed. R is an input argument if FACT = 'F';
217*> otherwise, R is an output argument. If FACT = 'F' and
218*> EQUED = 'R' or 'B', each element of R must be positive.
219*> \endverbatim
220*>
221*> \param[in,out] C
222*> \verbatim
223*> C is DOUBLE PRECISION array, dimension (N)
224*> The column scale factors for A. If EQUED = 'C' or 'B', A is
225*> multiplied on the right by diag(C); if EQUED = 'N' or 'R', C
226*> is not accessed. C is an input argument if FACT = 'F';
227*> otherwise, C is an output argument. If FACT = 'F' and
228*> EQUED = 'C' or 'B', each element of C must be positive.
229*> \endverbatim
230*>
231*> \param[in,out] B
232*> \verbatim
233*> B is DOUBLE PRECISION array, dimension (LDB,NRHS)
234*> On entry, the N-by-NRHS right hand side matrix B.
235*> On exit,
236*> if EQUED = 'N', B is not modified;
237*> if TRANS = 'N' and EQUED = 'R' or 'B', B is overwritten by
238*> diag(R)*B;
239*> if TRANS = 'T' or 'C' and EQUED = 'C' or 'B', B is
240*> overwritten by diag(C)*B.
241*> \endverbatim
242*>
243*> \param[in] LDB
244*> \verbatim
245*> LDB is INTEGER
246*> The leading dimension of the array B. LDB >= max(1,N).
247*> \endverbatim
248*>
249*> \param[out] X
250*> \verbatim
251*> X is DOUBLE PRECISION array, dimension (LDX,NRHS)
252*> If INFO = 0 or INFO = N+1, the N-by-NRHS solution matrix X
253*> to the original system of equations. Note that A and B are
254*> modified on exit if EQUED .ne. 'N', and the solution to the
255*> equilibrated system is inv(diag(C))*X if TRANS = 'N' and
256*> EQUED = 'C' or 'B', or inv(diag(R))*X if TRANS = 'T' or 'C'
257*> and EQUED = 'R' or 'B'.
258*> \endverbatim
259*>
260*> \param[in] LDX
261*> \verbatim
262*> LDX is INTEGER
263*> The leading dimension of the array X. LDX >= max(1,N).
264*> \endverbatim
265*>
266*> \param[out] RCOND
267*> \verbatim
268*> RCOND is DOUBLE PRECISION
269*> The estimate of the reciprocal condition number of the matrix
270*> A after equilibration (if done). If RCOND is less than the
271*> machine precision (in particular, if RCOND = 0), the matrix
272*> is singular to working precision. This condition is
273*> indicated by a return code of INFO > 0.
274*> \endverbatim
275*>
276*> \param[out] FERR
277*> \verbatim
278*> FERR is DOUBLE PRECISION array, dimension (NRHS)
279*> The estimated forward error bound for each solution vector
280*> X(j) (the j-th column of the solution matrix X).
281*> If XTRUE is the true solution corresponding to X(j), FERR(j)
282*> is an estimated upper bound for the magnitude of the largest
283*> element in (X(j) - XTRUE) divided by the magnitude of the
284*> largest element in X(j). The estimate is as reliable as
285*> the estimate for RCOND, and is almost always a slight
286*> overestimate of the true error.
287*> \endverbatim
288*>
289*> \param[out] BERR
290*> \verbatim
291*> BERR is DOUBLE PRECISION array, dimension (NRHS)
292*> The componentwise relative backward error of each solution
293*> vector X(j) (i.e., the smallest relative change in
294*> any element of A or B that makes X(j) an exact solution).
295*> \endverbatim
296*>
297*> \param[out] WORK
298*> \verbatim
299*> WORK is DOUBLE PRECISION array, dimension (MAX(1,4*N))
300*> On exit, WORK(1) contains the reciprocal pivot growth
301*> factor norm(A)/norm(U). The "max absolute element" norm is
302*> used. If WORK(1) is much less than 1, then the stability
303*> of the LU factorization of the (equilibrated) matrix A
304*> could be poor. This also means that the solution X, condition
305*> estimator RCOND, and forward error bound FERR could be
306*> unreliable. If factorization fails with 0<INFO<=N, then
307*> WORK(1) contains the reciprocal pivot growth factor for the
308*> leading INFO columns of A.
309*> \endverbatim
310*>
311*> \param[out] IWORK
312*> \verbatim
313*> IWORK is INTEGER array, dimension (N)
314*> \endverbatim
315*>
316*> \param[out] INFO
317*> \verbatim
318*> INFO is INTEGER
319*> = 0: successful exit
320*> < 0: if INFO = -i, the i-th argument had an illegal value
321*> > 0: if INFO = i, and i is
322*> <= N: U(i,i) is exactly zero. The factorization has
323*> been completed, but the factor U is exactly
324*> singular, so the solution and error bounds
325*> could not be computed. RCOND = 0 is returned.
326*> = N+1: U is nonsingular, but RCOND is less than machine
327*> precision, meaning that the matrix is singular
328*> to working precision. Nevertheless, the
329*> solution and error bounds are computed because
330*> there are a number of situations where the
331*> computed solution can be more accurate than the
332*> value of RCOND would suggest.
333*> \endverbatim
334*
335* Authors:
336* ========
337*
338*> \author Univ. of Tennessee
339*> \author Univ. of California Berkeley
340*> \author Univ. of Colorado Denver
341*> \author NAG Ltd.
342*
343*> \ingroup gesvx
344*
345* =====================================================================
346 SUBROUTINE dgesvx( FACT, TRANS, N, NRHS, A, LDA, AF, LDAF, IPIV,
347 $ EQUED, R, C, B, LDB, X, LDX, RCOND, FERR, BERR,
348 $ WORK, IWORK, INFO )
349*
350* -- LAPACK driver routine --
351* -- LAPACK is a software package provided by Univ. of Tennessee, --
352* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
353*
354* .. Scalar Arguments ..
355 CHARACTER EQUED, FACT, TRANS
356 INTEGER INFO, LDA, LDAF, LDB, LDX, N, NRHS
357 DOUBLE PRECISION RCOND
358* ..
359* .. Array Arguments ..
360 INTEGER IPIV( * ), IWORK( * )
361 DOUBLE PRECISION A( LDA, * ), AF( LDAF, * ), B( LDB, * ),
362 $ berr( * ), c( * ), ferr( * ), r( * ),
363 $ work( * ), x( ldx, * )
364* ..
365*
366* =====================================================================
367*
368* .. Parameters ..
369 DOUBLE PRECISION ZERO, ONE
370 PARAMETER ( ZERO = 0.0d+0, one = 1.0d+0 )
371* ..
372* .. Local Scalars ..
373 LOGICAL COLEQU, EQUIL, NOFACT, NOTRAN, ROWEQU
374 CHARACTER NORM
375 INTEGER I, INFEQU, J
376 DOUBLE PRECISION AMAX, ANORM, BIGNUM, COLCND, RCMAX, RCMIN,
377 $ rowcnd, rpvgrw, smlnum
378* ..
379* .. External Functions ..
380 LOGICAL LSAME
381 DOUBLE PRECISION DLAMCH, DLANGE, DLANTR
382 EXTERNAL lsame, dlamch, dlange, dlantr
383* ..
384* .. External Subroutines ..
385 EXTERNAL dgecon, dgeequ, dgerfs, dgetrf, dgetrs, dlacpy,
386 $ dlaqge, xerbla
387* ..
388* .. Intrinsic Functions ..
389 INTRINSIC max, min
390* ..
391* .. Executable Statements ..
392*
393 info = 0
394 nofact = lsame( fact, 'N' )
395 equil = lsame( fact, 'E' )
396 notran = lsame( trans, 'N' )
397 IF( nofact .OR. equil ) THEN
398 equed = 'N'
399 rowequ = .false.
400 colequ = .false.
401 ELSE
402 rowequ = lsame( equed, 'R' ) .OR. lsame( equed, 'B' )
403 colequ = lsame( equed, 'C' ) .OR. lsame( equed, 'B' )
404 smlnum = dlamch( 'Safe minimum' )
405 bignum = one / smlnum
406 END IF
407*
408* Test the input parameters.
409*
410 IF( .NOT.nofact .AND. .NOT.equil .AND. .NOT.lsame( fact, 'F' ) )
411 $ THEN
412 info = -1
413 ELSE IF( .NOT.notran .AND. .NOT.lsame( trans, 'T' ) .AND. .NOT.
414 $ lsame( trans, 'C' ) ) THEN
415 info = -2
416 ELSE IF( n.LT.0 ) THEN
417 info = -3
418 ELSE IF( nrhs.LT.0 ) THEN
419 info = -4
420 ELSE IF( lda.LT.max( 1, n ) ) THEN
421 info = -6
422 ELSE IF( ldaf.LT.max( 1, n ) ) THEN
423 info = -8
424 ELSE IF( lsame( fact, 'F' ) .AND. .NOT.
425 $ ( rowequ .OR. colequ .OR. lsame( equed, 'N' ) ) ) THEN
426 info = -10
427 ELSE
428 IF( rowequ ) THEN
429 rcmin = bignum
430 rcmax = zero
431 DO 10 j = 1, n
432 rcmin = min( rcmin, r( j ) )
433 rcmax = max( rcmax, r( j ) )
434 10 CONTINUE
435 IF( rcmin.LE.zero ) THEN
436 info = -11
437 ELSE IF( n.GT.0 ) THEN
438 rowcnd = max( rcmin, smlnum ) / min( rcmax, bignum )
439 ELSE
440 rowcnd = one
441 END IF
442 END IF
443 IF( colequ .AND. info.EQ.0 ) THEN
444 rcmin = bignum
445 rcmax = zero
446 DO 20 j = 1, n
447 rcmin = min( rcmin, c( j ) )
448 rcmax = max( rcmax, c( j ) )
449 20 CONTINUE
450 IF( rcmin.LE.zero ) THEN
451 info = -12
452 ELSE IF( n.GT.0 ) THEN
453 colcnd = max( rcmin, smlnum ) / min( rcmax, bignum )
454 ELSE
455 colcnd = one
456 END IF
457 END IF
458 IF( info.EQ.0 ) THEN
459 IF( ldb.LT.max( 1, n ) ) THEN
460 info = -14
461 ELSE IF( ldx.LT.max( 1, n ) ) THEN
462 info = -16
463 END IF
464 END IF
465 END IF
466*
467 IF( info.NE.0 ) THEN
468 CALL xerbla( 'DGESVX', -info )
469 RETURN
470 END IF
471*
472 IF( equil ) THEN
473*
474* Compute row and column scalings to equilibrate the matrix A.
475*
476 CALL dgeequ( n, n, a, lda, r, c, rowcnd, colcnd, amax, infequ )
477 IF( infequ.EQ.0 ) THEN
478*
479* Equilibrate the matrix.
480*
481 CALL dlaqge( n, n, a, lda, r, c, rowcnd, colcnd, amax,
482 $ equed )
483 rowequ = lsame( equed, 'R' ) .OR. lsame( equed, 'B' )
484 colequ = lsame( equed, 'C' ) .OR. lsame( equed, 'B' )
485 END IF
486 END IF
487*
488* Scale the right hand side.
489*
490 IF( notran ) THEN
491 IF( rowequ ) THEN
492 DO 40 j = 1, nrhs
493 DO 30 i = 1, n
494 b( i, j ) = r( i )*b( i, j )
495 30 CONTINUE
496 40 CONTINUE
497 END IF
498 ELSE IF( colequ ) THEN
499 DO 60 j = 1, nrhs
500 DO 50 i = 1, n
501 b( i, j ) = c( i )*b( i, j )
502 50 CONTINUE
503 60 CONTINUE
504 END IF
505*
506 IF( nofact .OR. equil ) THEN
507*
508* Compute the LU factorization of A.
509*
510 CALL dlacpy( 'Full', n, n, a, lda, af, ldaf )
511 CALL dgetrf( n, n, af, ldaf, ipiv, info )
512*
513* Return if INFO is non-zero.
514*
515 IF( info.GT.0 ) THEN
516*
517* Compute the reciprocal pivot growth factor of the
518* leading rank-deficient INFO columns of A.
519*
520 rpvgrw = dlantr( 'M', 'U', 'N', info, info, af, ldaf,
521 $ work )
522 IF( rpvgrw.EQ.zero ) THEN
523 rpvgrw = one
524 ELSE
525 rpvgrw = dlange( 'M', n, info, a, lda, work ) / rpvgrw
526 END IF
527 work( 1 ) = rpvgrw
528 rcond = zero
529 RETURN
530 END IF
531 END IF
532*
533* Compute the norm of the matrix A and the
534* reciprocal pivot growth factor RPVGRW.
535*
536 IF( notran ) THEN
537 norm = '1'
538 ELSE
539 norm = 'I'
540 END IF
541 anorm = dlange( norm, n, n, a, lda, work )
542 rpvgrw = dlantr( 'M', 'U', 'N', n, n, af, ldaf, work )
543 IF( rpvgrw.EQ.zero ) THEN
544 rpvgrw = one
545 ELSE
546 rpvgrw = dlange( 'M', n, n, a, lda, work ) / rpvgrw
547 END IF
548*
549* Compute the reciprocal of the condition number of A.
550*
551 CALL dgecon( norm, n, af, ldaf, anorm, rcond, work, iwork, info )
552*
553* Compute the solution matrix X.
554*
555 CALL dlacpy( 'Full', n, nrhs, b, ldb, x, ldx )
556 CALL dgetrs( trans, n, nrhs, af, ldaf, ipiv, x, ldx, info )
557*
558* Use iterative refinement to improve the computed solution and
559* compute error bounds and backward error estimates for it.
560*
561 CALL dgerfs( trans, n, nrhs, a, lda, af, ldaf, ipiv, b, ldb, x,
562 $ ldx, ferr, berr, work, iwork, info )
563*
564* Transform the solution matrix X to a solution of the original
565* system.
566*
567 IF( notran ) THEN
568 IF( colequ ) THEN
569 DO 80 j = 1, nrhs
570 DO 70 i = 1, n
571 x( i, j ) = c( i )*x( i, j )
572 70 CONTINUE
573 80 CONTINUE
574 DO 90 j = 1, nrhs
575 ferr( j ) = ferr( j ) / colcnd
576 90 CONTINUE
577 END IF
578 ELSE IF( rowequ ) THEN
579 DO 110 j = 1, nrhs
580 DO 100 i = 1, n
581 x( i, j ) = r( i )*x( i, j )
582 100 CONTINUE
583 110 CONTINUE
584 DO 120 j = 1, nrhs
585 ferr( j ) = ferr( j ) / rowcnd
586 120 CONTINUE
587 END IF
588*
589 work( 1 ) = rpvgrw
590*
591* Set INFO = N+1 if the matrix is singular to working precision.
592*
593 IF( rcond.LT.dlamch( 'Epsilon' ) )
594 $ info = n + 1
595 RETURN
596*
597* End of DGESVX
598*
599 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dgecon(norm, n, a, lda, anorm, rcond, work, iwork, info)
DGECON
Definition dgecon.f:132
subroutine dgeequ(m, n, a, lda, r, c, rowcnd, colcnd, amax, info)
DGEEQU
Definition dgeequ.f:139
subroutine dgerfs(trans, n, nrhs, a, lda, af, ldaf, ipiv, b, ldb, x, ldx, ferr, berr, work, iwork, info)
DGERFS
Definition dgerfs.f:185
subroutine dgesvx(fact, trans, n, nrhs, a, lda, af, ldaf, ipiv, equed, r, c, b, ldb, x, ldx, rcond, ferr, berr, work, iwork, info)
DGESVX computes the solution to system of linear equations A * X = B for GE matrices
Definition dgesvx.f:349
subroutine dgetrf(m, n, a, lda, ipiv, info)
DGETRF
Definition dgetrf.f:108
subroutine dgetrs(trans, n, nrhs, a, lda, ipiv, b, ldb, info)
DGETRS
Definition dgetrs.f:121
subroutine dlacpy(uplo, m, n, a, lda, b, ldb)
DLACPY copies all or part of one two-dimensional array to another.
Definition dlacpy.f:103
subroutine dlaqge(m, n, a, lda, r, c, rowcnd, colcnd, amax, equed)
DLAQGE scales a general rectangular matrix, using row and column scaling factors computed by sgeequ.
Definition dlaqge.f:142