A = [2 0 4 3; -4 5 -7 -10; 1 15 2 -4.5; -2 0 2 -13];
[L,U] = lufact(A)
L =
1.0000e+00 0 0 0 -2.0000e+00 1.0000e+00 0 0 5.0000e-01 3.0000e+00 1.0000e+00 0 -1.0000e+00 0 -2.0000e+00 1.0000e+00
U =
2 0 4 3 0 5 1 -4 0 0 -3 6 0 0 0 2
LtimesU = L*U
LtimesU =
2.0000e+00 0 4.0000e+00 3.0000e+00 -4.0000e+00 5.0000e+00 -7.0000e+00 -1.0000e+01 1.0000e+00 1.5000e+01 2.0000e+00 -4.5000e+00 -2.0000e+00 0 2.0000e+00 -1.3000e+01
Because MATLAB doesn't show all the digits by default, it's best to compare two quantities by taking their difference.
A - LtimesU
ans =
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
(Usually we can expect ``zero'' only up to machine precision. However, all the exact numbers in this example are also floating-point numbers.)
To solve a linear system, we no longer need the matrix .
b = [4;9;29;40];
z = forwardsub(L,b);
x = backsub(U,z)
x =
-3 1 4 -2
b - A*x
ans =
0 0 0 0