3 * =========== DOCUMENTATION ===========
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
9 *> Download CPFTRS + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cpftrs.f">
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cpftrs.f">
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cpftrs.f">
21 * SUBROUTINE CPFTRS( TRANSR, UPLO, N, NRHS, A, B, LDB, INFO )
23 * .. Scalar Arguments ..
24 * CHARACTER TRANSR, UPLO
25 * INTEGER INFO, LDB, N, NRHS
27 * .. Array Arguments ..
28 * COMPLEX A( 0: * ), B( LDB, * )
37 *> CPFTRS solves a system of linear equations A*X = B with a Hermitian
38 *> positive definite matrix A using the Cholesky factorization
39 *> A = U**H*U or A = L*L**H computed by CPFTRF.
47 *> TRANSR is CHARACTER*1
48 *> = 'N': The Normal TRANSR of RFP A is stored;
49 *> = 'C': The Conjugate-transpose TRANSR of RFP A is stored.
54 *> UPLO is CHARACTER*1
55 *> = 'U': Upper triangle of RFP A is stored;
56 *> = 'L': Lower triangle of RFP A is stored.
62 *> The order of the matrix A. N >= 0.
68 *> The number of right hand sides, i.e., the number of columns
69 *> of the matrix B. NRHS >= 0.
74 *> A is COMPLEX array, dimension ( N*(N+1)/2 );
75 *> The triangular factor U or L from the Cholesky factorization
76 *> of RFP A = U**H*U or RFP A = L*L**H, as computed by CPFTRF.
77 *> See note below for more details about RFP A.
82 *> B is COMPLEX array, dimension (LDB,NRHS)
83 *> On entry, the right hand side matrix B.
84 *> On exit, the solution matrix X.
90 *> The leading dimension of the array B. LDB >= max(1,N).
96 *> = 0: successful exit
97 *> < 0: if INFO = -i, the i-th argument had an illegal value
103 *> \author Univ. of Tennessee
104 *> \author Univ. of California Berkeley
105 *> \author Univ. of Colorado Denver
108 *> \date November 2011
110 *> \ingroup complexOTHERcomputational
112 *> \par Further Details:
113 * =====================
117 *> We first consider Standard Packed Format when N is even.
118 *> We give an example where N = 6.
120 *> AP is Upper AP is Lower
122 *> 00 01 02 03 04 05 00
123 *> 11 12 13 14 15 10 11
124 *> 22 23 24 25 20 21 22
125 *> 33 34 35 30 31 32 33
126 *> 44 45 40 41 42 43 44
127 *> 55 50 51 52 53 54 55
130 *> Let TRANSR = 'N'. RFP holds AP as follows:
131 *> For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
132 *> three columns of AP upper. The lower triangle A(4:6,0:2) consists of
133 *> conjugate-transpose of the first three columns of AP upper.
134 *> For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
135 *> three columns of AP lower. The upper triangle A(0:2,0:2) consists of
136 *> conjugate-transpose of the last three columns of AP lower.
137 *> To denote conjugate we place -- above the element. This covers the
138 *> case N even and TRANSR = 'N'.
157 *> Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
158 *> transpose of RFP A above. One therefore gets:
163 *> -- -- -- -- -- -- -- -- -- --
164 *> 03 13 23 33 00 01 02 33 00 10 20 30 40 50
165 *> -- -- -- -- -- -- -- -- -- --
166 *> 04 14 24 34 44 11 12 43 44 11 21 31 41 51
167 *> -- -- -- -- -- -- -- -- -- --
168 *> 05 15 25 35 45 55 22 53 54 55 22 32 42 52
171 *> We next consider Standard Packed Format when N is odd.
172 *> We give an example where N = 5.
174 *> AP is Upper AP is Lower
183 *> Let TRANSR = 'N'. RFP holds AP as follows:
184 *> For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
185 *> three columns of AP upper. The lower triangle A(3:4,0:1) consists of
186 *> conjugate-transpose of the first two columns of AP upper.
187 *> For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
188 *> three columns of AP lower. The upper triangle A(0:1,1:2) consists of
189 *> conjugate-transpose of the last two columns of AP lower.
190 *> To denote conjugate we place -- above the element. This covers the
191 *> case N odd and TRANSR = 'N'.
206 *> Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
207 *> transpose of RFP A above. One therefore gets:
212 *> -- -- -- -- -- -- -- -- --
213 *> 02 12 22 00 01 00 10 20 30 40 50
214 *> -- -- -- -- -- -- -- -- --
215 *> 03 13 23 33 11 33 11 21 31 41 51
216 *> -- -- -- -- -- -- -- -- --
217 *> 04 14 24 34 44 43 44 22 32 42 52
220 * =====================================================================
221 SUBROUTINE CPFTRS( TRANSR, UPLO, N, NRHS, A, B, LDB, INFO )
223 * -- LAPACK computational routine (version 3.4.0) --
224 * -- LAPACK is a software package provided by Univ. of Tennessee, --
225 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
228 * .. Scalar Arguments ..
229 CHARACTER TRANSR, UPLO
230 INTEGER INFO, LDB, N, NRHS
232 * .. Array Arguments ..
233 COMPLEX A( 0: * ), B( LDB, * )
236 * =====================================================================
240 PARAMETER ( CONE = ( 1.0E+0, 0.0E+0 ) )
242 * .. Local Scalars ..
243 LOGICAL LOWER, NORMALTRANSR
245 * .. External Functions ..
249 * .. External Subroutines ..
250 EXTERNAL XERBLA, CTFSM
252 * .. Intrinsic Functions ..
255 * .. Executable Statements ..
257 * Test the input parameters.
260 NORMALTRANSR = LSAME( TRANSR, 'N' )
261 LOWER = LSAME( UPLO, 'L' )
262 IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'C' ) ) THEN
264 ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
266 ELSE IF( N.LT.0 ) THEN
268 ELSE IF( NRHS.LT.0 ) THEN
270 ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
274 CALL XERBLA( 'CPFTRS', -INFO )
278 * Quick return if possible
280 IF( N.EQ.0 .OR. NRHS.EQ.0 )
283 * start execution: there are two triangular solves
286 CALL CTFSM( TRANSR, 'L', UPLO, 'N', 'N', N, NRHS, CONE, A, B,
288 CALL CTFSM( TRANSR, 'L', UPLO, 'C', 'N', N, NRHS, CONE, A, B,
291 CALL CTFSM( TRANSR, 'L', UPLO, 'C', 'N', N, NRHS, CONE, A, B,
293 CALL CTFSM( TRANSR, 'L', UPLO, 'N', 'N', N, NRHS, CONE, A, B,