Question

write a Matlab function file to solve system Ax=b by using the output of the function lufac2a your function should have inputs f=matrix return from lufac2a, piv=array return by lufac2a and b=right hand side of your system.the only output for your system should be x

guideline

1.use the column access for the matrix entries

2. do not create any other matrix in your function-get your data directly from the matrix passed into your function

3.do not use Matlab command designed for solving system

function file LUFAC2A GIVEN BELOW

function [f, rp, flag] = lufac2a(a)

% The purpose of this function is to apply Gaussian elimination with

% partial pivoting to the input matrix a . (This follows the LINPACK

% algorithm except uses elementary Gaussian transformations from Matrix

% Computations). This function returns:

%

% f - matrix containing the information about the L and U matrices in the

% factorization PA=LU

%

% rp - array containing information about the row interchanges used in the

% elimination process

%

% flag - error flag (set to 0 if a is invertible, and set to k>0 if a

% nonzero pivot could not be found for column k)

%

% The calling sequence is [f, rp, flag] = lufac2(a)

[m,n] = size(a);

if m ~= n

disp('The matrix must be a square matrix.')

return

end

f = a;

rp = zeros(n,1);

for j=1:n-1

[mx,p] = max(abs(f(j:n,j)));

if mx == 0

flag = j;

return

end

p = p + j - 1;

rp(j) = p;

if p~=j

temp = f(j,j:n);

f(j,j:n) = f(p,j:n);

f(p,j:n) = temp;

end

i=(j+1):n;

f(i,j) = f(i,j)/f(j,j);

f(i,i) = f(i,i) - f(i,j)*f(j,i);

end

if f(n,n)==0

flag = n;

else

flag = 0;

end

return

Answer 1

a = input('Enter Input matrix a = \n')
[m,n] = size(a);

if m ~= n

disp('The matrix must be a square matrix.')

return

end

f = a;

rp = zeros(n,1);

for j=1:n-1

[mx,p] = max(abs(f(j:n,j)));

if mx == 0

flag = j;

return

end

p = p + j - 1;

rp(j) = p;

if p~=j

temp = f(j,j:n);

f(j,j:n) = f(p,j:n);

f(p,j:n) = temp;

end

i=(j+1):n;

f(i,j) = f(i,j)/f(j,j);

f(i,i) = f(i,i) - f(i,j)*f(j,i);

end

if f(n,n)==0

flag = n;

else

flag = 0;

end

return

%%%_____________________________________________

%% This code is LU

a = input('enter the a matrix a')
[m,n] = size(a);
L = eye(n,n);
%conditon chak...
if (m ~= n)
error('a is not a square matrix')
end
%partial_pivoting..
for i = 1:n-1
for j = i+1:n
[dummy,pos] = max(abs(a(i:n,i)));
if dummy~=0
     pos=pos+(i-1)                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                              
a([pos,i],:) = a([i,pos],:);
else
     break
end
%lover triangular matrix
multi = a(j,i)/a(i,i);
a(j,:) = a(j,:) - multi*a(i,:);
L(j,i) = multi ;
end
end
%disp('the upper tringular matrix');
disp(a);
disp(L)