Computing Assignment { Gaussian Elimination on a large random matrix To complete this assignment, I suggest you download the GERandom.m from Canvas (Computing Assignment page). In class you have seen the nite-precision computations with Gaussian Elimination, implemented via Matlab's backslash command, leads to small errors for small matrices but possibly larger errors for larger matrices. The purpose of this assignment is to quantify the growth of this error for matrices with random entries. Let A be a random N N matrix with entries in [1; 1] , x = 1 N (1; 1; : : : ; 1)T be an N-vector (with unit length measured in the L1-norm) and b = Ax be the right-hand side vector. Let y = (yi) 2 RN be the result of solving the system Ax = b in nite precision using the backslash command i.e. y = Anb. To measure the error between x and y, we let = max 1iN jxi yij be the maximum componentwise dierence between the two vectors. Since A is a matrix with random values, we need to run this calculation a number of times with dierent realizations of A in order to get a reasonable value for . Let M be the number of trials and suppose that for the kth trial the error is (k). We dene the mean error as follows: EN = 1 M (1) + (2) + ::: + (M) The goal of this assignment is to investigate the size of the matrix N = N at which the mean error for Gaussian Elimination is EN 1. In other words, the point at which round-o error in Gaussian elimination is of the same magnitude as the vector x. In practice, your computer will not have the processing power to nd N exactly. Instead, you should extrapolate your data. Find EN for reasonable values of N, make a plot of log(N) versus log(EN) and then perform a suitable extrapolation. Your conclusions should be explained in a one-page report. Your report must include the following: A plot of N versus EN. 1 Justication for the values of N and M you chose. Explanation of how you do the extrapolation. An estimation of the number N. The following are useful Matlab built-in functions: A=2*rand(N,N)-1; a random N N matrix with entries in [1; 1] loglog(Ns, Err); loglog plot for N versus error p=polyfit(log(Ns), log(Err), 1); returns the coecients for the straight line y = p1x + p2 which ts to the loglog data

GeRandom.m file:

% Gaussian Elimination for a random matrix % Description: Computes the mean solution error over Ntr trials for the % system Ax=b where A is a random NxN matrix and x is an N-vector (with unit length measured in the L1-norm). % Plot the errors % File name: GERandom.m % YOU NEED MODIFY THIS CODE FOR COMPUTING ASSIGNMENT 3

N=10; % Matrix size Ntr=100; % Number of trials

errs=zeros(Ntr,1); % Vector of errors x=ones(N,1)/N; % exact solution vector which is an N-vector (with unit length measured in the L1-norm) for i=1:Ntr A=2*rand(N,N)-1; % Construct a random NxN matrix with entries in [-1, 1] b=A*x; % Compute the right-hand side vector y=A\b; % Solve the linear system errs(i)=max(abs(y-x)); % Compute the error end

% Compute the mean and standard deviation of the error mean_err = mean(errs) sdev_err = sqrt(var(errs))

% Plot the errors semilogy(1:Ntr,errs) title(['log(error) for ' num2str(Ntr) ' solves with a ' num2str(N) 'x' num2str(N) ' matrix'],'fontsize',14) xlabel('number','fontsize',12) ylabel('log(error)','fontsize',12)