フォトニック結晶のバンドダイアグラムの数値計算
これまでフォトニック結晶のバンドダイアグラムの計算には、「ＰＨバンド」を使用してきましたが、今回、FEMでバンドダイアグラムの計算をしてみましたので結果を掲載します。
FEMでバンドダイアグラムを計算する方法としては、導波路の固有値問題と同じく、界Φ(x, y)を直接解く方法と、Φ(x, y) = φ(x, y) exp(- jβx x) exp(- jβy y) と置く方法があるようです。ここでは、Φ(x, y)を直接解く方法で計算してみました。
正方格子
誘電体ロッド(半径 r = 0.18a、屈折率 n = 3.4)の正方格子のTEモード(電界がz方向成分のみ、多くの文献ではTMモードと記載されています。)のバンド構造を計算します。
バンドダイアグラムの計算結果
「ＰＨバンド」で計算した結果
FEMで計算した場合、フォトニックバンドギャップ(PBG)は、a/λ = 0.3092～0.4514となりました。一方、「PHバンド」で計算した場合0.317～0.459でした。多少の差異はありますがおよそ一致していると思います。
三角形格子
エアホール型三角形格子のバンドダイアグラムを計算しました。エアホールの半径r = 0.30a、基板の実効屈折率 n = 2.76です。
図面(a)
第一ブリュアンゾーンのΓ点を(0, 0)、K点を(4π/(3a), 0)、Ｍ点を(π/a, √3π/(3a))として計算しました。
バンドダイアグラムの計算結果
「PHバンド」の結果
こちらも問題なく解けているようです。
上記図面の場合、エアホールの半径が大きくなると計算できなくなるので周期をずらした場合も計算してみました。
図面(b)
図面(c)
三角形格子の第一ブリュアンゾーン(BZ)は六角形なので解析領域を六角形領域に取るのがいいと思いますが、境界条件を課す辺が多くなり少し面倒なので今回は計算しませんでしたが、またの機会に計算してみたいと思います。
さて、導波路の固有値問題では長方形領域をとって計算しましたがバンドダイアグラムの計算の場合、どうなるか試してみました。
図面(d)
バンドダイアグラムの計算結果
これまでの計算にはなかったバンドが出現しました。なんだろうと思い調べたところ、下記の文献で説明されている現象のようです。つまり、図面(d)の解析領域は、第一BZを１つ以上含んでおり、その為折りたたんだバンド("folded" band)が混入しているようです。文献のFig.2が分かりやすかったです。
(1) Bartlomiej Salski
"The unfolding of bandgap diagrams of hexagonal photonic crystals computed with FDTD"
http://www.jpier.org/PIERM/pierm27/03.12081313.pdf
Progress In Electromagnetics Research M, Vol. 27, 27 - 39, 2012
文献にはこのような折りたたんだバンドを元に戻す方法が記されていますが、それなりに手間がかかりそうです。
一方、長方形領域として計算したいだけなら図面(d)の半分の領域にすることで計算できるようです。下記文献(2)に載っていました。
(2) Ashutosh and P. K. Jain
"FDTD analysis of the dispersion characteristics of the metal pbg structures"
http://jpier.org/PIERB/pierb39/04.11120601.pdf
Progress In Electromagnetics Research B, Vol. 39, 71 - 88, 2012
比較的簡単な構造なのでFEMで試してみました。
周期境界条件は、左右の境界はこれまでと同じくΦ(a, y) = exp(-jβx a) Φ(0, y)を適用しますが、上下の境界は、文献(2)のFig.4のように辺を２つに分割して別々の境界条件を適用します。つまり、
上の境界の右半分 = exp( - jβx a/2) exp( - jβy √3a/2) ×下の境界の左半分
上の境界の左半分 = exp( + jβx a/2) exp( - jβy √3a/2) ×下の境界の右半分
という境界条件を適用します。
バンドダイアグラムの計算結果
どうやらこの解析領域だと不要なバンドの混入は発生していないようです。
第一BZを過不足なく解析領域に取り込むことが正しいバンドダイアグラムの計算には必要なようです。
【追記(2013-04-10)】
六角形領域で計算してみました。六角形領域の場合、右半分の境界上の点の界は左半分の境界上の点の界で表すことができます。また、左半分の境界の頂点の内、２つの頂点だけが独立で残りの頂点は従属変数となります。
左上の境界をB１、左の境界をB2、左下の境界をB3、右下の境界をB4、右の境界をB5、右上の境界をB6とし、頂点を上から反時計回りにv1(上)～v6(右上)とすると、
ΦB4 = a1 ΦB1
ΦB5 = a2 ΦB2
ΦB6 = a3 ΦB3
Φv3 = (1/a3) Φv1
Φv4 = a1 Φv2
Φv5 = a1 Φv1
Φv6 = a2 Φv2
ここで、a:格子定数として、
a1 = exp( - jβx(a/2))exp( +jβy(3a/(2√3)) )
a2 = exp( - jβxa)
a3 = exp( - jβx(a/2))exp( -jβy(3a/(2√3)) )
です。領域は一辺がa/√3 の正六角形となります。
バンドダイアグラム計算結果
他の領域分割の場合とほぼ同じ結果を得ることができました。
