Example: LU Factorization of a Sparse Matrix

The LU Factorization of the sparse 6 \times 6 matrix

 A=\begin{pmatrix} 10.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 \\ 0.0 & 10.0 & -3.0 & -1.0 & 0.0 & 0.0 \\ 0.0 & 0.0 & 15.0 & 0.0 & 0.0 & 0.0 \\ -2.0 & 0.0 & 0.0 & 10.0 & -1.0 & 0.0 \\ -1.0 & 0.0 & 0.0 & -5.0 & 1.0 & -3.0 \\ -1.0 & -2.0 & 0.0 & 0.0 & 0.0 & 6.0 \end{pmatrix}
is computed. The sparse coordinate form for A is given by a series of row, column, value triplets:

row column value
0 0 10.0
1 1 10.0
1 2 -3.0
1 3 -1.0
2 2 15.0
3 0 -2.0
3 3 10.0
3 4 -1.0
4 0 -1.0
4 3 -5.0
4 4 1.0
4 5 -3.0
5 0 -1.0
5 1 -2.0
5 5 6.0

Let

y^T = (1.0, 2.0, 3.0, 4.0, 5.0, 6.0)
so that b_1:=Ay = {(10.0, 7.0, 45.0, 33.0, -34.0, 31.0)}^T, and
b_2:=A^Ty = {(-9.0, 8.0, 39.0, 13.0, 1.0, 21.0)}^T
.

The LU factorization of A is used to solve the sparse linear systems Ax=b_1 and A^Tx=b_2 with iterative refinement. The reciprocal pivot growth factor and the reciprocal condition number are also computed.

import com.imsl.math.*;

public class SuperLUEx1
{
	public static void main(String args[]) throws Exception
	{
		int m;
		SuperLU slu;
		double conditionNumber, recip_pivot_growth;
		double[] sol = null;
		double Ferr, Berr;

		double[] b1 = { 10.0, 7.0, 45.0, 33.0, -34.0, 31.0 };
		double[] b2 = { -9.0, 8.0, 39.0, 13.0, 1.0, 21.0 };

		// Initialize matrix A.
		m = 6;
		SparseMatrix a = new SparseMatrix(m,m);

		a.set(0, 0, 10.0);
		a.set(1, 1, 10.0);
		a.set(1, 2, -3.0);
		a.set(1, 3, -1.0);
		a.set(2, 2, 15.0);
		a.set(3, 0, -2.0);
		a.set(3, 3, 10.0);
		a.set(3, 4, -1.0);
		a.set(4, 0, -1.0);
		a.set(4, 3, -5.0);
		a.set(4, 4, 1.0);
		a.set(4, 5, -3.0);
		a.set(5, 0, -1.0);
		a.set(5, 1, -2.0);
		a.set(5, 5, 6.0);

		// Compute the sparse LU factorization of a

		slu = new SuperLU(a);

		slu.setEquilibrate(false);
		slu.setColumnPermutationMethod(SuperLU.NATURAL_ORDERING);
		slu.setPivotGrowth(true);

		//slu.setPerformanceTuningParameters(1,2);
		//slu.setPerformanceTuningParameters(2,1);

		// Set option of iterative refinement
		slu.setIterativeRefinement(true);

		// Solve sparse system A*x = b1
		System.out.println();
		System.out.println("Solve sparse System Ax=b1");
		System.out.println("=========================");
		System.out.println();

		sol = slu.solve(b1);
		new PrintMatrix("Solution").print(sol);

		Ferr = slu.getForwardErrorBound();
		Berr = slu.getRelativeBackwardError();

		System.out.println();
		System.out.println("Forward error bound: "+Ferr);
		System.out.println("Relative backward error: "+ Berr);
		System.out.println();
		System.out.println();


		// Solve sparse system A^T*x = b2
		System.out.println();
		System.out.println("Solve sparse System A^Tx=b2");
		System.out.println("===========================");
		System.out.println();

		sol = slu.solveTranspose(b2);
		new PrintMatrix("Solution").print(sol);

		Ferr = slu.getForwardErrorBound();
		Berr = slu.getRelativeBackwardError();

		System.out.println();
		System.out.println("Forward error bound: "+Ferr);
		System.out.println("Relative backward error: "+ Berr);
		System.out.println();
		System.out.println();

		// Compute reciprocal pivot growth factor and condition number

		recip_pivot_growth = slu.getReciprocalPivotGrowthFactor();
		conditionNumber = slu.getConditionNumber();

		System.out.println("Pivot growth factor and condition number");
		System.out.println("========================================");
		System.out.println();
		
		System.out.println("Reciprocal pivot growth factor: "+ recip_pivot_growth);
		System.out.println("Reciprocal condition number: "+ conditionNumber);
		System.out.println();
   	}
}

Output


Solve sparse System Ax=b1
=========================

Solution
   0  
0  1  
1  2  
2  3  
3  4  
4  5  
5  6  


Forward error bound: 2.7448942589780064E-14
Relative backward error: 0.0



Solve sparse System A^Tx=b2
===========================

Solution
   0  
0  1  
1  2  
2  3  
3  4  
4  5  
5  6  


Forward error bound: 2.413791011361905E-15
Relative backward error: 0.0


Pivot growth factor and condition number
========================================

Reciprocal pivot growth factor: 1.0
Reciprocal condition number: 0.02445109780439122

Link to Java source.