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));

數值解：

解析解：

兩種解的差的絕對值：