1*> \brief \b CLA_GBRCOND_C computes the infinity norm condition number of op(A)*inv(diag(c)) for general banded matrices.
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download CLA_GBRCOND_C + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cla_gbrcond_c.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cla_gbrcond_c.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cla_gbrcond_c.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18*  Definition:
19*  ===========
20*
21*       REAL FUNCTION CLA_GBRCOND_C( TRANS, N, KL, KU, AB, LDAB, AFB,
22*                                    LDAFB, IPIV, C, CAPPLY, INFO, WORK,
23*                                    RWORK )
24*
25*       .. Scalar Arguments ..
26*       CHARACTER          TRANS
27*       LOGICAL            CAPPLY
28*       INTEGER            N, KL, KU, KD, KE, LDAB, LDAFB, INFO
29*       ..
30*       .. Array Arguments ..
31*       INTEGER            IPIV( * )
32*       COMPLEX            AB( LDAB, * ), AFB( LDAFB, * ), WORK( * )
33*       REAL               C( * ), RWORK( * )
34*       ..
35*
36*
37*> \par Purpose:
38*  =============
39*>
40*> \verbatim
41*>
42*>    CLA_GBRCOND_C Computes the infinity norm condition number of
43*>    op(A) * inv(diag(C)) where C is a REAL vector.
44*> \endverbatim
45*
46*  Arguments:
47*  ==========
48*
49*> \param[in] TRANS
50*> \verbatim
51*>          TRANS is CHARACTER*1
52*>     Specifies the form of the system of equations:
53*>       = 'N':  A * X = B     (No transpose)
54*>       = 'T':  A**T * X = B  (Transpose)
55*>       = 'C':  A**H * X = B  (Conjugate Transpose = Transpose)
56*> \endverbatim
57*>
58*> \param[in] N
59*> \verbatim
60*>          N is INTEGER
61*>     The number of linear equations, i.e., the order of the
62*>     matrix A.  N >= 0.
63*> \endverbatim
64*>
65*> \param[in] KL
66*> \verbatim
67*>          KL is INTEGER
68*>     The number of subdiagonals within the band of A.  KL >= 0.
69*> \endverbatim
70*>
71*> \param[in] KU
72*> \verbatim
73*>          KU is INTEGER
74*>     The number of superdiagonals within the band of A.  KU >= 0.
75*> \endverbatim
76*>
77*> \param[in] AB
78*> \verbatim
79*>          AB is COMPLEX array, dimension (LDAB,N)
80*>     On entry, the matrix A in band storage, in rows 1 to KL+KU+1.
81*>     The j-th column of A is stored in the j-th column of the
82*>     array AB as follows:
83*>     AB(KU+1+i-j,j) = A(i,j) for max(1,j-KU)<=i<=min(N,j+kl)
84*> \endverbatim
85*>
86*> \param[in] LDAB
87*> \verbatim
88*>          LDAB is INTEGER
89*>     The leading dimension of the array AB.  LDAB >= KL+KU+1.
90*> \endverbatim
91*>
92*> \param[in] AFB
93*> \verbatim
94*>          AFB is COMPLEX array, dimension (LDAFB,N)
95*>     Details of the LU factorization of the band matrix A, as
96*>     computed by CGBTRF.  U is stored as an upper triangular
97*>     band matrix with KL+KU superdiagonals in rows 1 to KL+KU+1,
98*>     and the multipliers used during the factorization are stored
99*>     in rows KL+KU+2 to 2*KL+KU+1.
100*> \endverbatim
101*>
102*> \param[in] LDAFB
103*> \verbatim
104*>          LDAFB is INTEGER
105*>     The leading dimension of the array AFB.  LDAFB >= 2*KL+KU+1.
106*> \endverbatim
107*>
108*> \param[in] IPIV
109*> \verbatim
110*>          IPIV is INTEGER array, dimension (N)
111*>     The pivot indices from the factorization A = P*L*U
112*>     as computed by CGBTRF; row i of the matrix was interchanged
113*>     with row IPIV(i).
114*> \endverbatim
115*>
116*> \param[in] C
117*> \verbatim
118*>          C is REAL array, dimension (N)
119*>     The vector C in the formula op(A) * inv(diag(C)).
120*> \endverbatim
121*>
122*> \param[in] CAPPLY
123*> \verbatim
124*>          CAPPLY is LOGICAL
125*>     If .TRUE. then access the vector C in the formula above.
126*> \endverbatim
127*>
128*> \param[out] INFO
129*> \verbatim
130*>          INFO is INTEGER
131*>       = 0:  Successful exit.
132*>     i > 0:  The ith argument is invalid.
133*> \endverbatim
134*>
135*> \param[in] WORK
136*> \verbatim
137*>          WORK is COMPLEX array, dimension (2*N).
138*>     Workspace.
139*> \endverbatim
140*>
141*> \param[in] RWORK
142*> \verbatim
143*>          RWORK is REAL array, dimension (N).
144*>     Workspace.
145*> \endverbatim
146*
147*  Authors:
148*  ========
149*
150*> \author Univ. of Tennessee
151*> \author Univ. of California Berkeley
152*> \author Univ. of Colorado Denver
153*> \author NAG Ltd.
154*
155*> \date September 2012
156*
157*> \ingroup complexGBcomputational
158*
159*  =====================================================================
160      REAL FUNCTION CLA_GBRCOND_C( TRANS, N, KL, KU, AB, LDAB, AFB,
161     $                             LDAFB, IPIV, C, CAPPLY, INFO, WORK,
162     $                             RWORK )
163*
164*  -- LAPACK computational routine (version 3.4.2) --
165*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
166*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
167*     September 2012
168*
169*     .. Scalar Arguments ..
170      CHARACTER          TRANS
171      LOGICAL            CAPPLY
172      INTEGER            N, KL, KU, KD, KE, LDAB, LDAFB, INFO
173*     ..
174*     .. Array Arguments ..
175      INTEGER            IPIV( * )
176      COMPLEX            AB( LDAB, * ), AFB( LDAFB, * ), WORK( * )
177      REAL               C( * ), RWORK( * )
178*     ..
179*
180*  =====================================================================
181*
182*     .. Local Scalars ..
183      LOGICAL            NOTRANS
184      INTEGER            KASE, I, J
185      REAL               AINVNM, ANORM, TMP
186      COMPLEX            ZDUM
187*     ..
188*     .. Local Arrays ..
189      INTEGER            ISAVE( 3 )
190*     ..
191*     .. External Functions ..
192      LOGICAL            LSAME
193      EXTERNAL           LSAME
194*     ..
195*     .. External Subroutines ..
196      EXTERNAL           CLACN2, CGBTRS, XERBLA
197*     ..
198*     .. Intrinsic Functions ..
199      INTRINSIC          ABS, MAX
200*     ..
201*     .. Statement Functions ..
202      REAL               CABS1
203*     ..
204*     .. Statement Function Definitions ..
205      CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
206*     ..
207*     .. Executable Statements ..
208      CLA_GBRCOND_C = 0.0E+0
209*
210      INFO = 0
211      NOTRANS = LSAME( TRANS, 'N' )
212      IF ( .NOT. NOTRANS .AND. .NOT. LSAME( TRANS, 'T' ) .AND. .NOT.
213     $     LSAME( TRANS, 'C' ) ) THEN
214         INFO = -1
215      ELSE IF( N.LT.0 ) THEN
216         INFO = -2
217      ELSE IF( KL.LT.0 .OR. KL.GT.N-1 ) THEN
218         INFO = -3
219      ELSE IF( KU.LT.0 .OR. KU.GT.N-1 ) THEN
220         INFO = -4
221      ELSE IF( LDAB.LT.KL+KU+1 ) THEN
222         INFO = -6
223      ELSE IF( LDAFB.LT.2*KL+KU+1 ) THEN
224         INFO = -8
225      END IF
226      IF( INFO.NE.0 ) THEN
227         CALL XERBLA( 'CLA_GBRCOND_C', -INFO )
228         RETURN
229      END IF
230*
231*     Compute norm of op(A)*op2(C).
232*
233      ANORM = 0.0E+0
234      KD = KU + 1
235      KE = KL + 1
236      IF ( NOTRANS ) THEN
237         DO I = 1, N
238            TMP = 0.0E+0
239            IF ( CAPPLY ) THEN
240               DO J = MAX( I-KL, 1 ), MIN( I+KU, N )
241                  TMP = TMP + CABS1( AB( KD+I-J, J ) ) / C( J )
242               END DO
243            ELSE
244               DO J = MAX( I-KL, 1 ), MIN( I+KU, N )
245                  TMP = TMP + CABS1( AB( KD+I-J, J ) )
246               END DO
247            END IF
248            RWORK( I ) = TMP
249            ANORM = MAX( ANORM, TMP )
250         END DO
251      ELSE
252         DO I = 1, N
253            TMP = 0.0E+0
254            IF ( CAPPLY ) THEN
255               DO J = MAX( I-KL, 1 ), MIN( I+KU, N )
256                  TMP = TMP + CABS1( AB( KE-I+J, I ) ) / C( J )
257               END DO
258            ELSE
259               DO J = MAX( I-KL, 1 ), MIN( I+KU, N )
260                  TMP = TMP + CABS1( AB( KE-I+J, I ) )
261               END DO
262            END IF
263            RWORK( I ) = TMP
264            ANORM = MAX( ANORM, TMP )
265         END DO
266      END IF
267*
268*     Quick return if possible.
269*
270      IF( N.EQ.0 ) THEN
271         CLA_GBRCOND_C = 1.0E+0
272         RETURN
273      ELSE IF( ANORM .EQ. 0.0E+0 ) THEN
274         RETURN
275      END IF
276*
277*     Estimate the norm of inv(op(A)).
278*
279      AINVNM = 0.0E+0
280*
281      KASE = 0
282   10 CONTINUE
283      CALL CLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE )
284      IF( KASE.NE.0 ) THEN
285         IF( KASE.EQ.2 ) THEN
286*
287*           Multiply by R.
288*
289            DO I = 1, N
290               WORK( I ) = WORK( I ) * RWORK( I )
291            END DO
292*
293            IF ( NOTRANS ) THEN
294               CALL CGBTRS( 'No transpose', N, KL, KU, 1, AFB, LDAFB,
295     $              IPIV, WORK, N, INFO )
296            ELSE
297               CALL CGBTRS( 'Conjugate transpose', N, KL, KU, 1, AFB,
298     $              LDAFB, IPIV, WORK, N, INFO )
299            ENDIF
300*
301*           Multiply by inv(C).
302*
303            IF ( CAPPLY ) THEN
304               DO I = 1, N
305                  WORK( I ) = WORK( I ) * C( I )
306               END DO
307            END IF
308         ELSE
309*
310*           Multiply by inv(C**H).
311*
312            IF ( CAPPLY ) THEN
313               DO I = 1, N
314                  WORK( I ) = WORK( I ) * C( I )
315               END DO
316            END IF
317*
318            IF ( NOTRANS ) THEN
319               CALL CGBTRS( 'Conjugate transpose', N, KL, KU, 1, AFB,
320     $              LDAFB, IPIV,  WORK, N, INFO )
321            ELSE
322               CALL CGBTRS( 'No transpose', N, KL, KU, 1, AFB, LDAFB,
323     $              IPIV, WORK, N, INFO )
324            END IF
325*
326*           Multiply by R.
327*
328            DO I = 1, N
329               WORK( I ) = WORK( I ) * RWORK( I )
330            END DO
331         END IF
332         GO TO 10
333      END IF
334*
335*     Compute the estimate of the reciprocal condition number.
336*
337      IF( AINVNM .NE. 0.0E+0 )
338     $   CLA_GBRCOND_C = 1.0E+0 / AINVNM
339*
340      RETURN
341*
342      END
343