matlab练习程序(差分法解一维热传导方程)
差分法计算一维热传导方程是计算偏微分方程数值解的一个经典例子。
热传导方程也是一种抛物型偏微分方程。
一维热传导方程如下:
该方程的解析解为:
通过对比解析解和数值解,我们能够知道数值解的是否正确。
下面根据微分写出差分形式:
整理得:
已知网格平面三条边的边界条件,根据上面递推公式,不断递推就能计算出每个网格的值。
matlab代码如下:
clear all;close all;clc; t = 0.03; %时间范围,计算到0.03秒 x = 1; %空间范围,0-1米 m = 320; %时间方向分320个格子 n = 64; %空间方向分64个格子 ht = t/(m-1); %时间步长dt hx = x/(n-1); %空间步长dx u = zeros(m,n); %设置边界条件 i=2:n-1; xx = (i-1)*x/(n-1); u(1,2:n-1) = sin(4*pi*xx); u(:,1) = 0; u(:,end) = 0; %根据推导的差分公式计算 for i=1:m-1 for j=2:n-1 u(i+1,j) = ht*(u(i,j+1)+u(i,j-1)-2*u(i,j))/hx^2 + u(i,j); end end %画出数值解 [x,t] = meshgrid(0:x/(n-1):1,0:0.03/(m-1):0.03); mesh(x,t,u) %画出解析解 u1 = exp(-(4*pi)^2*t).*sin(4*pi*x); figure; mesh(x,t,u1); %数值解与解析解的差 figure; mesh(abs(u-u1));
数值解:
解析解:
两种解的差的绝对值: