周期结构指定频率的本征分析方法

摘要

本发明涉及周期结构指定频率的本征分析方法，包括步骤：A.建立考虑导体和介质损耗周期结构内的电磁场边值问题，通过有限元法的标准变分原理得到电磁场边值问题的泛函方程。B.采用四面体网格剖分求解域，并保证周期边界主面和从面上的网格匹配。C.选择基函数，将电场强度矢量在所有网格内用基函数展开，并运用里兹方法得到展开系数的矩阵方程。D.将展开系数的矩阵方程写成关于传播常数的线性广义本征方程。E.给定一个频率，求解线性广义本征方程，根据传播常数，得到与给定频率相对应的衰减常数和相位常数，通过后处理，得到互作用阻抗。F.重复步骤E，得到不同频率对应的衰减常数，相位常数和互作用阳抗。

说明书

技术领域

本发明属于三维电磁场数值求解的技术领域，具体涉及一种周期结构指定频率的本征分 析方法。

背景技术

周期结构在微波管中应用非常广泛，包括波导截面的周期性变化，波导周期加载膜片， 周期填充介质等。在微波管中，普遍采用周期结构作为器件的高频电路，形成电子注与高频 场相互作用进行能量交换以实现微波振荡或放大的场所。周期结构的高频特性直接影响器件 的工作频率、频带宽度、换能效率和输出功率，以及其他一系列整管性能。高精度地获得周 期结构的高频特性有着极其重要的意义。

周期结构的高频特性包括色散特性、阻抗特性与衰减特性。其中，衰减特性主要由周期 结构中有限电导率导体的欧姆损耗与非理想介质的介质损耗引起，由衰减常数来描述。在高 效率微波管的设计中，准确计算由有限电导率导体和非理想介质引入的衰减特性至关重要。

目前，利用各种计算电磁学方法对周期结构进行三维建模，采用指定相移的本征分析方 法求出其高频特性，已经成为一种常见而有效的方法。在指定相移本征分析方法中，首先在 周期结构周期边界上指定一个相移，然后求解数值方法建立的本征方程获得本征频率与本征 场分布，得出色散特性，通过后续处理求出阻抗特性与衰减特性。重复这一过程，即可得到 周期结构的高频特性。

与不考虑损耗(假定导体为理想导体，介质为理想介质)周期结构的指定相移本征分析 方法相比，考虑损耗(有限电导率导体与非理想介质)周期结构的指定相移本征分析方法会 复杂很多，精度和效率也会低很多。以有限元方法为例，首先，由于周期边界的存在，引入 考虑介质损耗的介电常数，或考虑导体损耗的表面阻抗边界，都会使周期结构有限元分析产 生的有限元矩阵从无损耗情况下的厄米对称矩阵变成复不对称矩阵，从而增加求解的时间和 内存。其次，在考虑损耗的情况下，导体的表面阻抗是频率的函数，采用指定相移的本征分 析方法需要求解一个关于频率的非线性广义本征问题，这比求解无损耗情况下的线性广义本 征问题复杂地多。最后，在考虑损耗的情况下，介质的损耗角正切随频率变化，且没有一个 确定的函数关系。因此在采用指定相移的本征分析方法时，只能采用预估的介质损耗角正切 对周期结构进行本征分析，预估介质的损耗角正切值的精确性将严重影响周期结构本征分析 的精度。

发明内容

本发明的目的是为了克服现有的指定相移本征分析方法在考虑导体和介质损耗的周期结 构高频特性分析中存在的问题，提出了周期结构指定频率的本征分析方法，该方法能够更精 确和高效的求解出考虑导体和介质损耗的周期结构的高频特性。

为了实现上述目的，本发明的技术方案是：周期结构指定频率的本征分析方法，包括以 下步骤：

A.建立考虑导体和介质损耗周期结构内的电磁场边值问题，通过有限元法的标准变分原 理得到电磁场边值问题的泛函方程。

B.采用四面体网格剖分求解域，并保证周期边界主面和从面上的网格匹配。

C.选择基函数，将电场强度矢量在所有网格内用基函数展开，并运用里兹方法得到展开 系数的矩阵方程。

D.将展开系数的矩阵方程写成关于传播常数的线性广义本征方程。

E.给定一个频率，求解线性广义本征方程，得到传播常数与电场展开系数。根据传播常 数，得到与给定频率相对应的衰减常数和相位常数，通过后处理，得到互作用阻抗。

F.重复步骤E，得到不同频率对应的衰减常数，相位常数和互作用阻抗，即可得到周期 结构的高频特性。

本发明的有益效果：利用本发明提出的周期结构指定频率的本征分析方法可以准确快速 求出考虑导体和介质损耗行周期结构的高频特性。这是因为本发明提出的指定频率的本征分 析方法和现有的指定相移本征分析方法相比有以下四点优势：1)指定频率的本征分析方法克 服了指定相移的本征分析方法在考虑导体和介质损耗周期结构高频特性分析中存在的求解非 线性广义本征值方程的问题。在考虑导体和介质损耗周期结构高频特性分析中，指定频率的 本征分析方法仅仅需要求解一个线性广义本征值方程；2)在指定频率的本征分析方法中，衰 减常数是直接通过求解一个线性广义本征值方程的本征值得到的，避免了指定相移本征分析 方法中通过复杂后处理得到衰减常数过程中产生的二次累计数值误差；3)指定频率的本征分 析方法能够精确地确定随频率变化的介质的相对介电常数，而在指定相移的本征分析方法中 随频率变化的介质相对介电常数只能采用一个预估的近似值，从而影响周期结构高频特性分 析的精度；4)采用指定频率的本征分析方法，可以直接得到指定频率下的高频特性，避免了 指定相移的本征分析方法首先在不同相移点求出高频特性，然后通过插值近似求出指定频率 高频特性的复杂过程。

附图说明

图1是本发明的主流程图。

具体实施方式

下面结合附图和具体实施例对本发明作进一步说明。

如图1所示，周期结构指定频率的本征分析方法，包括以下步骤：

A.建立考虑导体和介质损耗周期结构内的电磁场边值问题，通过有限元法的标准变分 原理得到电磁场边值问题的泛函方程。

首先，根据麦克斯韦方程组与周期结构的边界与导体属性得到考虑导体和介质损耗时 周期结构内电磁场的边值问题，如下：

$\left(\begin{array}{cc}▿\times {\mu }_r^{-1}▿\times E-k_0^2{ϵ}_{r}E=0& \mathrm{in\Omega }\\ \hat{n}\times (E_s\times \hat{n})=\hat{n}\times (E_m\times \hat{n})e^{-(\alpha +\mathrm{j\beta })L}& \mathrm{on}{\Gamma }_{\mathrm{PBC}}\\ \hat{n}\times {\mu }_r^{-1}▿\times E=({\mathrm{jk}}_0{\eta }_0/Z_s)\hat{n}\times (E\times \hat{n})& \mathrm{on}{\Gamma }_{\mathrm{SIBC}}\end{array}\right)---\left(1\right)$

公式(1)中第一个式子为频域矢量波动方程，是周期结构有限元仿真中的主方程；其 中，Ω为周期结构的仿真区域空间范围，即为公式(1)的求解域，是矢性偏微分算子符 号，μ_r为求解域Ω中介质的相对磁导率，E为求解域Ω的电场强度矢量，k₀为自由空间 波数，ε_r为求解域Ω中介质的相对介电常数，在考虑介质损耗时，ε_r＝ε′_r(1-jtanδ)， 为复数，-jε′_rtanδ为ε_r的虚部，tanδ是用以描述介质损耗的损耗角正切，通常随频率 升高而增加，ε′_r为ε_r的实部，通常为不随频率变化的常数。

公式(1)中第二个式子为准周期边界条件，其中，Γ_PBC表示准周期边界，由主面及从 面组成；为边界的外法向单位矢量；E_m和E_s分别表示周期边界主面和从面上的电场；j 为虚数单位符号；α和β分别为衰减常数和相位常数，α+jβ为传播常数；L为周期长 度。在考虑导体和介质损耗的情况下，准周期边界条件的物理意义是在周期边界从面上的 电磁场和主面上的电磁场，除了相差一个复数相位系数e^-jβL外，还存在的幅度的衰减系数 e^-αL；

公式(1)中第三个式子为导体的阻抗边界条件，其中，Γ_SIBC表示阻抗边界；η₀为自 由空间波阻抗；Z_s为良导体的表面阻抗，满足：

$Z_s=(1+j)\sqrt{\frac{\mathrm{\pi f\mu }}{\sigma }}---\left(2\right)$

公式(2)中f为频率，σ为导体的有限电导率，显然，导体的表面阻抗Z_s是频率的 函数。由于周期结构中的导体一般都是良导体，因此在边值问题中强加表面阻抗边界条件 可以更加精确地计算导体损耗对周期结构高频特性的影响。

周期边界Γ_PBC和阻抗边界Γ_SIBC组成了求解域Ω的外边界；

从考虑导体和介质损耗时周期结构内电磁场的边值问题，即从公式(1)出发，通过有 限元法的标准变分原理得到电磁场边值问题的泛函方程F(E)，即公式(3)：

$F\left(E\right)=\frac{1}{2}{\int }_{\Omega }[{(▿\times E)}^{*}\frac{1}{{\mu }_r}(▿\times E)-k_0^2{ϵ}_r{E}^{*}·E]\mathrm{d\Omega }$

$+\frac{{\mathrm{jk}}_0{\eta }_0}{z_s}{\int }_{\Gamma }[\hat{n}\times (E^{*}\times \hat{n})]·[\hat{n}\times (E\times \hat{n})]\mathrm{d\Gamma }---\left(3\right)$

公式(3)中，上标*表示对物理量取共轭。使泛函方程即公式(3)取极小值，并且满足公 式(1)中的第二个方程的电场函数E即为周期结构内电磁场边值问题即公式(1)的解，dΓ 表示二维面积分的微元，dΩ表示三维体积分的微元。

B.采用四面体网格剖分求解域，并保证周期边界主面和从面上的网格匹配。剖分后的 求解域被人为分割为多个三维四面体网格，每个网格内的电场由多个基函数的线性组合来 近似表达。由于四面体网格剖分是有限元方法中的一种公知过程，因此本步骤不再详细描 述。

C.选择基函数，将电场强度矢量在所有网格内用基函数展开，并运用里兹方法得到展 开系数的矩阵方程。

选择合适的基函数N_m，将公式(3)和公式(1)中第二个方程中的电场E在所有网格内 用基函数N_m展开，即

$E=\underset{m}{\Sigma }{x_m}^t{N}_m^t---\left(4\right)$

公式(4)中，下标m取值从0到M_T，M_T为网格中所有基函数的个数。上标t将电场展 开系数x_m和基函数N_m按照所在区域进行区分。t∈{I，M，S}，I表示基函数N_i在除周期边 界外的计算区域内，M表示基函数N_m在周期边界主面上，S表示基函数N_m在周期边界从 面上。

将所有网格内的电场E用基函数N_m展开后，代入泛函方程，即公式(3)，运用里兹方 法，并利用周期边界条件将展开系数用替代，得到展开系数满足的矩阵方程如下：

$\left(\begin{array}{cc}{K}^{\mathrm{II}}\left(f\right)+\sqrt{f}{C}^{\mathrm{II}}& K^{\mathrm{IM}}\left(f\right)+K^{\mathrm{IS}}\left(f\right)e^{-\Phi }\\ K^{\mathrm{MI}}\left(f\right)+K^{\mathrm{SI}}\left(f\right)e^{\Phi }& K^{\mathrm{MM}}\left(f\right)+K^{\mathrm{SS}}\left(f\right)+K^{\mathrm{SM}}\left(f\right)e^{\Phi }+K^{\mathrm{MS}}\left(f\right)e^{-\Phi }\end{array}\right)·\left(\begin{array}{c}{x}^I\\ x^M\end{array}\right)=0---\left(5\right)$

公式(5)中，$x=\left(\begin{array}{c}{x}^I\\ x^M\end{array}\right)$为电场展开系数x_m组成的列向量，其中x^I是由周期边界外的网 格内电场展开系数x_m组成的列向量，x^M是周期边界主面上各个网格内电场展开系数x_m组 成的列向量，K表示计算过程中的中间矩阵，其上标从t∈{I，M，S}中任意选取两种常量。 Φ＝(α+jβ)L，α和β分别为衰减常数和相位常数，L为周期长度。分块矩阵C^II中的元 素满足，

$C_{\mathrm{mn}}^{\mathrm{II}}=\frac{2j\sqrt{\pi {\mu }_0\sigma }}{(1+j)}{\int }_{\Gamma }[\hat{n}\times (N_m^I\times \hat{n})]·[\hat{n}\times (N_n^I\times \hat{n})]\mathrm{d\Gamma }---\left(6\right)$

公式(6)中，∫_Γ(·)dΓ表示在有限电导率导体表面的面积分。如果基函数或不在 有限电导率导体的表面，则式中下标m和n用于区分不同的基函数。

其他各分块矩阵可统一描述为：

K^xy(f)＝S^xy-f²·T^xy    (7)

公式(7)中，f为频率；S^xy与T^xy为分块矩阵；x与y为上标，用以区分基函数所在的区域， x∈{I，M，S}，y∈{I，M，S}，I，M，S的含义与前同。S^xy与T^xy中各元素分别满足公式(8) 与公式(9)。

$S_{\mathrm{mn}}^{\mathrm{xy}}={\int }_{\Omega }(▿\times N_m^x)·\frac{1}{{\mu }_r}(▿\times N_n^y)\mathrm{d\Omega }---\left(8\right)$

$T_{\mathrm{mn}}^{\mathrm{xy}}={4\pi }^2{\mu }_0{ϵ}_0{\int }_{\Omega }{N}_m^x·{ϵ}_r{N}_n^y\mathrm{d\Omega }---\left(9\right)$

公式(8)与公式(9)中，∫_Ω(·)dΩ表示在求解域内的体积。

D.将展开系数的矩阵方程写成关于传播常数(α+jβ)的线性广义本征方程。

假定在所有四面体网格中不存在跨越周期边界主面和从面的四面体，那么在公式(5)中的 分块矩阵K^SM(f)＝0，K^MS(f)＝0。然后，对公式(5)中矩阵的第二行乘以e^-Φ，将公式(5) 改写成如下形式

公式(10)是以-e^-Φ为本征值的标准线性广义本征方程。公式(10)中 表示将矩阵$\left(\begin{array}{cc}{K}^{\mathrm{II}}\left(f\right)+\sqrt{f}{C}^{\mathrm{II}}& K^{\mathrm{IM}}\left(f\right)\\ K^{\mathrm{SI}}\left(f\right)& 0\end{array}\right)$简记为A，表示将矩阵 $\left(\begin{array}{cc}0& K^{\mathrm{IS}}\left(f\right)\\ K^{\mathrm{MI}}\left(f\right)& K^{\mathrm{MM}}\left(f\right)+K^{\mathrm{SS}}\left(f\right)\end{array}\right)$简记为B，为后续使用。

E.给定一个频率，求解线性广义本征方程，得到传播常数(α+jβ)与电场展开系 数x。根据传播常数，得到与给定频率相对应的衰减常数和相位常数，通过后处理，得到 互作用阻抗。

给定频率f，系数矩阵A与B为确定的大型稀疏奇异矩阵，可以采用位移求逆隐式重启 Arnoldi算法求解本征方程。首先给定一个预估本征值将标准线性广义本征方程(10)转 化成如下形式：

公式(11)中，$x=\left(\begin{array}{c}{x}^I\\ x^M\end{array}\right)$为电场展开系数x_m组成的列向量。λ为待求本征值-e^-Φ。只 要预估本征值不精确等于λ，那么矩阵为非奇异矩阵，采用隐式重启Arnoldi算 法即可快速求解方程(11)，得到本征值与电场展开系数x。

由于为给定的预估本征值，由的值可以得到感兴趣的本征值λ，即-e^-Φ。由 Φ＝(α+jβ)L，得到与给定频率相对应的衰减常数α和相位常数β

$\alpha =-\frac{1}{L}\mathrm{Re}\left(\ln (-\lambda )\right)---\left(12\right)$

$\beta =-\frac{1}{L}\mathrm{Im}\left(\ln (-\lambda )\right)---\left(13\right)$

公式(12)与公式(13)中，Re(·)与Im(·)分别表示取实部和虚部运算。

根据得到的电场展开系数x_m，结合基函数，由公式(4)可以得到求解域内的电场分 布，进而可以由互作用阻抗的定义得到与指定频率f相对应的互作用阻抗，由电场分布获 得互作用阻抗的过程为本领域的公知过程，因此不再详细描述。

F.重复步骤E，得到不同频率对应的衰减常数α，相位常数β和互作用阻抗，即可得 到周期结构的高频特性。因为衰减特性、色散特性、和互作用阻抗特性是不同频率对应的 衰减常数α，相位常数β和互作用阻抗的变化趋势来描述的，因此可分别得到周期结构的衰 减特性、色散特性、和互作用阻抗特性。

本领域的普通技术人员将会意识到，这里所述的实施例是为了帮助读者理解本发明的 原理，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。本领域的普通 技术人员可以根据本发明公开的这些技术启示做出各种不脱离本发明实质的其它各种具体 变形和组合，这些变形和组合仍然在本发明的保护范围内。