第九章_热传导方程的差分解法_郑大昉

u(x,0) = ϕ(x),
| x |< +∞
(9-12)
此处,并未给定介质的边界限制 此处 并未给定介质的边界限制. 并未给定介质的边界限制
初、边值混合问题 一方面,给定初始条件: 一方面 给定初始条件 给定初始条件
u(x,0) = ϕ(x),
第一类: u(0, t) = g1(t) 第一类
0≤ x ≤l
∂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)
时间的一阶差商近似为 近似为: 另, 对时间的一阶差商近似为
uN, j,k = uN−1, j,k
j =1 2,..., M −1 , k =1 2,... ,
(9-37)
(2)恒温条件(恒温热源 恒温条件 恒温热源 恒温条件 恒温热源) 低温热源为 令高温热源温度为 1 , 低温热源为 0. 高温热源温度为
u0, j,k =1 , ui,0,k = ui,M ,k = 0
λ, l, h,α,T;
二维情形及差分方法
内部无热源均匀介质中的二维热传导方程为: 内部无热源均匀介质中的二维热传导方程为 无热源均匀介质中的二维热传导方程为
∂u ∂2u ∂2u = λ( 2 + 2 ) ∂t ∂x ∂y
(0 < x < l, 0 < y < s, 0 < t < T)
初始条件为 初始条件为:
(9-15)
第三类: 第三类
∂u(0, t) ∂x − λ1(t)u(0, t) = g1(t) ∂u(l, t) − λ (t)u(l, t) = g (t) 2 2 ∂x
0 ≤ t ≤T
(9-16)
其中, 均为给定函数 给定函数. 其中 g1(t), g1(t), λ (t), λ2 (t) 均为给定函数 1 差分方法 令:
对上述二阶偏微分方程如何做数值计算 对上述二阶偏微分方程如何做数值计算? 数值计算
K 其中: 其中 λ = cρ
一维情形及差分方法
一维情形下,方程 写为: 一维情形下 方程 (9-10) 写为
∂u ∂2u =λ , λ > 0,0 < t ≤ T 2 ∂t ∂x
(9-11)
求解上述方程需要给定初始条件和边界条件 即定解条件. 求解上述方程需要给定初始条件和边界条件,即定解条件 初始条件 初值问题 时刻空间各点的温度分布 温度分布, 给定 t=0 时刻空间各点的温度分布 即:
k =1 2,... ,
x = Nh 0< y < M , h ∂uN, j,k 即: =0 ∂x
b.
j = 1 2,..., M −1 ,
(9-33)
k =1 2,... ,
(9-32), (9-33) 给出 给出:
u1, j,k −u0, j,k −u h
=0
(9-34)
k =1 2,... ,
边界条件(以一类为例): 边界条件 以一类为例 为例
(9-24)
k u0,k = g1(kτ ), uN,k = g2 (kτ ), = 0,1,..., M −1
此外, 计算值与精确解误差为 误差为: 此外 设计算值与精确解误差为 若取: 若取 则, 当
增加时 减少或不增加. k 增加时, εi,k将减少或不增加
最后,其计算步骤可列举如下: 最后 其计算步骤可列举如下 可列举如下 1. 给定 给定:
XN αh2 , M = YM , τ = , kmax = T 2. 计算 N = 计算: h τ λ h 3. 计算 计算: x = ih, y = jh, t = kτ
几何意义如图示 如图示: 的几何意义如图示 其计算可按如下简单步骤进行: 其计算可按如下简单步骤进行
1. 给定 给定:
l αh2 M = T ; 2. 计算 N = , τ = 计算: , h λ τ
3. 计算初值 ui,0 = ϕ(ih); 计算初值: 计算边界值: 计算边界值 u0,k = g1(kτ ), uN,k = g2 (kτ ); 4. 计算 ui,k +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+1 −ui,k ∂u = τ ∂t i,k
(9-20), (9-21) 代入 (9-11) 得:
(9-21)
ui,k +1 = αui+1,k + (1− 2α)ui,k +αui−1,k
(i =1 2,..., N −1. ,
其中: 其中
(9-22)
k = 0,1 M −1.) ,...,
∆x ≡ h, ∆t ≡τ , 空间标号为 空间标号为:
∂ui,k ∂x
+
时间标号为: i , 时间标号为 k
则一阶向前差商可近似为 一阶向前差商可近似为: 可近似为
=
ui+1,k − ui,k h
(9-17)
一阶向后差商可近似为 一阶向后差商可近似为: 可近似为
∂ui,k ∂x
−
=
ui,k − ui−1,k h
j =1 2,..., M1 −1 M2 +1 M −1 , , ,..
uN, j,k − uN−1, j,k h
j =1 2,..., M −1 ,
=0
k =1 2,... ,
(9-35)
或:
u0, j ,k = u1, j ,k
k =1 2,... ,
(9-36)
j =1 2,..., M1 −1 M2 +1 M −1 , , ,..
i
λ, h,α,T, XN ,YM
j
k
4. 计算初值及边界值 计算初值 初值及 5. 计算 计算:
T (i = 0,1 N, j = 0,1 M, k = 0,1 ) ,..., ,..., ,..., τ
T M = . τ
α=
上述
τλ
h2
l N = , h
(9-23)
[]
意为取整. 意为取整
显然, 递推求解还需要定解条件: 求解还需要定解条件 显然 对 (9-22) 做递推求解还需要定解条件
, 初始条件: 初始条件 ui,0 = ϕ(ih), i = 1 2,..., N −1
依能量守恒,即 依能量守恒 即: 得:
t2
Q2 = Q + Q3 1
(9-6) (9-7) (9-8)
∂u ∫t1 dt∫∫∫[cρ ∂t −∇⋅ (K∇u) − F]dV = 0 (V ) ∂u 从而: 从而 cρ = ∇⋅ (K∇u) + F(x, y, z, t) ∂t
可得相应的三维齐次方程 若无热源,可得相应的三维齐次方程 无热源 可得相应的三维齐次方程: 或:
(9-31)
其几何意义如图示 几何意义如图示: 如图示 边界条件 (1)绝热条件(绝热壁 绝热条件 绝热壁) 绝热条件 绝热壁 有关区域: 有关区域 a.
x =0
∂u0, j,k ∂x
0 < y < M1h, M2h < y < M h
即:
=0
(9-32)
j =1 2,..., M1 −1 M2 +1 M −1 , , ,..
第九章 热传导方程的差分解法
热传导及二阶偏微分方程
处的温度 温度为 若 t 时刻点 (x, y, z) 处的温度为 u(x, y, z, t) ，则 横截面为 传导的热量可表示为： 热量可表示为 时间内通过该处横截面 时间内通过该处横截面为 ∆S 传导的热量可表示为： 其中: K(x, y, z, t) > 0 ---称介质的 热传导系数. 其中 ---称介质的 热传导系数. ∂u ---为温度梯度沿截面法向分量 沿截面法向分量 为温度梯度沿截面法向分量. ∂n 若考虑一闭合 闭合曲面 （ 内的热量变化． 若考虑一闭合曲面 S 其体积为 V ）内的热量变化． 由 (9-1), 在
1 α= 2 ≤ h 2
τλ
εi,k
(9-25)
(9-26)
稳定差分格式. 此时 (9-22) 称稳定差分格式 数值计算时,应选参数满足(9-26) . 数值计算时 应选参数满足 满足
递推解 (9-22) 即:
ui,k+1 = αui+1,k + (1− 2α)ui,k +αui−1,k
(i =1 2,..., N −1. = 0,1 M −1.) , k ,...,
Q = ∫ dt ∫∫∫ ∇⋅[K(x, y, z, t)∇u]dV 1
t1 (V )ຫໍສະໝຸດ t2(9-3)哈密顿算子. 其中 ∇为哈密顿算子 温度随时间变化引起热量的消耗为 引起热量的消耗 另, V 内温度随时间变化引起热量的消耗为:
Q2 = ∫
其中
t2
t1
c和 ρ
∂u dt ∫∫∫ cρ dV ∂t (V )
(9-4)
分别为介质的比热容及密度 分别为介质的比热容及密度. 比热容