Gaussian elimination¶
We create a 4-by-4 linear system with the matrix
A = [
2 0 4 3
-4 5 -7 -10
1 15 2 -4.5
-2 0 2 -13
];
and with the right-hand side
b = [ 4, 9, 29, 40 ];
We define an augmented matrix by tacking \(b\) on the end as a new column.
S = [A b]
4×5 Array{Float64,2}:
2.0 0.0 4.0 3.0 4.0
-4.0 5.0 -7.0 -10.0 9.0
1.0 15.0 2.0 -4.5 29.0
-2.0 0.0 2.0 -13.0 40.0
The goal is to introduce zeros into the lower triangle of this matrix. By using only elementary row operations, we ensure that the matrix \(S\) always represents a linear system that is equivalent to the original. We proceed from left to right and top to bottom. The first step is to put a zero in the (2,1) location using a multiple of row 1:
@show mult21 = S[2,1]/S[1,1];
S[2,:] -= mult21*S[1,:]; # -= means "subtract the following from"
S
mult21 = S[2, 1] / S[1, 1] = -2.0
4×5 Array{Float64,2}:
2.0 0.0 4.0 3.0 4.0
0.0 5.0 1.0 -4.0 17.0
1.0 15.0 2.0 -4.5 29.0
-2.0 0.0 2.0 -13.0 40.0
We repeat the process for the (3,1) and (4,1) entries.
@show mult31 = S[3,1]/S[1,1];
S[3,:] -= mult31*S[1,:];
@show mult41 = S[4,1]/S[1,1];
S[4,:] -= mult41*S[1,:];
S
mult31 = S[3, 1] / S[1, 1] = 0.5
mult41 = S[4, 1] / S[1, 1] = -1.0
4×5 Array{Float64,2}:
2.0 0.0 4.0 3.0 4.0
0.0 5.0 1.0 -4.0 17.0
0.0 15.0 0.0 -6.0 27.0
0.0 0.0 6.0 -10.0 44.0
The first column has the zero structure we want. To avoid interfering with that, we no longer add multiples of row 1 to anything. Instead, to handle column 2, we use multiples of row 2. We’ll also exploit the highly repetitive nature of the operations to write them as a loop.
for i = 3:4
mult = S[i,2]/S[2,2]
S[i,:] -= mult*S[2,:]
end
S
4×5 Array{Float64,2}:
2.0 0.0 4.0 3.0 4.0
0.0 5.0 1.0 -4.0 17.0
0.0 0.0 -3.0 6.0 -24.0
0.0 0.0 6.0 -10.0 44.0
We finish out the triangularization with a zero in the (4,3) place. It’s a little silly to use a loop for just one iteration, but the point is to establish a pattern.
for i = 4
mult = S[i,3]/S[3,3]
S[i,:] -= mult*S[3,:]
end
S
4×5 Array{Float64,2}:
2.0 0.0 4.0 3.0 4.0
0.0 5.0 1.0 -4.0 17.0
0.0 0.0 -3.0 6.0 -24.0
0.0 0.0 0.0 2.0 -4.0
Recall that \(S\) is an augmented matrix: it represents the system \(Ux=z\), where
U = S[:,1:4]
4×4 Array{Float64,2}:
2.0 0.0 4.0 3.0
0.0 5.0 1.0 -4.0
0.0 0.0 -3.0 6.0
0.0 0.0 0.0 2.0
z = S[:,5]
4-element Array{Float64,1}:
4.0
17.0
-24.0
-4.0
The solutions to this system are identical to those of the original system, but this one can be solved by backward substitution.
using FundamentalsNumericalComputation
x = FNC.backsub(U,z)
4-element Array{Float64,1}:
-3.0
1.0
4.0
-2.0
b - A*x
4-element Array{Float64,1}:
0.0
0.0
0.0
0.0