物理現象を記述する数式の世界

偏微分方程式の数学解析

　わたしたちの暮らしを支える多くの技術は、自然界で起こる様々な現象を人類が理解することによって発展していきました。科学の加速的な発展は、複雑な物理現象を理解しやすい数式という形で明瞭化し、定量的な予測を可能にしたことによるところが大きいと考えられます。多くの場合、この物理現象を記述する数式は、高校数学で学ぶ微分を用いた微分方程式 “differential equation” という形の数式で表現されます。そのため、多くの物理現象が微分方程式のもつ数学的側面から理解されてきました。数理学府博士後期課程の藤井さんは、空気や水など流体の運動を記述するナビエ-ストークス方程式 “Navier-Stokes equations” について研究しています。微分方程式の魅力や最近得られた研究成果について藤井さんに解説していただきます。

藤井幹大（⼤学院数理学府）
構成：中島涼輔（大学院理学研究院）

微分方程式とは

　微分方程式 “differential equation” とは、微分を含んだ関数の関係式のことです。その関係式に出てくる関数を具体的に求めることを「微分方程式を解く」と言います。微分方程式は物理学の発展によって現れた概念なので、まずは物理的な例で考えてみましょう。

　雨粒は非常に高い位置から落ちてきますが、空気抵抗があるおかげでゆっくりとしたスピードで落ちてくることは皆さんも経験的に納得できるかと存じます。それでは雨粒はどのくらいの速度で落ちてくるのかということについて考えてみましょう (図1)。雨粒の質量を \(m\)、重力加速度を \(g\) として、空気抵抗が雨粒の速度 \(v\) に比例して \(-kv\) (\(k\) は正定数) と表されるとします。このとき、時刻 \(t=0\) で上空から鉛直下向きに落下をはじめた雨粒の \(t\) 秒後の速度 \(v(t)\) はニュートンの運動方程式\[m\frac{\mathrm{d}}{\mathrm{d}t}v(t)\,=\,mg\,-\,kv(t)\tag{1}\]をみたします^[1]。(1) は時刻 \(t\) についての関数 \(v(t)\) の満たすべき微分方程式です。つまり、微分方程式とは未知の物理的な現象を関数で表したとき、物理法則から導かれるその関数の満たすべき関係式として現れるのです。したがって、微分方程式を解くことができれば知りたかった物理現象を解明することにつながります。

<dfn class="fig">図1</dfn>：<span class="qrinews-figure-title">雨粒の落下</span>　速度に比例した空気抵抗を受けながら、重力によって鉛直下向きに落下する。 — 図1：雨粒の落下　速度に比例した空気抵抗を受けながら、重力によって鉛直下向きに落下する。

虫眼鏡をもったきゅうりくん@右下 — 図1：雨粒の落下　速度に比例した空気抵抗を受けながら、重力によって鉛直下向きに落下する。

　それでは微分方程式 (1) の解を求めましょう。(1) の両辺を \(m\) で割って、右辺の \(v(t)\) を左辺に移項することにより、

\[\frac{\mathrm{d}}{\mathrm{d}t}v(t)\,+\,\frac{k}{m}v(t)\,=\,g\tag{2}\]

となります。(2) の両辺に \(\mathrm{e}^{\frac{k}{m}t}\) をかけると

\[\mathrm{e}^{\frac{k}{m}t}\frac{\mathrm{d}}{\mathrm{d}t}v(t)\,+\,\frac{k}{m}\mathrm{e}^{\frac{k}{m}t}\,v(t)\,=\,g\,\mathrm{e}^{\frac{k}{m}t}\tag{3}\]

となります。ここで積の微分公式 \((fg)'=fg'+f'g\) を利用すると (3) の左辺は

\[\mathrm{e}^{\frac{k}{m}t}\frac{\mathrm{d}}{\mathrm{d}t}v(t)\,+\,\frac{k}{m}\mathrm{e}^{\frac{k}{m}t}\,v(t)\,=\,\frac{\mathrm{d}}{\mathrm{d}t}\left(\mathrm{e}^{\frac{k}{m}t}\,v(t)\right)\]

と書き直せるので (3) は

\[\frac{\mathrm{d}}{\mathrm{d}t}\left(\mathrm{e}^{\frac{k}{m}t}\,v(t)\right)\,=\,g\,\mathrm{e}^{\frac{k}{m}t}\tag{4}\]

と変形できます。よって両辺を時間区間 \([0,t]\) で積分することにより

\[\mathrm{e}^{\frac{k}{m}t}\,v(t)\,-\,v(0)\,=\,\int_0^tg\,\mathrm{e}^{\frac{k}{m}\tau}\mathrm{d}\tau\,=\,\frac{mg}{k}\left(\mathrm{e}^{\frac{k}{m}t}-1\right)\]

となり、

\[v(t)\,=\,\mathrm{e}^{-\frac{k}{m}t}\,v(0)\,+\,\frac{mg}{k}\left(1-\mathrm{e}^{-\frac{k}{m}t}\right)\]

を得ます。このままでは右辺に \(v(0)\) が残ってしまい、\(v(t)\) が決定できません。これは物理的に考えると初期時刻における速度が定まらなければ時刻 \(t\) での速度が定まらないという極めて自然なことに対応しています。そこで微分方程式では初期時刻における値である初期値を指定した上で方程式を解くことが多いです。雨粒は時刻 \(t=0\) で落下を開始したので \(v(0)=0\) と指定することが自然です。よって、求めたかった雨粒の速度を表す関数は

\[v(t)\,=\,\frac{mg}{k}\left(1-\mathrm{e}^{-\frac{k}{m}t}\right)\]

です。この \(v(t)\) の式から

\[\lim_{t\rightarrow\infty}v(t)\,=\,\frac{mg}{k}\]

となり、雨粒はどれだけの時間落ち続けても速度がどんどん大きくなることはなく、\(mg/k\) という速度に収束することがわかります。(これを終端速度と言います。)

偏微分方程式の例

　先の雨粒の問題では考える対象が一点の動きであったため、その点の動きは時刻 \(t\) のみの関数として表現できました。ところが流体の運動や熱の伝導といった現象は時刻 \(t\) と位置 \(x\) の多変数関数 “multivariable function” として表されることになります。多変数関数に対する微分方程式では偏微分 “partial derivative” が現れるため、特に偏微分方程式 “partial differential equation” と呼ばれます。

　ここでは物理現象を記述する偏微分方程式の例として「熱方程式 “heat equation” 」をご紹介します。長い一本のまっすぐな針金を実数直線 \(\mathbb{R}\) と考えて時刻 \(t\)、位置 \(x\) における針金の温度を \(u(t,x)\) とすると、温度分布 \(u(t,x)\) は適当な仮定の元で以下の偏微分方程式に従います (図2)^[2]^[3]：

\[\begin{align}\begin{cases}\dfrac{\partial}{\partial t}u(t,x)\,-\,\dfrac{\partial^2}{\partial x^2} u(t,x)\,=\,0\,,\qquad & t>0\,,\,\,x\in \mathbb{R}\,,\\[5pt]u(0,x)\,=\,u_0(x)\,, \qquad & x\in \mathbb{R}\,.\end{cases}\end{align}\tag{5}\]

ここで、\(u_0(x)\) は初期時刻における温度分布を表す与えられた関数です。

<dfn class="fig">図2</dfn>：<span class="qrinews-figure-title">針金の温度分布のイメージ</span>　針金を部分的に加熱し、初期時刻 \(t=0\) で加熱を止めた時の温度分布の時間変化 (冷めていく様子) は (5) で表される。図の一部は、<a href="https://www.irasutoya.com" class="link-to-external-page" target="_blank"><cite class="article" lang="ja">いらすとや</cite></a>より引用。 — 図2：針金の温度分布のイメージ　針金を部分的に加熱し、初期時刻 \(t=0\) で加熱を止めた時の温度分布の時間変化 (冷めていく様子) は (5) で表される。図の一部は、いらすとやより引用。

　このような多変数関数を解とする偏微分方程式で解を求めるためにはフーリエ変換 “Fourier transform” と呼ばれる道具が有効にはたらきます。関数 \(f(x)\) のフーリエ変換とフーリエ逆変換はそれぞれ

\[\begin{align}\mathcal{F}[f](\xi)\,&=\,\widehat{f}(\xi)\,=\,\int_{-\infty}^{\infty}\mathrm{e}^{-\mathrm{i}x\xi}f(x)\mathrm{d}x\,,\\\mathcal{F}^{-1}[f](x)\,&=\,\frac{1}{2\pi}\int_{-\infty}^{\infty}\mathrm{e}^{\mathrm{i}x\xi}f(\xi)\mathrm{d}\xi\,,\end{align}\]

と定義されます^[4]^[5]。関数 \(f(x)\) のフーリエ変換 \(\widehat{f}(\xi)\) をフーリエ逆変換すると、元の関数 \(f(x)\) を復元できます。フーリエ変換の利点の一つとして「微分を掛け算に変換できる」ことが挙げられます。具体的には次が成立します：

\[\mathcal{F}\left[\frac{\mathrm{d}}{\mathrm{d}x}f\right](\xi)\,=\,\widehat{\frac{\mathrm{d}}{\mathrm{d}x}f}(\xi)\,=\,\mathrm{i}\xi \widehat{f}(\xi)\]

つまり、微分をすることはフーリエ変換後の世界において \(\mathrm{i}\xi\) をかけるという単純な操作となるのです。これに着目して熱方程式 (5) を空間変数 \(x\) に関してフーリエ変換すると

\[\begin{align}\begin{cases}\dfrac{\partial}{\partial t}\widehat{u}(t,\xi)\,+\,\xi^2\widehat{u}(t,\xi)\,=\,0\,,\qquad & t>0\,,\,\,\xi\in \mathbb{R}\,,\\[5pt]\widehat{u}(0,\xi)\,=\,\widehat{u_0}(\xi)\,, \qquad & \xi\in \mathbb{R}\end{cases}\end{align}\tag{6}\]

となり、方程式 (6) は \(t\) のみの微分方程式となります。したがってこの方程式は雨粒の微分方程式のときと同様に解くことができて^[6]、

\[\widehat{u}(t,\xi)=\mathrm{e}^{-t\xi^2}\widehat{u_0}(\xi)\tag{7}\]

となります。よって (7) にフーリエ逆変換を作用させることで

\[\begin{align}u(t,x)&\,=\,\mathcal{F}^{-1}[e^{-t\xi^2}\widehat{u_0}(\xi)](x)\\&\,=\,\frac{1}{2\pi}\int_{-\infty}^{\infty}\mathrm{e}^{\mathrm{i}x\xi}e^{-t\xi^2}\int_{-\infty}^{\infty}\mathrm{e}^{-\mathrm{i}y\xi}u_0(y)\mathrm{d}y\,\mathrm{d}\xi\\&\,=\,\frac{1}{2\pi}\int_{-\infty}^{\infty}\left(\int_{-\infty}^{\infty}\mathrm{e}^{\mathrm{i}(x-y)\xi}\mathrm{e}^{-t\xi^2}\mathrm{d}\xi\right)u_0(y)\mathrm{d}y\end{align}\tag{8}\]

となります。さらに、(8) の括弧内の積分を具体的に計算する^[7]ことで

\[u(t,x)\,=\,\int_{-\infty}^{\infty}(4\pi t)^{-\frac{1}{2}}\mathrm{e}^{-\frac{(x-y)^2}{4t}}u_0(y)\mathrm{d}y\tag{9}\]

という表示を得ます。

　このような解の表示を得ることで解の持つ様々な性質に気付くことができます。例えば、\(\mathrm{e}^{-\frac{(x-y)^2}{4t}}\) は \(x\) について何回でも微分できる滑らかな関数ですから、与えられた初期時刻 \(t=0\) における温度分布 \(u_0(x)\) が微分できないようなギザギザな関数であったとしても、公式 (9) により \(t>0\) になった瞬間に解 \(u(t,x)\) は \(x\) について滑らかな関数になります(図3)。また、十分時間が経ったとき解 \(u(t,x)\) はどのように振る舞うかを考えてみましょう。\(\mathrm{e}^{-\frac{x^2}{4t}}\leqslant1\) ですから

\[|u(t,x)|\,\leqslant\,(4\pi t)^{-\frac{1}{2}}\int_{-\infty}^{\infty}|u_0(y)|\mathrm{d}y\,=\,Ct^{-\frac{1}{2}}\qquad(C\mbox{は正定数})\]

となるので解 \(u(t,x)\) は \(t\to\infty\) において少なくとも \(t^{-\frac{1}{2}}\) 程度のスピードで \(0\) に収束することがわかります。

図3：熱方程式の解の一例　初期時刻 \(t=0\) では \(x=\pm1\) で不連続な分布をしているが、\(t>0\) では滑らかな分布になる。また、\(u(t,x)\) の最大値である \(u(t,0)\) は時間が経つとともに減衰している。グラフは、プログラミング言語の Python を用いて (9) を数値的に計算した。(再生ボタンを押すと動画が再生します。)

　上記の議論は多次元の場合も成り立ちます。詳しく述べると、自然数 \(d=1,2,3\) に対して \(d\) 次元熱方程式は

\[\begin{align}\begin{cases}\dfrac{\partial}{\partial t}u(t,\boldsymbol{x})\,-\,\Delta u(t,\boldsymbol{x})\,=\,0\,,\quad & t>0\,,\,\, \boldsymbol{x}=(x_1,...,x_d)\in \mathbb{R}^d\,,\\[5pt]u(0,\boldsymbol{x})\,=\,u_0(\boldsymbol{x})\,, \quad & \boldsymbol{x}=(x_1,...,x_d)\in \mathbb{R}^d\end{cases}\end{align}\]

で与えられ^[8]、ここで \(\Delta=\displaystyle \sum_{j=1}^d \dfrac{\partial^2}{\partial x_j^2}\) は \(d\) 次元ラプラシアンです^[9]。このとき解は

\[u(t,\boldsymbol{x})\,=\,\int_{\mathbb{R}^d}(4\pi t)^{-\frac{d}{2}}\mathrm{e}^{-\frac{|\boldsymbol{x}-\boldsymbol{y}|^2}{4t}}u_0(\boldsymbol{y})\mathrm{d}\boldsymbol{y}\]

となり^[10]、さらに解の時間減衰評価

\[|u(t,\boldsymbol{x})|\leqslant Ct^{-\frac{d}{2}}\tag{10}\]

が成り立ちます。

非線形偏微分方程式と私の研究

　前節の熱方程式は未知関数の微分の線型結合で表された「線形 “linear” 」偏微分方程式でした。この線型性のおかげでフーリエ変換との相性が良く、解をある程度明示的に表すことができました。一方で、未知関数同士の積が方程式中に含まれるような「非線形 “nonlinear” 」偏微分方程式では解を明示的に表すことは難しいです。それでも、解の具体的表示を使わない様々な解析手法が開発されており、解の具体的な表示がわからなくても時刻無限大での減衰評価などを考察することは可能です。

　ここでは少しマニアックになりますが私が最近の研究で得た結果について簡単にご紹介いたします。私の研究対象は、地球上の大気や海洋の流れ等の地球流体をモデル化する際に現れる、以下の異方的粘性項を持つ \(3\) 次元ナビエ-ストークス方程式 “Navier-Stokes equations” です：

\[\begin{align}\begin{cases}\dfrac{\partial}{\partial t}\boldsymbol{u}\,-\,\Delta_{\rm h} \boldsymbol{u}\,+\,(\boldsymbol{u} \cdot \nabla)\boldsymbol{u}\,+\,\nabla p\,=\,0\,, \quad & t>0\,,\,\,\boldsymbol{x}\in \mathbb{R}^3\,,\\[5pt]\nabla \cdot \boldsymbol{u}\,=\,0\,, \quad & t\geqslant 0\,,\,\, \boldsymbol{x}\in \mathbb{R}^3\,,\\[5pt]\boldsymbol{u}(0,\boldsymbol{x})=\boldsymbol{u}_0(\boldsymbol{x})\,, \quad & \boldsymbol{x}\in \mathbb{R}^3\,.\end{cases}\end{align}\tag{11}\]

ここで \(\boldsymbol{u}=(u_1(t,\boldsymbol{x}), u_2(t,\boldsymbol{x}), u_3(t,\boldsymbol{x}))\) は流体の流速を表す未知のベクトル場で、\(p=p(t,\boldsymbol{x})\) は流体の圧力を表す未知のスカラー場です。また、\(\nabla=\left(\dfrac{\partial}{\partial x_1}, \dfrac{\partial}{\partial x_2},\dfrac{\partial}{\partial x_3}\right)\) は \(3\) 次元の勾配 “gradient” であり、\(\Delta_{\rm h}= \dfrac{\partial^2}{\partial x_1^2}+\dfrac{\partial^2}{\partial x_2^2}\) は水平方向の \(2\) 次元ラプラシアンを表します。方程式 (11) の特徴は、流体の粘性を表す項 \(-\Delta_{\rm h}\boldsymbol{u}\) が鉛直方向の変数 \(x_3\) に関する微分を含まないため等方的^[11]でない ( \(=\) 異方的) という点です。

　数学的に解の性質を解析する際、粘性項に現れるラプラシアンは非常に重要な因子の一つです。ところが (11) には \(x_3\) 方向にラプラシアンの効果がはたらかないので、粘性項が \(3\) 次元ラプラシアンで表現される通常の \(3\) 次元ナビエ-ストークス方程式より数学的取り扱いが難しくなります。一方で方程式の持つ異方性により通常のナビエ-ストークス方程式では起こらない面白い現象を捉えられることがあります。

　方程式 (11) では非線形項 \((\boldsymbol{u} \cdot \nabla)\boldsymbol{u}\) の複雑さゆえ、式 (9) のように解を具体的に表示することはできませんが、小さい解^[12]が存在することはすでに証明されています^[13]。小さい解を考えると、非線形項 \((\boldsymbol{u} \cdot \nabla)\boldsymbol{u}\) は小さい関数同士の積となりますから、非線型項の効果は、線形項 \((\partial/\partial t)\boldsymbol{u}- \Delta_{\rm h} \boldsymbol{u}\) に比べて相対的に大きさが小さくなるので無視できるだろうという考察ができます。したがって非線形問題の解も線形の解から少しだけ「ずれた」ものとみなすことができるため解の存在や性質を議論しやすいのです。私の研究では (11) の小さい解 \(\boldsymbol{u}(t,\boldsymbol{x})\) に関する次の時間減衰評価を得ました：

\[\begin{align}|(u_1(t,\boldsymbol{x}),u_2(t,\boldsymbol{x}))|&\,\leqslant\,Ct^{-1}\,,\\|u_3(t,\boldsymbol{x})|&\,\leqslant\,Ct^{-\frac{3}{2}}\,.\end{align}\]

ここで、\(C\) は \(t\) や \(\boldsymbol{x}\) に依存しない正定数です。方程式の線形主部は \((\partial/\partial t)\boldsymbol{u}- \Delta_{\rm h} \boldsymbol{u}\) という \(2\) 次元の熱方程式的な振る舞いをすることから \(u_1(t,\boldsymbol{x})\)と \(u_2(t,\boldsymbol{x})\) の時間減衰評価が \(t^{-1}\) で与えられることは (10) で \(d=2\) としたものに対応しており極めて自然な結果です。ところが \(u_3(t,\boldsymbol{x})\) の時間減衰率は \(t^{-\frac{3}{2}}\) であり、これは (10) において \(d=3\) としたものに対応しており \(3\) 次元の熱方程式の解のように振る舞うことを意味しています。鉛直方向の速度場が水平方向の速度場より早く減衰することは、地球流体の持つ成層効果の観点から考えても自然な結果です。本研究ではさらに詳しく解の時刻無限大における挙動を解析しておりますが、ここでは詳細を述べられません。詳しくは私のプレプリント (arXiv:2108.11940)をご覧ください。

おわりに

　微分方程式は学部の授業で習うと、様々なテクニックを習い計算練習をすることが多いですが、非線形偏微分方程式の多くは直接計算によって解を求めることはできません。それは方程式の持つ構造が複雑であることを意味しており、裏を返すと方程式の解が多様な性質を内包しているということです。

　非線型偏微分方程式の解が持つ多様な性質は各方程式によって大きく異なることが多く、それらを明らかにするために用いられる解析手法も各方程式ごとに異なることが多いです。そのため非線形偏微分方程式の研究では、各論的な研究が盛んであり^[14]「一般化してより広くに応用の効く統一的な理論を構築する」という数学の思想に完全に沿うものではないかもしれません。しかし、自分の興味のある方程式固有の性質を発見する研究は、方程式と対話しながらその個性を知ろうとする楽しみがあり、人間同士のコミュニケーションと似た魅力を感じます。本記事をお読み頂いた皆さんの中で偏微分方程式論の研究に興味をお持ちいただけた方がいらっしゃれば幸いです。

研究こぼれ話

私は極限を扱う微分積分等の数学に魅力を感じ数学科に入学しました。その気持ちを持ち続けて今も勉強・研究をしています。

フーリエ解析の勉強はその理論の美しさに感嘆し、私の知的好奇心を掻き立ててくれるものです。また、研究ではうまくいかないことも多いですが、研究対象としている偏微分方程式の持つ固有の性質に気づけたときや、計算が綺麗にまとまった時などは何物にも替え難い喜びがあります。

今後も興味の赴くままにたくさん研究をしていきたいと思っています。

Note:

[1] 加速度 \(a\) は速度の変化の割合ですから時間微分で与えられます。したがって (1) は高校物理で習う言葉で書くと \(ma=mg-kv\) となります。
[2] 微分記号において通常 “\(\mathrm{d}\)” と書くべきところを “\(\partial\)” と書いていますが、これは偏微分というものです。例えば \(\frac{\partial}{\partial t}u(t,x)\) は \(t,x\) に関する \(2\) 変数関数 \(u(t,x)\) を (\(x\) を定数と思って) \(t\) について微分していることを意味します。同様に \(\frac{\partial^2}{\partial x^2}u(t,x)\) は \(x\) についての \(2\) 階導関数であることを意味します。
[3] \(\mathbb{R}\) は実数全体の集合を表します。
[4] 端点が \(-\infty\) や \(\infty\) の積分は広義積分と呼ばれ、
\[\qquad\int_{-\infty}^{\infty}f(x)\mathrm{d}x\,=\,\lim_{L\to -\infty}\int_{L}^{0} f(x)\mathrm{d}x\,+\,\lim_{R\to\infty}\int_{0}^{R} f(x)\mathrm{d}x\]
と定義されます。
[5] \(\mathrm{i}=\sqrt{-1}\) は虚数単位を表します。
[6] (6) は (2) で \(k/m=\xi^2, g=0\) としたものと同じです。
[7] 複素解析学のコーシーの積分定理を用いて計算します。
[8] 太文字 \(\boldsymbol{x}\) はベクトルであることを意味します。高校では矢印記号を用いていたもの ( \(\vec{x}\) ) を大学数学では太文字としてベクトルを表します。また、\(\mathbb{R}^d\) は \(d\) 次元空間で \(d=1\) のときは数直線、\(d=2\) のときは平面、そして \(d=3\) のときは (\(3\) 次元) 空間を表します。
[9] ラプラシアン \(\Delta=\displaystyle \sum_{j=1}^d \dfrac{\partial^2}{\partial x_j^2}\) は \(1\) 次元の際の空間変数の \(2\) 階微分 \(\dfrac{\mathrm{d}^2}{\mathrm{d} x^2}\) の多次元版の一つです。
[10] \(|\boldsymbol{x}|\) は \(\boldsymbol{x}\) のユークリッドノルム (高校数学の言葉ではベクトルの大きさ) \(|\boldsymbol{x}|=\sqrt{x_1^2+...+x_d^2}\) を表します。
[11] 等方的とは、どの方向の空間変数も等しく同じ役割を持っていることです。
[12] 小さい解とは、初期値 \(\boldsymbol{u}_0(\boldsymbol{x})\) のある意味での大きさが十分小さいときの対応する方程式の解のことを意味します。ある意味でよい初期値を選択した場合には、初期値の小ささに応じて時間発展した解 \(\boldsymbol{u}(t,\boldsymbol{x})\) の大きさも小さくなります。
[13] Chemin et al. (2000)
[14] 偏微分方程式を抽象化してより一般的に議論する方向の研究もあります。

より詳しく知りたい方は・・・

タイトル: Large time behavior of solutions to the 3D anisotropic Navier-Stokes equation
著者: Mikihiro Fujii
掲載誌: arXiv:2108.11940 (2021)

タイトル: Fluids with anisotropic viscosity
著者: Jean-Yves Chemin, Benoît Desjardins, Isabelle Gallagher, Emmanuel Grenier
掲載誌: Mathematical Modelling and Numerical Analysis 34, 2:315–335 (2000)

キーワード: 偏微分方程式、流体力学、ナビエ-ストークス方程式、熱方程式

九州大学理学研究院理学府理学部