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