Annotation of rpl/lapack/lapack/dla_porpvgrw.f, revision 1.4
1.1 bertrand 1: DOUBLE PRECISION FUNCTION DLA_PORPVGRW( UPLO, NCOLS, A, LDA, AF,
2: $ LDAF, WORK )
3: *
4: * -- LAPACK routine (version 3.2.2) --
5: * -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
6: * -- Jason Riedy of Univ. of California Berkeley. --
7: * -- June 2010 --
8: *
9: * -- LAPACK is a software package provided by Univ. of Tennessee, --
10: * -- Univ. of California Berkeley and NAG Ltd. --
11: *
12: IMPLICIT NONE
13: * ..
14: * .. Scalar Arguments ..
15: CHARACTER*1 UPLO
16: INTEGER NCOLS, LDA, LDAF
17: * ..
18: * .. Array Arguments ..
19: DOUBLE PRECISION A( LDA, * ), AF( LDAF, * ), WORK( * )
20: * ..
21: *
22: * Purpose
23: * =======
24: *
25: * DLA_PORPVGRW computes the reciprocal pivot growth factor
26: * norm(A)/norm(U). The "max absolute element" norm is used. If this is
27: * much less than 1, the stability of the LU factorization of the
28: * (equilibrated) matrix A could be poor. This also means that the
29: * solution X, estimated condition numbers, and error bounds could be
30: * unreliable.
31: *
32: * Arguments
33: * =========
34: *
35: * UPLO (input) CHARACTER*1
36: * = 'U': Upper triangle of A is stored;
37: * = 'L': Lower triangle of A is stored.
38: *
39: * NCOLS (input) INTEGER
40: * The number of columns of the matrix A. NCOLS >= 0.
41: *
42: * A (input) DOUBLE PRECISION array, dimension (LDA,N)
43: * On entry, the N-by-N matrix A.
44: *
45: * LDA (input) INTEGER
46: * The leading dimension of the array A. LDA >= max(1,N).
47: *
48: * AF (input) DOUBLE PRECISION array, dimension (LDAF,N)
49: * The triangular factor U or L from the Cholesky factorization
50: * A = U**T*U or A = L*L**T, as computed by DPOTRF.
51: *
52: * LDAF (input) INTEGER
53: * The leading dimension of the array AF. LDAF >= max(1,N).
54: *
55: * WORK (input) DOUBLE PRECISION array, dimension (2*N)
56: *
57: * =====================================================================
58: *
59: * .. Local Scalars ..
60: INTEGER I, J
61: DOUBLE PRECISION AMAX, UMAX, RPVGRW
62: LOGICAL UPPER
63: * ..
64: * .. Intrinsic Functions ..
65: INTRINSIC ABS, MAX, MIN
66: * ..
67: * .. External Functions ..
68: EXTERNAL LSAME, DLASET
69: LOGICAL LSAME
70: * ..
71: * .. Executable Statements ..
72: *
73: UPPER = LSAME( 'Upper', UPLO )
74: *
75: * DPOTRF will have factored only the NCOLSxNCOLS leading minor, so
76: * we restrict the growth search to that minor and use only the first
77: * 2*NCOLS workspace entries.
78: *
79: RPVGRW = 1.0D+0
80: DO I = 1, 2*NCOLS
81: WORK( I ) = 0.0D+0
82: END DO
83: *
84: * Find the max magnitude entry of each column.
85: *
86: IF ( UPPER ) THEN
87: DO J = 1, NCOLS
88: DO I = 1, J
89: WORK( NCOLS+J ) =
90: $ MAX( ABS( A( I, J ) ), WORK( NCOLS+J ) )
91: END DO
92: END DO
93: ELSE
94: DO J = 1, NCOLS
95: DO I = J, NCOLS
96: WORK( NCOLS+J ) =
97: $ MAX( ABS( A( I, J ) ), WORK( NCOLS+J ) )
98: END DO
99: END DO
100: END IF
101: *
102: * Now find the max magnitude entry of each column of the factor in
103: * AF. No pivoting, so no permutations.
104: *
105: IF ( LSAME( 'Upper', UPLO ) ) THEN
106: DO J = 1, NCOLS
107: DO I = 1, J
108: WORK( J ) = MAX( ABS( AF( I, J ) ), WORK( J ) )
109: END DO
110: END DO
111: ELSE
112: DO J = 1, NCOLS
113: DO I = J, NCOLS
114: WORK( J ) = MAX( ABS( AF( I, J ) ), WORK( J ) )
115: END DO
116: END DO
117: END IF
118: *
119: * Compute the *inverse* of the max element growth factor. Dividing
120: * by zero would imply the largest entry of the factor's column is
121: * zero. Than can happen when either the column of A is zero or
122: * massive pivots made the factor underflow to zero. Neither counts
123: * as growth in itself, so simply ignore terms with zero
124: * denominators.
125: *
126: IF ( LSAME( 'Upper', UPLO ) ) THEN
127: DO I = 1, NCOLS
128: UMAX = WORK( I )
129: AMAX = WORK( NCOLS+I )
130: IF ( UMAX /= 0.0D+0 ) THEN
131: RPVGRW = MIN( AMAX / UMAX, RPVGRW )
132: END IF
133: END DO
134: ELSE
135: DO I = 1, NCOLS
136: UMAX = WORK( I )
137: AMAX = WORK( NCOLS+I )
138: IF ( UMAX /= 0.0D+0 ) THEN
139: RPVGRW = MIN( AMAX / UMAX, RPVGRW )
140: END IF
141: END DO
142: END IF
143:
144: DLA_PORPVGRW = RPVGRW
145: END
CVSweb interface <joel.bertrand@systella.fr>