Gaussian Elimination in Matlab - Programmers Heaven

Howdy, Stranger!

It looks like you're new here. If you want to get involved, click one of these buttons!

Categories

Welcome to the new platform of Programmer's Heaven! We apologize for the inconvenience caused, if you visited us from a broken link of the previous version. The main reason to move to a new platform is to provide more effective and collaborative experience to you all. Please feel free to experience the new platform and use its exciting features. Contact us for any issue that you need to get clarified. We are more than happy to help you.

Gaussian Elimination in Matlab

veronicak5678veronicak5678 Posts: 4Member
I am trying to write a Matlab program to perform Gaussian Elimination with Partial Pivoting. I getting bad results, but I can't see where my code deviates from the algorithm I'm using. Here is a link to the book with the steps I'm trying to follow:

[link=http://books.google.com/books?id=wmcL0y2avuUC&pg=PA368&lpg=PA368&dq=gaussian+elimination+algorithm+6.2+numerical+analysis&source=bl&ots=1esN2ZdIPZ&sig=jQkt8_9lM8CloVzy3vWz8vVxoR4&hl=en&ei=RvssTd_UD42usAOZ9ImbBw&sa=X&oi=book_result&ct=result&resnum=1&ved=0CBcQ6AEwAA#v=onepage&q&f=false]http://books.google.com/books?id=wmcL0y2avuUC&pg=PA368&lpg=PA368&dq=gaussian+elimination+algorithm+6.2+numerical+analysis&source=bl&ots=1esN2ZdIPZ&sig=jQkt8_9lM8CloVzy3vWz8vVxoR4&hl=en&ei=RvssTd_UD42usAOZ9ImbBw&sa=X&oi=book_result&ct=result&resnum=1&ved=0CBcQ6AEwAA#v=onepage&q&f=false[/link]

This is algorithm 6.2 which begins on page 362. Here is my code:
[code]clear
A= [1 -1 2 -1 ; 2 -2 3 -3; 1 1 1 0; 1 -1 4 3 ]
[~,q]=size(A);
n=q-1;
x=zeros(n,1);
m=[n,q];
p=0;
nrow=ones(n,1);
ncopy = [1, q];
%-------------------------------Step 1
for i=1:n
nrow(i,1)=i;
end
%------------------------------Step 2
for i=1:n-1
%-----------------------------Step 3
maxval = 0;
index=0;
for j =i:n
if abs(A(nrow(j),i))>maxval
maxval=abs(A(nrow(j),i));
index = j;
end
end
p=index;
%-----------------------------Step 4
if A(nrow(p),i) ==0
display 'No!'
end
%-----------------------------Step 5
if nrow(i,:)~=nrow(p,:)
ncopy = nrow(i,:);
nrow(i,:)=nrow(p,:);
nrow(p,:)= ncopy;
end
%------------------------------Step 6
for j=i+1:n
%------------------------------Step 7
m(nrow(j),i) =A(nrow(j),i)/A(nrow(i),i);
%------------------------------Step 8
A(nrow(j),:) =A(nrow(j),:)- (m(nrow(j),i)*A(nrow(i),:));
end
end
%------------------------------Step 9
if A(nrow(n),n)==0
display 'No!'
end
%------------------------------Step 10
x(n,1)=A(nrow(n),q)/A(nrow(n),n);
%------------------------------Step 11
for i = n-1:-1:1
sig=0;
for j=i+1:n
sig =sig + A(nrow(i),j)*x(j,1);
end
x(i,1) = (A(nrow(i),q) - sig)/A(nrow(i),i);
end
%-----------------------------Step 10
display ' A is'
disp (A)
display ' x is'
disp (x)[/code]
Sign In or Register to comment.