音波の方程式

POINT

  • 音波を流体力学の基礎方程式から導出する.
  • 摂動の1次までを考えると,波動方程式が得られる.

音波が波動方程式に従うことは,流体力学の方程式において摂動論を適用し,その1次の寄与を取り出すことで示されます.摂動の2次以上を考えると,もはや線形方程式では表すことができなくなります.
【関連記事】

記法

均質媒質中において,音波が存在しない平衡状態を「添字$0$」で表す.例えば,音波が存在するときの密度,圧力は
\begin{aligned}
\rho(\boldsymbol{x},t) &= \rho_0 + \rho^\prime(\boldsymbol{x},t) \\
p(\boldsymbol{x},t) &= p_0 + p^\prime(\boldsymbol{x},t)
\end{aligned}
と表される.

ここで,「均質媒質」とは$\rho_0,p_0$が位置依存性を持たないことを意味し,「平衡状態」とは$\rho_0,p_0$が時間依存性を持たないことを意味する.

音速

状態方程式$p = g(\rho)$が成り立つ場合(バロトロピー性)には,音速を
\begin{aligned}
c(\rho)=\sqrt{\frac{\mathrm{d} p}{\mathrm{d} \rho}(\rho)}
\end{aligned}
として
\begin{aligned}
\frac{\partial p}{\partial t}
&=\frac{\mathrm{d} p}{\mathrm{d} \rho} \frac{\partial \rho}{\partial t}
=c^2 \frac{\partial \rho}{\partial t} \\
\boldsymbol{\nabla} p
&=\frac{\mathrm{d} p}{\mathrm{d} \rho} \boldsymbol{\nabla} \rho
=c^2 \boldsymbol{\nabla} \rho
\end{aligned}
などと計算ができる.この関係式は,後で用いる.

特に,音波が存在しない(音波振幅が無限小)の極限での音速を

\begin{aligned}
c_0 = c(\rho_0) = \sqrt{\frac{\mathrm{d}p}{\mathrm{d}\rho}(\rho_0)}
\end{aligned}
と表すことにする.これが「音速」の意味を持つことは,後で導出する波動方程式の表式で,$c_0$が波動の速度を表すパラメータとなっていることからわかる.

音波の方程式

マッハ数$\varepsilon$により摂動展開を行い,流体力学の方程式(運動方程式,連続の方程式)の各次数($\varepsilon^n$の項)を取り出すことにより満たすべき方程式が得られる.特に,$\varepsilon$の1次の項は線形方程式となる.

摂動展開

音響マッハ数を
\begin{aligned}
\varepsilon = \frac{|\boldsymbol{u}|}{c_0}
\end{aligned}
と定義し,密度と圧力,粒子速度が
\begin{aligned}
\rho &= \rho_0 + \varepsilon \rho_1 + \varepsilon^2 \rho_2 +\cdots\\
p &= p_0 + \varepsilon p_1 + \varepsilon^2 p_2 +\cdots\\
\boldsymbol{u} &= \cancel{\boldsymbol{u}_0} + \varepsilon\boldsymbol{u}_1 + \varepsilon^2 \boldsymbol{u}_2 +\cdots
\end{aligned}
と摂動展開できるとする.このとき,$\boldsymbol{u}_0$が存在しないのは,上で音響マッハ数を$\boldsymbol{u}=O(\varepsilon)$の形で定義したためである.

運動方程式

物質微分は
\begin{aligned}
\frac{\mathrm{D}\boldsymbol{u}}{\mathrm{D} t}
&=\frac{\partial \boldsymbol{u}}{\partial t} + (\boldsymbol{u}\cdot \boldsymbol{\nabla})\boldsymbol{u} \\
&=\varepsilon \frac{\partial \boldsymbol{u}_1}{\partial t}
+\varepsilon^2
\biggl[\frac{\partial \boldsymbol{u}_2}{\partial t}
+ (\boldsymbol{u}_1\cdot \boldsymbol{\nabla})\boldsymbol{u}_1 \biggr]
+O(\varepsilon^3)
\end{aligned}
となるから,運動方程式
\begin{aligned}
\frac{\mathrm{D}\boldsymbol{u}}{\mathrm{D} t}
=-\frac{1}{\rho}\boldsymbol{\nabla} p
\end{aligned}
\begin{aligned}
\begin{cases}
O(\varepsilon^0)&: \displaystyle \,
0=-\boldsymbol{\nabla} p_0 \\[5pt]
O(\varepsilon)&: \displaystyle \,
\rho_0 \frac{\partial \boldsymbol{u}_1}{\partial t}= - \boldsymbol{\nabla} p_1 \\[5pt]
O(\varepsilon^2)&: \displaystyle \,
\rho_0 \biggl[\frac{\partial \boldsymbol{u}_2}{\partial t}
+ (\boldsymbol{u}_1\cdot \boldsymbol{\nabla})\boldsymbol{u}_1 \biggr]
+ \rho_1\frac{\partial \boldsymbol{u}_1}{\partial t}
= - \boldsymbol{\nabla} p_2 \\
& \qquad \qquad \qquad \vdots
\end{cases}
\end{aligned}
となる.

連続の方程式

同様に連続の方程式
\begin{aligned}
\frac{\partial \rho}{\partial t}
+\mathrm{div\,} (\rho \boldsymbol{u})
=0
\end{aligned}
\begin{aligned}
\begin{cases}
O(\varepsilon^0)&: \displaystyle \,
\frac{\partial \rho_0}{\partial t} = 0 \\[5pt]
O(\varepsilon)&: \displaystyle \,
\frac{\partial \rho_1}{\partial t} +\mathrm{div\,} (\rho_0 \boldsymbol{u}_1) = 0 \\[5pt]
O(\varepsilon^2)&: \displaystyle \,
\frac{\partial \rho_2}{\partial t} +\mathrm{div\,} (\rho_0 \boldsymbol{u}_2 + \rho_1 \boldsymbol{u}_1) = 0 \\
& \qquad \qquad \qquad \vdots
\end{cases}
\end{aligned}
となる.特に,$O(\varepsilon^0)$の方程式は,$\rho_0$が平衡状態を表す(時間に依存しない)ことと整合している.

状態方程式

状態方程式$p = g(\rho)$が成り立つ場合(バロトロピー性),
\begin{aligned}
p &=p_0 + \varepsilon p_1 + \varepsilon^2 p_2 +\cdots\\
g(\rho) &= g(\rho_0) + g^\prime (\rho_0) (\rho-\rho_0) \\
&\qquad\qquad + \frac{1}{2} g^{\prime\prime} (\rho_0) (\rho-\rho_0)^2
+ \cdots
\end{aligned}
であるから$\rho-\rho_0 = \varepsilon \rho_1 + \varepsilon^2 \rho_2 +\cdots$に注意すれば
\begin{aligned}
\begin{cases}
O(\varepsilon^0)&: \displaystyle \,
p_0 = g(\rho_0) \\[5pt]
O(\varepsilon)&: \displaystyle \,
p_1 = \underbrace{\frac{\mathrm{d}p}{\mathrm{d} \rho}(\rho_0)}_{=c_0^2} \cdot \rho_1\\[5pt]
O(\varepsilon^2)&: \displaystyle \,
p_2 = \underbrace{\frac{\mathrm{d}p}{\mathrm{d} \rho}(\rho_0)}_{=c_0^2} \cdot \rho_2
+\frac{1}{2} \frac{\mathrm{d}^2 p}{\mathrm{d} \rho^2}(\rho_0)
\cdot \underbrace{\rho_1^2}_{\mathrlap{=(p_1/c_0^2)^2}} \\[5pt]
& \qquad \qquad \qquad \vdots
\end{cases}
\end{aligned}
となる.

0次

摂動の0次の方程式
摂動の0次の量が満たすべき方程式は
\begin{aligned}
&\boldsymbol{\nabla} \rho_0=\boldsymbol{\nabla} p_0 = 0 \\
&\frac{\partial \rho_0}{\partial t} = \frac{\partial p_0}{\partial t} = 0 \\
&p_0 = g(\rho_0)
\end{aligned}
となる.これは,「均質媒質における平衡状態の密度,圧力がそれぞれ$\rho_0$,$p_0$であること(空間的,時間的な依存性を持たないこと)」及び「密度が$\rho_0$のときに圧力が$p_0$となること」の仮定と整合している.
【解説】
$O(\varepsilon^0)$の運動方程式,連続の方程式,状態方程式はそれぞれ
\begin{aligned}
\begin{cases}
\,\displaystyle
0=-\boldsymbol{\nabla} p_0 \\
\,\displaystyle
\frac{\partial \rho_0}{\partial t} = 0\\
\,\displaystyle
p_0 = g(\rho_0)
\end{cases}
\end{aligned}
である.

第3式から$p_0(\boldsymbol{x},t) = g\bigl(\rho_0(\boldsymbol{x},t)\bigr)$なので,第1式より

\begin{aligned}
\boldsymbol{\nabla} \rho_0
&=\frac{1}{g^\prime\bigl(\rho_0(\boldsymbol{x},t)\bigr)} \boldsymbol{\nabla} p_0 \\
&=0
\end{aligned}
である.同様に第2式から
\begin{aligned}
\frac{\partial p_0}{\partial t}
&=g^\prime\bigl(\rho_0(\boldsymbol{x},t)\bigr) \frac{\partial \rho_0}{\partial t} \\
&=0
\end{aligned}
が成り立つ.//

1次

摂動の1次の方程式
$\mathrm{rot\,} \boldsymbol{u}_1 = 0$であるとき,速度,圧力,密度は速度ポテンシャル$\phi_1(\boldsymbol{x},t )$を用いてそれぞれ
\begin{aligned}
\begin{cases}
\, \boldsymbol{u}_1(\boldsymbol{x},t ) = \mathrm{grad\,} \phi_1(\boldsymbol{x},t ) \\[5pt]
\, p_1(\boldsymbol{x},t ) = -\rho_0 \dfrac{\partial \phi_1}{\partial t}(\boldsymbol{x},t ) \\[5pt]
\, \rho_1(\boldsymbol{x},t ) = -\dfrac{\rho_0}{c_0^2} \dfrac{\partial \phi_1}{\partial t}(\boldsymbol{x},t )
\end{cases}
\end{aligned}
と表すことができる.また,速度ポテンシャルは波動方程式
\begin{aligned}
\biggl(\frac{1}{c_0^2}\frac{\partial^2}{\partial t^2} - \Delta\biggr)
\phi_1(\boldsymbol{x},t )
=0
\end{aligned}
を満たす.これに$\boldsymbol{\nabla},\partial/\partial t$を作用させることで
\begin{aligned}
\begin{cases}
\,\displaystyle
\biggl(\frac{1}{c_0^2}\frac{\partial^2}{\partial t^2} - \Delta\biggr)
\boldsymbol{u}_1(\boldsymbol{x},t )
=0 \\
\,\displaystyle
\biggl(\frac{1}{c_0^2}\frac{\partial^2}{\partial t^2} - \Delta\biggr)
p_1(\boldsymbol{x},t )
=0 \\
\,\displaystyle
\biggl(\frac{1}{c_0^2}\frac{\partial^2}{\partial t^2} - \Delta\biggr)
\rho_1(\boldsymbol{x},t )
=0
\end{cases}
\end{aligned}
が成り立つこともわかる.
$\mathrm{rot\,} \boldsymbol{u}_1 \neq 0$の場合は
\begin{aligned}
\begin{cases}
\,\displaystyle
\biggl(\frac{1}{c_0^2} \frac{\partial^2}{\partial t^2} -\Delta - \mathrm{rot\,}\mathrm{rot\,} \biggr) \boldsymbol{u}_1
= 0 \\
\,\displaystyle
\biggl(\frac{1}{c_0^2}\frac{\partial^2}{\partial t^2} - \Delta\biggr)
p_1(\boldsymbol{x},t )
=0\\
\,\displaystyle
\biggl(\frac{1}{c_0^2}\frac{\partial^2}{\partial t^2} - \Delta\biggr)
\rho_1(\boldsymbol{x},t )
=0
\end{cases}
\end{aligned}
が成り立つ.

【解説】
$O(\varepsilon^1)$の運動方程式,連続の方程式,状態方程式はそれぞれ

\begin{aligned}
\begin{cases}
\,\displaystyle
\rho_0 \frac{\partial \boldsymbol{u}_1}{\partial t}= - \boldsymbol{\nabla} p_1 \\
\,\displaystyle
\frac{\partial \rho_1}{\partial t} +\mathrm{div\,} (\rho_0 \boldsymbol{u}_1) = 0 \\
\,\displaystyle
p_1 = c_0^2 \rho_1
\end{cases}
\end{aligned}
である.

第2式は,第3式から

\begin{aligned}
\frac{1}{c_0^2}\frac{\partial p_1}{\partial t} +\rho_0 \mathrm{div\,} \boldsymbol{u}_1 = 0
\end{aligned}
となる(0次の結果$\boldsymbol{\nabla}\rho_0=0$を用いた).この式に$\partial/\partial t$を作用させ,第1式を用いれば
\begin{aligned}
\biggl( \frac{1}{c_0^2}\frac{\partial^2 }{\partial t^2} -\Delta \biggr) p_1 = 0
\end{aligned}
が得られる.よって第3式から
\begin{aligned}
\biggl( \frac{1}{c_0^2}\frac{\partial^2 }{\partial t^2} -\Delta \biggr) \rho_1 = 0
\end{aligned}
が成り立つこともわかる.

同様に,$\boldsymbol{\nabla}$を作用させ,第1式を用いれば

\begin{aligned}
-\frac{\rho_0}{c_0^2}\frac{\partial^2 \boldsymbol{u}_1}{\partial t^2}
+\rho_0 \underbrace{\mathrm{grad\,} \mathrm{div\,} \boldsymbol{u}_1}
_{=\Delta \boldsymbol{u}_1 + \mathrm{rot\,}\mathrm{rot\,} \boldsymbol{u}_1 }
= 0
\end{aligned}
である(ベクトル解析の公式を用いた)から,
\begin{aligned}
\biggl(\frac{1}{c_0^2} \frac{\partial^2}{\partial t^2} -\Delta - \mathrm{rot\,}\mathrm{rot\,} \biggr) \boldsymbol{u}_1
= 0
\end{aligned}
が導かれる.特に,$\mathrm{rot\,}\boldsymbol{u}_1(\boldsymbol{x},t) = 0$の場合は波動方程式
\begin{aligned}
\biggl( \frac{1}{c_0^2}\frac{\partial^2 }{\partial t^2} -\Delta \biggr) \boldsymbol{u}_1 = 0
\end{aligned}
となる.



以上を見通しよく計算するには,速度ポテンシャルを用いるとよい.

第1式の$\mathrm{rot\,}$をとると

\begin{aligned}
\frac{\partial }{\partial t} \mathrm{rot\,} \boldsymbol{u}_1 = 0
\end{aligned}
となる.

もし,$\mathrm{rot\,} \boldsymbol{u}_1(\boldsymbol{x},t=0) = 0$であれば,時間によらず$\mathrm{rot\,}\boldsymbol{u}_1(\boldsymbol{x},t) = 0$となる.このとき,速度ポテンシャル$ \phi_1(\boldsymbol{x},t )$を使って

\begin{aligned}
\boldsymbol{u}_1(\boldsymbol{x},t )
= \mathrm{grad\,} \phi_1(\boldsymbol{x},t )
% + \mathrm{rot\,}\boldsymbol{A}(\boldsymbol{x},t )
\end{aligned}
と書ける.

このとき,第1式は

\begin{aligned}
\boldsymbol{\nabla}
\biggl( \rho_0 \frac{\partial \phi_1}{\partial t} + p_1 \biggr)
=0
\end{aligned}
となる.よって,任意関数$f(t)$を用いて
\begin{aligned}
\rho_0 \frac{\partial \phi_1}{\partial t} + p_1
=f(t)
\end{aligned}
と表せる.$\phi_1$を$\displaystyle \phi_1-\int^t f(s)\,\mathrm{d}s$と取り直せば(他の式の表式を変えてしまうことなく)
\begin{aligned}
\rho_0 \frac{\partial \phi_1}{\partial t} + p_1
=0
\end{aligned}
と表せる.

また,

\begin{aligned}
\frac{\partial^2 \phi_1}{\partial t^2}
&=-\frac{1}{\rho_0}\frac{\partial p_1}{\partial t} \\
&=-\frac{c_0^2}{\rho_0}
\underbrace{\frac{\partial \rho_1}{\partial t}}_{\mathrlap{=-\mathrm{div\,} (\rho_0 \boldsymbol{u}_1)}} \\
&=c_0^2 \Delta \phi_1
\end{aligned}
より波動方程式
\begin{aligned}
\biggl(\frac{1}{c_0^2}\frac{\partial^2}{\partial t^2} - \Delta\biggr)
\phi_1(\boldsymbol{x},t )
=0
\end{aligned}
が得られる.//

2次

摂動の2次の方程式
「運動エネルギー密度」と「ポテンシャルエネルギー密度」をそれぞれ
\begin{aligned}
K=\frac{1}{2} \rho_0|\boldsymbol{u}_1|^2,\quad
U=\frac{p_1^2}{2\rho_0 c_0^2}
\end{aligned}
と表すとき,
\begin{aligned}
&\rho_0\frac{\partial \boldsymbol{u}_2}{\partial t}
=-\boldsymbol{\nabla} (p_2 + \mathcal{L}) ,\quad(\mathcal{L}=K-U)\\
&\frac{\partial p_2}{\partial t} + \rho_0 c_0^2 \mathrm{div\,} \boldsymbol{u}_2 \\
&=\frac{\partial}{\partial t}
\biggl[\biggl(1+\frac{\rho_0}{c_0^2}\frac{\mathrm{d}^2 p}{\mathrm{d} \rho^2}(\rho_0)\biggr)U+K\biggr]
\end{aligned}


1次の場合と同様に速度ポテンシャル$\phi_2$を導入すれば
\begin{aligned}
\begin{cases}
\,\displaystyle
\boldsymbol{u}_2 = \mathrm{grad\,}\phi_2(\boldsymbol{x}, t) \\
\,\displaystyle
p_2 = -\rho_0 \frac{\partial \phi_2}{\partial t}(\boldsymbol{x}, t) - \mathcal{L} \\
\,\displaystyle
\mathcal{L} = \frac{1}{2} \rho_0 | \boldsymbol{\nabla} \phi_1 |^2
- \frac{1}{2} \frac{\rho_0}{c_0^2} \biggl(\frac{\partial \phi_1}{\partial t}\biggr)^2
\end{cases}
\end{aligned}
と表すことができる.
【解説】
$O(\varepsilon^2)$の運動方程式,連続の方程式,状態方程式はそれぞれ
\begin{aligned}
\begin{cases}
\,\displaystyle
\rho_0 \biggl[\frac{\partial \boldsymbol{u}_2}{\partial t}
+ \underbrace{(\boldsymbol{u}_1\cdot \boldsymbol{\nabla})\boldsymbol{u}_1}_{=\boldsymbol{\nabla} K/\rho_0} \biggr]
+ \underbrace{\rho_1\frac{\partial \boldsymbol{u}_1}{\partial t}}_{=-\boldsymbol{\nabla} U}
= - \boldsymbol{\nabla} p_2\\
\,\displaystyle
\frac{\partial \rho_2}{\partial t}
+\underbrace{\mathrm{div\,} (\rho_0 \boldsymbol{u}_2 + \rho_1 \boldsymbol{u}_1)}
_{=\rho_0\mathrm{div\,}\boldsymbol{u}_2 -\frac{1}{c_0^2}\frac{\partial}{\partial t} (K+U) }
= 0 \\
\,\displaystyle
p_2 = c_0^2 \rho_2
+\overbrace{\frac{1}{2} \frac{\mathrm{d}^2 p}{\mathrm{d} \rho^2}(\rho_0)
\biggl(\frac{p_1}{c_0^2}\biggr)^2}^{\mathrlap{=\frac{\rho_0}{c_0^2}\frac{\mathrm{d}^2 p}{\mathrm{d} \rho^2}(\rho_0) \cdot U}}
\end{cases}
\end{aligned}
である.これらから,簡単な計算で上にまとめた摂動の2次の方程式が導ける.

最後に,計算の詳細を補足する.
まず,ベクトル解析の公式

\begin{aligned}
\boldsymbol{\nabla} (\boldsymbol{A}\cdot\boldsymbol{B})
&=(\boldsymbol{A}\cdot\boldsymbol{\nabla})\boldsymbol{B}
+(\boldsymbol{B}\cdot\boldsymbol{\nabla})\boldsymbol{A} \\
&\qquad
+\boldsymbol{A}\times\mathrm{rot\,}\boldsymbol{B}
+\boldsymbol{B}\times\mathrm{rot\,}\boldsymbol{A}
\end{aligned}
から
\begin{aligned}
\rho_0(\boldsymbol{u}_1\cdot \boldsymbol{\nabla})\boldsymbol{u}_1
=\boldsymbol{\nabla}\biggl(\frac{1}{2} \rho_0 |\boldsymbol{u}_1|^2\biggr)
=\boldsymbol{\nabla} K
\end{aligned}
と変形できる.

また,1次の結果を使えば

\begin{aligned}
\rho_1\frac{\partial \boldsymbol{u}_1}{\partial t}
&=-\frac{p_1}{c_0^2} \frac{\boldsymbol{\nabla}p_1}{\rho_0}
=-\boldsymbol{\nabla} U \\
c_0^2 \cdot \mathrm{div\,} ( \rho_1 \boldsymbol{u}_1)
&= \mathrm{div\,} (p_1 \boldsymbol{u}_1) \\
&= \boldsymbol{\nabla}p_1 \cdot \boldsymbol{u}_1 + p_1 \mathrm{div\,} \boldsymbol{u}_1 \\
&=-\rho_0 \frac{\partial \boldsymbol{u}_1}{\partial t}\cdot \boldsymbol{u}_1
- p_1 \frac{1}{\rho_0 c_0^2}\frac{\partial p_1}{\partial t} \\
&=-\frac{\partial}{\partial t} (K+U)
\end{aligned}
がわかる.//


参考文献

[1]音と音波 (基礎物理学選書 (4))
[2]非線形音響―基礎と応用 (音響テクノロジーシリーズ):第2章の内容,特に「$\S 2.1.2$ 線形な波動方程式」「$\S 2.1.6$ 微小ひずみ」を参考にしました.より進んだ非線形現象の内容も議論している書籍です.
[3]振動・波動 (基礎物理学選書 (8)):$\S 6.8$「音波」
[4]流体力学 (前編) (物理学選書 (14)):$\S 13.$「運動方程式の第1積分」(2) 渦なしの流れ.