Sun Performance Library Reference Manual: 1

Equilibration Scale Factors for a Symmetric Positive Definite Matrix

The subroutines described in this section compute equilibration scale factors for a real symmetric (or Hermitian) positive definite matrix A. Equilibration is intended to reduce the condition number of A with respect to the 2-norm. SCALE contains scale factors chosen so that the scaled matrix B with elements B(i,j) = SCALE(i) × A(i,j) × SCALE(j) has ones on the diagonal. This choice of SCALE puts the condition number of B within a factor N of the smallest possible condition number over all possible diagonal scalings. Note that these subroutines only compute the scale factors for A. Additional code must be written to actually equilibrate A.

Calling Sequence


CALL DPOEQU	(N, DA, LDA, DSCALE, DSCOND, DAMAX, INFO)
CALL SPOEQU	(N, SA, LDA, SSCALE, SSCOND, SAMAX, INFO)
CALL ZPOEQU	(N, ZA, LDA, DSCALE, DSCOND, DAMAX, INFO)
CALL CPOEQU	(N, CA, LDA, SSCALE, SSCOND, SAMAX, INFO)

void dpoequ	(int n, double da, int lda, double dscale, double dscond, double damax, int *info)
void spoequ	(int n, float sa, int lda, float sscale, float sscond, float samax, int *info)
void zpoequ	(int n, doublecomplex za, int lda, double dscale, double dscond, double damax, int *info)
void cpoequ	(int n, complex ca, int lda, float sscale, float sscond, float samax, int *info)

Arguments

N

Order of the matrix A. N 0.

xA

Matrix A; only the diagonal elements of A are referenced.

LDA

Leading dimension of the array A as specified in a dimension or type statement. LDA max(1, N).

xSCALE

On exit, the N scale factors for A.

xSCOND

On exit, the ratio of the smallest SCALE(i) to the largest SCALE(i). If SCOND 0.1 and AMAX is not close to overflow or underflow, it is not worth scaling by SCALE.

xAMAX

On exit, the absolute value of the largest element of A. If AMAX is very close to overflow or underflow, A should be scaled regardless of the value of SCOND.

INFO

On exit:

INFO = 0

Subroutine completed normally.

INFO < 0

The ith argument, where i = |INFO|, had an illegal value.

INFO > 0

The ith diagonal entry, where i = INFO, is non-positive.

N	Order of the matrix A. N 0.
xA	Matrix A; only the diagonal elements of A are referenced.
LDA	Leading dimension of the array A as specified in a dimension or type statement. LDA max(1, N).
xSCALE	On exit, the N scale factors for A.
xSCOND	On exit, the ratio of the smallest SCALE(i) to the largest SCALE(i). If SCOND 0.1 and AMAX is not close to overflow or underflow, it is not worth scaling by SCALE.
xAMAX	On exit, the absolute value of the largest element of A. If AMAX is very close to overflow or underflow, A should be scaled regardless of the value of SCOND.
INFO	On exit:
	INFO = 0	Subroutine completed normally.
	INFO < 0	The ith argument, where i = \|INFO\|, had an illegal value.
	INFO > 0	The ith diagonal entry, where i = INFO, is non-positive.

Sample Program



PROGRAM TEST
IMPLICIT NONE
C
INTEGER LDA, LDSCAL, N
PARAMETER (N = 4)
PARAMETER (LDA = N)
PARAMETER (LDSCAL = N)
C
DOUBLE PRECISION A(LDA,N), AMAX, SCALE(N), SCOND
INTEGER ICOL, INFO, IROW
C
EXTERNAL DPOEQU
INTRINSIC ABS
C
C Initialize the array A to store in symmetric form the 4x4
C symmetric positive definite matrix A shown below.
C
C 4096 -2 0 0
C A = -2 256 -2 0
C 0 -2 16 -2
C 0 0 -2 1
C
DATA A / 4.096D3, 38D8, -2.0D0, 2.56D2, 28D8, 0.0D0,
$ -2.0D0, 1.6D1, 8D8, 0.0D0, 0.0D0, -2.0D0, 1.0D0 /
C
C Print the initial values of the arrays.
C
PRINT 1000
DO 100, IROW = 1, N
PRINT 1010, (A(ICOL,IROW), ICOL = 1, IROW - 1),
$ (A(IROW,ICOL), ICOL = IROW, N)
100 CONTINUE
PRINT 1020
PRINT 1010, ((A(IROW,ICOL), ICOL = 1, N), IROW = 1, LDA)
C
C Compute the scale factors to use to equilibrate A.
C
CALL DPOEQU (N, A, LDA, SCALE, SCOND, AMAX, INFO)
IF (INFO .NE. 0) THEN
PRINT 1030, INFO
STOP 1
END IF
C
C Apply the scale factors to A then print the result.
C
DO 210, ICOL = 1, N
DO 200, IROW = 1, ICOL
A(IROW,ICOL) = A(IROW,ICOL) * SCALE(IROW) * SCALE(ICOL)
200 CONTINUE
210 CONTINUE
PRINT 1040
PRINT 1050, SCALE
PRINT 1060
DO 220, IROW = 1, N
PRINT 1010, (A(ICOL,IROW), ICOL = 1, IROW - 1),
$ (A(IROW,ICOL), ICOL = IROW, N)
220 CONTINUE
C
1000 FORMAT (1X, 'A in full form:')
1010 FORMAT (4(3X, F10.5))
1020 FORMAT (/1X, 'A in symmetric form: (* in unused elements)')
1030 FORMAT (1X, 'Error computing scale factors for A, INFO =',
$ 1X, I5)
1040 FORMAT (/1X, 'Scale factors:')
1050 FORMAT (1X, F8.5)
1060 FORMAT (/1X, 'Equilibrated A:')
C
END

Sample Output



A in full form:
4096.00000 -2.00000 0.00000 0.00000
-2.00000 256.00000 -2.00000 0.00000
0.00000 -2.00000 16.00000 -2.00000
0.00000 0.00000 -2.00000 1.00000

A in symmetric form: (* in unused elements)
4096.00000 -2.00000 0.00000 0.00000
********** 256.00000 -2.00000 0.00000
******** ******** 16.00000 -2.00000
******** ****** ******** 1.00000

Scale factors:
0.01562
0.06250
0.25000
1.00000

Equilibrated A:
1.00000 -0.00195 0.00000 0.00000
-0.00195 1.00000 -0.03125 0.00000
0.00000 -0.03125 1.00000 -0.50000
0.00000 0.00000 -0.50000 1.00000