% glu.m: general LU factorization without pivoting. clear; format short; % Step 0: Assign the matrix A and initialize two vectors d and w. % The components of the vector d are prescribed (nonzero) % values of the diagonal entries. The vector w is a logical % array: if w(k) = 1, then L(k,k) = d(k); if w(k) = -1, then % U(k,k) = d(k); if w(k) = 0, then L(k,k) = U(k,k). The last % case corresponds to the Cholesky factorization. n = 4; A = zeros(n,n); for i = 1:n for j = 1:n A(i,j) = 1/(i+j-1); end end d = zeros(n,1); for i = 1:n d(i) = 1; end w = zeros(n,1); for i = 1:n w(i) = 1; % this gives the Doolittle factorization % w(i) = -1; % this gives the Crout factorization % w(i) = 0; % this gives the Cholesky factorization end % Step 1: Basic Gaussian Elimination. L = zeros(n,n); % initialize the matrix L U = zeros(n,n); % initialize the matrix U for k = 1:n if w(k) == 1 L(k,k) = d(k); U(k,k) = A(k,k); for s = 1:(k-1) U(k,k) = U(k,k)-L(k,s)*U(s,k); end U(k,k) = U(k,k)/L(k,k); for j = (k+1):n U(k,j) = A(k,j); for s = 1:(k-1) U(k,j) = U(k,j)-L(k,s)*U(s,j); end U(k,j) = U(k,j)/L(k,k); end for i = (k+1):n L(i,k) = A(i,k); for s = 1:(k-1) L(i,k) = L(i,k)-L(i,s)*U(s,k); end L(i,k) = L(i,k)/U(k,k); end elseif w(k) == -1 % (fill out this part for the Crout factorization) else % (fill out this part for the Cholesky factorization) end end L % output L U % output U