001:       SUBROUTINE CGESC2( N, A, LDA, RHS, IPIV, JPIV, SCALE )
002: *
003: *  -- LAPACK auxiliary routine (version 3.2) --
004: *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
005: *     November 2006
006: *
007: *     .. Scalar Arguments ..
008:       INTEGER            LDA, N
009:       REAL               SCALE
010: *     ..
011: *     .. Array Arguments ..
012:       INTEGER            IPIV( * ), JPIV( * )
013:       COMPLEX            A( LDA, * ), RHS( * )
014: *     ..
015: *
016: *  Purpose
017: *  =======
018: *
019: *  CGESC2 solves a system of linear equations
020: *
021: *            A * X = scale* RHS
022: *
023: *  with a general N-by-N matrix A using the LU factorization with
024: *  complete pivoting computed by CGETC2.
025: *
026: *
027: *  Arguments
028: *  =========
029: *
030: *  N       (input) INTEGER
031: *          The number of columns of the matrix A.
032: *
033: *  A       (input) COMPLEX array, dimension (LDA, N)
034: *          On entry, the  LU part of the factorization of the n-by-n
035: *          matrix A computed by CGETC2:  A = P * L * U * Q
036: *
037: *  LDA     (input) INTEGER
038: *          The leading dimension of the array A.  LDA >= max(1, N).
039: *
040: *  RHS     (input/output) COMPLEX array, dimension N.
041: *          On entry, the right hand side vector b.
042: *          On exit, the solution vector X.
043: *
044: *  IPIV    (input) INTEGER array, dimension (N).
045: *          The pivot indices; for 1 <= i <= N, row i of the
046: *          matrix has been interchanged with row IPIV(i).
047: *
048: *  JPIV    (input) INTEGER array, dimension (N).
049: *          The pivot indices; for 1 <= j <= N, column j of the
050: *          matrix has been interchanged with column JPIV(j).
051: *
052: *  SCALE    (output) REAL
053: *           On exit, SCALE contains the scale factor. SCALE is chosen
054: *           0 <= SCALE <= 1 to prevent owerflow in the solution.
055: *
056: *  Further Details
057: *  ===============
058: *
059: *  Based on contributions by
060: *     Bo Kagstrom and Peter Poromaa, Department of Computing Science,
061: *     Umea University, S-901 87 Umea, Sweden.
062: *
063: *  =====================================================================
064: *
065: *     .. Parameters ..
066:       REAL               ZERO, ONE, TWO
067:       PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0, TWO = 2.0E+0 )
068: *     ..
069: *     .. Local Scalars ..
070:       INTEGER            I, J
071:       REAL               BIGNUM, EPS, SMLNUM
072:       COMPLEX            TEMP
073: *     ..
074: *     .. External Subroutines ..
075:       EXTERNAL           CLASWP, CSCAL, SLABAD
076: *     ..
077: *     .. External Functions ..
078:       INTEGER            ICAMAX
079:       REAL               SLAMCH
080:       EXTERNAL           ICAMAX, SLAMCH
081: *     ..
082: *     .. Intrinsic Functions ..
083:       INTRINSIC          ABS, CMPLX, REAL
084: *     ..
085: *     .. Executable Statements ..
086: *
087: *     Set constant to control overflow
088: *
089:       EPS = SLAMCH( 'P' )
090:       SMLNUM = SLAMCH( 'S' ) / EPS
091:       BIGNUM = ONE / SMLNUM
092:       CALL SLABAD( SMLNUM, BIGNUM )
093: *
094: *     Apply permutations IPIV to RHS
095: *
096:       CALL CLASWP( 1, RHS, LDA, 1, N-1, IPIV, 1 )
097: *
098: *     Solve for L part
099: *
100:       DO 20 I = 1, N - 1
101:          DO 10 J = I + 1, N
102:             RHS( J ) = RHS( J ) - A( J, I )*RHS( I )
103:    10    CONTINUE
104:    20 CONTINUE
105: *
106: *     Solve for U part
107: *
108:       SCALE = ONE
109: *
110: *     Check for scaling
111: *
112:       I = ICAMAX( N, RHS, 1 )
113:       IF( TWO*SMLNUM*ABS( RHS( I ) ).GT.ABS( A( N, N ) ) ) THEN
114:          TEMP = CMPLX( ONE / TWO, ZERO ) / ABS( RHS( I ) )
115:          CALL CSCAL( N, TEMP, RHS( 1 ), 1 )
116:          SCALE = SCALE*REAL( TEMP )
117:       END IF
118:       DO 40 I = N, 1, -1
119:          TEMP = CMPLX( ONE, ZERO ) / A( I, I )
120:          RHS( I ) = RHS( I )*TEMP
121:          DO 30 J = I + 1, N
122:             RHS( I ) = RHS( I ) - RHS( J )*( A( I, J )*TEMP )
123:    30    CONTINUE
124:    40 CONTINUE
125: *
126: *     Apply permutations JPIV to the solution (RHS)
127: *
128:       CALL CLASWP( 1, RHS, LDA, 1, N-1, JPIV, -1 )
129:       RETURN
130: *
131: *     End of CGESC2
132: *
133:       END
134: