Program does not show the matrix
6 views (last 30 days)
Show older comments
This program is made to approximate the solution to the Poisson equation. I have no errors with the equations, but it won't display the matrix and proceeds with displaying "End". Please help
a=0;
b=2;
c=0;
d=1;
m=5;
n=6;
N= 100;
TOL = 10^-4;
syms ga(y) gb(y) gc(x) gd(x) f(x,y);
f(x,y)=x*exp(y);
ga(y)=0;
gb(y)=2*exp(y);
gc(x)=x;
gd(x) =x*exp(1);
h=(b-a)/n;
k=(d-c)/m;
for i = 1:n-1
xi(i)=a+i*h;
end
for j = 1:m-1
yj(j)=c+j*k;
end
for i=1:n-1
for j=1:m-1
wij(i,j)=0;
end
end
lambda=h^2/k^2;
mu=2*(1+lambda);
l=1;
while l <= N
z=(-(h^2)*f(xi(1),yj(m-1))+ga(yj(m-1))+lambda*gd(xi(1))+lambda*wij(1,m-2)+wij(2,m-1))/mu;
NORM=abs(z-wij(1,m-1));
wij(1,m-1)=z;
for i=2:n-2
z = (-h^2*f(xi(i),yj(m-1))+lambda*gd(xi(i))+wij(i-1,m-1)+wij(i+1,m-1)+lambda*wij(i,m-2))/mu;
if abs(wij(i,m-1)-z)>NORM
NORM = abs(wij(i,m-1)-z);
end
wij(i,m-1) = z;
end
z=(-h^2*f(xi(n-1),yj(m-1))+gb(yj(m-1))+lambda*gd(xi(n-1))+wij(n-2,m-1)+lambda*wij(n-1,m-2))/mu;
if abs (wij(n-1,m-1)-z) > NORM
NORM = abs (wij(n-1,m-1)-z);
end
wij(n-1,m-1)=z;
for j=m-2:-1:2
z=(-(h^2)*f(xi(1),yj(j))+ga(yj(j))+lambda*wij(1,j+1)+lambda*wij(1,j-1)+wij(2,j))/mu;
if abs(wij(1,j)-z)>NORM
NORM = abs(wij(1,j)-z);
end
wij(i,j)=z;
for i=2:n-2
z=(-(h^2)*f(xi(i),yj(j))+wij(i-1,j)+lambda*wij(i,j+1)+wij(i+1,j)+lambda*wij(i,j-1))/mu;
if abs (wij(i,j)-z)>NORM
NORM=abs(wij(i,j)-z);
end
wij(i,j)=z;
end
z=(-h^2*f(xi(n-1),yj(j))+gb(yj(j))+wij(n-2,j)+lambda*wij(n-1,j+1)+lambda*wij(n-1,j-1))/mu;
if abs(wij(n-1,j)-z)>NORM
NORM=abs(wij(n-1,j)-z);
end
wij(n-1,j)=z;
end
z = (-h^2*f(xi(1),yj(1))+ga(yj(1))+lambda*gc(xi(1))+lambda*wij(1,2)+wij(2,1))/mu;
if abs(wij(1,1)-z) > NORM
NORM = abs(wij(1,1)-z);
end
wij(1,1) = z;
for i = 2:n-2
z = (-h^2*f(xi(i),yj(1))+lambda*gc(xi(i))+wij(i-1,1)+lambda*wij(i,2)+wij(i+1,1))/mu;
if abs(wij(i,1)-z)>NORM
NORM = abs(wij(i,1)-z);
end
wij(i,1) = z;
end
z=(-h^2*f(xi(n-1),yj(1))+gb(yj(1))+lambda*gc(xi(n-1))+wij(n-2,1)+lambda*wij(n-1,2))/mu;
if abs (wij(n-1,1)-z) > NORM
NORM= abs(wij(n-1,1)-z) ;
end
wij(n-1,1)=z;
if NORM<= TOL
for i=1:n-1
for j=1:m-1
MATRIX(i)=[xi(i) yj(j) wij(i,j)];
end
end
display(MATRIX)
end
l=l+1;
end
display("End")
2 Comments
Rik
on 29 Jan 2021
@LightFury Yeji it is very rude to edit away major parts of your question after receiving an answer. Please don't do it again.
Answers (1)
Walter Roberson
on 27 Jan 2021
a=0;
b=2;
c=0;
d=1;
m=5;
n=6;
N= 100;
TOL = 10^-4;
syms ga(y) gb(y) gc(x) gd(x) f(x,y);
f(x,y)=x*exp(y);
ga(y)=0;
gb(y)=2*exp(y);
gc(x)=x;
gd(x) =x*exp(1);
h=(b-a)/n;
k=(d-c)/m;
for i = 1:n-1
xi(i)=a+i*h;
end
for j = 1:m-1
yj(j)=c+j*k;
end
for i=1:n-1
for j=1:m-1
wij(i,j)=0;
end
end
lambda=h^2/k^2;
mu=2*(1+lambda);
l=1;
found = 0;
while l <= N
if found > 0
fprintf('solution was already found at iteration #%d, why did we start iteration %#d\n', found, l);
end
z=(-(h^2)*f(xi(1),yj(m-1))+ga(yj(m-1))+lambda*gd(xi(1))+lambda*wij(1,m-2)+wij(2,m-1))/mu;
NORM=abs(z-wij(1,m-1));
wij(1,m-1)=z;
for i=2:n-2
z = (-h^2*f(xi(i),yj(m-1))+lambda*gd(xi(i))+wij(i-1,m-1)+wij(i+1,m-1)+lambda*wij(i,m-2))/mu;
if abs(wij(i,m-1)-z)>NORM
NORM = abs(wij(i,m-1)-z);
end
wij(i,m-1) = z;
end
z=(-h^2*f(xi(n-1),yj(m-1))+gb(yj(m-1))+lambda*gd(xi(n-1))+wij(n-2,m-1)+lambda*wij(n-1,m-2))/mu;
if abs (wij(n-1,m-1)-z) > NORM
NORM = abs (wij(n-1,m-1)-z);
end
wij(n-1,m-1)=z;
for j=m-2:-1:2
z=(-(h^2)*f(xi(1),yj(j))+ga(yj(j))+lambda*wij(1,j+1)+lambda*wij(1,j-1)+wij(2,j))/mu;
if abs(wij(1,j)-z)>NORM
NORM = abs(wij(1,j)-z);
end
wij(i,j)=z;
for i=2:n-2
z=(-(h^2)*f(xi(i),yj(j))+wij(i-1,j)+lambda*wij(i,j+1)+wij(i+1,j)+lambda*wij(i,j-1))/mu;
if abs (wij(i,j)-z)>NORM
NORM=abs(wij(i,j)-z);
end
wij(i,j)=z;
end
z=(-h^2*f(xi(n-1),yj(j))+gb(yj(j))+wij(n-2,j)+lambda*wij(n-1,j+1)+lambda*wij(n-1,j-1))/mu;
if abs(wij(n-1,j)-z)>NORM
NORM=abs(wij(n-1,j)-z);
end
wij(n-1,j)=z;
end
z = (-h^2*f(xi(1),yj(1))+ga(yj(1))+lambda*gc(xi(1))+lambda*wij(1,2)+wij(2,1))/mu;
if abs(wij(1,1)-z) > NORM
NORM = abs(wij(1,1)-z);
end
wij(1,1) = z;
for i = 2:n-2
z = (-h^2*f(xi(i),yj(1))+lambda*gc(xi(i))+wij(i-1,1)+lambda*wij(i,2)+wij(i+1,1))/mu;
if abs(wij(i,1)-z)>NORM
NORM = abs(wij(i,1)-z);
end
wij(i,1) = z;
end
z=(-h^2*f(xi(n-1),yj(1))+gb(yj(1))+lambda*gc(xi(n-1))+wij(n-2,1)+lambda*wij(n-1,2))/mu;
if abs (wij(n-1,1)-z) > NORM
NORM= abs(wij(n-1,1)-z) ;
end
wij(n-1,1)=z;
if NORM<= TOL
for i=1:n-1
for j=1:m-1
MATRIX(i)=[xi(i) yj(j) wij(i,j)];
end
end
display(MATRIX)
found = l;
end
l=l+1;
end
display("End")
if found == 0
fprintf('Solution was never found, ran out of iterations\n');
else
fprintf('Solution was found at iteration #%d\n', found);
end
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!