LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
 All Classes Files Functions Variables Typedefs Macros
dgesc2.f
Go to the documentation of this file.
1 *> \brief \b DGESC2 solves a system of linear equations using the LU factorization with complete pivoting computed by sgetc2.
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download DGESC2 + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgesc2.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgesc2.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgesc2.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE DGESC2( N, A, LDA, RHS, IPIV, JPIV, SCALE )
22 *
23 * .. Scalar Arguments ..
24 * INTEGER LDA, N
25 * DOUBLE PRECISION SCALE
26 * ..
27 * .. Array Arguments ..
28 * INTEGER IPIV( * ), JPIV( * )
29 * DOUBLE PRECISION A( LDA, * ), RHS( * )
30 * ..
31 *
32 *
33 *> \par Purpose:
34 * =============
35 *>
36 *> \verbatim
37 *>
38 *> DGESC2 solves a system of linear equations
39 *>
40 *> A * X = scale* RHS
41 *>
42 *> with a general N-by-N matrix A using the LU factorization with
43 *> complete pivoting computed by DGETC2.
44 *> \endverbatim
45 *
46 * Arguments:
47 * ==========
48 *
49 *> \param[in] N
50 *> \verbatim
51 *> N is INTEGER
52 *> The order of the matrix A.
53 *> \endverbatim
54 *>
55 *> \param[in] A
56 *> \verbatim
57 *> A is DOUBLE PRECISION array, dimension (LDA,N)
58 *> On entry, the LU part of the factorization of the n-by-n
59 *> matrix A computed by DGETC2: A = P * L * U * Q
60 *> \endverbatim
61 *>
62 *> \param[in] LDA
63 *> \verbatim
64 *> LDA is INTEGER
65 *> The leading dimension of the array A. LDA >= max(1, N).
66 *> \endverbatim
67 *>
68 *> \param[in,out] RHS
69 *> \verbatim
70 *> RHS is DOUBLE PRECISION array, dimension (N).
71 *> On entry, the right hand side vector b.
72 *> On exit, the solution vector X.
73 *> \endverbatim
74 *>
75 *> \param[in] IPIV
76 *> \verbatim
77 *> IPIV is INTEGER array, dimension (N).
78 *> The pivot indices; for 1 <= i <= N, row i of the
79 *> matrix has been interchanged with row IPIV(i).
80 *> \endverbatim
81 *>
82 *> \param[in] JPIV
83 *> \verbatim
84 *> JPIV is INTEGER array, dimension (N).
85 *> The pivot indices; for 1 <= j <= N, column j of the
86 *> matrix has been interchanged with column JPIV(j).
87 *> \endverbatim
88 *>
89 *> \param[out] SCALE
90 *> \verbatim
91 *> SCALE is DOUBLE PRECISION
92 *> On exit, SCALE contains the scale factor. SCALE is chosen
93 *> 0 <= SCALE <= 1 to prevent owerflow in the solution.
94 *> \endverbatim
95 *
96 * Authors:
97 * ========
98 *
99 *> \author Univ. of Tennessee
100 *> \author Univ. of California Berkeley
101 *> \author Univ. of Colorado Denver
102 *> \author NAG Ltd.
103 *
104 *> \date September 2012
105 *
106 *> \ingroup doubleGEauxiliary
107 *
108 *> \par Contributors:
109 * ==================
110 *>
111 *> Bo Kagstrom and Peter Poromaa, Department of Computing Science,
112 *> Umea University, S-901 87 Umea, Sweden.
113 *
114 * =====================================================================
115  SUBROUTINE dgesc2( N, A, LDA, RHS, IPIV, JPIV, SCALE )
116 *
117 * -- LAPACK auxiliary routine (version 3.4.2) --
118 * -- LAPACK is a software package provided by Univ. of Tennessee, --
119 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
120 * September 2012
121 *
122 * .. Scalar Arguments ..
123  INTEGER lda, n
124  DOUBLE PRECISION scale
125 * ..
126 * .. Array Arguments ..
127  INTEGER ipiv( * ), jpiv( * )
128  DOUBLE PRECISION a( lda, * ), rhs( * )
129 * ..
130 *
131 * =====================================================================
132 *
133 * .. Parameters ..
134  DOUBLE PRECISION one, two
135  parameter( one = 1.0d+0, two = 2.0d+0 )
136 * ..
137 * .. Local Scalars ..
138  INTEGER i, j
139  DOUBLE PRECISION bignum, eps, smlnum, temp
140 * ..
141 * .. External Subroutines ..
142  EXTERNAL dlaswp, dscal
143 * ..
144 * .. External Functions ..
145  INTEGER idamax
146  DOUBLE PRECISION dlamch
147  EXTERNAL idamax, dlamch
148 * ..
149 * .. Intrinsic Functions ..
150  INTRINSIC abs
151 * ..
152 * .. Executable Statements ..
153 *
154 * Set constant to control owerflow
155 *
156  eps = dlamch( 'P' )
157  smlnum = dlamch( 'S' ) / eps
158  bignum = one / smlnum
159  CALL dlabad( smlnum, bignum )
160 *
161 * Apply permutations IPIV to RHS
162 *
163  CALL dlaswp( 1, rhs, lda, 1, n-1, ipiv, 1 )
164 *
165 * Solve for L part
166 *
167  DO 20 i = 1, n - 1
168  DO 10 j = i + 1, n
169  rhs( j ) = rhs( j ) - a( j, i )*rhs( i )
170  10 CONTINUE
171  20 CONTINUE
172 *
173 * Solve for U part
174 *
175  scale = one
176 *
177 * Check for scaling
178 *
179  i = idamax( n, rhs, 1 )
180  IF( two*smlnum*abs( rhs( i ) ).GT.abs( a( n, n ) ) ) THEN
181  temp = ( one / two ) / abs( rhs( i ) )
182  CALL dscal( n, temp, rhs( 1 ), 1 )
183  scale = scale*temp
184  END IF
185 *
186  DO 40 i = n, 1, -1
187  temp = one / a( i, i )
188  rhs( i ) = rhs( i )*temp
189  DO 30 j = i + 1, n
190  rhs( i ) = rhs( i ) - rhs( j )*( a( i, j )*temp )
191  30 CONTINUE
192  40 CONTINUE
193 *
194 * Apply permutations JPIV to the solution (RHS)
195 *
196  CALL dlaswp( 1, rhs, lda, 1, n-1, jpiv, -1 )
197  RETURN
198 *
199 * End of DGESC2
200 *
201  END