The matrix we had to solve was:

1 2 1 4 13 0 -4 2 -5 2 0 0 -5 -7.5 -35 0 0 0 -9 -18First we normalise the upper-triangle matrix, by simply dividing each row such that the leading coefficient is one:

1 2 1 4 13 0 1 -0.5 1.2 -0.5 0 0 1 1.5 7 0 0 0 1 2(this simplifies the back-substitution, but we can skip/combine this step with the back-substitution)

For back-substitution we work our way backwards from the bottom of the matrix to the top, progressively eliminating each variable. As with gaussian elimination we select a pivot row, and subtract that from the rows above it. First, we start with the last row, and subtract 1.5 times that row from the row above.

1 2 1 4 13 0 1 -0.5 1.2 -0.5 0 0 1 0 4 <-- subtract pivot row * 1.5 0 0 0 1 2 <-- pivotSimilarly, we continue on to the second row, subtracting 1.2 times, and the top row, subtracting four times.

1 2 1 0 5 <-- subtract pivot row * 4 0 1 -0.5 0 -3 <-- subtract pivot row * 1.2 0 0 1 0 4 0 0 0 1 2 <-- pivotAgain, we repeat the process for the third column:

1 2 0 0 1 0 1 0 0 -1 0 0 1 0 4 <-- pivot 0 0 0 1 2And finally, the second column:

1 0 0 0 3 0 1 0 0 -1 <-- pivot 0 0 1 0 4 0 0 0 1 2Now we have our solution to the system of equations from our original gaussian elimination problem.

a = 3, b = -1, c = 4 and d = 2.

In words/pseudo-code, the process is:

- Pivot through all the rows, starting from the bottom to the top
- For each row above the pivot, calculate how many times we need to subtract the pivot row from this row.
- For each element in the row, subtract the corresponding element from the pivot row, multiplied by the value above.

for (int p=n-1;p>0;p--) { //pivot backwards through all the rows for (int r=p-1;r>=0;r--) { //for each row above the pivot float multiple = mat[r][p] / mat[p][p]; //how many multiples of the pivot row do we need (to subtract)? for (int c=p-1;c<m;c++) { mat[r][c] = mat[r][c] - mat[p][c]*multiple; //subtract the pivot row element (multiple times) } } }(complete code here)

This process can be applied to find the inverse of a general matrix. Beginning with any matrix we want to invert, we augment it with the identity matrix. For example:

2 4 -2 1 0 0 4 9 -3 0 1 0 -2 -3 7 0 0 1Now we can apply gaussian elimination to generate:

2 4 -2 1 0 0 0 1 1 -2 1 0 0 0 4 3 -1 1The normalise the upper triangle to get:

1 2 -1 0.5 0 0 0 1 1 -2 1 0 0 0 1 0.75 -0.25 0.25And finally, back-substitution to get our solved inverse:

1 0 0 6.8 -2.8 0.75 0 1 0 -2.8 1.2 -0.25 0 0 1 0.75 -0.25 0.25In this entire discussion I have left out ill-conditioned and singular matrices, but I'll leave modifying the code for that as an exercise for the reader.