Numerically solving a diffusion equation with a piecewise initial condition
显示 更早的评论
Consider the diffusion equation given by 
with initial conditions 
and boundary conditions 
.
with initial conditions and boundary conditions .I wish to numerically solve this equation by the finite difference formula 
where 
with ∆x = 0.1 and ∆t = 0.01. The exact solution is given by 
where with ∆x = 0.1 and ∆t = 0.01. The exact solution is given by Here is the code that I have worked out so far.
clear all
D = 1/4;
dx = 0.1;
dt = 0.01;
x = 0:dx:1;
t = 0:dt:1;
Nx = length(x);
Nt = length(t);
u = zeros(Nx, Nt);
u(x<=1/2,1) = x(x<=1/2);
u(x>1/2,1) = 1 - x(x>1/2);
u(1,:) = 0;
u(Nx,:) = 0;
for j = 1:Nt-1
for i = 2:Nx-1
u(i,j+1) = D*(dt/dx^2)*(u(i+1,j) + u(i-1,j)) + (1-2*D*(dt/dx^2))*(u(i,j));
end
end
% Compute the analytical solution V(x,t)
k = 1:100;
V = zeros(length(x), length(t));
for j = 1:Nt
for i = 1:Nx
V(i,j) = sum((4./(k.*pi).^2).*sin(k.*pi/2).* sin(k.*pi.*x(i)).*exp(-D.*t(j).*(k.*pi).^2));
end
end
% Plot the numerical and exact solutions
figure;
[X,T] = meshgrid(x,t);
surf(X,T,u');
xlabel('x');
ylabel('t');
zlabel('u(x,t)');
title('Numerical solution');
figure;
surf(X,T,V');
xlabel('x');
ylabel('t');
zlabel('V(x,t)');
title('Exact solution');
However, there seems to be something wrong with the code. Any support would be greatly appreciated.
3 个评论
Bill Greene
2023-4-25
Bas123
2023-4-25
Alexander
2023-12-2
Hi! This was really helpful code for a HW assignment I'm working on. My graphs are a little wacky as well, I am using different differencing schemes, but I think the problem may lie in your parameters. I think by finding the correct dx and dt you may be able to generate smoother graphs.
However I am also pretty new to this and not sure exactly how to find the best spatial and time step sizes. Good luck!
采纳的回答
I changed
V(i,j) = sum((4./(k.*pi).^2).*sin(k.*pi/2).* sin(k.*pi.*x(j)).*exp(-0.5.*t(i).*(k.*pi).^2));
to
V(i,j) = sum((4./(k.*pi).^2).*sin(k.*pi/2).* sin(k.*pi.*x(i)).*exp(-D.*t(j).*(k.*pi).^2));
in your code and the results look at least similar.
You should check the results in detail, i.e. plots of u over x for some times t or plots of u over t for some positions x. A surface plot looks fine, but is not suited to find errors or unprecise results.
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Electromagnetics 的更多信息
Translated by 
