CN115392079A

CN115392079A - 一种基于完全矢量有限元的均匀波导模式计算方法

Info

Publication number: CN115392079A
Application number: CN202210992120.9A
Authority: CN
Inventors: 吴春华; 赵华鹏; 胡俊
Original assignee: University of Electronic Science and Technology of China
Current assignee: University of Electronic Science and Technology of China
Priority date: 2022-08-17
Filing date: 2022-08-17
Publication date: 2022-11-25

Abstract

本发明公开了一种基于完全矢量有限元的均匀波导模式计算方法，属于计算电磁学领域。本发明先将均匀波导横截面剖分成若干三角形离散单元取得各三角形单元顶点坐标；再利用顶点坐标计算出三角形单元各边对应的矢量基函数；以三角形各边的切向电场或切向磁场作为横向场插值参量，将其分别与各边对应的矢量基函数相乘后再求和得到插值函数。与此同时，将矢量波动方程转化为电场和磁场的方程，并对得到的电场和磁场方程进行变分泛函处理；最后将插值函数代入变分泛函后的电场和磁场方程中得到B矩阵为正定矩阵的均匀波导广义本征方程。以此求解出均匀波导的截止波数和截止频率。本发明计算耗时短、精度高等优点，适用于求解均匀波导的各模式。

Description

一种基于完全矢量有限元的均匀波导模式计算方法

技术领域

本发明涉及计算电磁学领域，具体涉及一种基于完全矢量有限元的均匀波导模式计算方法。

背景技术

波导不仅是信息的载体和能量的传输手段，还是各类通信和雷达系统组件的框架。截止波数k_c作为波导的重要参数之一，其值取决于波导系统的横截面和边界条件，决定着电磁波是否能在波导中传播，所以很多微波部件的优化分析往往也是以波导本征问题的求解为基础，也是计算电磁学最为基础和重要的组成部分。

波导本征问题的核心就是确定亥姆霍兹方程中的k_c，如今有很多方法可以求解波导本征问题，如：矩量法、有限差分法、有限元法、边界元法等。其中，有限元法应用最为广泛，又分为矢量有限元法和节点有限元法。而矢量有限元法因其克服了节点有限元法有时会出现“非物理解”、强加边界条件不便和处理媒质边缘困难等问题，更具有优势且更常见。

矢量有限元法中横向场以矢量基函数作为插值函数，以棱边切向电场作为插值参量，而纵向场以标量基函数作为插值函数，以节点纵向电场作为插值参量。通过将矢量波动方程等效为变分问题，然后离散单元再次变分后得到广义本征方程形式：Ax＝λBx，其中矩阵A和矩阵B是稀疏的，可由离散单元坐标值计算得到。通过求解广义本征方程得到λ的值，就可以求出截止波长k_c，但由此得到的矩阵B是非正定的，无法直接利用隐式重启Arnoldi(IRAM)法直接求解。如果用常规算法则计算时间太长且使用内存巨大，也可以使用稀疏矩阵快速LU分解技术来改进IRAM算法，不过其算法过程复杂且不直观。

可见，传统的矢量有限元法同时采用矢量基函数和标量基函数会导致B矩阵不正定，无法直接利用现有的算法或如ARPACK之类的软件直接求解，常规算法存在计算时间太长且使用内存巨大等诸多问题，有必要对其进行改进。

发明内容

本发明的目的在于：针对上述传统的矢量有限元方法的改进需求，本发明提出了一种基于完全矢量有限元的均匀波导模式计算方法，将波导中的TE波和TM波分开计算，利用TE波无纵向电场分量和TM波无纵向磁场分量的特点，横向场(电场或磁场)依然采用矢量基函数作为插值函数，以棱边切向场(电场或磁场)作为插值参量，不使用纵向场插值参量，同样最后可以得到广义本征方程，且矩阵B是正定的，能直接利用如IRAM等现有的算法直接求解，其方法更为简单，所需内存更少，精度更高，更适合用于计算均匀波导的截止波数。

为实现上述目的，本发明采用如下技术方案：

一种基于完全矢量有限元的均匀波导模式计算方法，包括以下步骤：

步骤1、将均匀波导横截面剖分成若干三角形离散单元，取得各三角形离散单元各顶点坐标；

步骤2、基于步骤1得到的顶点坐标计算得到三角形离散单元各边对应的矢量基函数；

步骤3、以三角形各边的切向电场或切向磁场作为横向场插值参量，纵向场不使用插值参量；将插值参量分别与步骤2得到的各边对应的矢量基函数相乘后再求和得到插值函数；

将矢量波动方程转化为电场和磁场的方程；基于变分原理，利用均匀波导TE波没有E_z分量和TM波没有H_z分量的特点对电场和磁场方程进行变分泛函处理；

步骤4、将插值函数代入变分泛函处理后的电场和磁场方程中，得到均匀波导的广义本征值方程；

步骤5、求解步骤4得到的广义本征值方程即可得到均匀波导的截止波数和截止频率。

进一步的，所述步骤2中三角形离散单元各边对应的矢量基函数的计算过程为：

步骤2.1、设三角形离散单元的面积为Δ，三个顶点按逆时针方向依次编号为1、2、3，p(x，y)为三角形内任意一点，该点与三角形各顶点连线将三角形划分为三个小三角形，并将各顶点对应变所在的三角形的面积分别表示为Δ₁、Δ₂、Δ₃；

步骤2.2、根据步骤2.1三角形离散单元面积Δ和任意一个小三角形面积Δ_i构建关于P点坐标的常函数L_i：

步骤2.3、基于步骤2.2构建的常函数L_i和i点对应边长度l_i构建三角形离散单元各边对应的矢量基函数

进一步的，所述步骤5中求解步骤4得到的广义本征方程采用了IRAM算法。

本发明提出的一种基于完全矢量有限元的均匀波导模式计算方法，利用均匀波导TE波没有E_z分量和TM波没有H_z分量的特点，只将横向场用棱边矢量基函数线性表示，从而可以分开计算TE波和TM波的截止频率，得到的广义本征方程的B矩阵为正定矩阵，利用已有的算法或软件即可直接求解广义本征方程。相较于传统矢量有限元方法，更简单直观且计算时间更短，求解过程中不需要纵向场的参与，相当于强行令纵向场为0，忽略了纵向场的影响，计算精度更高，内存使用更少，更适用于求解均匀波导的各模式。

附图说明

图1为三角形棱边元示意图；

图2为实施例中圆-矩形同轴波导结构示意图。

具体实施方式

下面结合附图和实施例对本发明作进一步的详细说明，所述是对本发明的解释而不是限定。

步骤1、将波导横截面剖分成若干三角形离散单元，取得各三角形离散单元各顶点坐标；

步骤2、基于步骤1得到的顶点坐标计算得到三角形离散单元各边对应的矢量基函数。

详细过程如下：

如图1所示，横向电场分量E_t或磁场分量H_t用三角形各边的切向电场E_ti或切向磁场H_ti作为插值参量，设三角形三个顶点按逆时针方向分别编号为1，2，3，面积为Δ。P(x,y)为其内的一点，与三个顶点的连线所形成的三个小三角形面积分别表示为Δ₁、Δ₂、Δ₃。结合步骤1得到的顶点坐标，得到四个三角形的面积表达式如下：

2.2、根据步骤2.1三角形离散单元面积Δ和任意一个小三角形面积Δ_i构建关于P点坐标的常函数L_i：

2.3、基于步骤2.2构建常函数L_i和顶点1对应边长长度l₁计算得到三角形离散单元顶点1对应边的矢量基函数：

其中l₁为三角形顶点1对应的边的长度。

设

为从顶点2到顶点3的单位向量，则有

和

由此在第一边上则有关系式

成立；

垂直于第二条边，且在第二条边上L₂为0，则在第二边上有关系式：

同理，第三边上也存在关系式：

由此可见，

具有作为对应于第一边的切向电场或切向磁场参量的基函数的特点。基于该特点即可得到对应于第二和第三边的基函数：

步骤3、计算插值函数，同时对矢量波动方程进行变分泛函处理。

本实施例仅采用三角形单元内的横向电场或横向磁场作为插值参量，为方便描述，三角形单元内的横向电场或横向磁场统一用P表示。将插值参量P分别与步骤2.3得到的各边对应的矢量基函数相乘后，再求和得到插值函数；插值函数的表达式为：

其中

为横向电场或磁场，

为步骤2中得到的矢量基函数，P_ti为沿着顶点i对应的边的电场或磁场分量。

对矢量波动方程进行变分泛函处理的过程为：

根据矢量波动方程得到电场和磁场的方程；基于变分原理，利用均匀波导TE波没有E_z分量和TM波没有H_z分量的特点对电场和磁场方程进行变分泛函处理。具体的：

在均匀填充波导的全波分析中，边值问题由矢量波动方程定义：

或

其中k₀为真空中电磁波的传播常数，μ_r为波导中填充介质的相对磁导率，ε_r为波导中填充介质的相对介电常数，由于实际波导的边界都是导电壁，所以对于电场来说，其边界条件为

对于磁场来说，其边界条件为

其中

为边界法向分量，根据变分原理，将该边值问题等效为下列变分问题：

或

其中

和

的函数表达式相同，将其统一用

来表示，得到表达式如下：

其中

为电场或磁场，*为共轭符号，k₀为真空中的波数，μ_r和ε_r分别为波导中填充介质的相对磁导率和相对介电常数。

根据电磁波在均匀波导中的传播特征，即TE波没有E_z分量，将TE波的电场改写为如下形式：

其中

为电场，

为横向电场，β为传播常数，将上式对z求偏微分，即可得到

将将三维微分算子表示为xy坐标系下的微分算子和对z偏微分的和

为xy坐标系下的微分算子，将其代入

得到：

上述泛函中只涉及k₀和β的平方项，不涉及它们的一次项。因此离散上述泛函可以得到由k₀确定传播常数β的广义本征值表达形式。

同理，TM波没有H_z分量，经过推导可得到下列泛函：

由此，得到了均匀波导本征模问题的等效泛函。

步骤4、将插值函数代入变分泛函处理后的电场和磁场方程中，得到均匀波导的广义本征值方程。具体的：

4.1、先将步骤2.4得到的插值函数

改写为行项量和列向量相乘的形式：

其中“T”表示转置，代入泛函

中得到：

其中M为离散的三角形单元总数，矩阵A^e和矩阵B^e的表达式如下：

矩阵A^e和B^e可通过解析计算或数值计算方法计算得到。

4.2、对插值参量进行统一编号得到：

然后求其变分后得到广义本征方程：

[A]{P_t}＝-β²[B]{P_t}

其中矩阵A和矩阵B都是稀疏对称矩阵，且矩阵B是正定矩阵，写成标准的广义本征值方程形式为：

[A]{P}＝λ[B]{P}

4.3、利用IRAM等现有方法直接求解广义本征方程[A]{P}＝λ[B]{P}；当

时，求解的是TM波，无需强加边界条件，直接求得的本征值λ就是-β²。当

时，求解的是TE波，需要考虑边界条件

也就是使边界上的切向电场为0，为此，我们先找到边界棱边，然后删除矩阵A和矩阵B中对应边界棱边编号的行和列，这样就保证了边界上的切向电场为0，求解强加边界条件后的广义本征方程，就能得到TE波的传播常数。如此，在给定k₀后，TE波和TM波的传播常数可分别通过直接利用IRAM等现有方法直接求解广义本征方程得到。最后通过公式可以求出截止波数和截止频率：

其中c为自由空间的光速。

为说明本发明方法的优点，下面给出具体的实施例：

本实施例以圆-矩形同轴波导为例进行建模，如图2所示，矩形长边a＝5.08cm,短边b＝2.54cm，内圆半径r₀＝0.635cm。使用本发明提供的方法计算此圆-矩形同轴波导的截止波数，剖分得到5302条棱边和3468个三角形单元，求解两个广义本征方程共用0.08s，计算结果与已发表文献的结果对比如下：

由此可见，本发明能够准确计算均匀波导各模式的截止波数，且计算精度高，计算时间短，有效解决了传统矢量有限元方法得到的广义本征方程中B矩阵不正定，导致普通方法求解时间长，占用内存大，而改进求解算法复杂不直观的问题。

以上所述，仅为本发明的具体实施方式，本说明书中所公开的任一特征，除非特别叙述，均可被其他等效或具有类似目的的替代特征加以替换；所公开的所有特征、或所有方法或过程中的步骤，除了互相排斥的特征和/或步骤以外，均可以任何方式组合；本领域的技术人员根据本发明技术方案的技术特征所做出的任何非本质的添加、替换，均属于本发明的保护范围。

Claims

1.一种基于完全矢量有限元的均匀波导模式计算方法，其特征在于：包括以下步骤：

根据矢量波动方程得到电场和磁场的方程；基于变分原理，利用均匀波导TE波没有E_z分量和TM波没有H_z分量的特点对电场和磁场方程进行变分泛函处理；

步骤5、求解步骤4得到的广义本征值方程即可得到截止波数和截止频率。

2.根据权利要求1所述的一种基于完全矢量有限元的均匀波导模式计算方法，其特征在于：所述步骤2中三角形离散单元各边对应的矢量基函数的计算过程为：

3.根据权利要求1所述的一种基于完全矢量有限元的均匀波导模式计算方法，其特征在于：所述步骤5中求解步骤4得到的广义本征方程采用了IRAM算法。