FDTD境界条件

1. はじめに

有限差分時間領域法（FDTD）シミュレーションでは、境界条件は計算領域の外側境界において電磁場がどのように振る舞うかを定義します。物理的には、波がシミュレーションセルの端に到達したときに何が起こるかを決定します。これはFDTDモデルの最も重要な側面の1つです。なぜなら、数値領域は常に有限である一方、物理問題はしばしば開放空間での伝搬、繰り返し構造、または閉じ込められた空洞を表現するからです。したがって、境界条件の役割は、切り詰められたメッシュの外側における場の正しい数学的延長を課すことです。

境界の選択が不適切であると、波はメッシュ内へ反射して戻り、物理解に干渉し、強い定在波アーティファクトを生じる可能性があります。最悪の場合、エネルギーがシミュレーション体積の内部に閉じ込められ、計算された場分布はもはや意図したデバイスを表さなくなります。これは特に、開放フォトニック構造、散乱問題、放射系、およびインパルス応答計算において重要であり、そこでは外向きの波が領域をきれいに離れる必要があります。逆に、周期構造では正しい挙動は吸収ではなく繰り返しであり、いくつかのテスト問題では境界で場を明示的に固定することが望ましい場合があります。したがって、境界条件の選択はモデルの物理と一致していなければなりません。

数学的には、Maxwellの回転方程式は有限のYee格子上で解かれます。ソースのない等方領域では、これらは次のように書けます

\[ \mu \frac{\partial \mathbf{H}}{\partial t} = - \nabla \times \mathbf{E}, \qquad \varepsilon \frac{\partial \mathbf{E}}{\partial t} + \sigma \mathbf{E} = \nabla \times \mathbf{H}. \]

内部の更新方程式はこれらの関係から直接導かれますが、最外層セルでは差分ステンシルが保存されたメッシュの外側に達します。境界条件は、更新を完了するために必要な欠けた情報を与えます。実際には、これはシミュレーションボックスの各外面に対して、接線方向の場を固定するのか、吸収するのか、周期的に折り返すのか、あるいは整合した吸収層を通じて減衰させるのかを決定する数値規則を割り当てることを意味します。

OghmaNanoでは、境界条件は Optical リボン上の Boundary Conditions ボタンから設定します（Figure ?? を参照）。対応するエディタ（Figure ?? を参照）では、FDTD領域の各外面に異なる条件を割り当てることができます。利用可能な境界タイプは Dirichlet、Mur ABC、Periodic、および PML です。

選択可能な6つの境界は、デカルトシミュレーションボックスの6つの面に対応します。OghmaNanoで使用されるデフォルトの表示構成では、\(y_0\) はシミュレーション領域の上面、\(y_1\) は下面、\(x_0\) は左面、\(x_1\) は右面、\(z_0\) はユーザーに最も近い面、\(z_1\) はユーザーから最も遠い面です。この方向マッピングを Figure ?? に示します。

境界条件エディタは、これら6つの面のそれぞれに独立して1つの規則を割り当てます。これにより、たとえば横方向に周期境界、伝搬方向に吸収境界を持つ混合領域や、自由空間への放射のためにすべての面にPMLを持つ領域を構成することが可能になります。PMLが選択されると、その面に対して追加の深さパラメータが利用可能になります。この深さはメッシュセル数で与えられ、吸収媒質のためにどれだけの外側FDTD層を確保するかを決定します。

2. 境界条件エディタ

境界条件エディタは、\(y_{\min}\)、\(y_{\max}\)、\(x_{\min}\)、\(x_{\max}\)、\(z_{\min}\)、および \(z_{\max}\) の各面に対して個別の制御を提供し、それぞれ Figure ?? に示された方向規約における \(y_0\)、\(y_1\)、\(x_0\)、\(x_1\)、\(z_0\)、および \(z_1\) に対応します。各面について、選択されたオプションは、格子外縁における接線方向の電場および磁場がどのように更新されるかを決定します。

3. Dirichlet境界条件

Dirichlet境界条件は、境界上の場を所定の値に固定します。最も一般的なFDTDの使用では、この値はゼロであるため、境界は

\[ \mathbf{E}_{\mathrm{tan}} = 0 \qquad \text{on the selected face,} \]

を満たします。したがって、この境界は電場に対して完全反射壁として振る舞い、外向きの波はすべてシミュレーション領域へ反射して戻ります。

3. Mur ABC境界条件

Mur吸収境界条件は、計算領域の端における外向き進行波の振る舞いを近似することにより、反射を減少させるように設計されています。その考え方は、境界法線方向に1方向波動方程式を満たすことで、格子を離れるエネルギーが反射して戻るのではなく外向きに伝搬し続けるようにすることです。

\(+x\) 方向に伝搬する一次元波に対して、正確な外向き波の関係は

\[ \frac{\partial u}{\partial x} + \frac{1}{c}\frac{\partial u}{\partial t}=0. \]

です。Murの一次吸収境界条件は、この方程式をFDTDメッシュ上で離散化することで得られます。左境界では、標準形は

\[ u_0^{n+1} = u_1^{n} + \frac{c\Delta t-\Delta x}{c\Delta t+\Delta x} \left( u_1^{n+1}-u_0^{n} \right), \]

であり、反対側にも同様の式があります。ここで \(u_i^n\) は空間インデックス \(i\) と時間インデックス \(n\) における場、\(\Delta x\) はメッシュ間隔、\(\Delta t\) はFDTD時間ステップです。

Mur吸収境界条件は計算コストが低く、追加メモリもほとんど必要としないため、軽量シミュレーションや探索的計算に適しています。これらは外向き波を近似しようとするものであり、エネルギーが外側境界で反射するのではなく、シミュレーション領域を離れることを可能にします。

Mur条件はあくまで近似にすぎないため、ある程度の反射は常に残ります。斜入射、広帯域パルス、および複雑な多次元散乱場では精度が低下します。したがって、Mur境界は、波が境界にほぼ法線方向から到達し、非常に高い精度が要求されない単純な伝搬問題に最も適しています。

高度な光学シミュレーション、放射構造、または波が複数角度から境界に入射する場合には、代わりに完全整合層（PML）を使用すべきです。PML境界は、Mur条件よりもはるかに効果的で、かつ著しく低い反射で外向き波を除去する吸収領域を導入します。

4. 周期境界条件

周期境界条件は、シミュレーション領域の一方の面を反対側の面と同一視し、一方の側から出た場が他方の側から再び入るようにします。これは、物理構造が空間的に無限に繰り返されており、計算セルがその繰り返しジオメトリの単位セルを表す場合に適しています。

最も単純な形では、周期性は

\[ \mathbf{E}(x+L_x,y,z,t)=\mathbf{E}(x,y,z,t), \qquad \mathbf{H}(x+L_x,y,z,t)=\mathbf{H}(x,y,z,t), \]

を課し、\(y\) または \(z\) 方向についても、それらの境界が周期として指定されている場合は同様です。数値的には、領域の一方の側のすぐ外側で必要となる場の値は反対側から取得されます。したがって、対向する面は整合した対として扱わなければなりません。物理的には、これは解が孤立した物体ではなく、シミュレーションされたセルの無限格子列に対応することを意味します。

周期境界は、回折格子、フォトニック結晶、メタマテリアルの単位セル、および周期的導波路ジオメトリのように、空間的に繰り返す構造に一般的に使用されます。単一の単位セルのみをモデル化することで、ソースとジオメトリが仮定された周期性と整合している限り、完全な繰り返し構造の振る舞いを表現しながら計算コストを大幅に削減できます。

周期境界では、シミュレーション領域の一方の側を離れた場は反対側から再び入ります。したがって、エネルギーは吸収されず、計算セル内を再循環します。そのため、周期境界は開放空間を表現するために使用すべきではありません。シミュレーション領域のすべての面に周期境界を適用すると、エネルギーは領域の外へ出ることができず、全電磁エネルギーが時間とともに蓄積する可能性があります。多くの実用的なFDTDシミュレーションでは、そのため横方向に周期境界を適用し、伝搬方向にはPMLのような吸収境界を使用します。

最も一般的な場合、周期場は隣接セル間に位相シフトを含むことがあります。これはBloch条件

\[ \mathbf{E}(x+L_x,y,z,t)=\mathbf{E}(x,y,z,t)e^{ik_xL_x}, \]

によって記述され、磁場についても同様の関係があります。ここで実装されている標準的な周期境界は、特別な場合 \(k_x=0\) に対応し、場がシミュレーション領域の対向面で同一に繰り返されることを意味します。

5. PML境界条件

完全整合層（PML）は、現代のFDTDシミュレーションで使用される標準的な高性能吸収境界です。その目的は、広い周波数範囲と広い角度範囲にわたって、最小限の反射で外向きの電磁波を除去することです。外縁で1ステップの境界式を課すのではなく、PMLは物理シミュレーション領域を取り囲む人工的な吸収媒質を導入します。この層に入る波は、界面で鋭いインピーダンス不整合を見ることなく減衰を受けます。

その本質的な考え方は、PMLが入口面で内部媒質に整合しているため、入射波がそれ自体で反射を生むような不連続に遭遇しないということです。層の内部では、場は指数関数的に減衰します。連続体の形では、これは複素座標伸張によって解釈でき、たとえば

\[ x \;\rightarrow\; \int_0^x s_x(\xi)\,d\xi, \qquad s_x(\xi)=\kappa_x(\xi)+\frac{\sigma_x(\xi)}{j\omega\varepsilon_0}, \]

のように表されます。他の方向についても対応する伸張係数があります。ここで \(\sigma_x\) は人工的な導電率プロファイル、\(\kappa_x\) は伸張座標PML形式で用いられるスケーリング係数です。その効果は、伝搬場が物理領域との界面でほぼ無反射のまま、層を通過する間に減衰することです。

OghmaNanoでは、ある面に対して PML を選択すると、シミュレーション領域のその側でこの吸収層が有効になります。対応する PML depth パラメータは、その層に割り当てるメッシュセル数を設定します。一般に、より厚いPMLは、場が減衰する距離が長くなるため、より強い減衰を与えます。PMLが薄すぎると、いくらかのエネルギーが外側の切り詰め境界に到達して反射する可能性があります。十分に厚ければ、外向きの場は外縁に達する前に無視できる値まで減衰します。

厚さ依存性を理解する便利な方法は、吸収媒質内の場成分が概ね

\[ E(d)\sim E(0)\exp(-\alpha d), \]

として減衰することに着目することです。ここで \(d\) はPML内を進んだ距離、\(\alpha\) は導電率グレーディングによって設定される有効減衰係数です。したがって、PMLセル数を増やすと、外側の切り詰め境界に到達する残留振幅が減少します。実用的には、PML深さは、層の背面からの反射が目的とする物理信号に比べて無視できる程度に十分大きくする必要があります。

PMLは通常、導波路放射、孤立物体からの散乱、アンテナ様放射、および過渡パルス伝搬を含む、開放フォトニックおよび電磁シミュレーションにおいて好まれる境界条件です。場が斜めに境界へ入射したり、広いスペクトル成分を含んだりする場合でも、単純な吸収境界条件よりはるかに頑健です。このため、波が計算ボックスから出ていくことが想定される多くのFDTDモデルで、PMLがデフォルトの選択肢となっています。

PMLを用いる場合でも、適切なモデリング実践は重要です。吸収層は、強いエバネッセント場や高共振場の領域に近づけすぎてはいけません。なぜなら、非伝搬の近接場が吸収媒質と相互作用して解を変化させる可能性があるからです。一般に、物理構造とPMLの間には、外向き波面がきれいに形成されるのに十分な自由空間を設けるべきです。このように用いることで、PMLはFDTDにおいて、開放かつ無反射の外部領域に最も近い近似を提供します。