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
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