距離が 2 倍なら時間は 4 倍？熱方程式と熱核でわかる熱の広がり

1 点に集めた熱は、時刻 $t$ には $\frac{1}{\sqrt{4\pi kt}}{e}^{-x^2/(4kt)}$ の形で広がっています。山の幅は $\sqrt{4kt}$ で、時間の平方根でしか伸びません。

距離を 2 倍にすると、届くまでの時間は 4 倍。銅の棒なら 1 cm に 1 秒ほど、10 cm には 90 秒ほどかかります^[4]。

この $\sqrt{t}$ が、熱方程式のほとんどすべてを決めています。

方程式が出てくるところ

棒の微小な区間 $[x,x+\Delta x]$ を考えます。ここに蓄えられた熱量は、比熱 $c$、密度 $\rho$、断面積 $1$ として $c\rho u\Delta x$。

蓄えの変化は、左端から入る量と右端から出る量の差です。熱の流れる量を $q(x,t)$ と書きます。

$$c\rho \frac{\partial u}{\partial t}\Delta x=q(x,t)-q(x+\Delta x,t)$$

右辺は $-\frac{\partial q}{\partial x}\Delta x$ とほぼ等しい。$\Delta x$ で割れば保存則が出ます。

$$c\rho \frac{\partial u}{\partial t}=-\frac{\partial q}{\partial x}$$

ここまでは熱の出入りを数えただけ。もう 1 本、$q$ と $u$ を結ぶ関係が要ります。

フーリエの法則

熱は温度の高いほうから低いほうへ、勾配に比例して流れます。比例定数を熱伝導率 $K$ とする。

$$q=-K\frac{\partial u}{\partial x}$$

マイナスは向きの約束です。温度が右上がり（$u_x>0$）なら熱は左へ流れる。

保存則に入れると、$u$ だけの式になります。

$$c\rho \frac{\partial u}{\partial t}=K\frac{{\partial }^{2}u}{\partial x^2}⟹\frac{\partial u}{\partial t}=k\frac{{\partial }^{2}u}{\partial x^2},k=\frac{K}{c\rho }$$

この $k$ を熱拡散率と呼びます^[2]。伝導率そのものではなく、熱容量で割った量が効く。

単位が答えを決める

$k$ の単位は ${\mathrm{m}}^2/s$ です。長さの 2 乗を時間で割った形。

この方程式に出てくる量は $x$、$t$、$k$ の 3 つだけなので、時間と長さを結ぶ組み合わせは $\frac{x^2}{k}$ しかありません。

つまり距離 $x$ まで熱が届く時間は $x^2$ に比例する。倍の距離なら 4 倍、10 倍の距離なら 100 倍の時間がかかります。

数値は熱物性のデータベースによります^[4]。銅と水で 760 倍の開きがある。金属の柄が熱くなるまでの時間と、鍋の水が温まるまでの時間の差はここから来ます。

熱核

初期の熱が原点の 1 点に集まっている場合を解きます。答えの形は、次元の話から予想できる。

長さのスケールが $\sqrt{kt}$ しかないので、解は $\eta =\frac{x}{\sqrt{4kt}}$ という無次元量の関数になるはずです。

熱の総量が変わらない条件と合わせて、次の形を仮定します。

$$u(x,t)=\frac{1}{\sqrt{t}}F(\eta ),\eta =\frac{x}{\sqrt{4kt}}$$

これを方程式に代入します^[1]。

相似解の計算

$\eta$ の微分は $\frac{\partial \eta }{\partial t}=-\frac{\eta }{2t}$、$\frac{\partial \eta }{\partial x}=\frac{1}{\sqrt{4kt}}$ です。

$$\begin{array}{rl}{u}_t& =-\frac{1}{2}{t}^{-3/2}(F+\eta F^{\prime })\\ k{u}_{xx}& =\frac{1}{4}{t}^{-3/2}{F}^{\prime \prime }\end{array}$$

等しいと置いて $t^{-3/2}$ を払うと、$\eta$ だけの常微分方程式になります。

$$F^{\prime \prime }+2\eta F^{\prime }+2F=0$$

左辺は $(F^{\prime }+2\eta F{)}^{\prime }$ とまとめられます。積の微分を展開すれば $F^{\prime \prime }+2F+2\eta F^{\prime }$ で、確かに一致する。

したがって $F^{\prime }+2\eta F$ は定数です。遠方で $F$ も $F^{\prime }$ も $0$ に近づくので、その定数は $0$。

$$F^{\prime }=-2\eta F⟹F=C{e}^{-{\eta }^2}$$

変数分離で解けました。総量を $1$ にそろえると $C=\frac{1}{\sqrt{4\pi k}}$ になります。

$$\Phi (x,t)=\frac{1}{\sqrt{4\pi kt}}\exp \left(-\frac{x^2}{4kt}\right)$$

これが熱核です^[1]。偏微分方程式が、相似変数を通すと 1 階の常微分方程式に落ちました。

面積は変わりません。熱の総量が保たれるので、広がったぶんだけ高さが落ちます。

一般の初期温度

熱核が点熱源の答えなので、初期温度 $f$ は点熱源の重ね合わせと見られます^[1]。

$$u(x,t)={\int }_{-\infty }^{\infty }\Phi (x-\xi ,t)f(\xi )d\xi$$

$\xi$ にあった熱が $\Phi$ の形で広がり、それを全部足した。方程式が線形だからできる操作です。

例：跳びのある初期温度

$x<0$ で温度 $1$、$x>0$ で温度 $0$ という初期状態をとります。積分は $\xi <0$ の範囲だけ残る。

$$u(x,t)={\int }_{-\infty }^0\frac{1}{\sqrt{4\pi kt}}{e}^{-(x-\xi {)}^2/(4kt)}d\xi$$

$s=\frac{x-\xi }{\sqrt{4kt}}$ と置き換えると、標準的なガウス積分になります。

$$u(x,t)=\frac{1}{2}\mathrm{erfc}\left(\frac{x}{\sqrt{4kt}}\right)$$

$\mathrm{erfc}$ は補誤差関数です。$x=0$ では $\mathrm{erfc}(0)=1$ なので、温度はいつでも $\frac{1}{2}$。

跳びの真ん中の値が、時間によらず保たれます。

熱は無限の速さで伝わる

熱核は $t>0$ のどの $x$ でも正です。指数関数は $0$ になりません。

原点だけを温めた瞬間に、$1$ km 先の温度も上がっている計算になる。値は $e^{-10^{12}}$ のような桁ですが、厳密には $0$ ではありません。

物理としては近似の限界です。分子の運動を連続の量でならしたモデルなので、極端に短い時間や距離では成り立ちません。

最大値原理

長方形の領域 $0\le x\le L$、$0\le t\le T$ で考えます。温度の最大値がどこに現れるか、という問いです。

答えは、下辺（$t=0$）と左右の辺（$x=0$ と $x=L$）のどこかです。この 3 辺を放物型境界と呼びます。

内部と上辺には最大値が来ない。あとから熱が湧き出すことはない、という当たり前の内容を式にしたものです。

証明

そのままだと、内部の最大点で $u_t=0$、$u_{xx}\le 0$ となって矛盾が出ません。等号があるためです。

そこで小さな $\epsilon >0$ をとり、$v=u+\epsilon x^2$ を考えます。この $v$ は方程式を少しだけ破る。

$$v_t-k{v}_{xx}=u_t-k{u}_{xx}-2k\epsilon =-2k\epsilon <0$$

$v$ の最大値が内部か上辺の点 $(x_0,t_0)$ にあったとします。$x$ について最大なので $v_{xx}\le 0$。

$t$ については、$t_0<T$ なら $v_t=0$、$t_0=T$ なら $v_t\ge 0$ です。どちらでも $v_t\ge 0$ になる。

すると $v_t-k{v}_{xx}\ge 0$ となり、さきほどの $-2k\epsilon <0$ と食い違います。

よって $v$ の最大値は放物型境界の上。$\epsilon \to 0$ とすれば、$u$ についても同じことが言えます。

一意性

同じ初期温度と同じ境界温度をもつ解が 2 つあったとして、差を $w$ と置きます。

$w$ は熱方程式を満たし、放物型境界の上で $0$ です。最大値原理から $w\le 0$。

$-w$ にも同じ議論を当てると $w\ge 0$ になります。したがって $w\equiv 0$ で、解は 1 つに決まる。

時間を逆に回せない

熱方程式は時間の向きを変えると別ものになります。$\sin nx$ の形の温度は $e^{-k{n}^{2}t}$ で減りますが、逆向きでは $e^{+k{n}^{2}t}$ で増える。

$n=100$ なら、$t$ を少し戻すだけで $e^{10^{4}kt}$ の桁になります。初期データのわずかな誤差が、桁違いに拡大される。

過去の温度分布を求める問題が難しいのはこのためです。前向きには滑らかになる一方で、後ろ向きには手がつけられません^[1]。

数値で解いてみる

$x$ を $\Delta x$ 刻み、$t$ を $\Delta t$ 刻みにして、2 階微分を差分に置き換えます^[3]。

$$u_j^{n+1}=u_j^n+r\left(u_{j+1}^n-2u_j^n+u_{j-1}^n\right),r=\frac{k\Delta t}{(\Delta x{)}^2}$$

前の時刻の 3 点から次の時刻を作るだけ。連立方程式を解く必要がありません。

ところが $r$ を大きくとると、計算が発散します。刻みを細かくして精度を上げようとしたときに、まず出会う壁です。

安定条件の導出

$u_j^n=g^n{e}^{\mathrm{i}j\theta }$ という形を入れて、1 ステップあたりの倍率 $g$ を調べます。

$$g=1+r\left(e^{\mathrm{i}\theta }-2+e^{-\mathrm{i}\theta }\right)=1+2r(\cos \theta -1)=1-4r{\sin }^2\frac{\theta }{2}$$

$|g|\le 1$ が全部の $\theta$ で成り立てばよい。${\sin }^2\frac{\theta }{2}$ は最大で $1$ なので、条件は次のとおりです。

$$1-4r\ge -1⟺r\le \frac{1}{2}$$

$r>\frac{1}{2}$ では $|g|>1$ の成分が出て、1 ステップごとに符号を変えながら増えていきます^[3]。

$\Delta x$ を半分にすると $\Delta t$ は 4 分の 1 にしなければならない。ここでも $x^2$ と $t$ の対応が顔を出します。

上下で違うのは $r$ の値だけです。$0.4$ と $0.55$ の差が、この結果を分けます。

確率の言葉に移す

熱核は、平均 $0$、分散 $2kt$ の正規分布の密度そのものです。

$$\Phi (x,t)=\frac{1}{\sqrt{2\pi \cdot 2kt}}\exp \left(-\frac{x^2}{2\cdot 2kt}\right)$$

分散が時間に比例するのは、独立な小さい移動を足し合わせた結果です。ランダムウォークの極限がブラウン運動で、その位置の分布が熱核になる。

標準偏差は $\sqrt{2kt}$ で時間の平方根。熱の広がり方と、酔歩が原点から離れる速さが同じ式に乗っています。

検算のしかた

熱核を得たら、方程式に入れ直すより先に、総量を確かめると早い。

$${\int }_{-\infty }^{\infty }\Phi (x,t)dx=1$$

これが $1$ でなければ、係数のどこかを間違えています。$\int e^{-a{x}^2}dx=\sqrt{\pi /a}$ に $a=\frac{1}{4kt}$ を入れて確かめる。

次元も見ます。$\Phi$ は温度の密度なので $\frac{1}{\text{長さ}}$ の次元。$\frac{1}{\sqrt{kt}}$ は ${\mathrm{m}}^{-1}$ で、確かに合っています。

保存則とフーリエの法則を並べれば方程式が出て、次元を見れば $\sqrt{t}$ が出る。熱方程式の性質は、ほとんどがこの 2 つから読めます。

出典