IvyFEM開発日誌（10）.NET 有限要素法 IvyFEM 0.0.0.21 - 時間領域FEMによるH面導波管直角コーナーベンドの散乱係数の計算(Givoli-Neta-Patlashenko's High order ABC)
IvyFEM.dll 0.0.0.21をリリースしました。
Givoli-Neta-Patlashenkoの高次ABCを用いた時間領域FEMに対応しました。
この方法は以前に実装してみたことがありました。PML（これも試した）のように追加領域を必要とせず解けるのでメモリの節約になります。
今回この方法をIvyFEMに実装しました。
〇直線導波管
励振面（下図の領域内線分上）から正弦変調ガウシアンパルスを励振したときを調べ、左右の端面で反射波が発生しないことを確かめました。
直線導波管（励振中）
4096ステップ完了時
時間領域の計算で得られたポート1付近の観測点のEzの時間変化ez1(t)、ポート2付近の観測点のEzの時間変化ez2(t)のグラフからは（入射波より少し遅れてくる）反射波は見られません。
またez1(t)、ez2(t)をFFTしてEz1(f) (= Inc(f))、Ez2(f)(=T(f))を計算し、S21 = T(f)/Inc(f)を計算した結果も合わせて示しています。得られた|S21|を見るとほぼ1.0になりました。
〇直角コーナーベンド
以前に計算した結果は次の記事です。
今回計算した結果は同等の結果が得られました。
励振中１
励振中２
励振終了
反射波、透過波到達
反射波、透過波通過中
反射波、透過波通過完了
4096ステップ完了
このベンドの時間領域の計算で得られたポート1付近の観測点のEzの時間変化ez1(t)、ポート2付近の観測点のEzの時間変化ez2(t)、および直線導波路のポート1付近の観測点のEzの時間変化ezInc(t)を用いると、
反射波の時間変化 r(t) = ez1(t) - ezInc(t)
透過波の時間変化 t(t) = ez2(t)
で、ezInc(t)、r(t)、t(t)をFFTすると複素振幅Inc(f)、R(f)、T(f)が得られます。そうすると、
S11(f) = R(f) / Inc(f)
S21(f) = T(f) / Inc(f)
から散乱係数が求められます。
|S11|、|S21|の計算結果を見るとW/λ = 1.5以下で少し劣化が見られますがまずまずの結果が得られました。
動作確認に使ったサンプルソースコードも公開しています。
