作者: ls_本性
专题简介：本文章将从零开始学习并完成2018年全国大学生数学建模的A题，由于本人也是新手第一次接触偏微分方程题目（可能会有很多错误，希望大家指出）
😊😊😊如果觉得文章不错或能帮助到你学习，可以点赞👍收藏📁评论📒+关注哦！👍👍
本文末尾有惊喜哦！
高温作业专用服装设计
背景知识
热传导方程的导出（控制方程）
热传导方程是基于傅里叶定律和能量守恒定律推导
傅里叶定律
热流密度（q）：单位时候内通过截面面积S的截面所传递的热量Q。
傅里叶定律：单位时候内通过截面面积S的截面所传递的热量Q正比于垂直该截面方向上的温度变化率.
温度： u ( x , y , z , t ) u(x,y,z,t) u(x,y,z,t)表示一点 ( x , y , z ) (x,y,z) (x,y,z)在t时刻的温度
对于x,y,z方向傅里叶定律：
q
x
=
−
λ
∂
u
∂
x
{q_x} = - \lambda \frac{{\partial u}}{{\partial x}}
qx=−λ∂x∂u
q
y
=
−
λ
∂
u
∂
y
{q_y} = - \lambda \frac{{\partial u}}{{\partial y}}
qy=−λ∂y∂u
q
z
=
−
λ
∂
u
∂
z
{q_z} = - \lambda \frac{{\partial u}}{{\partial z}}
qz=−λ∂z∂u
对于整体：
q
→
=
−
λ
∇
\overrightarrow q = - \lambda \nabla
q=−λ∇
其中
λ
\lambda
λ表示热传导率,
∇
\nabla
∇为散度
我们假定：
- 该物体三维各项同性的均值介质（ λ \lambda λ都相同）
- 温度变化不剧烈（温度不会影响 λ \lambda λ）
微元法推导方程
我们选取物题的极小一部分进行分析。
沿x方向流过的单位时间热量：
(
q
x
∣
x
−
q
x
∣
x
+
Δ
x
)
Δ
y
Δ
ｚ
({{q}_{x}}{{|}_{x}}-{{q}_{x}}{{|}_{x+\Delta x}})\Delta y\Deltaｚ
(qx∣x−qx∣x+Δx)ΔyΔｚ
该式子表示单位时间在
x
x
x方向流入的热量减去在
x
+
Δ
x
x+\Delta x
x+Δx流出的热量
对该公式变形转化成偏导格式：
(
q
x
∣
x
−
q
x
∣
x
+
Δ
x
)
Δ
y
Δ
z
Δ
x
Δ
x
\frac{({{q}_{x}}{{|}_{x}}-{{q}_{x}}{{|}_{x+\Delta x}})\Delta y\Delta z\Delta x}{\Delta x}
Δx(qx∣x−qx∣x+Δx)ΔyΔzΔx
⟶
\longrightarrow
⟶
−
∂
q
x
∂
x
Δ
x
Δ
y
Δ
z
-\frac{\partial {{q}_{x}}}{\partial x}\Delta x\Delta y\Delta z
−∂x∂qxΔxΔyΔz
对于y,z方向同理，得到流过单位时间整个物体的热量：
−
(
∂
q
x
∂
x
+
∂
q
y
∂
y
+
∂
q
z
∂
z
)
Δ
x
Δ
y
Δ
z
=
−
∇
⋅
q
→
Δ
x
Δ
y
Δ
z
-(\frac{\partial {{q}_{x}}}{\partial x}+\frac{\partial {{q}_{y}}}{\partial y}+\frac{\partial {{q}_{z}}}{\partial z})\Delta x\Delta y\Delta z=-\nabla \cdot \overset{\to }{\mathop{q}}\,\Delta x\Delta y\Delta z
−(∂x∂qx+∂y∂qy+∂z∂qz)ΔxΔyΔz=−∇⋅q→ΔxΔyΔz
其中
∇
⋅
q
→
\nabla \cdot \overset{\to }{\mathop{q}}\,
∇⋅q→表示
q
→
\overset{\to }{\mathop{q}}\,
q→的散度。
定义 F ( x , y , z , t ) F(x,y,z,t) F(x,y,z,t)为单位时间，单位体积，物体自身的产热量。
由能量守恒定律得：
−
∇
⋅
q
→
Δ
x
Δ
y
Δ
z
Δ
t
+
F
Δ
x
Δ
y
Δ
z
Δ
t
=
c
ρ
Δ
x
Δ
y
Δ
z
Δ
u
-\nabla \cdot \overset{\to }{\mathop{q}}\,\Delta x\Delta y\Delta z\Delta t+F\Delta x\Delta y\Delta z\Delta t=c\rho \Delta x\Delta y\Delta z\Delta u
−∇⋅q→ΔxΔyΔzΔt+FΔxΔyΔzΔt=cρΔxΔyΔzΔu
其中
c
c
c为比热容，
ρ
\rho
ρ为密度,
Δ
u
\Delta u
Δu为温度改变量
对等式右边的式子乘以
Δ
t
\Delta t
Δt除以
Δ
t
\Delta t
Δt，然后对整个式子约分得：
−
∇
⋅
q
→
+
F
=
c
ρ
∂
u
∂
t
-\nabla \cdot \overset{\to }{\mathop{q}}\,+F=c\rho \frac{\partial u}{\partial t}
−∇⋅q→+F=cρ∂t∂u
其傅里叶定律得到的方程代入得：
λ
∇
2
+
F
=
c
ρ
∂
u
∂
t
\lambda {{\nabla }^{2}}+F=c\rho \frac{\partial u}{\partial t}
λ∇2+F=cρ∂t∂u
其中
∇
2
{{\nabla }^{2}}
∇2表示拉普拉多算子
最后对该等式变形得到最终的热传导控制方程：
∂
u
∂
t
−
k
∇
2
=
f
\frac{\partial u}{\partial t}-k{{\nabla }^{2}}=\text{f}
∂t∂u−k∇2=f
其中
k
=
λ
c
ρ
k=\frac{\lambda}{c\rho}
k=cρλ表示为扩散系数,
f
=
F
c
ρ
f=\frac{F}{c\rho}
f=cρF
一般来说若物理内部不会产生，则该方程退化为：
k
∇
2
=
∂
u
∂
t
k{{\nabla }^{2}}=\frac{\partial u}{\partial t}
k∇2=∂t∂u
直观理解为，温度的改变量与其内部某点位置的四面八方都有关（这里可能解释的有些模糊，后面会对一维热传导直观理解用差分进一步的解释）
定解条件
对于偏微分方程只知道其控制方程的话该方程存在无穷多解，要得到其唯一解，还需要对其加上定解条件。
初始条件
相比常微分方程其初始条件是一个具体的数值而言，我们需要注意偏微分方程其初始条件应该是一个函数。
eg: 方程 ∂ u ( x , y ) ∂ 2 x = 0 \frac{\partial u(x,y)}{{{\partial }^{2}}x}=0 ∂2x∂u(x,y)=0其通解之一为： u = c 1 ( y ) + x c 2 ( y ) u={{c}_{1}}(y)+x{{c}_{2}}(y) u=c1(y)+xc2(y)，其中 c 1 ( y ) {c}_{1}(y) c1(y), c 2 ( y ) {c}_{2}(y) c2(y)为两个关于 y y y的函数。
边界条件
- 第一类边界条件又称狄利克雷条件，规定了所研究物理量在边界上的数值，对于热传导方程而言其物理意义为：规定了物体边界上的温度值。
u ( x , y , z , t ) ∣ x 0 , y 0 , z 0 = f ( x 0 , y 0 , z 0 , t ) u(x,y,z,t){{|}_{{{x}_{0}},{{y}_{0}},{{z}_{0}}}}=f({{x}_{0}},{{y}_{0}},{{z}_{0,}}t) u(x,y,z,t)∣x0,y0,z0=f(x0,y0,z0,t)其中 ( x 0 , y 0 , z 0 ) ({x}_{0},{y}_{0},{z}_{0}) (x0,y0,z0)是边界点的坐标， f f f是时间 t t t的已知函数。 - 第二类边界条件又称诺依曼条件，规定了所研究物理量在边界外法线方向上方向导数的数值，对于热传导方程而言其物理意义为：规定了物体边界上的热流密度(单位时间单位面积的热流量)。
谈到热流密度我们就应该想到傅里叶定律，在傅里叶定律中热流密度的方向是由高温流向低温，而并非物理量的法向方向。
定义温度变化的梯度方向为： ∂ u ∂ n = ∇ n → \frac{\partial u}{\partial n}=\nabla \overrightarrow{n} ∂n∂u=∇n
其中
n
→
\overrightarrow{n}
n为物体边界的法向向量，
∇
\nabla
∇为散度
在结合傅里叶定律，可得：
−
λ
∂
u
∂
n
∣
x
0
,
y
0
,
z
0
=
ψ
(
x
0
,
y
0
,
z
0
,
t
)
-\lambda\frac{\partial u}{\partial n}{{|}_{{{x}_{0}},{{y}_{0}},{{z}_{0}}}}=\psi ({{x}_{0}},{{y}_{0}},{{z}_{0,}}t)
−λ∂n∂u∣x0,y0,z0=ψ(x0,y0,z0,t)
其中
(
x
0
,
y
0
,
z
0
)
({x}_{0},{y}_{0},{z}_{0})
(x0,y0,z0)是边界点的坐标，
f
f
f是时间
t
t
t的已知函数。
- 第三类边界条件又称混合边界条件，规定了所研究物理量及其外法向导数的线性组合在边界上的数值。对于热传导方程第三类边界条件又称牛顿冷却定律，其表示：单位时间单位面积的散热量与物体边界温度和外界温度只差成正比。
− λ ∂ u ∂ n ∣ x 0 , y 0 , z 0 = h ( u ∣ x 0 , y 0 , z 0 − u 外界 ) -\lambda\frac{\partial u}{\partial n}{{|}_{{{x}_{0}},{{y}_{0}},{{z}_{0}}}}=h(u{{|}_{{{x}_{0}},{{y}_{0}},{{z}_{0}}}}-{{u}_{外界}}) −λ∂n∂u∣x0,y0,z0=h(u∣x0,y0,z0−u外界)
其中 h h h为转化系数
模型建立
我们假设防护服每一层：三维各项同性的均值介质。那么我们便可以把这个三维问题抽象成一维问题来解决。
这里给出其一维简化图，由于四层材料，厚度不同，故我们需要对每二层进行一一讨论，这里我们在将其简化为二层材料，来建立二层耦合介质温度分布模型，再将其推广到四层。
二层耦合介质温度分布模型
控制方程
对于介质a,b，其满足热传导的一般方程：
耦合条件（a,b接触的边界条件）
- 由于a,b介质边界处于同一处，故在稳定情况下（温度没有发生急剧变化时）认为二介质边界温度相同：
u 1 ∣ d 1 = u 2 ∣ d 1 {{u}_{1}}{{|}_{d}}_{_{1}}={{u}_{2}}{{|}_{d}}_{_{1}} u1∣d1=u2∣d1 - 边界的热流密度（单位时间单位面积热流量）相同:
λ 1 ∂ u ∂ n 12 ∣ d 1 = λ 2 ∂ u ∂ n 21 ∣ d 1 {{\lambda }_{1}}\frac{\partial u}{\partial {{n}_{12}}}{{|}_{d1}}={{\lambda }_{2}}\frac{\partial u}{\partial {{n}_{21}}}{{|}_{d1}} λ1∂n12∂u∣d1=λ2∂n21∂u∣d1
其中 n 12 {n}_{12} n12表示由介质a的法向方向（由介质a指向外部）， n 21 {n}_{21} n21表示由介质b的法向方向（由外部指向介质b）
初始条件
由于人穿着防护服，故防护服在
t
=
0
t=0
t=0时刻的温度近似等于人体温度,即：
{
u
1
(
x
,
0
)
=
37
x
∈
[
0
,
d
1
]
u
2
(
x
,
0
)
=
37
x
∈
[
d
1
,
d
2
]
\left\{ \begin{align} & {{u}_{1}}(x,0)=37 x\in [0,{{d}_{1}}] \nonumber \\ & {{u}_{2}}(x,0)=37 x\in [{{d}_{1}},{{d}_{2}}] \nonumber \end{align} \right.
{u1(x,0)=37 x∈[0,d1]u2(x,0)=37 x∈[d1,d2]
边界条件
-
左边界条件：左侧边界直接与外界接触，故其温度与外界恒定温度相同:
u 1 ( 0 , t ) = 75 {{u}_{1}}(0,t)=75 u1(0,t)=75 t ∈ [ 0 , T ] t\in [0,T] t∈[0,T] -
右边界条件：介质b右侧边界温度高于人体温度，故与人皮肤发生热量交换，满足牛顿冷却定律
− λ 2 ∂ u ∂ n 22 = h ( u 2 − u 人 ) -{{\lambda }_{2}}\frac{\partial u}{\partial {{n}_{22}}}=h({{u}_{2}}-{{u}_{人}}) −λ2∂n22∂u=h(u2−u人)
由于此处空气对流很弱，故不考虑热对流。 -
h h h的确定：查文献可知h的范围是 [ 5 , 25 ] [5,25] [5,25]。然后采用变步长法去枚举最佳的转换系数
整体模型建立
控制方程
{
∂
u
i
∂
t
=
a
i
2
∂
2
u
i
∂
x
2
x
∈
⋃
i
=
1
4
[
L
i
−
1
,
L
i
]
,
L
0
=
0
a
i
2
=
λ
i
c
i
ρ
i
i
=
1
,
2
,
3
,
4
\left\{ \begin{align} & \frac{\partial {{u}_{i}}}{\partial t}={{a}_{i}}^{2}\frac{{{\partial }^{2}}{{u}_{i}}}{\partial {{x}^{2}}}x\in \bigcup\limits_{i=1}^{4}{[{{L}_{i-1}},{{L}_{i}}],{{L}_{0}}=0}\nonumber \\ & {{a}_{i}}^{2}=\frac{{{\lambda }_{i}}}{{{c}_{i}}{{\rho }_{i}}} i=1,2,3,4 \nonumber\\ \end{align} \right.
⎩
⎨
⎧∂t∂ui=ai2∂x2∂2uix∈i=1⋃4[Li−1,Li],L0=0ai2=ciρiλi i=1,2,3,4
对四个不同介质建立不同的控制方程
耦合条件
- 对一二，二三，三四边界建立耦合方程：
{ λ i ∂ u i ∂ x ∣ Γ i = λ i + 1 ∂ u i + 1 ∂ x ∣ Γ i + 1 i = 1 , 2 , 3 u i ∣ Γ i = u i + 1 ∣ Γ i + 1 i = 1 , 2 , 3 \left\{ \begin{align} & {{\lambda }_{i}}\frac{\partial {{u}_{i}}}{\partial x}{{|}_{\Gamma i}}={{\lambda }_{i+1}}\frac{\partial {{u}_{i+1}}}{\partial x}{{|}_{{{\Gamma }_{i+1}}}} i=1,2,3\nonumber \\ & {{u}_{i}}{{|}_{\Gamma i}}={{u}_{i+1}}{{|}_{\Gamma i+1}} i=1,2,3\nonumber \\ \end{align} \right. ⎩ ⎨ ⎧λi∂x∂ui∣Γi=λi+1∂x∂ui+1∣Γi+1 i=1,2,3ui∣Γi=ui+1∣Γi+1 i=1,2,3
此处为什么是对 x x x偏导而不是对 n n n偏导，这里我们把边界处近似的当成一条指线，且我们在这里只考虑标量不考虑矢量，故将对 n n n的偏导退化为对 x x x的偏导
初值条件
u
i
(
x
,
0
)
=
37
i
=
1
,
2
,
3
,
4
{{u}_{i}}(x,0)=37 i=1,2,3,4
ui(x,0)=37 i=1,2,3,4
防护服初始温度37度
边界条件
-
左边界条件：左侧边界直接与外界接触，故其温度与外界恒定温度相同:
u 1 ( 0 , t ) = 75 {{u}_{1}}(0,t)=75 u1(0,t)=75 t ∈ [ 0 , T ] t\in [0,T] t∈[0,T] -
右边界条件：介质b右侧边界温度高于人体温度，故与人皮肤发生热量交换，满足牛顿冷却定律
− λ 2 ∂ u ∂ n 22 = h ( u 4 ( L 4 , t ) − u 人 ( L 4 , t ) ) -{{\lambda }_{2}}\frac{\partial u}{\partial {{n}_{22}}}=h({{u}_{4}}({L}_{4},t)-{{u}_{人}({L}_{4},t)}) −λ2∂n22∂u=h(u4(L4,t)−u人(L4,t))
1193
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
