private static double[] gaussianElimination(double[][] A, double[] b) {
		int n = b.length;
		// set up permutation
		int[] perm = new int[n];
		for (int i = 0; i < n; ++i)
			perm[i] = i;
		// convert to upper triangular form
		for (int i = 0; i < n; ++i) {
			// pivot
			for (int j = i + 1; j < n; ++j)
				if (Math.abs(A[perm[i]][i]) < Math.abs(A[perm[j]][i])) {
					int t = perm[j];
					perm[j] = perm[i];
					perm[i] = t;
				}
			// zero out column
			for (int j = i + 1; j < n; ++j) {
				double coeff = -A[perm[j]][i] / A[perm[i]][i];
				A[perm[j]][i] = 0;
				for (int k = i + 1; k < n; ++k)
					A[perm[j]][k] += coeff * A[perm[i]][k];
				b[perm[j]] += coeff * b[perm[i]];
			}
		}
		// run backwards substitution
		for (int i = n - 1; i >= 0; --i) {
			for (int j = i + 1; j < n; ++j)
				b[perm[i]] -= b[perm[j]] * A[perm[i]][j];
			b[perm[i]] /= A[perm[i]][i];
		}
		// un-permute
		double[] x = new double[n];
		for (int i = 0; i < n; ++i)
			x[i] = b[perm[i]];
		return x;
	}