how to solve 1D heat equation by Crank-Nicolson method
115 ビュー (過去 30 日間)
古いコメントを表示
I need to solve a 1D heat equation
by Crank-Nicolson method
. The tempeture on both ends of the interval is given as the fixed value u(0,t)=2, u(L,t)=0.5. I solve the equation through the below code, but the result is wrong. Attached figures are the correct result. I don't know why? Could you please anyone offer me a hand? Thanks a lot.
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/272299/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/272300/image.png)
clear all;
L=1;
dx=0.1;
imax=L/dx+1;
X=linspace(0,L,imax);
% inital value
f0 = 2-1.5.*X+sin(pi.*X);
dt=0.05;
nstep=15;
D=1.44;
alpha=D*dt/(dx)^2;
e = ones(imax,1);
B = [-alpha*e 2*(1+alpha)*e -alpha*e];
Lx = spdiags(B,[-1 0 1],imax,imax);
B = [alpha*e 2*(1-alpha)*e alpha*e];
Rx = spdiags(B,[-1 0 1],imax,imax);
%%CN method:
u=zeros(imax,nstep+1);
u(:,1)=(f0)';
for n=2:nstep+1
u(:,n)=Lx\(Rx*u(:,n-1));
% boundary value
u(1,n)=2;
u(end,n)=0.5;
end
% display contour plot
t=[0:nstep];
[Xg,tg] = meshgrid(X,t);
figure;
contourf(Xg,tg,u.');
4 件のコメント
reyhdh saleh
2022 年 1 月 1 日
%This program is meant to test out the Crank-Nicolson scheme using a simple
%nonlinear diffusion scheme.
n=5000;m=30;
t=linspace(0,20,n);% set time and distance scales
x=linspace(0,1,m);
dx=x(2)-x(1);dt=t(2)-t(1); %Increments
s=dt/(2*dx^2);%Useful for the solution.
u=zeros(n,m); %set up solution matrix
p=s*ones(1,m-1); q=-(1+2*s)*ones(1,m);
Z=diag(p,1)+diag(p,-1)+diag(q);
A=Z(2:m-1,2:m-1);
%Add in intial condition:
u(1,:)=exp(-5*(x-0.5).^2);
v=zeros(m-2,1);
%Solve the system
for i=2:n-1
%Construct the RHS for the solution
for j=2:m-1
v(j-1)=s*u(i-1,j+1)-(2*s+1)*u(i-1,j)-s*u(i-1,j-1);
end
%Solve for the new time step
w=A\v;
u(i,2:m-1)=w;
u(i,1)=u(i,2);
u(i,end)=u(i,end-1);
end
I have problem ....
Undefined function or variable 'Crank'
what does it mean
Torsten
2022 年 1 月 1 日
We don't know since the string "Crank" (most probably for Crank-Nicholson) does not appear in the code you posted.
You can run the code from above just as it is as a script in MATLAB.
採用された回答
その他の回答 (1 件)
参考
カテゴリ
Help Center および File Exchange で Ordinary Differential Equations についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!