  function Feb1Lec()
% Relative error when solving Ax = b via \, i.e., Gaussian Elimination
% Results will be displayed to full machin precision, i.e., 16 digits.
format long
clc
% We will use the Pascal matrices, e.g.,
A = Pascal(5)
disp('A = L*L'' where')
L = abs(Pascal(5,1))
pause
for n = 4:4:20
    % Set up the "true" linear system..
    A = Pascal(n);
    b = A*ones(n,1);
    % Note that A\b = ones(n,1) in exact arithmetic.
    % Let's perturb and b with roundoff-level noise
    delta_A = eps*norm(A,1)*(-1+2*rand(n,n));
    delta_b = eps*norm(b,1)*(-1+2*rand(n,1));
    % Solve the perturbed system and compute the condition number
    x = (A+delta_A)\(b+delta_b);
    kappa = cond(A);
    clc
    fprintf('n = %2d\n\ncond(A,1) = %10.3e\n\nx =\n',n,kappa)
    fprintf('  %20.15f\n',x)
    pause
end
