高斯消元法 - Java实现

高斯消元法

原理（CMU - PSP）：

1. I = 1

2. LET P=A[K] , I = max { |A[J, I]| , I <= J <= N}

3. IF P = 0 THEN EXIT

4. IF K > I THEN FOR J = I TO N + 1, SWAP A[I, J] AND A[K, J]

5. 5.a FOR J = I + 1 TO N + 1 LET A[I, J] = A[I, J] / A[I, I]

    5.b SET A[I, I]  = 1

6. FOR L = I + 1 to N, multiply row I by -A[L, I]   and add to row L

7. I = I + 1. IF I < N THEN GO TO 2

8. EXIT to back-substitution

/*cb*/public class GaussElimination {
	private double[][] matrix;
	
	/*mb*/public double[] gaussElimination(double[][] m) {
		this.matrix = m;

		final int N = matrix.length;
		final double[] ret = new double[N];
				
		int i = 0;
		int j = 0;
		int k;
		while(i < N - 1) {
			k = findMax(matrix, i);
			
			if(k > i) {
				for(j = i; j < N + 1; j ++) {
					swap(i, j, k);
				}
			}
						
			for(j = i + 1; j < N + 1; j ++) {
				matrix[i][j] /= matrix[i][i];
			}
			matrix[i][i] = 1.0;
						
			for(int l = i + 1; l < N; l ++) {
				double factor = -matrix[l][i];
				for(int c = 0; c < matrix[0].length; c ++) {
					matrix[l][c] += matrix[i][c] * factor;
				}
			}
			
			i ++;
		}
		 
		for(int s = N - 1; s >= 0; s --) {
			if(s == N - 1) {
				ret[s] = matrix[s][matrix[s].length - 1] / matrix[s][matrix[s].length - 2];
			} else {
				ret[s] = matrix[s][matrix[s].length - 1];
				int t = matrix[s].length - 2;
				for(; t > s; t --) {
					ret[s] -= ret[t] * matrix[s][t];
				}
				ret[s] /= matrix[s][t];
			}
		}
		
		return ret;
	}/*me*/
	
	/*mb*/private int findMax(double[][] matrix, int position) {
		int index = 0;
		
		double max = 0.0;
		for(int i = position; i < matrix.length; i ++) {
			if(Math.abs(matrix[i][position]) > max) {
				max = matrix[i][position];
				index = i;
			}
		}
		
		return index;
	}/*me*/
	
	/*mb*/private void swap(int i, int j, int k) {
		double temp = matrix[i][j];
		matrix[i][j] = matrix[k][j];
		matrix[k][j] = temp;
	}/*me*/
	
}/*ce*/