MATLAB Answers

Xianjie
0

This calculate error is why?

Asked by Xianjie
on 10 Nov 2012
When I calculate some integration, normally I use the for loop the code is as following:
M=10;N=10;P=4;
for m1=-P:M
for n1=-P:N
s=(m1+P)*(N+1+P)+n1+P+1;
for m=-P:M
for n=-P:N
t=(m+P)*(N+1+P)+n+P+1;
A2(s,t)=intefai1(a,m1,m);
end
end
end
end
for m1=-P:M
for n1=-P:N
s=(m1+P)*(N+1+P)+n1+P+1;
for m=-P:M
for n=-P:N
t=(m+P)*(N+1+P)+n+P+1;
B2(s,t)=intefai1(b,n1,n);
end
end
end
end
function y = intefai1(a,m1,m)
% fai(m1,x)*fai(m,x)
if m1 >= 0
if m >= 0
y = alpha1(a,m1,m);
else
y = -gam(a,m1,abs(m));
end
else
if m >= 0
y = -gam(a,m,abs(m1));
else
y = betass(a,abs(m),abs(m1));
end
end
function y=gam(a,m,n)
% cos(m*pi*x/a)*sin(n*pi*x/a)
if m==0
if n==m
y=0;
else
y=-a*(-1+(-1)^n)/(n*pi);
end
else
if n==m
y=0;
else
y=a*n*(-1+(-1)^(m+n))/(pi*(m^2-n^2));
end
end
function y=betass(a,m,n)
% sin(m*pi*x/a)*sin(n*pi*x/a)
y=a/2*dell(m,n)*(1-dell(0,m))*(1-dell(0,n));
function y=alpha1(a,m,n)
% cos(m*pi*x/a)*cos(n*pi*x/a)
y=a*(dell(0,m)*dell(0,n)+1/2*dell(m,n)*(1-dell(0,m)*dell(0,n)));
When I use logical operation like the following
function y=betass(a,m,n)
% sin(m*pi*x/a)*sin(n*pi*x/a)
y = a/2*(m==n).*(~(m==0)).*(~(n==0));
function y=alpha1(a,m,n)
y=a*((~m&~n)+1/2*(m==n).*(m|n));
function y=gam(a,m,n)
% cos(m*pi*x/a)*sin(n*pi*x/a)
y = a*n.*((-1).^(m+n)-1)./(m.^2-n.^2)/pi;
y(isnan(y)) = 0;
M=10;N=10;P=4;
[xx,yy]=ndgrid(-P:M,-P:N);
cc = -P:M;
dd = kron(cc',ones(N+P+1,1));
mm1 = repmat(dd,1,(N+P+1)*(M+P+1));
nn1 = repmat(xx,N+P+1,M+P+1);
ee = kron(cc,ones(1,N+P+1));
mm = repmat(ee,(M+P+1)*(N+P+1),1);
nn = repmat(yy,(N+P+1),(M+P+1));
A1=intefai1(a,mm1,mm);
B1=intefai1(b,nn1,nn);
Like this, the A1, B1, A2 and B2 are not the same fully. Some element is different. The difference is 6E-19. Why only some element is different, when mm=0, mm1=5,mm1=11, the element is different. Thank you for your sugestion.

  2 Comments

Please properly format your code so it is readable.
Ok. Thank you. Now I format my code. Can you give me some sugestions? Thank you for your concern.

Sign in to comment.

Tags

2 Answers

Answer by Walter Roberson
on 10 Nov 2012

  2 Comments

Thank you for your advice. But only some function can encounter this phenomenon.
Any time you compute the same quantity in two different ways, the problem can arise.

Sign in to comment.


Answer by Dr. Seis
on 12 Nov 2012

Instead of determining whether things are equal... determine whether the difference between two numbers are within a tolerance level. For example, instead of
if a==b
% do something
end
try
if abs(a-b) < 1e-10
% do something
end

  1 Comment

Incidentally, if my suggestion for trying to vectorize you code did not help in your previous question, feel free to ask follow-up questions. That is how we learn more about what you are looking for and how you learn more about the hints we suggest.
For example, I would have benefited from knowing the information you provided here in this question for answering your previous question.

Sign in to comment.