Poisson equation with finite difference method

AI Thread Summary
The discussion centers on solving the Poisson equation using the finite difference method, specifically regarding the behavior of maximum temperature as mesh size changes. It is noted that the maximum temperature should converge to a constant value as the mesh size decreases, indicating that the solution is stabilizing. Participants discuss the implementation of the finite difference equation and share code snippets, with suggestions to include print statements for debugging. A key point raised is that plotting maximum temperature against 1/n² should yield a straight line, helping to identify the true maximum temperature as the mesh size approaches zero. The conversation emphasizes the importance of verifying numerical stability and convergence in finite difference methods.
range.rover
Messages
15
Reaction score
0
Hi guys , i am solving this equation by Finite difference method.
(dt2/dx2 + dt2/dy2 )= -Q(x,y)
i have developed a program on this to calculate the maximum temperature, when i change the mesh size the maximum temperature is also changing,

Should the maximum temperature change with mesh size (or) it should be constant ?

an anybody help me out on this...please
 
Physics news on Phys.org
range.rover said:
Hi guys , i am solving this equation by Finite difference method.
(dt2/dx2 + dt2/dy2 )= -Q(x,y)
i have developed a program on this to calculate the maximum temperature, when i change the mesh size the maximum temperature is also changing,

Should the maximum temperature change with mesh size (or) it should be constant ?

an anybody help me out on this...please
When the mesh size is reduced to a small enough value, the temperatures should approach final constant values. Show us your finite difference equation.

Chet
 
Thanx for your reply Chestermiller,
B(i,j) is the right hand side matrix,
I just applied central defference to the equation and took out the value of point(i,j)

Tnew(i,j) = ((T((i+1,j) +T(i-1,j) +T(i,j+1)+T(i,j-1) )*k)+(B(i,j)*〖dx〗^2)))/(4*k)
 
range.rover said:
Thanx for your reply Chestermiller,
B(i,j) is the right hand side matrix,
I just applied central defference to the equation and took out the value of point(i,j)

Tnew(i,j) = ((T((i+1,j) +T(i-1,j) +T(i,j+1)+T(i,j-1) )*k)+(B(i,j)*〖dx〗^2)))/(4*k)
The term in parenthesis involving the T's should also be divide by 4k.

Chet
 
yes, you are true the whole term is divided by 4k
my above equation is like

T = (((a+b+c+d)*k)+(B(i,j)*dx^2))/(4*k)

a = T((i+1,j)
b =T(i-1,j)
c =T(i,j+1)
d =T(i,j-1)
B(i,j) is right hand side where dx is space between the nodes.
 
range.rover said:
yes, you are true the whole term is divided by 4k
my above equation is like

T = (((a+b+c+d)*k)+(B(i,j)*dx^2))/(4*k)

a = T((i+1,j)
b =T(i-1,j)
c =T(i,j+1)
d =T(i,j-1)
B(i,j) is right hand side where dx is space between the nodes.
Yes, better. Does this resolve your problem now?

Chet
 
NO it has'nt solved my problem. if you can be kind enough i am posting my code have a look at it%%%%%% STRICTLY ONLY FOR ODD NUMBERS %%%%%%%%%

clear all
close all
clc
n = 153;% n grids and has n - 2 interior points per dimension
x = linspace(0,1,n);% line in x direction of 0 to 1 and divided into n parts
dx = x(2)-x(1); % distance between grids
y = x;
dy = dx;
k=1; %thermal conductivity
Q=8000 ; % heat source
m=(n-2)^2;

%Right hand side matrix
B=zeros(m,1);
if (mod(m,2)==0)
B(m/2,1)=Q;
else
B((m+1)/2,1)=Q;
end

%% Main Matrix

%Initialize a heat source matrix
A = zeros(n);
A(2:n-1,2:n-1) = reshape(B,n-2,n-2);
A = (dx^2*A)/(4*k);

Tsol = A;
%boundary conditions
Tsol(1,1:n) = 0; %TOP
Tsol(n,1:n) = 0; %BOTTOM
Tsol(1:n,1) = 0; %LEFT
Tsol(1:n,n) = 0; %RIGHT

tic
error = 1;
z = 0;
i = 2:n-1;
j = 2:n-1;

while error > 0.000001
z = z+1;
Torg = Tsol;
Tsol(i,j) = abs(((Torg(i+1,j)+Torg(i-1,j)+Torg(i,j+1)+Torg(i,j-1))*k))/(4*k);
Tsol = Tsol + A;
error = max(max(abs(Torg-Tsol)));
end
toc
figure,surf(x,y,Tsol),
figure, contour(x,y,Tsol),
title('Temperature (Steady State)'),xlabel('x'),ylabel('y'),colorbar
figure,pcolor(x,y,Tsol),shading interp,
title('Temperature (Steady State)'),xlabel('x'),ylabel('y'),colorbar
 
Hi Range.Rover,

I looked over your code briefly, and it looks like you have the right idea. But I'm not going to debug your code for you. That's up to you. I would suggest that you start out by putting lots of print statements into the code, and seeing what the different variables are doing. You can also go through the calculation by hand, according to what the code dictates.

What is the symptom of your issue? Is the maximum temperature changing as you change n? If you plot the maximum temperature from successive solutions as a function of 1/n2, the plot should approach a straight line, and the true maximum temperature should be the intercept at 1/n2 approaching zero.

Chet
 
Back
Top