滤波对角化方法的简单探讨及应用

傅霞　刘楠

摘要：介绍了滤波对角化计算本征模式的方法，其中一个重要的过程是求解广义逆矩阵。通過matlab编程实现了谐振腔电磁场本征值和本征模式的计算，并对比了采用奇异值分解和高斯消元法求解广义逆矩阵对本征值计算结果的影响。

关键词：滤波对角化；本征值；本征模式；奇异值分解；高斯消元法

中图分类号：TB

文献标识码：A

doi：10.19311/j.cnki.16723198.2017.09.098

0引言

滤波对角化算法（Filter Diagonalization Method，简称FDM）是在1990年由Neuhauser提出的用于计算量子力学本征值的算法，对求解各种系统本征值得到较好的效果。随着FDM算法在应用中的长时间探索与研究，其求解过程简化得到重大发展，并且在事件分析、收敛性、可靠性等方面都有很大改善，FDM算法在物理、化学、生物科学等领域得到迅速应用。

谐振腔本征值和本征模式问题是电磁领域的一个重要问题，多数谐振腔的优化都要以此为基础，获得本征值和本征模式后也利于进行工作模式分析。将FDM方法应用于电磁波本征模式的计算，其基本思想是：首先通过给定激励源，采用时域有限差分（FDTD）求解电磁场，由于腔体本身的谐振而产生进行滤波作用，然后利于矩阵对角化方法处理电磁场时域结果，求得特定结构中电磁波的本征值及其对应的本征模式。该方法是一种混合方法，在弥补了解析方法和数值方法的弱点的同时发挥它们各自的长处。

1滤波

对于存在本征振荡模式的系统，其参量具有exp（-iωmt）形式，其本征频率为ωm，该系统可以用状态矢量s描述，则s满足的方程的一般形式为

ids（t）dt=Ls（t）（1）

其差分离散形式为

s（t+Δt）-s（t）Δt=-iLs（t）（2）

如果给定一个0时刻的初始值s（0），通过式（2），可以得到s（Δt），然后可以再由s（Δt）得到s（2Δt），一直计算下去，在没有外加干扰的条件下，式（3）经过若干时间步长Δt的演变，系统中只会存在本征模式的线性组合，也就是说通过方程（2）实现了滤波过程，s将演变为本征模式的线性组合。从初始时刻开始，经过一个特定时间段T后，s（T）、s（T+Δt）、s（T+2Δt）…都会是共振模式的线性组合。前面t在0～T之间的推进过程即为滤波过程。

一般只关心某个频段的本征模式，可以通过使用具有频域特性的激励来区分希望的和不希望的模式，即通过激励，促使某些频率范围内的共振模式形成。只要在式（2）右端增加一个激励gf（t）即可

s（t+Δt）-s（t）Δt=-iHs（t）+gf（t）（3）

其中g=g（x，y，z），为空间位置的函数，gf（t）将激励频率在f（t）内且具有与g（x，y，z）相同形式的本征模式。在初始时间段内，经过充分的激励后，再采用方程（2）滤波，则经过较短的时间后，将使s成为希望存在的频率范围内的共振模式的线性组合。

2对角化

对角化部分是使用小规模的线性代数方法在一定的范围内区分相近的本征值。对于满足等式Lvm=λvm的vm称为H的本征模式，λ称为本征值。找出L个向量s1，…，sL（L必须大于本征模式数M），sl是所有本征模式的线性组合。

采用FDTD法计算出的电场E或B构造R和S矩阵，再根据对角化方法采用Matlab编程实现上述算法，其中关键的一个步骤就是求解S的广义逆矩阵S，其计算方法有多种，包括采用奇异值分解法（SVD）、迭代法。本文分别SVD和高斯消元法两种方法实现广义逆的求解。在Matlab中广义逆直接可以采用调用pinv函数实现，其本质是SVD算法。

另外，也可以根据下式求解S：

S=ST（SST）-1（12）

但是由于先求SST的逆矩阵，数值误差较大，在Matlab中通常转换为采用以下公式求解

S=ST/（SST）（13）

其中“/”为Matlab的右除，其本质是采用高斯消元法。

采用频率在c（23.6m-1）～c（25.3m-1）之间的激励，对1m×1.01m二维矩形腔体进行FDTD求解，实现滤波过程。

采用的网格剖分为500×500。取P和L分别为100和20，采用SVD和高斯消元法两种方法计算了本征值，计算结果及误差如表1所示。可以看出采用高斯消元法，对角化过程的误差更小，其计算广义逆矩阵的时间为1.7ms，而采用SVD计算广义逆矩阵的时间为0.9ms，相比滤波过程的计算时间（与网格数量有关，一般需要几十分钟），它们之间计算时间的差别可以忽略。同时可以看到两种计算得到的本征值与理论值误差基本一致，这是由于本例中采用FDTD进行滤波过程的数值误差更大。

4结论

本文介绍了采用滤波对角化计算本征模式的方法，通过matlab编程实现了谐振腔电磁场本征值和本征模式的计算，通过对比奇异值分解法和高斯消元法，发现高斯消元法用于对角化更加准确。

参考文献

[1]Neuhauser D.Time‐dependent reactive scattering in the presence of narrow resonances：Avoiding long propagation times[J].Journal of Chemical Physics，1991，95（7）：49274932.

[2]方剑委.基于滤波对角化方法提高傅立叶变换质谱数据质量[D].长沙：国防科学技术大学，2013.

[3]曾亚文.微波器件谐振腔特征模式分析的算法研究[D].北京：电子科技大学，2008.

[4]Werner G R，Cary J R.Extracting degenerate modes and frequencies from time-domain simulations with filter-diagonalization[J].Journal of Computational Physics，2008，227（10）：52005214.

[5]郭文彬.奇异值分解及其在广义逆理论中的应用[D].上海：华东师范大学，2004.

[6]施吉林，史晓发.广义逆矩阵的两个迭代算法[J].大连理工大学学报，1986，（S1）：3745.