## Equilibrate Matrices

The equilibration of a matrix is used to transform that matrix to a matrix with a lower condition number. This is usually recommended for using linear equation solvers in which pivoting is used, e.g. a Gaussian elimination routine in which pivoting is performed, a Doolittle routine using pivoting, a Crout routine using pivoting, etc.Given an equation

**A x**=

**B**where

**A**is an

*n × n*matrix and

**B**and

**x**are

*n*×1 column vectors, this routine calculates

*n×n*non-singular diagonal matrices

**R**and

**C**such that the max norm of the rows and columns of the

*n×n*matrix

**R A C**is bounded by 1. By choosing the diagonal elements of both

**R**and

**C**to be powers of 2, the significands of the corresponding elements of

**A**and

**R A C**agree and the significands of the corresponding elements of

**B**and

**R B**agree. The system

**A x**=

**B**, then is transformed to

**R A C**(

**C**) =

^{-1}x**R B**. This transformed system is then solved for

**y**=

**C**by using any of the linear equation solvers, in which pivoting is used. After solving for

^{-1}x**y**=

**C**, the value of

^{-1}x**x**is then solved by

**x**=

**C y**.

### Function List

- int Equilibrate_Matrix( double *A, int nrows, int ncols, double r[ ], double c[ ] )

Given the matrix equation**A x**=**B**, the routine Equilibrate_Matrix( ) determines the diagonal matrix**R**with diagonal given by the array*r[ ]*and the diagonal matrix**C**with diagonal given by the array*c[ ]*such that the product**R A C**equilibrates the matrix*A*. The array*r*should be dimensioned at least by*nrows*in the calling routine and the array*c*should be dimensioned at least by*ncols*in the calling routine. Although Equilibrate_Matrix( ) is primarily used to equilibrate an ill-conditioned square matrix, the routine is written to equilibrate an arbitrary*nrows×ncols*matrix. The original matrix*A*is replaced with the equilibrated matrix. The function*Equilibrate_Matrix( )*returns 0 for success or 1 if an underflow is detected which indicates a loss of significance.

- int Equilibrate_Right_Hand_Side( double B[ ], double r[ ], int n )

After the matrix**A**is equilibrated, the function Equilibrate_Right_Hand_Side( ) multiplies the right hand side,**B**, of the linear equation**A x**=**B**by the diagonal matrix**R**. I.e. after calling Equilibrate_Matrix(A,n,n,r,c) and then calling Equilibrate_Right_Hand_Side(B,r,n), the original equation**A x**=**B**is transformed to the system (**R A C**)**C**^{-1}**x**=**R B**. The input array*r*to Equilibrate_Right_Hand_Side( ) is the same array*r*calculated by Equilibrate_Matrix. The original column vector**B**is replaced with the product**RB**. The function Equilibrate_Right_Hand_Side( ) returns 0 for success, 1 if an underflow is detected, 2 if an overflow is detected or 3 if both an underflow and an overflow are detected. If an overflow is detected, then the process is halted immediately and the results are unusable, the column vector**B**has been destroyed.

- int Unequilibrate_Solution( double x[ ], double c[ ], int n )

After equilibrating the matrix**A**and the right hand side,**B**, of the matrix equation**A x**=**B**and then solving the resulting matrix equation, Unequilibrate_Solution( ) transforms the solution of the equilibrated system to a solution of the original system. The input array*c*to Unequilibrate_Solution( ) is the same array*c*calculated by Equilibrate_Matrix( ). The input column vector**x**is the solution to the equilibrated system of linear equations as calculated by a linear equation solver. The original column vector**x**is replaced with the product**Cx**. The function Unequilibrate_Solution( ) returns 0 for success, 1 if an underflow is detected, 2 if an overflow is detected or 3 if both an underflow and an overflow are detected. If an overflow is detected, then the results are unusable and the column vector**x**has been destroyed.

- int Equilibrate_CMatrix( double complex *A, int nrows, int ncols, double r[ ], double c[ ] )

Given the matrix equation**A x**=**B**, the routine Equilibrate_CMatrix determines the diagonal matrix**R**with diagonal given by the array*r[ ]*and the diagonal matrix**C**with diagonal given by the array*c[ ]*such that the product**R A C**equilibrates the matrix*A*. The array*r*should be dimensioned at least by*nrows*in the calling routine and the array*c*should be dimensioned at least by*ncols*in the calling routine. Although Equilibrate_CMatrix is primarily used to equilibrate an ill-conditioned square matrix, the routine is written to equilibrate an arbitrary*nrows×ncols*matrix. The original matrix*A*is replaced with the equilibrated matrix. The function*Equilibrate_CMatrix( )*returns 0 for success, 1 if an underflow is detected which indicates a loss of significance, 2 if an overflow is detected or 3 if both an underflow and an overflow are detected. If an overflow is detected, then the process is halted immediately and the results are unusable, the matrix**A**has been destroyed.

- int Equilibrate_CRight_Hand_Side( double complex B[ ], double r[ ], int n )

After the matrix**A**is equilibrated, the function Equilibrate_CRight_Hand_Side multiplies the right hand side,**B**, of the linear equation**A x**=**B**by the diagonal matrix**R**. I.e. after calling Equilibrate_CMatrix(A,n,n,r,c) and then calling Equilibrate_CRight_Hand_Side(B,r,n), the original equation**A x**=**B**is transformed to the system (**R A C**)**C**^{-1}**x**=**R B**. The input array*r*to Equilibrate_CRight_Hand_Side is the same array*r*calculated by Equilibrate_CMatrix. The original column vector**B**is replaced with the product**RB**. The function Equilibrate_CRight_Hand_Side( ) returns 0 for success, 1 if an underflow is detected, 2 if an overflow is detected or 3 if both an underflow and an overflow are detected. If an overflow is detected, then the process is halted immediately and the results are unusable, the column vector**B**has been destroyed.

- int Unequilibrate_CSolution( double complex x[ ], double c[ ], int n )

After equilibrating the matrix**A**and the right hand side,**B**, of the matrix equation**A x**=**B**and then solving the resulting matrix equation, Unequilibrate_CSolution( ) transforms the solution of the equilibrated system to a solution of the original system. The input array*c*to Unequilibrate_CSolution( ) is the same array*c*calculated by Equilibrate_CMatrix( ). The input column vector**x**is the solution to the equilibrated system of linear equations as calculated by a linear equation solver. The original column vector**x**is replaced with the product**Cx**. The function Unequilibrate_CSolution( ) returns 0 for success, 1 if an underflow is detected, 2 if an overflow is detected or 3 if both an underflow and an overflow are detected. If an overflow is detected, then the results are unusable and the column vector**x**has been destroyed.

*C* Source Code

- The file, equilibrate_matrix.c, contains the routines Equilibrate_Matrix( ), Equilibrate_Right_Hand_Side( ), and Unequilibrate_Solution( ).

- The file, equilibrate_cmatrix.c, contains the routines Equilibrate_CMatrix( ), Equilibrate_CRight_Hand_Side( ), and Unequilibrate_CSolution( ).

*C* Test Code, Test Results, and Build Shell Script

- The file, testequilibratematrix.c, contains a test program of Equilibrate_Matrix( ), Equilibrate_Right_Hand_Side( ), and Unequilibrate_Solution( ). This program requires the file equilibrate_matrix.c listed above.

- The file, EquilibrateMatrixTests.txt, contains the results of the test program testequilibratematrix.c.

- The file, testequilibratematrix.sh, contains the shell script used to compile, link, and execute the test program testequilibratematrix.c.

- The file, testequilibratecmatrix.c, contains a test program of Equilibrate_CMatrix( ), Equilibrate_CRight_Hand_Side( ), and Unequilibrate_CSolution( ). This program requires the file equilibrate_cmatrix.c listed above.

- The file, EquilibrateCMatrixTests.txt, contains the results of the test program testequilibratecmatrix.c.

- The file, testequilibratecmatrix.sh, contains the shell script used to compile, link, and execute the test program testequilibratecmatrix.c.