熱傳導方程式 (或稱熱方程式 )是一個重要的偏微分方程式 ,它描述一個區域內的溫度 如何隨時間 變化。
一維熱方程式圖解(觀看動畫版 )
熱傳導在三維的各向同性 介質裡的傳播可用以下方程式表達:
∂
u
∂
t
=
d
i
v
(
U
u
)
=
k
(
∂
2
u
∂
x
2
+
∂
2
u
∂
y
2
+
∂
2
u
∂
z
2
)
=
k
(
u
x
x
+
u
y
y
+
u
z
z
)
{\displaystyle {\partial u \over \partial t}=\mathrm {div} (Uu)=k\left({\partial ^{2}u \over \partial x^{2}}+{\partial ^{2}u \over \partial y^{2}}+{\partial ^{2}u \over \partial z^{2}}\right)=k(u_{xx}+u_{yy}+u_{zz})\quad }
其中:
u =u (t , x , y , z )表溫度,它是時間變數t與空間變數(x,y,z)的函數。
∂
u
{\displaystyle {\partial u}}
/
∂
t
{\displaystyle {\partial t}}
是空間中一點的溫度對時間的變化率。
u
x
x
{\displaystyle u_{xx}}
,
u
y
y
{\displaystyle u_{yy}}
與
u
z
z
{\displaystyle u_{zz}}
溫度對三個空間座標軸的二次導數。
k 是熱擴散率 ,決定於材料的熱導率 、密度 與熱容 。
熱方程式是傅立葉冷卻律的一個推論(詳見條目熱傳導 )。
如果考慮的介質不是整個空間,則為了得到方程式的唯一解,必須指定u的邊界條件 。如果介質是整個空間,為了得到唯一性,必須假定解的增長速度有個指數型的上界,此假定吻合實驗結果。
熱方程式的解具有將初始溫度平滑化的特質,這代表熱從高溫處向低溫處傳播。一般而言,許多不同的初始狀態會趨向同一個穩態(熱平衡 )。因此我們很難從現存的熱分佈反解初始狀態,即使對極短的時間間隔也一樣。
熱方程式也是拋物線偏微分方程式 最簡單的例子。
利用拉普拉斯算子 ,熱方程式可推廣為下述形式
u
t
=
k
Δ
u
,
{\displaystyle u_{t}=k\Delta u,\quad }
其中的
Δ
{\displaystyle \Delta }
是對空間變數的拉普拉斯算子。
熱方程式支配熱傳導及其它擴散過程,諸如粒子擴散或神經細胞的動作電位 。
就技術上來說,熱方程式違背狹義相對論 ,因為它的解表達一個「擾動可以在瞬間傳播至空間各處」的情況。擾動在前方光錐 外的影響通常可忽略不計,但是若要為熱傳導推出一個合理的速度,則須轉而考慮一個雙曲線型偏微分方程式。
在理想狀態下一根棍子的熱傳導,配上均勻的邊界條件。
以下解法首先由約瑟夫·傅立葉 在他於1822年出版的著作Théorie analytique de la chaleur (中譯:解析熱學)給出。先考慮只有一個空間變數的熱方程式,這可以當作棍子的熱傳導之模型。方程式如下:
(
1
)
u
t
=
k
u
x
x
{\displaystyle (1)\ u_{t}=ku_{xx}\quad }
其中u = u (t , x )是t 和x 的雙變數函數。
x 是空間變數,所以x ∈ [0,L ],其中L 表示棍子長度。
t 是時間變數,所以t ≥ 0。
假設下述初始條件
(
2
)
u
(
0
,
x
)
=
f
(
x
)
∀
x
∈
[
0
,
L
]
{\displaystyle (2)\ u(0,x)=f(x)\quad \forall x\in [0,L]\quad }
其中函數f 是給定的。再配合下述邊界條件
(
3
)
u
(
t
,
0
)
=
0
=
u
(
t
,
L
)
∀
t
>
0
{\displaystyle (3)\ u(t,0)=0=u(t,L)\quad \forall t>0\quad }
.
讓我們試著找一個非恆等於零的解,使之滿足邊界條件(3)並具備以下形式:
(
4
)
u
(
t
,
x
)
=
X
(
x
)
T
(
t
)
.
{\displaystyle (4)\ u(t,x)=X(x)T(t).\quad }
這套技術稱作分離變數法 。現在將u 代回方程式(1),
T
′
(
t
)
k
T
(
t
)
=
X
″
(
x
)
X
(
x
)
.
{\displaystyle {\frac {T'(t)}{kT(t)}}={\frac {X''(x)}{X(x)}}.\quad }
由於等式右邊只依賴x ,而左邊只依賴t ,兩邊都等於某個常數− λ,於是:
(
5
)
T
′
(
t
)
=
−
λ
k
T
(
t
)
{\displaystyle (5)\ T'(t)=-\lambda kT(t)\quad }
(
6
)
X
″
(
x
)
=
−
λ
X
(
x
)
.
{\displaystyle (6)\ X''(x)=-\lambda X(x).\quad }
以下將證明(6)沒有λ ≤ 0的解:
假設λ < 0,則存在實數B 、C 使得
X
(
x
)
=
B
e
−
λ
x
+
C
e
−
−
λ
x
{\displaystyle X(x)=Be^{{\sqrt {-\lambda }}\,x}+Ce^{-{\sqrt {-\lambda }}\,x}}
。
從(3)得到
X
(
0
)
=
0
=
X
(
L
)
.
{\displaystyle X(0)=0=X(L).\quad }
於是有B = 0 = C ,這蘊含u 恆等於零。
假設λ = 0,則存在實數B 、C 使得
X
(
x
)
=
B
x
+
C
.
{\displaystyle X(x)=Bx+C.\quad }
仿上述辦法可從等式(3)推出u 恆等於零。
因此必然有λ > 0,此時存在實數A 、B 、C 使得
T
(
t
)
=
A
e
−
λ
k
t
{\displaystyle T(t)=Ae^{-\lambda kt}\quad }
X
(
x
)
=
B
sin
(
λ
x
)
+
C
cos
(
λ
x
)
{\displaystyle X(x)=B\sin({\sqrt {\lambda }}\,x)+C\cos({\sqrt {\lambda }}\,x)}
。
從等式(3)可知C = 0,因此存在正整數n 使得
λ
=
n
π
L
{\displaystyle {\sqrt {\lambda }}=n{\frac {\pi }{L}}}
。
由此得到熱方程式形如(4)的解。
一般而言,滿足(1)與(3)的解相加後仍是滿足(1)與(3)的解。事實上可以證明滿足(1)、(2)、(3)的解由下述公式給出:
u
(
t
,
x
)
=
∑
n
=
1
+
∞
D
n
(
sin
n
π
x
L
)
e
−
n
2
π
2
k
t
L
2
{\displaystyle u(t,x)=\sum _{n=1}^{+\infty }D_{n}\left(\sin {\frac {n\pi x}{L}}\right)e^{-{\frac {n^{2}\pi ^{2}kt}{L^{2}}}}}
其中
D
n
=
2
L
∫
0
L
f
(
x
)
sin
n
π
x
L
d
x
{\displaystyle D_{n}={\frac {2}{L}}\int _{0}^{L}f(x)\sin {\frac {n\pi x}{L}}\,dx}
。
上面採用的方法可以推廣到許多不同方程式。想法是:在適當的函數空間 上,算子
u
↦
u
x
x
{\displaystyle u\mapsto u_{xx}}
可以用它的特徵向量 表示。這就自然地導向線性自伴算子 的譜理論 。
考慮線性算子 Δ u = u x x ,以下函數序列
e
n
(
x
)
=
2
L
sin
n
π
x
L
{\displaystyle e_{n}(x)={\sqrt {\frac {2}{L}}}\sin {\frac {n\pi x}{L}}}
(n ≥ 1)
是Δ的特徵向量。誠然:
Δ
e
n
=
−
n
2
π
2
L
2
e
n
{\displaystyle \Delta e_{n}=-{\frac {n^{2}\pi ^{2}}{L^{2}}}e_{n}}
。
此外,任何滿足邊界條件f (0)=f (L )=0的Δ的特徵向量都是某個e n 。令L2 (0, L )表 [0, L ]上全體平方可積函數 的向量空間 。這些函數e n 構成L2 (0, L )的一組正交歸一基 。更明白地說:
⟨
e
n
,
e
m
⟩
=
∫
0
L
e
n
(
x
)
e
m
(
x
)
d
x
=
{
0
n
≠
m
1
m
=
n
{\displaystyle \langle e_{n},e_{m}\rangle =\int _{0}^{L}e_{n}(x)e_{m}(x)dx=\left\{{\begin{matrix}0&n\neq m\\1&m=n\end{matrix}}\right.}
最後,序列{e n }n ∈ N 張出L2 (0, L )的一個稠密的線性子空間 。這就表明我們實際上已將算子Δ 對角化 。
一般而言,熱傳導的研究奠基於以下幾個原理。首先注意到熱流是能量 流的一種形式,因此可以談論單位時間內流進空間中一塊區域的熱量。
單位時間內流入區域 V 的熱量由一個依賴於時間的量q t (V )給出。假設q 有個密度 Q(t,x) ,於是
q
t
(
V
)
=
∫
V
Q
(
t
,
x
)
d
x
{\displaystyle q_{t}(V)=\int _{V}Q(t,x)\,dx\quad }
熱流是個依賴於時間的向量函數H (x ),其刻劃如下:單位時間內流經一個面積為dS 而單位法向量為n 的無窮小曲面元素的熱量是
H
(
x
)
⋅
n
(
x
)
d
S
{\displaystyle \mathbf {H} (x)\cdot \mathbf {n} (x)\,dS}
因此單位時間內進入V 的熱流量也由以下的面積分給出
q
t
(
V
)
=
−
∫
∂
V
H
(
x
)
⋅
n
(
x
)
d
S
{\displaystyle q_{t}(V)=-\int _{\partial V}\mathbf {H} (x)\cdot \mathbf {n} (x)\,dS}
其中n (x)是在x 點的向外單位法向量。
H
(
x
)
=
−
A
(
x
)
⋅
∇
u
(
x
)
{\displaystyle \mathbf {H} (x)=-\mathbf {A} (x)\cdot \nabla u(x)}
其中A (x )是個3×3實對稱正定矩陣 。
利用格林定理 可將之前的面積分轉成一個體積分
q
t
(
V
)
=
−
∫
∂
V
H
(
x
)
⋅
n
(
x
)
d
S
{\displaystyle q_{t}(V)=-\int _{\partial V}\mathbf {H} (x)\cdot \mathbf {n} (x)\,dS}
=
∫
∂
V
A
(
x
)
⋅
∇
u
(
x
)
⋅
n
(
x
)
d
S
{\displaystyle =\int _{\partial V}\mathbf {A} (x)\cdot \nabla u(x)\cdot \mathbf {n} (x)\,dS}
=
∫
V
∑
i
,
j
∂
x
i
a
i
j
(
x
)
∂
x
j
u
(
t
,
x
)
d
x
{\displaystyle =\int _{V}\sum _{i,j}\partial _{x_{i}}a_{ij}(x)\partial _{x_{j}}u(t,x)\,dx}
溫度在x 點對時間的改變率與流進x 點所在的無窮小區域的熱量成正比,此比例常數與時間無關,而可能與空間有關,寫作κ(x)。
∂
t
u
(
t
,
x
)
=
κ
(
x
)
Q
(
t
,
x
)
{\displaystyle \partial _{t}u(t,x)=\kappa (x)Q(t,x)\,}
將以上所有等式合併,便獲得支配熱流的一般公式。
∂
t
u
(
t
,
x
)
=
κ
(
x
)
∑
i
,
j
∂
x
i
a
i
j
(
x
)
∂
x
j
u
(
t
,
x
)
{\displaystyle \partial _{t}u(t,x)=\kappa (x)\sum _{i,j}\partial _{x_{i}}a_{ij}(x)\partial _{x_{j}}u(t,x)}
註記 :
係數κ(x )是該材料在x 點的密度 和比熱 的積的倒數。
在等方向性介質的情況,矩陣A 只是個純量,等於材料的導熱率。
在非等向的情況,A 不一定是純量,我們鮮少能明確寫出熱方程式的解。然而通常可考慮相應的抽象柯西問題 ,證明它是適定 的,並(或)導出若干定性結果(諸如初始值保持正性、無窮傳播速度、收斂至平衡態或一些平滑化性質)。這些論證通常有賴於單參數半群 理論:舉例來說,如果A 是個對稱矩陣,那麼由
A
u
(
x
)
:=
∑
i
,
j
∂
x
i
a
i
j
(
x
)
∂
x
j
u
(
x
)
{\displaystyle Au(x):=\sum _{i,j}\partial _{x_{i}}a_{ij}(x)\partial _{x_{j}}u(x)}
定義的橢圓算子 是自伴而且耗散的,因此由譜定理 導出它生成一個單參數半群。
在粒子擴散的模型中,我們考慮的方程式涉及
在大量粒子集體擴散的情況:粒子的體積濃度 ,記作c 。或者
在單一粒子的情況:單一粒子對位置的機率密度函數 ,記作P 。
不同情況下的方程式:
c
t
=
D
Δ
c
,
{\displaystyle c_{t}=D\Delta c,\quad }
或者
P
t
=
D
Δ
P
.
{\displaystyle P_{t}=D\Delta P.\quad }
c 與P 都是位置與時間的函數。D 是擴散係數,它控制擴散速度,通常以公尺/秒為單位。
如果擴散係數D 依賴於濃度c (或第二種情況下的機率密度P ),則我們得到非線性擴散方程式 。
單一粒子在粒子擴散方程式下的隨機軌跡是個布朗運動 。
如果一個粒子在時間
t
=
0
{\displaystyle t=0}
時置於
R
→
=
0
→
{\displaystyle {\vec {R}}={\vec {0}}}
,則相應的機率密度函數具有以下形式:
P
(
R
→
,
t
)
=
G
(
R
→
,
t
)
=
1
(
4
π
D
t
)
3
/
2
e
−
R
→
2
4
D
t
{\displaystyle P({\vec {R}},t)=G({\vec {R}},t)={\frac {1}{(4\pi Dt)^{3/2}}}e^{-{\frac {{\vec {R}}^{2}}{4Dt}}}}
它與機率密度函數的各分量
R
x
{\displaystyle R_{x}}
、
R
y
{\displaystyle R_{y}}
和
R
z
{\displaystyle R_{z}}
的關係是:
P
(
R
→
,
t
)
=
1
(
4
π
D
t
)
3
/
2
e
−
R
x
2
+
R
y
2
+
R
z
2
4
D
t
=
P
(
R
x
,
t
)
P
(
R
y
,
t
)
P
(
R
z
,
t
)
{\displaystyle P({\vec {R}},t)={\frac {1}{(4\pi Dt)^{3/2}}}e^{-{\frac {R_{x}^{2}+R_{y}^{2}+R_{z}^{2}}{4Dt}}}=P(R_{x},t)P(R_{y},t)P(R_{z},t)}
隨機變數
R
x
,
R
y
,
R
z
{\displaystyle R_{x},R_{y},R_{z}}
服從平均數為0、變異數為
2
D
t
{\displaystyle 2\,D\,t}
的正態分佈 。在三維的情形,隨機向量
R
→
{\displaystyle {\vec {R}}}
服從平均數為
0
→
{\displaystyle {\vec {0}}}
、變異數為
6
D
t
{\displaystyle 6\,D\,t}
的正態分佈。
在t=0 時,上述
P
(
R
→
,
t
)
{\displaystyle P({\vec {R}},t)}
的表示式帶有奇異點。對應於粒子處在原點之初始條件,其機率密度函數是在原點的狄拉克δ函數 ,記為
δ
(
R
→
)
{\displaystyle \delta ({\vec {R}})}
(三維的推廣是
δ
(
R
→
)
=
δ
(
R
x
)
δ
(
R
y
)
δ
(
R
z
)
{\displaystyle \delta ({\vec {R}})=\delta (R_{x})\delta (R_{y})\delta (R_{z})}
);擴散方程式對此初始值的解也稱作格林函數 。
粒子擴散方程式 首先由Adolf Fick於1855年導得。
格林函數是擴散方程式在粒子位置已知時的解(數學家稱之為擴散方程式的基本解 )。當粒子初始位置在原點
0
→
{\displaystyle {\vec {0}}}
時,相應的格林函數記作
G
(
R
→
,
t
)
{\displaystyle G({\vec {R}},t)}
(t>0 );根據擴散方程式對平移的對稱性,對一般的已知初始位置
R
→
0
{\displaystyle {\vec {R}}^{0}}
,相應的格林函數是
G
(
R
→
−
R
→
0
,
t
)
{\displaystyle G({\vec {R}}-{\vec {R}}^{0},t)}
。
對於一般的初始條件,擴散方程式的解可以透過積分分解為一族格林函數的疊加 。
舉例來說,設t=0 時有一大群粒子,根據濃度分佈的初始值
c
(
R
→
,
0
)
{\displaystyle c({\vec {R}},0)}
分佈於空間中。擴散方程式的解將告訴我們濃度分佈如何隨時間演化。
跟任何(廣義)函數一樣,濃度分佈的初始值可以透過積分表為狄拉克δ函數的疊加:
c
(
R
→
,
t
=
0
)
=
∫
c
(
R
→
0
,
t
=
0
)
δ
(
R
→
−
R
→
0
)
d
R
x
0
d
R
y
0
d
R
z
0
{\displaystyle c({\vec {R}},t=0)=\int c({\vec {R}}^{0},t=0)\delta ({\vec {R}}-{\vec {R}}^{0})dR_{x}^{0}\,dR_{y}^{0}\,dR_{z}^{0}}
擴散方程式是線性的,因此在之後的任一時刻t ,濃度分佈變為:
c
(
R
→
,
t
)
=
∫
c
(
R
→
0
,
t
=
0
)
G
(
R
→
−
R
→
0
,
t
)
d
R
x
0
d
R
y
0
d
R
z
0
{\displaystyle c({\vec {R}},t)=\int c({\vec {R}}^{0},t=0)G({\vec {R}}-{\vec {R}}^{0},t)dR_{x}^{0}\,dR_{y}^{0}\,dR_{z}^{0}}
在粒子擴散的情形,我們可以將狄拉克δ函數對應的初始條件理解為粒子落在一個已知位置。一般而言,任何擴散過程的解都有這種表法,包括熱傳導或動量 的擴散;後者關係到流體的黏性 現象。
熱傳導方程式在許多隨機過程 的數學模型 中出現,諸如布萊克-斯科爾斯模型 與奧恩斯坦-烏倫貝克過程 。在金融數學 ,1973年發表的布萊克-斯科爾斯模型 作為期權 定價模型,當中的差分方程式 可以轉成熱傳導方程式,並從此導出較簡單的解 。許多簡單期權的延伸模型沒有解析解 ,因此必須以數值方法 計算模型給出的定價。熱傳導方程式可以用克蘭克-尼科爾森方法 有效地求數值解,此方法也可用於許多無解析解的模型(詳見文獻Wilmott,1995)。
熱傳導方程式及其非線性的推廣型式也被應用於圖像分析 。
量子力學 中的薛丁格方程式 雖然有類似熱傳導方程式的數學公式 (但時間參數為
i
{\displaystyle i}
是虛數單位 ),本質卻不是擴散問題 ,解的定性行為也完全不同。
熱傳導方程式在流形 上的推廣 是處理阿蒂亞-辛格指標定理 的主要工具之一,由此也導向熱傳導方程式在黎曼幾何 中有許多應用。