CP090-计算物理热传导方程的差分解法
差分方法
figure
h=plot(x,uu(:,1),'linewidth',5);
set(h,'EraseMode','xor')
axis([0,1,0,0.25]);
fork=2:200
uu(2:99,2)=(1-2*c)*uu(2:99,1)+c*(uu(3:100,1)+uu(1:98,1))-b*dt/dx*(uu(3:100,1)-uu(2:99,1));
plot(u(1,:))
subplot(2,1,2)
plot(u(end,:))
差分方程所得的数值解的图形如图4所示，其中（a）是开始状态，(b)是最后状态。
(a) 初始状态
(b) 最后状态
【程序】
N=500;dx=0.01;dt=0.000001;
c=50*dt/dx/dx;
A=500;b=5;
x=linspace(0,1,100)';
方程设置是parobolic型，系数取为 。
解题的时间范围为 ，初始条件是 。
为了有足够的精度，将初始化的网格作了两次细分。而作图的选项为Contour和Animation。
作为对比，可以更改初始条件为 ,即 。
资料来源：数学物理方程与Matlab可视化.
(a) 整体图
(b) 上图：初始状态，(c)下图：最后状态
图3 解析解的图形
【程序】：
a2=50;b=5;
[x,t]=meshgrid(0:0.01:1,0:0.000001:0.0005);
Anfun=inline('2*(x-0.5).^2.*exp(5*x./2./50).*sin(n*pi*x)','x','n');
6第六讲 典型模型方程-热传导方的差分格式程
ADE Methods
同一时间步，同时左右扫描
6、Hopscotch Method
Comments
2、Richardaon’s Method: CTCS
3、Simple Implicit(Laasonen) Method
算子表示：
其中：
未知，每一个时间步长需要求解三对角方程组
放大因子：
4、Crank-Nicolson Method: famous
1 N 1 1 uN uN j 1 2u j j 1
where
G 1
3-D: ADI Methods
假设
∆t左移，右端合并即为C-N格式
其中 rx x 2 , ry y 2
at
at
3-D: ADI Methods
4、Splitting or Fractional-step Methods
5、ADE Methods
1-D:
或
Simple Explicit Method
差分格式的放大因子与精确解的放大因子比较
精度高
Simple Explicit Method
Simple Explicit Method: Example
Simple Explicit Method: Example
没有相位误差，幅值误差：1.88%
2-D: Crank-Nicolson Scheme
通常采用迭代方 法，需要比三对 角更多的计算资 源
3、2-D: ADI Methods
time
1/ 2 n 1 / 2 1/ 2 uin , uin 1, j , ui , j 1, j
1 n 1 n 1 uin, j 1 , ui , j , ui , j 1
第九章_热传导方程的差分解法_郑大昉
类似地,其偏微分用差分近似为: 类似地 其偏微分用差分近似为 近似为
∂ui, j,k ui, j,k+1 − ui, j,k = ∂t τ 2 ∂ ui, j,k ui+1, j,k − 2ui, j,k + ui−1, j,k = 2 ∂x h2 ∂2ui, j,k ui, j+1,k − 2ui, j,k + ui, j−1,k = 2 ∂y h2
∂ui,k ∂x ∂ui,k − ∂x + h
(9-18)
二阶中心差商可近似为 二阶中心差商可近似为: 可近似为
∂2ui,k ∂x
即:
2
=
−
(9-19)
ui+1,k − 2ui,k + ui−1,k ∂2u = 2 2 ∂x i,k h
(9-20)
时间的一阶差商近似为 近似为: 另, 对时间的一阶差商近似为
(9-27)
u(x, y,0) = ϕ(x, y)
(9-28)
其边界条件留待后面给出 边界条件留待后面给出. 留待后面给出
差分方法 仍设空间步长 h 仍设空间步长: 空间步长 时间步长: 时间步长 空间为: 网格. 空间为 N× M 网格
τ
则:
Nh = l,
M =s h
t = kτ , k = 0,1 2,... , x = ih, i = 0,1 N ,..., y = jh, j = 0,1,..., M
∆t
∂u ∆Q = −K(x, y, z, t)∆t∆S ∂n
(9-1)
t1
t2 t1
和
Q =∫ 1
∂u dt ∫∫ K(x, y, z, t) dS ∂n (S)
热传导方程的差分格式汇总
热传导方程的差分格式汇总1.显式差分格式：显式差分格式是最简单的一种方法，通过将导热方程时间和空间上的导数进行近似，引入差分算子，将方程转化为差分格式。
其中最常见的差分格式有：a. 前向差分法（Forward Difference Method）：利用当前节点和其相邻节点的温度值进行计算。
例如，在一维离散情况下，可以使用公式：u(i,j+1)=u(i,j)+α(u(i+1,j)-2u(i,j)+u(i-1,j))b. 后向差分法（Backward Difference Method）：利用当前节点和其相邻节点的温度值进行计算。
例如，在一维离散情况下，可以使用公式：u(i,j+1)=u(i,j)+α(u(i+1,j+1)-2u(i,j+1)+u(i-1,j+1))c. 中心差分法（Central Difference Method）：利用当前和其相邻节点的温度值进行计算。
例如，在一维离散情况下，可以使用公式：u(i,j+1)=u(i,j)+α(u(i+1,j)-2u(i,j)+u(i-1,j))+β(u(i+1,j)-u(i-1,j))其中α和β是时间和空间步长的比例因子。
2.隐式差分格式：显式差分格式具有较大的稳定性限制。
为了克服这个问题，可以使用隐式差分格式，其中使用下一个时间步长的温度值来求解当前时间步长。
常见的隐式差分格式有：a. C-N差分法（Crank-Nicolson Method）：利用前后两个时间步长的温度值进行计算。
例如，在一维离散情况下，可以使用公式：u(i,j+1)=u(i,j)+0.5α(u(i+1,j+1)-2u(i,j+1)+u(i-1,j+1))+0.5α(u(i+1,j)-2u(i,j)+u(i-1,j))b. 力学模拟法（Finite Element Method）：将空间离散化后，通过引入有限元方法，将热传导问题转化为线性方程组，再通过求解线性方程组得到温度分布。
一维热传导方程的差分法
一维热传导方程的差分法一维热传导方程描述了一个物体内部热的传递规律。
这个方程可用于解决各种问题，如材料的温度分布、传热速率等。
对于一维热传导方程，可以通过差分法来求解。
差分法是一种数值求解法，通过将原方程离散化成差分形式，将导数转化为有限差分，从而得到差分方程组。
通过求解差分方程组就可以得到离散点上的数值解。
关于一维热传导方程的差分法，以下是具体步骤。
1. 确定精度和空间网格数在差分法中，需要首先确定精度和空间离散化的步长。
通常情况下，精度越高，计算量越大，但是结果也越接近真实情况。
空间网格数越多，计算量也会越大，但是离散化的结果也越接近真实情况。
因此，需要在计算效率和结果准确性之间做出权衡。
2. 离散化热传导方程将一维热传导方程离散化，得到差分方程组。
通过 Taylor 展开，将导数转化为有限差分的形式，得到如下式子：$$ \frac{T_{i+1}-2T_{i}+T_{i-1}}{\Deltax^{2}}=\frac{\partial^{2}T}{\partial x^{2}}|_{x=i\Delta x,t}=\frac{1}{\alpha}\frac{\partial T}{\partial t}|_{x=i\Delta x,t} $$其中，$T_i$ 表示在 $x=i\Delta x$ 处的温度值，$\Delta x$ 表示空间分割步长，$\frac{1}{\alpha}$ 表示材料的热扩散系数。
3. 构建差分方程组通过对差分方程组进行简单的变形，得到一个带有时间变化的差分方程组：其中，$n$ 表示时间步长，$\Delta t$ 表示时间离散化步长。
4. 初始条件和边界条件为了有效地求解差分方程组，我们需要知道初始条件和给定的边界条件。
在一维热传导方程中，初始条件是物体最初的温度分布，而边界条件通常包括物体边界的温度和热流量。
5. 使用迭代算法求解差分方程组通过使用迭代算法（如欧拉法、隐式迭代法、迭代加速法等），可以求解差分方程组的数值解。
差分方法
解题的时间范围为 ，初始条件是 。
为了有足够的精度，将初始化的网格作了两次细分。而作图的选项为Contour和Animation。
作为对比，可以更改初始条件为 ,即 。
资料来源：数学物理方程与Matlab可视化.
end
mesh(X,Y,un);
axis([-1 1 -1 1 -0.4 0]);
pause(0.1)
end
figure(2)
wn=0;
fork=1:N
wnn=2*(U0-u0)*(-1)^k.*sin(k.*pi.*RR).*exp(-k^2*pi^2*a2*TT)./(pi*k.*RR);
wn=wn+wnn;
，即
稳定且非振荡的条件为
截断误差为
另一种格式为
即
该式称为隐式格式。对任何步长都是恒稳定的。在 上取值的唯一限制是，要将截断误差保持在合理的程度上从而节约计算时间。
截断误差为
。
二、一维热传导方问题
2．1 无限长细杆的热传导
无限长细杆的热传导的定解问题是
利用Fourier变换求得问题的解是
其中取初始温度分布如下：
end
2．3输运问题（非齐次方程）
定解问题是
解析解为
其中
我们分别用解析解与差分方程的数值解画图。
计算中取 ，这个时间很短，因为这个过程的时间很短。解析解的图形如图3所示，其中图（a）是以 为变量所画的表面图，从图中可以看出变化的全貌，图(b)是初始状态,图(c)是最后的状态。解析解在初始状态所画出的图形与差分方程的解有一定的偏差。
为了作图，取 ，在级数求和中只取了前面10项。
【程序】
热传导方程三层并行差分格式初始条件的计算
１００１－２４６Ｘ （ ２０１１ ）０４－０４８８－０５热传导方程三层并行差分格式初始条件的计算左风丽１崔霞２袁光伟２（１．北京应用物理与计算数学研究所高性能计算中心，北京１０００８８；２．北京应用物理与计算数学研究所计算物理重点实验室，北京１０００８８）摘要：给出二维热传导问题的三层差分格式初始条件的一种显式计算方法，对于由此形成的内边界预估校正三层并行差分算法，证明稳定性和收敛性定理．并行数值试验表明，方法稳定，且与通常采用隐式格式计算初始条件的方法相比，易于程序实现；与已有的扰动算法相比，能大幅度减小误差．三层差分格式；初始条件；稳定显式计算方法 Ｏ２４１．８２；Ｏ２４６Ａ２０１０－０８－１２２０１１－０２－１８国防基础科研项目（ Ｂ１５２０１１００１１）、中国工程物理研究院科学技术基金（２０１０Ａ０２０２０１０）、计算物理实验室基金、国家自然科学基金（６０９７３１５１，６１０３３００９）、８６３课题（２００９ＡＡ０１Ａ１３４，２０１０ＡＡ０１２３０３）及９７３ 课题（２０１１ＣＢ３０９７０２）资助项目左风丽（１９７１ －），女，河南，副研究员，硕士，主要从事大规模科学与工程并行计算研究．第４期４９０＠＠［ １ ］ Ｔｈｏｍａｓ Ｊ Ｗ． Ｎｕｍｅｒｉｃａｌ ｐａｒｔｉａｌ ｄｉｆｆｅｒｅｎｔｉａｌ ｅｑｕａｔｉｏｎｓ ｆｉｎｉｔｅ ｄｉｆｆｅｒｅｎｃｅ ｍｅｔｈｏｄｓ ［ Ｍ］． Ｓｐｒｉｎｇｅｒ－Ｖｅｒｌａｇ， １９９７．＠＠［２］ Ｒｉｃｈｔｍｙｅｒ Ｒ Ｄ，Ｍｏｒｔｏｎ Ｋ Ｗ．初值问题的差分方法［Ｍ］．袁国兴，杜明笙，王汉强，译．广州：中山大学出版社，１９９２．＠＠［３］李德元，陈光南．抛物型方程差分方法引论［Ｍ］．北京：科学出版社，１９９５．＠＠［４］陆金甫，关治．偏微分方程数值解法［Ｍ］．第二版．北京：清华大学出版社，２００４．＠＠［ ５ ］ Ｙｕａｎ Ｇ Ｗ， Ｚｕｏ Ｆ Ｌ． Ｐａｒａｌｌｅｌ ｄｉｆｆｅｒｅｎｃｅ ｓｃｈｅｍｅｓ ｆｏｒ ｈｅａｔ ｃｏｎｄｕｃｔｉｏｎ ｅｑｕａｔｉｏｎ ［ Ｊ ］． Ｉｎｔｅｒｎａｔｉｏｎａｌ Ｊｏｕｒｎａｌ ｏｆ Ｃｏｍｐｕｔｅｒ Ｍａｔｈｅｍａｔｉｃｓ， ２００３， ８０ （ ８ ） ： ９９５ － ９９９．＠＠［ ６ ］ Ｇüｎｔｅｒ Ｓ， Ｌａｃｋｎｅｒ Ｋ． Ａ ｍｉｘｅｄ ｉｍｐｌｉｃｉｔ－ｅｘｐｌｉｃｉｔ ｆｉｎｉｔｅ ｄｉｆｆｅｒｅｎｃｅ ｓｃｈｅｍｅ ｆｏｒ ｈｅａｔ ｔｒａｎｓｐｏｒｔ ｉｎ ｍａｇｎｅｔｉｚｅｄ ｐｌａｓｍａｓ ［ Ｊ ］． Ｊｏｕｒｎａｌ ｏｆ Ｃｏｍｐｕｔａｔｉｏｎａｌ Ｐｈｙｓｉｃｓ， ２００９， ２２８ ： ２８２ － ２９３．＠＠［ ７ ］ Ｙｕａｎ Ｇ Ｗ， Ｈａｎｇ Ｘ Ｄ， Ｓｈｅｎｇ Ｚ Ｑ． Ｐａｒａｌｌｅｌ ｄｉｆｆｅｒｅｎｃｅ ｓｃｈｅｍｅｓ ｗｉｔｈ ｉｎｔｅｒｆａｃｅ ｅｘｔｒａｐｏｌａｔｉｏｎ ｔｅｒｍｓ ｆｏｒ ｑｕａｓｉ－ｌｉｎｅａｒ ｐａｒａｂｏｌｉｃ ｓｙｓｔｅｍｓ ［Ｊ］． Ｓｃｉｅｎｃｅ ｉｎ Ｃｈｉｎａ Ｓｅｒｉｅｓ Ａ： Ｍａｔｈｅｍａｔｉｃｓ， ２００７， ５０（２） ： ２５３ －２７５．＠＠［ ８ ］ Ｓｈｅｎｇ Ｚ Ｑ， Ｙｕａｎ Ｇ Ｗ， Ｈａｎｇ Ｘ Ｄ． Ｕｎｃｏｎｄｉｔｉｏｎａｌ ｓｔａｂｉｌｉｔｙ ｏｆ ｐａｒａｌｌｅｌ ｄｉｆｆｅｒｅｎｃｅ ｓｃｈｅｍｅｓ ｗｉｔｈ ｓｅｃｏｎｄ ｏｒｄｅｒ ａｃｃｕｒａｃｙ ｆｏｒ ｐａｒａｂｏｌｉｃ ｅｑｕａｔｉｏｎ ［ Ｊ ］． Ａｐｐｌｉｅｄ Ｍａｔｈｅｍａｔｉｃｓ ａｎｄ Ｃｏｍｐｕｔａｔｉｏｎ， ２００７， １８４ ： １０１５ － １０３１．＠＠［ ９ ］ Ｙｕａｎ Ｇ Ｗ， Ｓｈｅｎｇ Ｚ Ｑ， Ｈａｎｇ Ｘ Ｄ． Ｔｈｅ ｕｎｃｏｎｄｉｔｉｏｎａｌ ｓｔａｂｉｌｉｔｙ ｏｆ ｐａｒａｌｌｅｌ ｄｉｆｆｅｒｅｎｃｅ ｓｃｈｅｍｅｓ ｗｉｔｈ ｓｅｃｏｎｄ ｏｒｄｅｒ ｃｏｎｖｅｒｇｅｎｃｅ ｆｏｒ ｎｏｎｌｉｎｅａｒ ｐａｒａｂｏｌｉｃ ｓｖｓｔｅｍ ［ Ｊ ］． Ｊｏｕｒｎａｌ ｏｆ Ｐａｒｔｉａｌ Ｄｉｆｆｅｒｅｎｔｉａｌ Ｅｑｕａｔｉｏｎｓ， ２００７， ２０ （ １ ） ： １ － ２０．＠＠［１０］许士良．计算机常用算法［Ｍ］．第二版．北京：清华大学出版社，１９９５．Ｉｎｉｔｉａｌｉｚａｔｉｏｎ Ｍｅｔｈｏｄ ｉｎ Ｔｈｒｅｅ－ｌａｙｅｒ Ｐａｒａｌｌｅｌ Ｄｉｆｆｅｒｅｎｃｅ Ｓｃｈｅｍｅ ｆｏｒ Ｈｅａｔ ＥｑｕａｔｉｏｎＺＵＯ ＦｅｎｇｌｉＣＵＩ Ｘｉａ ＹＵＡＮ Ｇｕａｎｇｗｅｉ。
2有限差分法及热传导数值计算PPT演示课件
t1
1 a11
(b1
a12t2
a13t3 )
t2
1 a 22
(b2
a 21t1 a 23t3 )
1 t3 a 33 (b3 a 31t1 a 32t2 )
•24
（2）假设一组解（迭代初场），记为： t1(0)、t2(0)并、t代3(0) 入迭代方程求得第一 次解
每次计算t1(1)均、t用2(1)、最t3(1新) 值代入。
(1) 平直边界上的节点
如图所示 边界节点 (m,n) 只能代表半个元体，若边界上有向 该元体传递的热流密度为q ，据能量守恒定律对该元体有：
tm1,n tm,n ytm,n1tm,n x
x
y 2
x
2
tm,n1 tm,n y
Φm,n
2xyyqw
0
Байду номын сангаас
xy tm ,n1 4 2 tm 1 ,n tm ,n 1 tm ,n 1 x2 Φ m ,n 2 x q w
非稳态项 的离散有三种不同的格式。如果将函数在节 点(n，i+1)对点(n，i)作泰勒展开，可有
•30
•31
由式(b)可得在点(n，i)处一阶导数的一种差分表示式 , 的向前差分:
类似地，将t在点(n，i-1)对点(n，i)作泰勒展开，可得 的向后差分的表达式：
如果将t在点(n，i+1)及(n，i-1)处的展开式相加，则可得 一阶导数的中心差分的表达式：
qw
y x
•16
(3) 内部角点
如图所示内部角点代表了 3/4 个元体，在同样的假设条 件下有
tm1,ntm,ny tm,n1tm,nx tm,n1tm,n x
x
热传导方程的差分格式讲解
热传导方程的左分格式—上机卖习报告二零一gg年五月一维抛物方程的初边值问题分别用向前差分格式、向后差分格式、六点对称格式，求解下列问题:du d2u”(兀0) = sin兀X、0 <x <1w(0,O = z/(l，O = 0, r >0在f = 0.05,0.1和0.2时刻的数值解，并与解析解u^t) = e-7：l sm(^x)进行比较。
1差分格式形式设空间步长h = l/N，时间步长r>0, T=M T,网比r = r/h2.(1)向前差分格式向前差分格式，即Z = /C\) ‘“；=0 =心),必=吆=0其中，丿= 1,2,…,N —1,R = 1,2,…,M—l. ^r^at/h2表示网比。
(1)式可改写成如下:M*+1 = + (i-2r)Uj + rw*_! + tfj此格式为显格式。
其矩阵表达式如下：Q-2r r)r l-2r(j、r 1一2广rl吐7、厂1一2、用丿加(2)向后差分格式(1)向后差分格式，即=0=久形)上：=WN =a其中j = 12・・\N_l,k = H,M_L (2)式可改写成- rw ：［： + (l+2r )叶' -中；；=0 + 叭此种差分格式被称为隐格式。
其矩阵表达式如下：rl + 2r -r( j \ I”-r l + 2r-r l + 2r -rW.V-1-r 1 + 2广丿MJ< UN >(3) 六点对称格式六点差分格式：喟-0 _ a加:-2喟+唸;唏- 2”； +吃,—T2L戸 戸 J眄=0产久XJM=H ；=O.将(3)式改写成-g 唸;+ (1 + 时-1 昭=g 略 + (1 - 诃 * * 咯 + /其矩阵表达式如下：(1 + r -r/2<l-r r/2 ) ( j\ -r/2 l + rr/2 1-rui-r/2 l + r -r/2r/2 1-r r/2X-r 1+2匚M丿r/2 l-2r ；E >2利用MATLAB 求解问题的过程对每种差分格式依次取N = 40., r=l/1600, r=l/3200, el/6400,用 MATLAB 求解并图形比较数值解与精确解，用表格列出不同剖分时的Z?误差。
热传导方程的推导与求解
热传导方程的推导与求解热传导方程是描述物体内部温度分布随时间变化的方程，常用于研究热传导过程和热能传递的问题。
在物理学和工程学中，热传导是一种重要的热传递方式，热传导方程的推导与求解对于理解热传导现象和解决实际问题具有重要意义。
热传导方程基于热传导定律，即热量在热传导过程中沿温度梯度方向从高温区传向低温区。
假设我们考虑一个一维热传导问题，研究物体中某一点的温度随时间的变化。
我们使用x轴表示物体的空间坐标，t表示时间。
首先，我们需要建立热传导方程的基本框架。
根据热传导定律，我们可以得到热传导方程的一般形式：∂T/∂t = α ∂²T/∂x²其中，T表示温度，t表示时间，α表示热扩散系数。
该方程说明了温度随时间和空间的变化率与热扩散系数α和温度梯度的平方成正比。
热扩散系数α反映了物体对热传导的难易程度，是与物体材料性质相关的参数。
根据热传导方程的一般形式，我们可以继续推导具体问题的热传导方程。
以一根长为L的均匀杆以及杆的初始温度分布T(x,0)为例，我们可以推导出热传导方程的初始和边界条件。
首先，我们考虑初始条件，即t=0时刻的温度分布。
假设杆的初始温度分布为T(x,0) = f(x)，其中f(x)是一个已知函数。
那么在t=0时刻，温度分布满足T(x,0) = f(x)。
其次，我们需要确定边界条件。
根据实际问题的不同特点，边界条件可以是温度的固定值或者温度梯度的固定值。
以杆的两端温度固定为T(0,t) = T0和T(L,t) = TL为例，我们可以得到边界条件。
有了初始条件和边界条件，我们可以开始求解热传导方程。
一种常用的方法是使用分离变量法。
假设温度分布可以表示为T(x,t) = X(x)T(t)，其中X(x)是与x有关的函数，T(t)是与t有关的函数。
将该形式的温度分布代入热传导方程，我们可以得到两个方程：X(x)T'(t) = αX''(x)T(t)将这两个方程变量分离，并将常数项记为-k²，我们可以得到两个独立的常微分方程：T'(t)/T(t) = αk²，X''(x)/X(x) = -k²分别求解这两个常微分方程，我们可以得到X(x)和T(t)的解。
热传导方程的差分解法
热传导方程的差分解法物理学中对热传导现象和扩散现象等物理过程的描述, 通常采用二阶偏微分方程, 统称为热传导方程.9.1. 热传导方程概述一般而言, 在介质内部传导的热量与传热时间、传热截面及温度梯度成正比. 设t 时刻, 点(),,x y z 处的温度为(),,,u x y z t , 则t ∆时间内通过S ∆横截面积传导的热量为(),,,uQ k x y z t t S n∆∆∆∂=-∂其中(),,,0k x y z t >, 是介质的热传导系数. un ∂∂是温度沿S ∆面的法向微商, 即温度梯度的法向分量. 为讨论热传导的规律, 设在介质中任取一小区域V , 其边界面S 为一封闭曲面. 现讨论自1t 至2t 时间内, 小体积V 内热量变化的情况. 首先, 通过包面S 传入V 的热量为()211,,,t t S u Q dt k x y z t ds n ∂=∂⎰⎰⎰ 由矢量积分定理可得()211,,,t t VQ dt k x y z t u dV =∇⋅∇⎡⎤⎣⎦⎰⎰⎰⎰ 其中∇是哈密顿算子.设介质的比热容为c , 密度为ρ, 则V 内温度变化所消耗的热量为212t t V u Q dt c dV tρ∂=∂⎰⎰⎰⎰设体积V 内部热源密度为(),,,F x y z t , 其物理意义是, t 时刻, 点(),,x y z 处, 单位体积热源在单位时间内产生的热量. 所有内部热源产生的热量为()213,,,t t VQ dt F x y z t dV =⎰⎰⎰⎰由能量守恒定律, 即213Q Q Q =+可得()2110t t Vu Q dt c k u F dV t ρ∂⎡⎤=-∇⋅∇-=⎢⎥∂⎣⎦⎰⎰⎰⎰因为体积和时间都是任取的, 所以有 ()u c k u F t ρ∂=∇⋅∇+∂ (9.1) 式(9.1)称为各向同性介质有热源的热传导方程, 也叫做三维非齐次热传导方程. 为简单起见, 设介质是均匀的, 即c 、ρ和k 都是常量. 再设体积V 内无热源, 即(),,,0F x y z t =, 则有u c k u t ρ∂=∆∂ (9.2) 式(9.2)称为各向同性介质无热源的热传导方程, 也叫做三维齐次热做传导方程. 其中∆是拉普拉斯算子. 式(9.2)也可表示为2222222u u u u a t xy z ⎛⎫∂∂∂∂=++ ⎪∂∂∂∂⎝⎭ (9.3)其中2k a c ρ=. 9.2. 一维热传导方程的差分解法各向同性介质中无热源的一维热传导方程为22220,0u u a a t T t x ∂∂=><≤∂∂ (9.4) 其中T 表明时间的有限范围. 要求解方程(9.4), 需要一定的初始条件和边界条件, 统称为定解条件.9.2.1 初值问题()(),0u x x x ϕ=<+∞ (9.5)即初始时刻空间各点的温度颁布函数.9.2.2 初、边值混合问题初始条件为()(),00u x x x l ϕ=≤≤ (9.6)0x =和x l =两端的边界条件有三种情况:第一类边界条件()()()()120,0,u t g t t u l t g t =⎧⎪≥⎨=⎪⎩(9.7) 第二类边界条件()()()()120,0,u t g t xt T u l t g t x∂⎧=⎪⎪∂≤≤⎨∂⎪=⎪∂⎩ (9.8) 其中()1g t 、()2g t 为给定函数.第三类边界条件()()()()()()()()11220,0,0,,u t t u t g t xt T u l t t u l t g t xλλ∂⎧-=⎪⎪∂≤≤⎨∂⎪-=⎪∂⎩ (9.9) 其中1λ、2λ、()1g t 、()2g t 为给定函数, 其中10λ≥, 20λ≥, 且不同时为零.用差分方法求解偏微分方程式(9.4), 首先要建立差分格式. 通常取空间步长和时间步长均为常量. 设空间步长为h , 时间步长为τ, 计算时的步序号空间用i 表示, 时间用k 表示.定义一阶向前商近似为1kk k i i i u u u xh ++∂-=∂一阶向后差商近似为1k k k i i i u u u xh--∂-=∂二阶中心差商作为二阶微商近似为21122,2k k k i i i i k u u u ux h +--+∂=∂ (9.10) 对时间的一阶差分近似为1,k k i i i k u u ut τ+-∂=∂ (9.11) 将(9.10)和式(9.11)代入(9.4), 并令22a h τα=(9.12)即可得一维热传导方程的差分格式为()111121,2,,10,1,,k k k k i i i i u u u u i N k Mααα++-=+-+=-= (9.13)其中,l T N M h τ⎡⎤⎡⎤==⎢⎥⎢⎥⎣⎦⎣⎦, “[]”表示取整.定解条件为()()()()()()001211,2,,11,11,,1i kk Nu i h i N u g k u g k k M ϕττ=-=+=-=-=+差分公式(9.13)为显式格式, 可由初始条件和边界条件逐次计算出任一时刻各点的温度. 习惯上把时刻计算的各点称为一层, 而计算则是一层一层进行的. 计算过程中层间各点的关系如图9.1所示. 从图中可直观地看出, 1k +时刻第i 个点的值是由k 时刻1i -, i 和1i +三点的值算出来.由于初始条件和边界条件的误差及其计算中的舍入误差, 用式(9.13)计算出的值并非该式的精确解k i u . 设计算值与其精确之间的误差为k i ε, 若当k 增加时, k i ε有减小的趋势, 或至少不增加, 则称其差分格式为稳定差分格式. 可以证明, 对于一维热传导方程, 差分格式(9.13)为稳定差分格式的充分条件是2212a h τα=≤(9.14) 差分格式(9.13)计算的具体步骤如下: 1. 给定2,,,,a l h T α2. 计算初始值: ,l T N M h τ⎡⎤⎡⎤==⎢⎥⎢⎥⎣⎦⎣⎦, 计算22h aατ=3. 计算初始值:()()011,2,,1i u i h i N ϕ=-=+ ;计算边界值:()()()()0121,11,,1k k N u g k u g k k M ττ=-=-=+ 4. 用差分格式(9.13)计算1k i u +. 泛定方程2201,0u ux t t x ∂∂=<<<∂∂初始条件()(),04101u x x x x =-≤≤边界条件()()0,001,0u t t u t =⎧⎪≤⎨=⎪⎩程序设计: clear %设置参数 lambda=1; alpha=1/6; L=1; h=0.01; T=0.6;tao=alpha*h^2/lambda; N=fix(L/h); M=fix(T/tao);%设置u 矩阵及x 的值 I=N+1; K=M+1; for i=1:I x(i)=(i-1)*h; endu(I,K)=zeros; %设置初始条件 u(:,1)=4.*x.*(1-x);%设置左端第一类边界条件 u(1,:)=0;%设置右端第一类边界条件 u(I,:)=0; %计算矩阵u for k=1:K for i=2:I-1u(i,k+1)=1/6*u(i+1,k)+2/3*u(i,k)+1/6*u(i-1,k); end end %u ； for k=1:1000:Kplot(x,u(:,k),'-k','LineWidth',2) hold on endx/cmT C Oaxis([0,1,0,1])xlabel ('\fontsize{14}\bfx/cm') ylabel ('\fontsize{14}\bfTC^O') grid on8.6 一维扩散方程的有限差分格式8.6.1 隐式六点差分格式(C —N 格式)以下介绍一维扩散方程或热传导方程的有限差分解法, 考虑一维扩散方程的定解问题()()()()()()()22max 2002111222,,0,0,0t u x t u x t a x l t t t x u x t u x k a u c a u b c x n ua ubc x l n ρ=⎧∂∂=≤≤<<⎪∂∂⎪=⎪⎛⎫⎪=⎨ ⎪∂⎝⎭+==⎪∂⎪⎪∂+==⎪∂⎩ (8.62) 取,x h t ∆∆τ==进行离散化, 如图8.12所示, 结点坐标为()()()()11,2,11,2,i kx i h i N t k k K τ=-=⎧⎪⎨=-=⎪⎩ (8.63) 结点处的函数为(),ki k i u x t u =. 在(),12i k +点, u t∂∂用中心差商,22ux∂∂用(),i k 和(),1i k +两点的中心差商的平均值代替, 则(8.62)中的偏微分方程变为()()()1111111121222k k k k k k k k i i i i i i i i u u u u u u u u hλτ+++++-+-⎡⎤-=-++-+⎣⎦(8.64) 引入212211,1,1a P P P h P Pτ==+=-, 将上式中的含()1k u +项移至等号左边, 将含()ku 项移至等号右边, 式(8.64)变为11111112122k k k k k k i i i i i i u Pu u u Pu u ++++-+--+-=++ (8.65) 上式表明由k 时的值可求得1k +时的u 值, 但要解联立方程组, 所以这种差分格式是隐式的. 整个方程涉及到六个点处的u 值, 所以称为隐式差分格式, 又称为Crank_Nicolson 格式, 简称C_N 格式, 误差为()()22O O h τ+, 是无条件稳定的.8.6.2 边界条件的差分格式由式(8.62)知, 一维扩散方程的边界条件为()()11122200u a u b c x n u a u b c x a n ∂⎧+==⎪⎪∂⎨∂⎪+==⎪∂⎩(8.66)在x 轴上设置两个虚格点0i =和1i N =+(见图8.13). 用中心差商代替.66)中的un∂∂, 则得()()1110212211222N N N b a u u u c hb a u u uc h +-⎧+-=⎪⎪⎨⎪+-=⎪⎩(8.67) 由式(8.67)解出011111222u hc b ha u b u =-+tk +k 图8.12和12222122N N N u hc b ha u b u +-=-+,代入1i =的式(8.65), 有()()11112111111122211111121111121111112121211111222222222242k k k k k k k k k k k k k k k k u Pu hc ha u u u P u hc ha u b u b ub Puha ub ub u b P u hc ha u b u ++++++++-+--+=++-+-++-=++-+整理得到()()1111111212111212k k k kb P ha u b u b P ha u b u hc +++-=-++ (8.68(a)) 同理, 代入i N =的式(8.65), 得到 ()()()()1111111211111222211122221211111222121212212222222222222222k k k k k k N N N N N N k k k k k k k kN N N N N N N N k k k k k k k N N N N N N N u Pu u u P u u hc b ha u b u Pu u hc ha u b u P u u hc ha u b u b Pu b u hc ha u u b P u b ++++-+-++++----++++----+-=++--++-=-+++--++-=-+++21k N u - 整理得()()11212122122222k k k kN N N N b u b P ha u b u b P ha u hc ++---++=+-+ (8.68(b))8.6.3 差分方程组及其求解把式(8.65)和式(8.68(a))和(8.68(b))结合起来, 构成差分方程组, 其形式为AU R = (8.69)其中, ()12,,N U u u u = 是未知量组成的矢量. 系数矩阵A 是三对角的, 而R 是由前一时刻的u 值组成的矢量()12,,N R R R R = . 该方程可利用MA TLAB 求解. 由式(8.65)和式(8.68(a))和(8.68(b))可知()()11211121121212222222k kk k ki i i i k k NN N R b P ha u b u hc R u P u u R b u b P ha u hc-+-⎧=-++⎪=++⎨⎪=+-+⎩ (8.70) 11111112213121121121b P ha b P P A P b b P ha +-⎛⎫ ⎪-- ⎪ ⎪-- ⎪= ⎪ ⎪ ⎪-- ⎪ ⎪-+⎝⎭ (8.70)8.6.4 计算实例研究细杆导热问题. 杆的初始温度是均匀的0u , 保持杆的一端的温度为不变的0u ,至于另一端则有强度为恒定的0q 的热流进入. (解析解见数理方法P214)杆上温度(),u x t 满足下列泛定方程和定解条件(数理方法P214)()20t xx u a u a k c ρ-== (8.71)00x x x lu u u q k ==⎧=⎪⎨=⎪⎩ (8.72) ()000t u u x l ==<< (8.73)边界条件不是齐次的, 首先要处理这个问题. 取一个既满足边界条件(8.72)又满足泛定方程(8.71)的函数(),v x t ,()00,q v x t u x k=+(8.74)计算程序: clear%设置边界条件参数 u0=0; q0k=10; D=1; a1=1.0; b1=0.0; a2=0.0; b2=1.0; c1=u0; c2=q0k;%设置u 矩阵及计算解方程系数 I=101; K=101; h=0.1; tao=0.1; P=tao*D/h^2; P1=1/P+1; P2=1/P-1;for i=1:I x(i)=(i-1)*h; end for k=1:K t=(k-1)*tao; endu(I,K)=zeros; %设置初始条件 u(:,1)=u0;%设置左端第一类边界条件 u(1,:)=u0; %设置系数矩阵A A(I,K)=zeros; A(1,1)=b1*P1+h*a1;x/cmu /u 0A(1,2)=-b1;A(I,K-1)=-b2;A(I,K)=b2*P1+h*a2;for i=2:K-1A(i,i)=2*P1;A(i,i-1)=-1;A(i,i+1)=-1;end%解方程for k=1:K-1R(1,1)=(b1*P2-h*a1)*u(1,k)+b1*u(2,k)+2*h*c1;R(I,1)=b2*u(I-1,k)+(b2*P2-h*a2)*u(I,k)+2*h*c2;for i=2:I-1R(i,1)=u(i-1,k)+2*P2*u(i,k)+u(i+1,k);endc=rank(A)==rank([A R]);u(:,k+1)=A\R;end%作图程序for k=10:20:100plot(x,u(:,k),'-k','LineWidth',2)hold onendaxis([0,10,0,35])xlabel ('\fontsize{14}\bfx/cm')ylabel ('\fontsize{14}\bfu/u_0')grid on%理论结果作图程序clearu0=0;q0k=10;I=101;h=0.1;D=1;for i=1:Ix(i)=(i-1)*h;endl=10;a=sqrt(D);for k=10:20:100t=0.1*k;u=0;for n=1:1000u=u+2*q0k*l/pi^2*(-1).^(n)./(n-1/2)^2*exp(-(n-1/2).^2*pi^2*a^2.*t/l^2).*sin((n-1/2).*pi.*x/l);endU=u+u0+q0k.*x;plot(x,U,':r','LineWidth',2)hold onendaxis([0,10,0,35])grid on例：clear%设置边界条件参数u0=0;q0k=10;D=1;a1=1.0;b1=0.0;a2=0.0;b2=1.0;c1=u0;c2= 0;%设置u矩阵及计算解方程系数I=101;K=101;h=0.1;tao=0.1;P=tao*D/h^2;P1=1/P+1;P2=1/P-1;for i=1:Ix(i)=(i-1)*h;endfor k=1:Kt=(k-1)*tao;endu(I,K)=zeros;%设置初始条件u(:,1)=-q0k.*x;%设置左端第一类边界条件u(1,:)=u0;%设置系数矩阵AA(I,K)=zeros;A(1,1)=b1*P1+h*a1;A(1,2)=-b1;A(I,K-1)=-b2;A(I,K)=b2*P1+h*a2;for i=2:K-1A(i,i)=2*P1;A(i,i-1)=-1;A(i,i+1)=-1;end%解方程for k=1:K-1R(1,1)=(b1*P2-h*a1)*u(1,k)+b1*u(2,k)+2*h*c1;R(I,1)=b2*u(I-1,k)+(b2*P2-h*a2)*u(I,k)+2*h*c2;for i=2:I-1R(i,1)=u(i-1,k)+2*P2*u(i,k)+u(i+1,k);endc=rank(A)==rank([A R]);u(:,k+1)=A\R;end%作图程序for k=1:10:101plot(x,u(:,k),'-k','LineWidth',2)hold onend%axis([0,10,0,35])xlabel ('\fontsize{14}\bfx/cm')ylabel ('\fontsize{14}\bfu/u_0')grid on%理论结果作图程序clearu0=0;q0k=10;I=101;h=0.1;D=1;for i=1:Ix(i)=(i-1)*h;endl=10;a=sqrt(D);for k=1:10:101t=0.1*(k-1);u=0;for n=1:10000u=u+2*q0k*l/pi^2*(-1).^(n)./(n-1/2)^2*exp(-(n-1/2).^2*pi^2*a^2.*t/l^2).*sin((n-1/2).*pi.*x/l);endU=u;plot(x,U,':r','LineWidth',2)hold onend%axis([0,10,0,35])grid onx/cmu /u 0热传导方程的混合问题泛定方程2201,0u u x t t x ∂∂=<<<∂∂初始条件 ()(),04101u x x x x =-≤≤边界条件()()0,001,0u t t u t =⎧⎪≤⎨=⎪⎩ 问题的数值解.clear%设置边界条件参数u0=0;D=1;a1=1.0;b1=0.0;a2=1.0;b2=0.0;c1=u0;c2=u0;%设置u 矩阵及计算解方程系数I=101;K=101;h=0.01;tao=0.01;P=tao*D/h^2;P1=1/P+1;P2=1/P-1;for i=1:Ix(i)=(i-1)*h;endfor k=1:Kt=(k-1)*tao;endu(I,K)=zeros;%设置初始条件u(:,1)=4.*x.*(1-x);%设置左端第一类边界条件u(1,:)=u0;%设置右端第一类边界条件u(101,:)=u0;x/cmu /u 0%设置系数矩阵AA(I,K)=zeros;A(1,1)=b1*P1+h*a1;A(1,2)=-b1;A(I,K-1)=-b2;A(I,K)=b2*P1+h*a2;for i=2:K-1A(i,i)=2*P1;A(i,i-1)=-1;A(i,i+1)=-1;end%解方程for k=1:K-1R(1,1)=(b1*P2-h*a1)*u(1,k)+b1*u(2,k)+2*h*c1;R(I,1)=b2*u(I-1,k)+(b2*P2-h*a2)*u(I,k)+2*h*c2;for i=2:I-1R(i,1)=u(i-1,k)+2*P2*u(i,k)+u(i+1,k);endc=rank(A)==rank([A R]);u(:,k+1)=A\R;end%作图程序for k=1:5:101plot(x,u(:,k),'-k','LineWidth',2)hold onend%axis([0,10,0,35])xlabel ('\fontsize{14}\bfx/cm')ylabel ('\fontsize{14}\bfu/u_0')grid on泛定方程2201,0u u x t t x ∂∂=<<<∂∂初始条件 ()()(),0sin 4101u x x x x =-≤≤边界条件 ()()0,001,0u t t u t =⎧⎪≤⎨=⎪⎩问题的数值解.clear%设置参数lambda=1;alpha=1/6;L=1;h=0.01;T=0.6;tao=alpha*h^2/lambda;N=fix(L/h);M=fix(T/tao);%设置u 矩阵及x 的值I=N+1;K=M+1;for i=1:Ix(i)=(i-1)*h;endu(I,K)=zeros;%设置初始条件u(:,1)=sin(4*pi.*x.*(1-x));%设置左端第一类边界条件u(1,:)=0;%设置右端第一类边界条件u(I,:)=0;%计算矩阵ufor k=1:Kfor i=2:I-1 u(i,k+1)=1/6*u(i+1,k)+2/3*u(i,k)+1/6*u(i-1,k); endendu;for k=1:100:1000plot(x,u(:,k),'-k','LineWidth',2)hold onend%axis([0,1,0,1])xlabel ('\fontsize{14}\bfx/cm') ylabel ('\fontsize{14}\bfTC^O') grid on。
一维热传导方程的差分法
一维热传导方程的差分法1. 引言1.1 介绍一维热传导方程的差分法一维热传导方程是描述物体内部温度分布随时间变化的数学模型。
差分法是一种常用的数值解法，通过将时间和空间进行离散化，将偏微分方程转化为差分方程，从而可以通过计算机进行数值求解。
在一维热传导方程的差分法中，我们通常将时间和空间分别进行离散化，将连续的温度变化转化为离散的温度值。
通过迭代计算，可以得到物体内各个离散点的温度随时间的变化情况。
差分法的优点在于可以较好地模拟物体内部温度分布的变化，同时可以较快地得到数值解，对于复杂的边界条件和非线性问题也有较好的适用性。
通过研究一维热传导方程的差分法，可以更好地理解物体内部温度分布的变化规律，为工程实践提供有效的数值模拟手段。
同时也可以探讨数值解法的稳定性和收敛性，为进一步的数值模拟研究提供参考。
通过不断改进差分法的算法和技术，可以更准确地预测物体内部温度变化，为工程设计和科学研究提供有力支持。
1.2 研究背景一维热传导方程是描述热量在一维空间内传递和分布的数学模型，广泛应用于工程领域和物理学中。
研究热传导方程的差分法是为了解决实际问题中复杂边界条件和非线性情况下的热传导问题，以及对传热过程进行数值模拟和分析。
在工程实践中，热传导问题经常出现在各种材料的传热过程中，例如石油钻井中地下油层的温度分布、金属材料的焊接过程中的温度控制等。
研究热传导方程的差分法可以帮助工程师们更好地理解热传导过程，优化工程设计，提高生产效率。
研究热传导方程的差分法还可以为其他科学领域提供理论支持和数值计算方法。
在地质学中用于模拟地热传导过程、在气象学中用于模拟大气环流等。
深入研究一维热传导方程的差分法对于推动科学研究和解决实际问题具有重要意义。
1.3 研究目的研究目的是通过对一维热传导方程的差分法进行深入分析和研究，探索其在实际工程和科学问题中的应用潜力。
具体来说，我们的研究目的包括以下几个方面：我们希望能够建立一种有效的数学模型，用以描述和解决一维热传导问题，为实际问题的数值模拟提供理论基础。
CP090-计算物理热传导方程的差分解法
t T
（9.6）
9.2 一维热传导方程的差分解法
边界条件： 3、 第三类边界条件：
u (, t ) x (t )u (, t ) g (t ) u (l , t ) (t )u (l , t ) g (t ) x
其中， (t )
9.1 热传导方程概述
四、三维非齐次热传导方程 由能量守恒定律，即
Q Q Q
可得：
或
u t dt V [c t ( Ku) F ]dV u c ( Ku ) F ( x, y, z, t ) t
t
（9.1）
式（9.1）称为各向同性介质有热源的热传导方程，也叫三维 非齐次热传导方程。
9.1 热传导方程概述
五、三维齐次热传导方程 当介质均匀（ c 、 和 K 为常数） 内无热源（ F ( x, y, z, t ) ）时： 、V
u c Ku, t
上式可表示为：
其中 x y z
u K u u u ( ) t c x y z
9.2 一维热传导方程的差分解法
例 9.1 求热传导方程混合问题：
u u t x u ( x,) x( x) u (, t ) , u (, t )
x , t x t
的数值解，取 N=10，h=0.1,计算到 k=36 为止。
9.2 一维热传导方程的差分解法
2、差分格式的稳定条件
h
3、具体步骤 （1）给定 , l , h, , T ； （2）计算 h / , N l / h, M T / ；
热传导方程德有限差分法
热传导方程德有限差分法在热传导领域，热传导方程是一个非常重要的数学模型。
而德有限差分法则是一种广泛使用的数值求解方法，用于解决热传导方程。
本文将介绍德有限差分法在热传导方程中的应用，包括方法的基本原理、求解过程及实际应用。
一、德有限差分法的基本原理德有限差分法是一种常用的数值程序，用于解决偏微分方程问题，尤其是热传导方程问题。
其基本思想是将二阶偏微分方程的差分替换为有限差分，再将有限差分数列的递推公式表示出来。
用这些公式代替偏微分方程中的导数，然后将其转化为一组线性方程组求解，从而得到数值解。
具体来说，偏微分方程可以通过一组一阶方程表示为：∂u/∂t = α ∂²u/∂x² (1)其中，u(x,t)表示物理量在时空域里的分布状态，α 表示热扩散系数，t表示时间。
热传导方程本质上是一个物理问题，而这里的关键在于如何求解其数值解。
德有限差分法的核心思想是将时间和空间分别分成大小相等的网格，将连续曲线上的点离散成一组点，从而转化为一个差分方程解析模型。
具体过程如下：1.选择网格网格的大小和数量；2.确定初始条件和边界条件；3.用有限差分逼近原方程；4.计算节点上的值；5.实现迭代算法。
二、对热传导方程应用德有限差分法当使用德有限差分法时，我们将网格分为水平和垂直方向，用i 和j分别表示各自的索引。
时间变量t用k表示。
由此可得，差分方程数列如下：uij(k + 1) = uij(k) + α(t/k)[ui-1,j(k) - 2ui,j(k) + ui+1,j(k) + ui,j-1(k) - 2ui,j(k) +ui,j+1(k)]这个式子表明，一个节点的表面积将受到其周围节点温度的影响，所以该节点的温度会发生变化。
通过迭代计算，我们可以得到数值解。
数值解可以通过散点图进行可视化展示，以便更好地理解结果，并作为之后评估模拟结果的基础。
三、热传导方程德有限差分法在实际应用中的举例在实际应用中，热传导方程德有限差分法可以用于解决多种问题。
一维热传导方程的差分法
一维热传导方程的差分法一维热传导方程是描述物体内部温度分布随时间演化的数学模型。
在工程领域，热传导方程经常被用来分析物体在不同热边界条件下的温度分布和热传导速度。
为了求解一维热传导方程，常常会采用差分法来进行数值计算。
差分法是一种利用差分逼近代替微分运算的数值方法，通过将空间和时间均匀划分为若干个小区间，用离散的点代表连续的物理量，在这些离散点上建立差分方程，最终得到一个离散的求解方程组。
通过求解这个方程组，可以得到不同时间和空间点上的温度分布。
一维热传导方程的一般形式可以写作：\[ \frac{\partial u}{\partial t} = \alpha \frac{\partial^2 u}{\partial x^2} \]\(u(x, t)\)为温度场变量，\(x\)为空间坐标，\(t\)为时间，\(\alpha\)为热传导系数。
为了利用差分法进行数值计算，首先需要将一维热传导方程离散化。
空间坐标可以划分为若干个网格点，记为\(x_i\)，时间可以划分为若干个时间步长，记为\(\Delta t\)。
通过差分法，可以用以下二阶中心差分逼近代替偏导数：将上述离散化的结果代入一维热传导方程，可以得到如下的差分方程：\(u_i^n\)表示在空间点\(x_i\)和时间点\(t_n\)处的温度值。
通过求解上述差分方程，可以得到物体在不同的时间和空间点的温度分布。
为了求解这个差分方程，可以采用显式差分法或者隐式差分法。
显式差分法是一种迭代数值计算的方法，通过某一时刻的温度值计算下一个时刻的温度值；隐式差分法是一种同时求解多个时刻温度值的方法，需要通过线性方程组的求解来得到下一个时间点的温度分布。
在实际工程中，差分法常常会遇到数值稳定性和收敛性的问题，需要谨慎选择时间步长和空间步长以保证数值计算的准确性。
还需要考虑边界条件和初始条件的选择，对于不同类型的热传导问题，需要考虑不同的求解策略。
差分法是求解一维热传导方程的一种重要的数值方法，通过将连续的偏导数转化为离散的差分方程进行数值计算，可以得到物体在不同时间和空间点的温度分布，为工程实践中的热传导问题提供重要的数值模拟手段。
一维热传导方程的差分法
一维热传导方程的差分法一维热传导方程是描述热传导现象的重要方程，它在材料科学、热科学以及工程领域中有着广泛的应用。
在实际工程中，为了求解一维热传导方程，常常会采用差分法这一数值计算方法。
本文将详细介绍一维热传导方程的差分法，并给出具体的数值计算步骤。
一维热传导方程可以表示为：∂u/∂t = α * ∂²u/∂x²u表示物体温度随时间和空间的变化，t为时间，x为空间坐标，α为热传导系数。
在一维情况下，我们只考虑温度随空间坐标x和时间t的变化，不考虑y和z方向的变化。
为了求解这个方程，我们需要给定初始条件和边界条件，通常我们会给定物体的初始温度分布以及物体边界的温度变化情况。
差分法是一种常用的数值计算方法，它将连续的变量用离散的形式来表示，并通过有限差分近似连续微分方程。
在一维热传导方程中，我们可以采用差分法来离散化空间和时间，然后通过迭代计算来求解温度的变化情况。
我们将空间和时间进行离散化。
假设我们将空间坐标x分成N个小段，时间t分成M 个小段，那么我们可以将空间坐标和时间分别表示为x_i = i*Δx和t_n = n*Δt，其中i = 0,1,2,...,N，n = 0,1,2,...,M，Δx和Δt分别为空间和时间的步长。
然后我们用u_i^n 来表示在空间坐标x_i和时间t_n处的温度。
接下来，我们用有限差分方法来离散化一维热传导方程。
我们可以采用中心差分法来逼近二阶空间导数：∂²u/∂x² ≈ (u_{i-1}^n - 2u_i^n + u_{i+1}^n) / Δx²通过对u_i^{n+1}进行求解，我们可以得到迭代更新方程：通过迭代计算，我们可以得到物体在空间和时间上的温度变化情况。
在实际工程中，我们通常会根据具体问题的要求来选择合适的空间步长和时间步长，并通过迭代计算来求解一维热传导方程。
九热传导方程的差分解法
u c ( K u ) F ( x , y , z , t ) t
u c ( K u ) K u t
若 F(x,y,z,t)0, c, , K, 为常数,则:
其中: 为拉普拉斯算子:
u u u u 2 2 2 x y z 所以热传导方程为: 2 2 2 u u u u ( 2 2 2) t x y z 其中: Kc.
u g ( t ) , k 0 , 1 , 2 , . . . , M 0 , k 1 k
u g () t k 0 , 1 , 2 , . . . , M N , k 2 k,
t
u0, M u M ,1 u0,1 u1,1 u0,0 u1,0
u N 1, M u N 1,1 uN 1,0
即:
gt ( ) , k 0 , 1 , 2 , . . . , M 1 k
u i,k u i 1 ,k
对时间一阶向前插商:
u i,k t
u i,k 1 u i,k
代入热传导方程:
u u u 2 u u ik , 1 ik , i 1 , k ik , i 1 , k 2 h 迭代公式:
u u ( 1 2 ) u u ; i , k 1 i 1 , k i , k i 1 , k
量. 取空间中的一个小区域 V, 其边界面 S 为一封闭曲面. 则 t1 到 t2 时刻通过包面 S 传入 V 的热量为:
u Q dt K ( x ,y , z , t ) ds 1 t S 1 n
t 2
由高斯公式:
Q [ K ( x ,y , z , t ) u ] dV 1 dt
《计算物理(本科)》[第9章]
![《计算物理(本科)》[第9章]](https://img.taocdn.com/s3/m/eed4d5196c175f0e7cd13764.png)
u i , j ,0 0 i 0,1,...,N; j 0,1,...M (初始条件, t 0时u 0) u0, j ,k u1, j ,k j 1,2,..., M 1 1, M 2 1,.., M 1, k 0,1,2,...
u N , j ,k u N 1, j ,k j 1,2,..., M 1, k 0,1,2,...
四、边界条件(如图所示的具体问题) y
第 九 章 热 传 导 方 程 的 数 值 解 法
M
绝热边界：粉红色部分绝热壁， {x=0,y∈(0,M1h).and.(M2h,Mh)}; {x=Nh,y∈(0,Mh)}
绝热位置应满足
u 0, j ,k x u N , j ,k x 0 0
M2
M1
合肥工业大学电子科学与应用物理学院
恒温热源边界：红色部分，即
第 九 章 热 传 导 方 程 的 数 值 解 法
{x=0, y ∈[M1h，M2h]} 归一化后高温源温度取“1”，即
u 0, j , k 1 j M 1 , M 1 1,..., M 2 1, M 2 ; k 0,1,2,...
M1 N x
Nh=l
则有
Mh=s
(时间t序号) (空间x序号) (空间y序号)
t=kτ k=0,1,2,… x=ih y=jh i=0,1,…,N j=0,1,…,M
合肥工业大学电子科学与应用物理学院
对节点(i，j)，在k时刻(即t=kτ)的差分式
第 九 章 热 传 导 方 程 的 数 值 解 法
1.给定λ，l，h，α，T ； 2.计算N=[l/h]，M=[T/τ ] ，τ=αh2/λ； 3.计算初始值： ui ,0 (ih ) , i 0,1,2,.......,N
热传导方程二阶并行区域分解差分算法
热传导方程二阶并行区域分解差分算法
田敏;羊丹平
【期刊名称】《山东大学学报：理学版》
【年(卷),期】2006(41)5
【摘要】提出了一类新的计算热传导方程数值解的并行差分算法.算法基于区域分解和子区域校正,在每个子区域上进行残量修正,各子域之间可以并行计算.证明了算法的收敛性,并且理论分析表明,在每一时间步,只需校正一次或两次,即可达到最优的收敛阶.数值试验表明了算法的有效性和优越性.
【总页数】9页(P12-19)
【关键词】区域分解;子区域校正;单位分解;中心差分格式;热传导方程
【作者】田敏;羊丹平
【作者单位】山东大学数学与系统科学学院
【正文语种】中文
【中图分类】O241.82
【相关文献】
1.热传导方程紧差分格式的区域分解算法 [J], 张红梅;岳素芳;许娟
2.热传导方程的一类区域分解差分算法 [J], 张红梅
3.关于热传导方程有限差分区域分解并行算法精度的注记 [J], 万正苏;方春华;张再云
4.热传导方程紧差分格式的区域分解算法 [J], 张红梅;岳素芳;许娟
5.热传导方程的一类新型重叠型并行区域分解有限差分算法 [J], 田敏;羊丹平因版权原因，仅展示原文概要，查看原文内容请购买。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
function u=fai(x) u=4*x*(1-x); function u=g1(x) u=0; function u=g2(x) u=0;
9.3 二维热传导方程的差分解法
一、二维热传导方程 各向同性介质中无热源的二维热传导方程为：
u u u K （ ）， , t T （9.9） t x x c x l , y s,
9.3 二维热传导方程的差分解法
例 9.2 求热传导方程混合问题：
u u u x , y , t y t x u ( x, y,) u ( x,, t ) u ( x,, t ) x , t u (, y, t ) y , t x u (, y, t ) y M , M y , t x u (, y, t ) M y M , t
9.2 一维热传导方程的差分解法
例 9.1 求热传导方程混合问题：
u u t x u ( x,) x( x) u (, t ) , u (, t )
x , t x t
的数值解，取 N=10，h=0.1,计算到 k=36 为止。
function main %热传导方程的差分解法 lda=1;l=1;h=0.05;alpha=0.5;tao=alpha* h^2/lda;T=tao*100;N=l/h;M=T/tao; for i=1:N+1 u(i,1)=fai((i-1)*h); end for k=1:M u(1,k)=g1(k*tao); u(N,k)=g2(k*tao); end for k=1:M for i=2:N u(i,k+1)=alpha*u(i+1,k)+(12*alpha)*u(i,k)+alpha*u(i-1,k); end plot([0:h:l],u(:,k+1)); hold on; pause(0.05); end
将差分格式（9.10）代入偏微分方程中得：
（9.10）
ui , j ,k ( )u i , j ,k (ui , j ,k u i , j ,k u i , j ,k ui , j ,k )
h
i ,,..., N , j ,,...,M , k ,,,...
k ,,,... i ,,,..., N j ,,,...,M
9.3 二维热传导方程的差分解法
对于节点 (i, j ) ，在 k 时刻有：
u i , j ,k u i , j ,k u i , j ,k t u i , j ,k u i , j ,k u i , j ,k u i , j ,k x h u i , j ,k u i , j ,k u i , j ,k u i , j ,k y h
第九章 热传导方程的差分解法
9.1 热传导方程概述
一、传入热量
t 时间内通过 S 横截面积传导的热量为： u Q K ( x, y, z , t )tS n u 其中，K ( x, y, z, t ) 是介质的热传导系数。 是温度沿 S 面的 n
法向微商，即温度梯度的法向分量。 通过包面 S 传入 V 的热量为：
9.1 热传导方程概述
五、三维齐次热传导方程 当介质均匀（ c 、 和 K 为常数） 内无热源（ F ( x, y, z, t ) ）时： 、V
u c Ku, t
上式可表示为：
其中 x y z
u K u u u ( ) t c x y z
4、 一维热传导方程的差分格式
u i ,k u i ,k
整理得：
u i ,k u i ,k u i ,k h
u i ,k u i ,k ( - )u i ,k u i ,k
h
i ,,...,N , k ,,,...,M
上式为三维齐次热传导方程。
（9.2）
9.2 一维热传导方程的差分解法
一、一维热传导方程 各向同性介质中无热源的一维热传导方程为：
u u t x
二、初始、边界条件 初始条件：
u( x,) ( x)
xl
（9.4）
3、具体步骤
h
（1）给定 , h, , T , XN, YM ； （2）计算 h / , N XN / h, M YM / h, K T / ；
（3）计算初始值： ui , j , (ih, jh) ； 计算边界值： u, j ,k g (k , jh),u N , j ,k g (k , jh) ； ui ,,k g (k , ih),ui , N ,k g (k , ih) ； （4）用差分格式计算 ui , j ,k 。
(ui 1, j ,k ui 1, j ,k ui , j 1,k ui , j 1,k ) i 0,1,..., N , j 0,1,...,M i 0,1,2,..., N j 1,2,...,M 1 1, M 2 1,...,M 1 j 0,1,...,M j M 1 , M 1 1,...,M 2
function chafen2 %二维热传导方程的差分解法 lda=1;l=1;s=2;h=0.05;alpha=0.25;tao=alpha*h^2/lda; T=tao*1000;N=l/h;M=s/h;K=T/tao;M1=0.4*M;M2=0.6*M; for i=1:N+1 for j=1:M+1 u(i,j,1)=0; for k=1:K u(i,1,k)=0; u(i,M,k)=0; if M1<=j & j<=M2 u(1,j,k)=1; end end end end for k=1:K-1 for j=2:M for i=2:N u(i,j,k+1)=(1-4*alpha)*u(i,j,k)+alpha*(u(i+1,j,k)+u(i1,j,k)+u(i,j+1,k)+u(i,j-1,k)); end u(N+1,j,k+1)=u(N,j,k+1); if (2<=j&j<=M1) || (M2<=j&j<=M) u(1,j,k+1)=u(2,j,k+1); else u(1,j,k+1)=1; end end surf(u(:,:,k+1)); pause(0.05); end
的数值解，取 N=20，M=20，h=0.05,计算到 k=100 为止。
9.3 二维热传导方程的差分解法
解：将初始边界混合问题转化为差分格式
ui , j ,k 1 (1 4 )ui , j ,k ui , j , 0 0 u i , 0 , k ui , M , k 0 u0, j ,k u1, j ,k u N 1, j , k u N , j , k u0, j ,k 1
u i ,k x
和
u i ,k t
u i ,k u i ,k
u i ,k u i ,k h
u i ,k x
3、 二阶中心差商
u i ,k ui ,k x
x
h
ui ,k ui ,k u i ,k h
9.2 一维热传导方程的差分解法
Q dt F ( x, y, z, t )dV
t V
t
其中， F ( x, y, z, t ) 为内部热源密度，表示单位时间单位体积所产生 的热量。 三、消耗热量 V 内消耗热量：
Q
t
t
其中， c 为介质的比热容， 为质量密度。
u dt c dV V t
（9.11）
9.3 二维热传导方程的差分解法
初始、边界条件的差分格式：
u i , j , (ih, jh) u , j ,k g (k , jh) u N , j ,k g (k , jh) u i , ,k g (k , ih) u i , N ,k g (k , ih)
Q dt
t
t
S
u K ( x, y, z , t ) ds n
由矢量分析（高斯散度定理）可得：
Q dt [K ( x, y, z, t )u]dV
t V
t
其中， 是哈密顿算子。
9.1 热传导方程概述
二、产生热量 V 内所有热源产生的热量：
9.1 热传导方程概述
四、三维非齐次热传导方程 由能量守恒定律，即
Q Q Q
可得：
或
u t dt V [c t ( Ku) F ]dV u c ( Ku ) F ( x, y, z, t ) t
t
（9.1）
式（9.1）称为各向同性介质有热源的热传导方程，也叫三维 非齐次热传导方程。
t T
（9.6）
9.2 一维热传导方程的差分解法
边界条件： 3、 第三类边界条件：
u (, t ) x (t )u (, t ) g (t ) u (l , t ) (t )u (l , t ) g (t ) x
其中， (t )
i ,,..., N , j ,,...,M k ,,,..., j ,,...,M k ,,...,M , i ,,..., N