綜合上面數據的結果後,我們於本小節比較 Lanczos 與 MLCIRR 的精度與 計算代價。其計算成本之數據如 Table 1,而精準度表現如 Figure 12,MLCIRR 的部分的參數為 N = 16, N2 = 8, M = 4, M2 = 4, L = 8, Γ = 25, ρ = 25, ρmin =
Table 1: Lanczso 與 MLCIRR 各階段的時間比較
程序 CGS 迭代 (次) CGS 時間 (秒) 迭代倍數 時間倍數
Lanczos 2914 6214 1 1
MLCIRRinitial 10626 21438 3.65 3.45 MLCIRRrefine2 13542 28661 4.64 4.61 MLCIRRrefine2 14572 30767 5.01 4.95 MLCIRRtotal 38740 80997 13.29 13.03
實際 8 線程 實際值 (理論值)
MLCIRRinitial 10602 3449 3.65 0.55(0.46) MLCIRRrefine1 13517 4642 4.64 0.75(0.57) MLCIRRrefine2 14538 4826 5.01 0.78 (0.61) MLCIRRtotal 38657 12973 13.29 2.08(1.62)
6 結論
我們逐項的探討並且實際測試了 CIRR 中各項參數的交互影響,進一步比 較了 CIRR 和 FEAST 的異同及優缺點,最後提出一個整合兩個方法的演算法 MLCIRR。此演算法的優勢在於採用 CIRR 的觀點,以參數 M 去創造虛擬的向 量空間,得到了相當於加速 FEAST 第一次收斂的效果。之後的迭代中,我們 重新劃分了路徑,使得同樣分割數的條件下,更逼近實際積分值,同時內部存 在的特徵對相較於原本的問題,顯得更加的分散,更加重要的是,這些特徵對 將落在不易遺漏的積分路徑旁,而非容易遺漏的中央區域。並且由於同樣擁有 M 參數的關係,我們進行迭代時的向量,將不像 FEAST 一直維持在估計量的 2 倍左右,而會逐漸與目標向量完全相等,達到完全不浪費計算量的情況。又因 為積分路徑急速的收縮,導致外迭代收斂速度亦會是 FEAST 的數倍。而此演 算法存在的風險為,當選取範圍內擁有大量特徵簇時,它並不能像 FEAST 慢 慢的修正,而會因為積分路徑的重新劃分和起始向量不足得原因,造成無法挽 回的特徵對遺漏,同樣因為積分路徑的重劃關係,當問題本身可以使用 LU 分 解時,此演算法並不能像 FEAST 反覆的利用已解過的 LU 分解。
關於 CIRR 的建議參數部分,我們在這邊下一個總結。
• N,M,L 的取捨
當不進行迭代時,相同的線性系統計算成本下,N 與 L 之間的取捨將會是
以 N 為重,原因在於若 N 不足將無法一次性的得到較好的殘差,並且同 時限制了 M 的大小,進一步使得殘差得精準度降低,而在有迭代的情形 下,增加 L 可以有效的確保特徵對不被遺漏。
• 長短軸比 r 的影響
當的積分路徑 (平移值) 遠離特徵值時,能有效的降低迭代次數,但同時也 使得無論是單一次的收斂速度或者是迭代時的收斂速度大幅下滑。並且此 參數在開拓投影空間的 M 和 L 越小時,對迭代時收斂速度的影響越明顯。
因此,除非僅是用以估計的情況下,否則盡量取小。
• 建議參數
對於一個大範圍且尚未預估特徵對數的區域,我們建議的初始測試參數 是,N = 16, M = 8, L = 16 這組參數一次容納約 50 組準確的特徵對,並 且由於 L = 16 在 32 線程以下,將不會因為平移值的不同,因解線性系統 時間不同而造成等帶的其情況,並且 N, M 亦都在理論上的最佳範圍內。
若是對於一以知的特徵對求更精確的解,則為 N = 8, M = 4, L = 1 其中 初始向量為已知的特徵向量。
而關於整體 CIRR 方面:
• CIRR 在大範圍下有著容易遺漏特徵對的問題
• CIRR 於良好的目標區域下 (線性系統迭代次數尚未大幅攀升, 內部特徵對 稀疏) 表現將會比 Lanczos 來的好
• MLCIRR 大幅度的加速了 FEAST,但於起始的圈選範圍和預估的 M, L 有更嚴苛的要求
• 沒有任何的方式能保證 CIRR 不遺漏特徵對,故 CIRR 最安全的使用方式 仍侷限於殘差的縮小
於參考文獻中 [2, 20] 中可以發現,無論是 T. Sakura 或者 E. Polizzi 的團隊,
於近期的研究中接有著墨於如何預估某特定範圍內的特徵對。因此,如何有效 的估計範圍內的特徵對仍是 CIRR 當今最重要的課題,而本文中第四節可以靠 著觀察平移值與迭代次數的關係,來避開大量的特徵簇存在,可以說是非常特 列的作法。若是問題本身是經由微分方程離散而來者,靠著較粗的切割計算較 小的矩陣書來估計特徵對的分佈也不失為一個好的方式。
References
[1] J. Asakura, T. Sakurai, H. Tadano, T. Ikegami and K. Kimura, A numerical method for nonlinear eigenvalue problems using contour integrals, JSIAM Letters, Vol.1,52-55 (2009).
[2] Y. Futamura, H. Tadano, and T. Sakurai, Parallel stochastic estimation method of eigenvalue distribution, JSIAM Letters, Vol. 2, 127-130, (2010).
[3] T.-M. Huang, H.-E. Hsieh, W.-W. Lin and W. Wang, Eigendecomposition of the discrete double-curl operator with application to fast eigensolver for three dimensional photonic crystals, SIAM J. Matrix Anal. Appl., Vol. 34, 369-391(2013)
[4] T.-M. Huang, H.-E. Hsieh, W.-W. Lin, and W. Wang.Fast lanczos eigenvalue solvers for band structures of three dimensional photonic Crystals with face-centered cubic lattice(preprint)
[5] T.-M. HUANG, Y.-L. Huang.W.-W. Lin, and W.-C. Wang, A Null space free Jacobi-Davidson iteration for maxwell’s operator(preprint)
[6] T.-M. Huang, W.-J. Chang, Y.-L. Huang, W.-W. Lin, W.-C. Wang and W. Wang, Preconditioning bandgap eigenvalue problems in three dimen-sional photonic crystals simulations, Journal of Computational Physics, Vol.
229,8684-8703(2010)
[7] T.-M. Huang, Y.-C. Kuo and W. Wang, Computing extremal eigenvalues for three-dimensional photonic crystals with wave vectors near the Brillouin zone center, J. Sci. Comput., Vol. 55, 529-551(2013)
[8] T. Ikegami and T. Sakurai, Contour integral eigensolver for non-Hermitian systems: a Rayleigh-Ritz-type approach, RANMEP2008Taiwanese J. Math., Vol. 14, pp. 825–837 (2010).
[9] T. Ikegami, T. Sakurai, U. Nagashima, A filter diagonalization for generalized eigenvalue problems based on the Sakurai-Sugiura projection method,CS-TR-08-13, Tsukuba (2008)
[10] A. Levin, D. Zhang, E. Polizzi, FEAST fundamental framework for electronic structure calculations: reformulation and solution of the muffintin problem computer physics communications, V183, I11, 2370-2375 (2012)
[11] E. Polizzi, Density-matrix-based algorithms for solving eingenvalue problems Phys. Rev. B., Vol. 79, 115112 (2009)
[12] E. Polizzi,A high-performance numerical library for solving eigenvalue prob-lems: FEAST solver user’s guide v2.1s,(2013)
[13] T. Sakurai, J. Asakura, H. Tadano, T. Ikegami and K. Kimura, A method for finding zeros of polynomial equations using a contour integral based eigen-solver, Proc. Symbolic Numeric Computations 2009, Kyoto. 143-147 (2009) [14] T. Sakurai, J. Asakura, H. Tadano and T. Ikegami, Error analysis for a matrix
pencil of Hankel matrices with perturbed complex moments, JSIAM Letters, Vol. 1,76-79 (2009).
[15] T. Sakurai, H. Sugiura, A projection method for generalized eigenvalue prob-lems, ISE-TR-02-189, Tsukuba,(2002)
[16] T. Sakurai, H. Tadano, CIRR: a Rayleigh-Ritz type method with contour integral for generalized eigenvalue problems, Proc. The First China-Japan-Korea Joint Conference on Numerical Mathematics, Vol. 36,745-757 (2007).
[17] T. Sakurai, H. Tadano, T. Ikegami, U. Nagashima, A parallel eigensolver using contour integration for generalized eigenvalue problems in molecular simulation, CS-TR-08-14, Tsukuba(2008)
[18] K. Senzaki, H. Tadano, T. Sakurai and Z. Bai, A method for profiling the distribution of eigenvalues using the AS method,RANMEP 2008 Taiwanese J. Math, Vol. 14,839-853 (2010).
[19] M. Suzuki, I. Suzuk, Bloch theorem and Energy band, Lecture Note on Solid State Physics(2006)
[20] P.-T. Tang, E. Polizzi,FEAST as a subspace iteration eigensolver Aacelerated by approximate spectral projection,SIMAX, (accepted)
[21] I. Yamazaki, T. Ikegami, H. Tadano, T. Sakurai, Performance comparison of parallel eigensolvers based on a contour integral method and a Lanczos method, Parallel Computing ,39,280–290.(2013)
[22] C. Yang, W.-G. Gao, Z.-J. Bai, X.-S. Li, L.-Q. Lee, P. Husbands. An algebraic substructuring method for large-scale eigenvalue calculation. SIAM J. SCI.
COMPUT. Vol. 27, No. 3. 873–892(2005)
[23] K.-S. Yee, Numerical solution of initial boundary value problems involving Maxwell equations in isotropic media,IEEE Trans. Antennas Propagat.AP-14(3), 302-307. (1966)