Improve the Stability by Pivoting

Below is a sample of the symbolic example we worked out in class.

>> % We illustrate with a very small example the dramatic effects of
>> % neglecting the pivoting.
>> format short e
>> A = rand(2,2)  % pick a random 2-by-2 matrix

A =

   9.5013e-01   6.0684e-01
   2.3114e-01   4.8598e-01

>> A(1,1) = eps   % eps is the machine precision

A =

   2.2204e-16   6.0684e-01
   2.3114e-01   4.8598e-01

>> % As we can see, A(1,1) is not zero, so we can theoretically move
>> % ahead with our Gaussian elimination:
>> M21=eye(2);    % multiplier matrix to eliminate A(2,1)
>> M21(2,1) = -A(2,1)/A(1,1)

M21 =

   1.0000e+00            0
  -1.0410e+15   1.0000e+00

>> U = M21*A

U =

   2.2204e-16   6.0684e-01
            0  -6.3170e+14

>> L = M21;
>> L(2,1) = -M21(2,1)

L =

   1.0000e+00            0
   1.0410e+15   1.0000e+00

>> L*U

ans =

   2.2204e-16   6.0684e-01
   2.3114e-01   5.0000e-01

>> A

A =

   2.2204e-16   6.0684e-01
   2.3114e-01   4.8598e-01

>> % Comparing A(2,2) with the (2,2) element of L*U we see a large
>> % relative error.
>> % Let us investigate the effect on solving a linear system.
>> b = rand(2,1)  % random right hand side vector

b =

   8.9130e-01
   7.6210e-01

>> % Now we use L and U to solve the system by back substitution:
>> y = L\b

y =

   8.9130e-01
  -9.2780e+14

>> x = U\y
Warning: Matrix is close to singular or badly scaled.
         Results may be inaccurate. RCOND = 3.515055e-31.

x =

            0
   1.4687e+00

>> % As you can see, MATLAB issues a warning: U is close to singular.
>> % Let us compare with the solution matlab gives us:
>> sol = A\b

sol =

   2.0901e-01
   1.4687e+00

>> % compute the residual
>> b - A*sol

ans =

   1.1102e-16
            0

>> b - A*x

ans =

            0
   4.8311e-02