從上一小節中，了解到了 CIRR 與 FEAST 兩個演算法的優缺點，在此我 們提出一個折衷的辦法來改善 CIRR 對大型區域內的特徵向量有所遺漏的問 題，卻又不花費有如 FEAST 般龐大的計算代價，由於我們並未搜尋到相關的論 文，所以在此暫且稱之為 MLCIRR(Muti-Level CIRR)。由上小節得知，CIRR 由於初始向量空間無法容納下所有的特徵向量，所以無法迭代，而 FEAST 又 因為龐大的計算代價和精度問題而讓人詬病。顯然得，若想要讓整個演算法得 迭代得以成立，我們必須將架構偏向 FEAST，但我們可以於一些地方做出改進。
• 將 CIRR 中的 SVD 步驟移除
CIRR 的問題在於 M 會使得內部得特徵向量經過 SVD 而消逝，而又因無 法迭代而不能取消 SVD。因此，若今日的架構可以迭代，則我們可以去除
SVD 的步驟，進而使得縝密性上升。如此一來就直接改善了 FEAST 第一 次的代價和精度，即使從這步開始進行 FEAST，都已經是有所改善。
• 隨著每次的迭代改變積分路徑
由於需要迭代的緣故，不得不把架構偏向 FEAST，但這不代表不能再度 加快其收斂。回憶積分路徑在是一個大的橢圓曲線，而我們所求的特徵對 又由於矩陣是漢米頓矩陣的緣故，所有的特徵值都在於實數軸上。若將大 橢圓在進行迭代時分化成小橢圓，而依然包覆著整個實數軸，且這個舉 動，比起直接再用大橢圓進行迭代，並不會增加任何的計算量，甚至會有 減少，因為僅將所對應該區域的的特徵向量成為該小橢圓中的初始向量，
而不在所求區域內之向量，自然就不再被進行計算。而注意到由於線性系 統的計算量是 N∗ L 所以比起 FEAST 的迭代，計算量顯然是下降的 (因 為 N 可以固定，而 L 只會減少)。更棒的事情是，比起 FEAST 的大橢圓 用非常少的 N 去切割，如今小橢圓一樣用 N 去切割，則相對於原本的積 分值逼近程度，是有非常大幅度的提升的，因此迭代收斂速度也會遠快於 原始的 FEAST，於我們的實際測試結果而言，最多 2 次就可以直接達到 解線性系統的精準度，而若以原始的 FEAST 則需要約 6 次左右的迭代。
• 視殘差決定迭代所需包覆的區域
關於如何選擇小橢圓，我們視第一輪 CIRR 的結果，來決定該如何選擇。
如果第一次的 CIRR，有部分的特徵對殘差過高 (> δ)，則我們建議使用 將整條實數軸覆蓋的方式，然後以每兩個特徵值為端點，並設置最小的距 離以避免特徵簇或重根問題。假設 l 為分割的小橢圓數，而該區域內有
ˆ
ml 個 CIRR 算出的特徵值，則給予此小橢圓 ˆml− 1 個隨機向量，並且以 2M ≤ N 的原則決定 M。
• 殘差過高時，以隨機向量為初始值
由於殘差過高的原因，此處的特徵向量本身通常並沒有太好的修正性值，
而且由於不能精準的知道其到底是在兩個相鄰小橢圓中的哪一個，導致我 們無法得知將其至入何者，為了避免這些煩惱，故我們採取的是隨機向 量。以數值結果來說，即使並沒有使用到特徵向量來修正，卻因為積分路 徑的大幅度縮小，以及內部的特徵對稀疏的原因，於此步驟就可以得到解 線性系統精度的殘差，而鮮少需要用到第三步，亦即是若殘差足夠小的情 況下所使用的第二種劃分法。
• 殘差足夠小時，以極小的區域包覆單根，並導入特徵向量為初始向量 如果已有足夠小的殘差 (ε)，則可以直接以每個特徵值為圓心，用非常小 的數值為半徑 (ρ = 10ε)，分別將每個特徵值包覆，以求直接達到計算目 標精度 (ϵ)，若有重根或特徵簇現象，則置入相對應之特徵向量，且不再 重新計算一次該區域。該注意的是，如果此廣義特徵值問題中，有重根，
亦或是非常鄰近之特徵簇，那就會使得我們的半徑幾乎區近於 0，雖然並 不會使整個程序無法運行，卻會造成大量的計算浪費。因此，我們需要再 附加兩個條件來避開此類問題。第一個條件是，我們設定一個最小的半徑 rmin，若兩個特徵值 λj 與 λj+1 的距離小於 rmin 則 λj+1 = λj+2 以此類推 直至 |λj − λj+l| > rmin。如此一來就解決了過度浪費計算量的問題。然而 若問題本身擁有特徵簇，則無論始任何的特徵問題解法都會對其感到棘 手，而感謝區塊形式的解法被推導完成，我們可以藉由增加初始向量的數 量 L，來達到開拓足夠的投影空間來容納其數量，雖然其也同時造成了數 倍的計算量，但其是無法被避免的。
以有限的核心數來說，比起 FEAST 對全部的特徵向量再做一次迭代，我們 此處可以僅對不足的區域做迭代，而若核心數足夠的情況下，亦可同時改善所 有的區域，故多了非常大的彈性，然而此處可能牽扯到技術層面上如何分配核 心的問題，但由於我們的數值測試是用 8 核心進行測試，故只要 8≤ N 就已經 可以達到全分配，所以並無處理此問題，但應該並非無法解決的課題。
此演算法最大的變動，除了將第一輪的 FEAST 以 CIRR 取代而節省計算成 本並提升精準度外，在於之後的迭代將積分路徑切為小圈使得精準度大幅提升，
並其代價與原始 FEAST 的迭代相較不增反減。然而，重新劃分積分路徑會使得 整個演算法與傳統的線性系統解法 LU 分解無法共用，因為若不改變積分路徑 和化分數，則整個 FEAST 中都僅需解同一組線性系統問題，而當在迭代時一 併改變積分路徑的話，則無法反覆的利用於第一輪解出的 LU 分解之結果，但 在矩陣問題都非常龐大的情況下，通常並不會使用 LU 分解來解線性系統問題，
因此這並不能算上非常大的缺憾。
3 馬克斯威爾方程與建構的矩陣
在本節中，將介紹如何生成並離散三維光子晶體的馬克斯威爾方程，並如何 倚靠快速傅立葉 (FFT) 變換為預條件 (preconditioner) 快速解出線性系統
Algorithm 3 MLCIRR method
1: 輸入:A, B ∈ Rn×n, V ∈ Rn×L, N, N2, M, M2, γ, ρ, ρmin, r, δ
2: 輸出:(ˆλ1, ˆx1), (ˆλ2, ˆx2)...(ˆλm, ˆxm)
3: 令 ωj = γ + ρ(exp(2πi(j+
1 2)
N ))1,≤ i ≤ N
4: 解 Yi = (ωiB− A)−1BV, 1≤ i ≤ N
5: 計算 ˆSk= N1 ∑N−1
j=0 (ωjρ−γ)k+1Yj, k = 0, 1, ..., M − 1
6: 由 ˆS 建構正交基 ˜Q
7: 使用瑞雷 -瑞茲程序投影 ˜A = ˜QTA ˜Q, ˜B = ˜QTB ˜Q
8: 計算投影矩陣書 ( ˜A, ˜B), 1≤ i ≤ ˆm 之特徵對 (ˆλi, ˆwi)
9: 還原並計算特徵向量 ˆxi = ˜Q ˆwi, 1≤ i ≤ ˆm
10: 計算殘差 εi
11: if ∃εi > δ then
12: i = 1, j = 1, k = 1
13: while i≤ m − 1 do
14: k = i, ρj = λ˜k−˜λ2i+1
15: while ρi < ρmin 且 i < m− 1 do
16: i = i + 1
17: ρj = λ˜k−˜λ2i+1
18: end while
19: γj = λ˜k+˜2λi+1,i = i + 1, Lj = i− k
20: end while
21: 以 γj, ρj 為中心和半徑，N2, M2, Lj 以 CIRR 解各區域 GEP。
22: end if
23: while εi > ϵ do
24: 令 γ = λi, ρi = 10εi, Vi = xi，以 FEAST 解該區域 GEP
25: end while
Figure 1: 左為簡易立方晶格之單一原求與圓柱垂直圓柱示意圖，右為面心立晶 格之單一圓球與橢圓柱示意圖
3.1 從馬克斯威爾方程到離散出的矩陣
光子晶體的馬克斯威爾方程可以寫成
▽ × H = iωE,
▽ × E = −iωu0H
▽ · εE = 0,
▽ · H = 0,
H 和 E 分別代表磁場和電場，ω 代表頻率，u0 代表磁場常數，ε 為材質的介電 系數。令 λ = u0ω2 並將磁場 H 以電場 E 的關係式帶入，則我們可以得到一個 特徵值問題如下
▽ × ▽ × E = λεE,
▽ · εE = 0, 取時諧波後又由於布洛赫定理 [19]，可得
E(x + aℓ) = ei2πk·aℓE(x)
其中 l = 1, 2, 3 分別為三個維度方向的序號，而 aℓ 為光子晶體的晶格移動的 方向向量，因此我們所要解的三維光子晶格的兩種標準型態即為當 al 的夾角
為 π2 時的簡易晶格 (SC,Simple Cube)，與夾腳為 π3 的面心立方晶格 (FCC,Face Centered Cubic)。顯然的 SC 的三個移動向量剛好分別為各方向的單位向量，而 FCC 的三個分量則分別為
其中若是簡易立方晶格則 J2 = In1, J3 = In1n2，若是面心立方晶格則 (positive difine)，因此可以使用 CIRR 對其求解。然而回憶起 CIRR 最大量的 計算代價，在於解大量 Yi = (ωiB− A)−1BV 的線性系統問題。雖然 LU 分解 在於 FEAST 時可以達到重複利用的效果，但由於矩陣大小過於龐大，實際上 連開始運算都會造成記憶體不足的問題，更加不用談論效率。而較傳統的迭代 法諸如 Jacobi，GS，SOR，SSOR 如老牛拖車，而 MINRES，GMRES 亦差強 人意。由於 (ωiB− A) 並非對稱正定，故無法使用 CG，但我們可以退而求其 次選擇廣義的 CG(GCG)，其中最廣為人知的是 CGS 和 BiCGstab。而要加速 GCG 的方法，是找一個良好的預處理矩陣 M，於 GCG 中插入解 M 的線性系