Reputation: 203
I am trying to write a program in Matlab to solve a system of linear equations using LU decomposition which employs the use of gauss elimination and hence a lot of arithmetic steps.
The answers are close to the right solution but the round-off errors are pretty high as compared to other languages like Python etc.
For instance, one of the solutions is exactly 3 but I get 2.9877.
I know in-built functions should be used for such trivial things since Matlab is a high-computing language but if I still wanted to do it with loops etc would I always get round-off errors or is there a way to reduce those while doing numerical calculations?
I am attaching the code but it's big and not worth reading. I have still attached it for the sake of completion. One can just note the use of many arithmetic operations which introduce a lot of round-off errors.
Are these round-off errors intrinsic to Matlab and unavoidable?
clc
clear
%No of equations is n
n=3;
%WRITING THE Coefficients
A(1,1)=3;
A(1,2)=-0.1;
A(1,3)=-0.2;
B(1)=7.85;
A(2,1)=0.1;
A(2,2)=7;
A(2,3)=-0.3;
B(2)=-19.3;
A(3,1)=0.3;
A(3,2)=-0.2;
A(3,3)=10;
B(3)=71.4;
%Forward Elimination
for i=1:n-1
for j=i+1:n
fact=A(j,i)/A(i,i);
A(j,i)=fact;
A(j,j:n)=A(j,j:n)-fact*A(i,j:n);
B(j)=B(j)-fact*B(i);
end
end
disp(A)
% Calculating d matrix
sum=0;
D(1)=B(1);
for i=2:n
for j=1:i-1
sum=sum+A(i,j)*B(j);
D(i)=B(i)-sum;
end
end
disp("D =")
disp(transpose(D))
%Back Substitution
X(n)=D(n)/A(n,n);
for z=n-1:-1:1
sum=0;
for w=z+1:n
sum=sum+A(z,w)*X(w);
end
X(z)=(D(z)-sum)/(A(z,z));
end
disp("Solution X is:")
disp(transpose(X))
Upvotes: 0
Views: 110
Reputation:
Never forget to distrust your coding.
If you comment out line 24 and replace 26 with A(j,i:n)=A(j,i:n)-fact*A(i,i:n);
you get the nice solution (long format)
3.000000000000000
-2.500000000000000
7.000000000000002
I am not saying that this is the best fix (it is not), but it clearly shows that roundoff is not guilty. The wrong solution remained close because the system is strongly diagonal-dominant.
Upvotes: 5