Lots of trailing whitespaces in the files of Syd. Cleaning this. No big deal.
[platform/upstream/lapack.git] / SRC / cpocon.f
1 *> \brief \b CPOCON
2 *
3 *  =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 *            http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download CPOCON + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cpocon.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cpocon.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cpocon.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 *  Definition:
19 *  ===========
20 *
21 *       SUBROUTINE CPOCON( UPLO, N, A, LDA, ANORM, RCOND, WORK, RWORK,
22 *                          INFO )
23 *
24 *       .. Scalar Arguments ..
25 *       CHARACTER          UPLO
26 *       INTEGER            INFO, LDA, N
27 *       REAL               ANORM, RCOND
28 *       ..
29 *       .. Array Arguments ..
30 *       REAL               RWORK( * )
31 *       COMPLEX            A( LDA, * ), WORK( * )
32 *       ..
33 *
34 *
35 *> \par Purpose:
36 *  =============
37 *>
38 *> \verbatim
39 *>
40 *> CPOCON estimates the reciprocal of the condition number (in the
41 *> 1-norm) of a complex Hermitian positive definite matrix using the
42 *> Cholesky factorization A = U**H*U or A = L*L**H computed by CPOTRF.
43 *>
44 *> An estimate is obtained for norm(inv(A)), and the reciprocal of the
45 *> condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))).
46 *> \endverbatim
47 *
48 *  Arguments:
49 *  ==========
50 *
51 *> \param[in] UPLO
52 *> \verbatim
53 *>          UPLO is CHARACTER*1
54 *>          = 'U':  Upper triangle of A is stored;
55 *>          = 'L':  Lower triangle of A is stored.
56 *> \endverbatim
57 *>
58 *> \param[in] N
59 *> \verbatim
60 *>          N is INTEGER
61 *>          The order of the matrix A.  N >= 0.
62 *> \endverbatim
63 *>
64 *> \param[in] A
65 *> \verbatim
66 *>          A is COMPLEX array, dimension (LDA,N)
67 *>          The triangular factor U or L from the Cholesky factorization
68 *>          A = U**H*U or A = L*L**H, as computed by CPOTRF.
69 *> \endverbatim
70 *>
71 *> \param[in] LDA
72 *> \verbatim
73 *>          LDA is INTEGER
74 *>          The leading dimension of the array A.  LDA >= max(1,N).
75 *> \endverbatim
76 *>
77 *> \param[in] ANORM
78 *> \verbatim
79 *>          ANORM is REAL
80 *>          The 1-norm (or infinity-norm) of the Hermitian matrix A.
81 *> \endverbatim
82 *>
83 *> \param[out] RCOND
84 *> \verbatim
85 *>          RCOND is REAL
86 *>          The reciprocal of the condition number of the matrix A,
87 *>          computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an
88 *>          estimate of the 1-norm of inv(A) computed in this routine.
89 *> \endverbatim
90 *>
91 *> \param[out] WORK
92 *> \verbatim
93 *>          WORK is COMPLEX array, dimension (2*N)
94 *> \endverbatim
95 *>
96 *> \param[out] RWORK
97 *> \verbatim
98 *>          RWORK is REAL array, dimension (N)
99 *> \endverbatim
100 *>
101 *> \param[out] INFO
102 *> \verbatim
103 *>          INFO is INTEGER
104 *>          = 0:  successful exit
105 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
106 *> \endverbatim
107 *
108 *  Authors:
109 *  ========
110 *
111 *> \author Univ. of Tennessee
112 *> \author Univ. of California Berkeley
113 *> \author Univ. of Colorado Denver
114 *> \author NAG Ltd.
115 *
116 *> \date November 2011
117 *
118 *> \ingroup complexPOcomputational
119 *
120 *  =====================================================================
121       SUBROUTINE CPOCON( UPLO, N, A, LDA, ANORM, RCOND, WORK, RWORK,
122      $                   INFO )
123 *
124 *  -- LAPACK computational routine (version 3.4.0) --
125 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
126 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
127 *     November 2011
128 *
129 *     .. Scalar Arguments ..
130       CHARACTER          UPLO
131       INTEGER            INFO, LDA, N
132       REAL               ANORM, RCOND
133 *     ..
134 *     .. Array Arguments ..
135       REAL               RWORK( * )
136       COMPLEX            A( LDA, * ), WORK( * )
137 *     ..
138 *
139 *  =====================================================================
140 *
141 *     .. Parameters ..
142       REAL               ONE, ZERO
143       PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
144 *     ..
145 *     .. Local Scalars ..
146       LOGICAL            UPPER
147       CHARACTER          NORMIN
148       INTEGER            IX, KASE
149       REAL               AINVNM, SCALE, SCALEL, SCALEU, SMLNUM
150       COMPLEX            ZDUM
151 *     ..
152 *     .. Local Arrays ..
153       INTEGER            ISAVE( 3 )
154 *     ..
155 *     .. External Functions ..
156       LOGICAL            LSAME
157       INTEGER            ICAMAX
158       REAL               SLAMCH
159       EXTERNAL           LSAME, ICAMAX, SLAMCH
160 *     ..
161 *     .. External Subroutines ..
162       EXTERNAL           CLACN2, CLATRS, CSRSCL, XERBLA
163 *     ..
164 *     .. Intrinsic Functions ..
165       INTRINSIC          ABS, AIMAG, MAX, REAL
166 *     ..
167 *     .. Statement Functions ..
168       REAL               CABS1
169 *     ..
170 *     .. Statement Function definitions ..
171       CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
172 *     ..
173 *     .. Executable Statements ..
174 *
175 *     Test the input parameters.
176 *
177       INFO = 0
178       UPPER = LSAME( UPLO, 'U' )
179       IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
180          INFO = -1
181       ELSE IF( N.LT.0 ) THEN
182          INFO = -2
183       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
184          INFO = -4
185       ELSE IF( ANORM.LT.ZERO ) THEN
186          INFO = -5
187       END IF
188       IF( INFO.NE.0 ) THEN
189          CALL XERBLA( 'CPOCON', -INFO )
190          RETURN
191       END IF
192 *
193 *     Quick return if possible
194 *
195       RCOND = ZERO
196       IF( N.EQ.0 ) THEN
197          RCOND = ONE
198          RETURN
199       ELSE IF( ANORM.EQ.ZERO ) THEN
200          RETURN
201       END IF
202 *
203       SMLNUM = SLAMCH( 'Safe minimum' )
204 *
205 *     Estimate the 1-norm of inv(A).
206 *
207       KASE = 0
208       NORMIN = 'N'
209    10 CONTINUE
210       CALL CLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
211       IF( KASE.NE.0 ) THEN
212          IF( UPPER ) THEN
213 *
214 *           Multiply by inv(U**H).
215 *
216             CALL CLATRS( 'Upper', 'Conjugate transpose', 'Non-unit',
217      $                   NORMIN, N, A, LDA, WORK, SCALEL, RWORK, INFO )
218             NORMIN = 'Y'
219 *
220 *           Multiply by inv(U).
221 *
222             CALL CLATRS( 'Upper', 'No transpose', 'Non-unit', NORMIN, N,
223      $                   A, LDA, WORK, SCALEU, RWORK, INFO )
224          ELSE
225 *
226 *           Multiply by inv(L).
227 *
228             CALL CLATRS( 'Lower', 'No transpose', 'Non-unit', NORMIN, N,
229      $                   A, LDA, WORK, SCALEL, RWORK, INFO )
230             NORMIN = 'Y'
231 *
232 *           Multiply by inv(L**H).
233 *
234             CALL CLATRS( 'Lower', 'Conjugate transpose', 'Non-unit',
235      $                   NORMIN, N, A, LDA, WORK, SCALEU, RWORK, INFO )
236          END IF
237 *
238 *        Multiply by 1/SCALE if doing so will not cause overflow.
239 *
240          SCALE = SCALEL*SCALEU
241          IF( SCALE.NE.ONE ) THEN
242             IX = ICAMAX( N, WORK, 1 )
243             IF( SCALE.LT.CABS1( WORK( IX ) )*SMLNUM .OR. SCALE.EQ.ZERO )
244      $         GO TO 20
245             CALL CSRSCL( N, SCALE, WORK, 1 )
246          END IF
247          GO TO 10
248       END IF
249 *
250 *     Compute the estimate of the reciprocal condition number.
251 *
252       IF( AINVNM.NE.ZERO )
253      $   RCOND = ( ONE / AINVNM ) / ANORM
254 *
255    20 CONTINUE
256       RETURN
257 *
258 *     End of CPOCON
259 *
260       END