function A = Aoneseps(n)
% A = Aoneseps(n)
%
% Calculates a matrix A of size nxn for which Gaussian
% elimination is unstable.


disp([' generate data A,x, b=Ax, for special A, with n= ',num2str(n)]);
A=-tril(ones(n));
A(1:n+1:end)=ones(n,1);
A(:,end)=ones(n,1);
xhat=randn(n,1);
b=A*xhat;

disp(' ')
disp('Solve by x = A\b:');
pause(.1)

tic, x = A\b; toc
res = norm(A*x-b);  disp(['residual ',num2str(res)]);
err = norm(x-xhat); disp(['error ',num2str(err)]);



disp(' ')
E=randn(n);
Aeps=A+   100*eps*  (  E/ (max(max(abs(E))))  );
disp('add small random perturbation to A --> Aeps')
disp('Solve by xeps = Aeps\b:');
pause(.1)

tic, xeps = Aeps\b; toc
res = norm(Aeps*xeps-b);  disp(['residual ||A*xeps-b|| ',num2str(res)]);
