## 概要

波動の伝播を数値的に計算する方法には[有限要素法(FDM)](../1628)、[$k$-space method](../1627)などがある。これらの方法は波動が無限に伝播することを仮定しているが、実際のシミュレーションでは有限の格子内で波動が伝播する。そのため、境界では波動が反射する問題が生じる。これを解決するためには、実際にシミュレーションしたい領域よりもはるかに広く格子を設定するか、境界条件をうまく設定して波動の反射が起こらないようにする。そのような境界条件を**吸収境界条件**<sup>absorbing boundary condition (ABC)</sup>という。

## 公式

次のような波動方程式を考えよう。

$$
\dfrac{\partial ^{2}}{\partial t^{2}} u(x,t) = c^{2} \dfrac{\partial ^{2}}{\partial x^{2}} u(x,t) \quad \text{for } (x, t) \in \mathbb{R} \times [0, \infty)
$$

これを有限な格子空間$[0, 1] \times [0, T]$でシミュレーションするとしよう。$x$軸空間を$M$等分すると$\Delta x = 1/M$、時間を$N$等分すると$\Delta t = T/N$である。波動方程式の解$u$の離散化を$U\_{i}^{n} = u(i \ast \Delta x, n \ast \Delta t)$という。このとき$i = 0, 1, \cdots, M$、$n = 0, 1, \cdots, N$である。吸収境界条件は次のようになる。

$$
U\_{0}^{n+1} = U\_{1}^{n} + \dfrac{\left(c \Delta t - \Delta x \right)}{\left(c \Delta t + \Delta x \right)} \left( U\_{1}^{n+1} - U\_{0}^{n} \right)
$$

$$
U\_{M}^{n+1} = U\_{M-1}^{n} + \dfrac{\left(c \Delta t - \Delta x \right)}{\left(c \Delta t + \Delta x \right)} \left( U\_{M-1}^{n+1} - U\_{M}^{n} \right)
$$

## 説明

1次元と2次元で吸収境界条件を適用した結果は次の通りである。

![](1631_1.webp#center)

![](1631_2.webp#center)

## 導出[^1] 
[^1]: Mur, Gerrit. *"Absorbing boundary conditions for the finite-difference approximation of the time-domain electromagnetic-field equations."* IEEE transactions on Electromagnetic Compatibility 4 (1981): 377-382.

波動方程式を次のように整理しよう。

$$
\begin{align*}
         && \dfrac{\partial ^{2}}{\partial t^{2}} u(x,t) &= c^{2} \dfrac{\partial ^{2}}{\partial x^{2}} u(x,t) \\\
\implies && \dfrac{\partial ^{2}}{\partial t^{2}} u(x,t) - c^{2} \dfrac{\partial ^{2}}{\partial x^{2}} u(x,t) &= 0 \\\
\implies && \left(\dfrac{\partial }{\partial t} - c \dfrac{\partial }{\partial x} \right) \left(\dfrac{\partial }{\partial t} + c \dfrac{\partial }{\partial x} \right) u(x,t) &= 0 \\\
\end{align*}
$$

このとき波動方程式の解$u(x,t)$は次の2つの方程式のいずれかを満たす。

$$
\left(\dfrac{\partial }{\partial t} - c \dfrac{\partial }{\partial x} \right) u(x,t) = 0 \tag{1}
$$

$$
\left(\dfrac{\partial }{\partial t} + c \dfrac{\partial }{\partial x} \right) u(x,t) = 0 \tag{2}
$$

$(1)$を満たす$u\_{+}(x, t) = f(x + ct)$は$-x$方向に移動する波動である。反対に$(2)$を満たす$u\_{-}(x, t) = g(x - ct)$は$+x$方向に移動する波動である。境界の$x = 0$の場所では$U$は左方向に行く波動だけが存在するはずなので、$U$は$(1)$を満たさなければならない。$x = 0$での$t$に対する微分係数を求めると、

$$
\dfrac{\partial u}{\partial t}|\_{x = 0} \approx \dfrac{U\_{0}^{n+1} - U\_{0}^{n}}{\Delta t} \approx \dfrac{\dfrac{U\_{0}^{n+1} + U\_{1}^{n+1}}{2} - \dfrac{U\_{0}^{n} + U\_{1}^{n}}{2}}{\Delta t}
$$

2つ目の近似式で関数値を隣接する2点の平均としたのは、より良い近似のためである[^2]。同様に$x$に対する微分係数は次のように求めることができる。
[^2]: https://eecs.wsu.edu/~schneidj/ufdtd/chap6.pdf

$$
c\dfrac{\partial u}{\partial x}|\_{x = 0} \approx c \dfrac{U\_{1}^{n} - U\_{0}^{n}}{\Delta x} \approx c\dfrac{\dfrac{U\_{1}^{n+1} + U\_{1}^{n}}{2} - \dfrac{U\_{0}^{n+1} + U\_{0}^{n}}{2}}{\Delta x}
$$

$U$は$x = 0$で$(1)$を満たさなければならないので、

$$
\dfrac{\dfrac{U\_{0}^{n+1} + U\_{1}^{n+1}}{2} - \dfrac{U\_{0}^{n} + U\_{1}^{n}}{2}}{\Delta t} - c\dfrac{\dfrac{U\_{1}^{n+1} + U\_{1}^{n}}{2} - \dfrac{U\_{0}^{n+1} + U\_{0}^{n}}{2}}{\Delta x} = 0
$$

$$
\implies \dfrac{ U\_{0}^{n+1} + U\_{1}^{n+1} - U\_{0}^{n} - U\_{1}^{n}}{2 \Delta t} - c\dfrac{ U\_{1}^{n+1} + U\_{1}^{n} - U\_{0}^{n+1} - U\_{0}^{n}}{2 \Delta x} = 0
$$

$$
\implies U\_{0}^{n+1} + U\_{1}^{n+1} - U\_{0}^{n} - U\_{1}^{n} - \dfrac{c\Delta t}{\Delta x} \left( U\_{1}^{n+1} + U\_{1}^{n} - U\_{0}^{n+1} - U\_{0}^{n} \right) = 0
$$

$$
\implies \left(\dfrac{c \Delta t}{\Delta x} + 1 \right) U\_{0}^{n+1} = \left(\dfrac{c \Delta t}{\Delta x} + 1 \right) U\_{1}^{n} - (U\_{1}^{n+1} - U\_{0}^{n}) + \dfrac{c\Delta t}{\Delta x} \left( U\_{1}^{n+1} - U\_{0}^{n} \right)
$$

$$
\implies \left(\dfrac{c \Delta t}{\Delta x} + 1 \right) U\_{0}^{n+1} = \left(\dfrac{c \Delta t}{\Delta x} + 1 \right) U\_{1}^{n} + \left(\dfrac{c \Delta t}{\Delta x} - 1 \right) \left( U\_{1}^{n+1} - U\_{0}^{n} \right)
$$

$$
\begin{align*}
\implies U\_{0}^{n+1} &= U\_{1}^{n} + \dfrac{\left(\dfrac{c \Delta t}{\Delta x} - 1 \right)}{\left(\dfrac{c \Delta t}{\Delta x} + 1 \right)} \left( U\_{1}^{n+1} - U\_{0}^{n} \right) \\\
&= U\_{1}^{n} + \dfrac{\left(c \Delta t - \Delta x \right)}{\left(c \Delta t + \Delta x \right)} \left( U\_{1}^{n+1} - U\_{0}^{n} \right)
\end{align*}
$$

同様に、$x=1$の場合の吸収境界条件は次の通りである。

$$
U\_{M}^{n+1} = U\_{M-1}^{n} + \dfrac{\left(c \Delta t - \Delta x \right)}{\left(c \Delta t + \Delta x \right)} \left( U\_{M-1}^{n+1} - U\_{M}^{n} \right)
$$

■

## コード

1次元波動方程式における吸収境界条件をシミュレートしたコードは次の通りである。