# Solving triangular systems¶

It’s easy to get just the lower triangular part of any matrix using the tril command.

using FundamentalsNumericalComputation

A = rand(1.:9.,5,5)
L = tril(A)

5×5 Array{Float64,2}:
8.0  0.0  0.0  0.0  0.0
4.0  4.0  0.0  0.0  0.0
7.0  3.0  5.0  0.0  0.0
5.0  8.0  2.0  3.0  0.0
4.0  8.0  2.0  7.0  8.0


We’ll set up and solve a linear system with this matrix.

b = ones(5)
x = FNC.forwardsub(L,b)

5-element Array{Float64,1}:
0.125
0.125
-0.05
-0.17499999999999996
0.10312499999999997


It’s not clear what the error in this answer is. However, the residual, while not zero, is comparable to $$\varepsilon_\text{mach}$$ in size.

b - L*x

5-element Array{Float64,1}:
0.0
0.0
0.0
0.0
0.0


Next we’ll engineer a problem to which we know the exact answer.

alpha = 0.3;
beta = 2.2;
U = diagm(0=>ones(5),1=>[-1,-1,-1,-1])
U[1,[4,5]] = [ alpha-beta, beta ]
U

5×5 Array{Float64,2}:
1.0  -1.0   0.0  -1.9   2.2
0.0   1.0  -1.0   0.0   0.0
0.0   0.0   1.0  -1.0   0.0
0.0   0.0   0.0   1.0  -1.0
0.0   0.0   0.0   0.0   1.0

x_exact = ones(5)
b = [alpha,0,0,0,1]

x = FNC.backsub(U,b)
err = x - x_exact

5-element Array{Float64,1}:
2.220446049250313e-16
0.0
0.0
0.0
0.0


Everything seems OK here. But another example, with a different value for $$\beta$$, is more troubling.

alpha = 0.3;
beta = 1e12;
U = diagm(0=>ones(5),1=>[-1,-1,-1,-1])
U[1,[4,5]] = [ alpha-beta, beta ]
b = [alpha,0,0,0,1]

x = FNC.backsub(U,b)
err = x - x_exact

5-element Array{Float64,1}:
-4.882812499995559e-5
0.0
0.0
0.0
0.0


It’s not so good to get four digits of accuracy after starting with sixteen! But the source of the error is not hard to track down. Solving for $$x_1$$ performs $$(\alpha-\beta)+\beta$$ in the first row. Since $$|\alpha|$$ is so much smaller than $$|\beta|$$, this a recipe for losing digits to subtractive cancellation.