有限元方法求解金属衬底介质片对平面波的反射问题

徐源，孙一月

(1.中国传媒大学 理学院，北京100024；2.中国人民武装警察部队工程大学 理学院，陕西 西安710086)

有限元方法求解金属衬底介质片对平面波的反射问题

徐源1，孙一月2

(1.中国传媒大学 理学院，北京100024；2.中国人民武装警察部队工程大学 理学院，陕西 西安710086)

任意平面波可以被分解为只有Z轴向电场的Ez-极化平面波和只有Z轴向磁场的Hz-极化平面波。我们需要研究的是金属衬底介质片对这两种波的反射能力。根据标量的赫尔姆霍兹方程分别得到相应的极化方程，在方程的边界上采用不同边界条件。最后我们通过采用有限元的思想来分析问题，通过Matlab语言实现有限元方法的过程来解决这个问题。文章的最后给出了两个数值实验来验证真解和数值解之间的误差关系保证了算法的正确性和时效性。

极化平面波；有限元；Matlab

1 问题的描述

在图1中，一束均匀的入射波在非均质介质片面上发生了反射。介质片的宽度为L，相对电导率为εr，相对的磁导率为μr都视为x的方程。周围的自由空间介质中εr=μr=1，接下来我们需要考虑的是介质片的反射能力。

众所周知任意平面波可以被分解为一个只有Z轴向电场的Ez-极化平面波和一个只有Z轴向磁场的Hz-极化平面波。因此我们可以充分考虑这两种波的案例。对于Ez-极化平面波例子，其入射波的表达式可以是：

(1)

图1 金属介质片平面波的反射图

E0为一常量，指的是入射场的大小，θ代表图1中的入射角度。显而易见为了满足在与x坐标轴相垂直交界面处场的连续性条件，那么总场必须要有一个共同的系数e-jk0ysinθ。根据标量的赫尔姆霍兹方程使得电场Ez满足：

(2)

对于边界条件我们采用的是齐次的狄利克雷边界条件

Ez|x=0=0

(3)

相似地，对于极化波，入射波可以表示为：

(4)

总的磁场满足标量的赫尔姆霍兹方程：

(5)

对应的纽曼边界条件为：

(6)

2 解的分析

此问题的解决第一步需要把介质片分成很多的薄层，定义为M层(图2)每一层的相对电导率εr和相对磁导率可以近似地看作为常值，用εrm和μrm其中(m=1，2，3，…，M)。因此关于方程2在第m层的解可以表示为：

Ezm=(Amejkxmx+Bme-jkxmx)e-jk0ysinθ

(7)

图2

我们强制Ez和Hy在两层之间的交界面处是连续的，指出第m层和第(m+1)层的时候，我们有：

(8)

这里的xm+1指的是交界面处的位置，Rm=Bm/Am，并且有：

(9)

由于边界条件式3，一旦有A1+B1=0成立，那么R1=-1。根据式8我们又可以计算出R2，R3等等。最后我们就可以计算出我们感兴趣的反射系数RM+1，我们用|RM+1|2来表示介质片反射能力的百分比。相类似的过程可以直接用在磁场的Hz极化波上面：

(10)

由边界条件式6，一旦有A1-B1=0成立，那么R1=1。根据式8我们又可以计算出R2，R3等等。我们就可以计算出我们感兴趣的Hz极化波反射系数RM+1。

3 有限元的解

首先，让我们用有限元的方法来解决这个问题。首先确定问题的求解区域。显然问题的求解区域是半无穷的(0≤x≤∞)，遗憾是有限元方法不能直接的用于这种半有界的求解域。于是我们需要适当的缩小下求解域，我们可以引进适当的人工边界条件。为了满足解的有效性，人工的边界和求解域直接应该足够的接近。所以针对此问题最可能的选择是令x=L，这里的L指的是板的厚度。

接下来要做的是在此边界上找一个适当的边界条件。需要指出的是在板的外部，总场可以表示为入射场和反射场的叠加。因此对于Ez极化波我们有：

Ez(x)=E0ejk0xcosθ+RE0e-jk0xcosθxgt;L

(11)

这里的R代表反射系数以及公共系数e-jk0xcosθ已经表示过了，对Ez进行对x求偏导数，我们有：

(12)

也可以表示为：

(13)

或者

(14)

此条件适用于在μr=1处介质片外表面的边界上。如果希望将边界放在介质片的内表面，那么可以用到式14的边界条件需要指出的是：

这是针对Ez和Hz连续边界条件的结果，因此我们不难发现：

(15)

作为边界条件仅仅在介质片的内表面是可行的。显然式13和式14是等价的，它们的选择并不会影响方程的解。相类似地针对极化波我们同样可以找到适当的边界条件强加在Hz如：

(16)

或

(17)

对于式2和式5的不同，我们分别给出了在x=0处的边界条件和在x=L处推导出的边界条件。我们就可以继续运用有限元方法来求解问题的数值解。我们首先是把介质片或者是区域(0，L)分成了一定数量的薄层，再次将M层每一层都看作是一个单元。因此我们就可以通过编写有限元方法的计算机语言可以高效地解决这个问题。一旦场在x=L处被解决了，那么方程11的反射系数也可以得到，或者更具体地对于Ez极化波的反射系数来说：

Hz极化波来说也是类似，就是把Ez用Hz来代替，E0用H0来代替。

4 数值实验

为了保证有限元解的有效性，我们令电导率εr=4+(2-0.1j)(1-x/L)2，相对应的磁导率μr=2-0.1j。对于介质片的厚度采取五倍于自由空间的波长。实验的结果已经展示在图3。采用50单元(51个节点)和100单元(101个节点)的有限元的数值解和真解的比较。正如我们所看到的节点的数量越多，数值解就越接近真解。

5 实验分析

通过Matlab语言实现的有限元算法来求解金属衬底介质板对平面波的反射问题，我们不难发现，随着求解单元数量的增加，数值解与真实解之间的误差会越来越小。这达到了我们预先的估计和目的，用有限元的方法来解决偏微分方程问题的时效性和准确性优势明显。根据完整的有限元的求解过程构编写程序语言是本文的最终目的。

(a)电场极化算例 (b)磁场极化算例图3 金属衬底介质片的反射系数εr=4+(2-0.1j)(1-x/L)2，μr=2-0.1j，L=5λ

[1]Harris，Edward G.Introduction to modern theoretical physics[M].V I，1975.

[2]Jianming Jin.The finite element method in electromagnetics[M].Canada：John Wiley amp; Sonsinc，2002：83-90.

[3]严登俊，黄学良，胡敏强，勾磊.有限元网格生成技术分析[J].微特电机，1999(01).

[4]傅为农，江建中.关于电磁场数值分析的若干认识[J].微特电机，1998(05).

[5]刘春太，杨晓东，申长雨，陈静波.任意平面区域的变尺寸有限元网格划分[J].计算力学学报，2000(01).

[6]Gerald Recktenwald.数值方法和MATLAB实现与应用[M].北京：机械工业出社，2004.

[7]Snyder A W，Love J D.光波导理论[M].北京：人民邮电出版社，1991.

(责任编辑：宋金宝)

FiniteElementMethodtoSolvePlane-waveReflectionbyaMetal-backedDielectricSlab

XU Yuan1，SUN Yi-yue2

(1.School of Science，Communication University of China，Beijing 100024；2.College of Science，Engineering University of the Chinese People’s Armed Police Force，Xi’an 710086)

It is well known that any plane wave can be decomposed into an Ez-polarized plane wave having only a z-compoment for the electric field and a Hz-polarized plane wave having only a z-compoment for the magnetic field，we need to find the power reflected by the slab for two polarization cases only.With observation，the scalar Helmholtz equation governing the functions，the boundary condition for equations is different.We consider the finite element solution to this problem and write Matlab code to solve the problem. Finally this dissertation gives the running result of this system and analyzes the accuracy of the algorithm and its stabilization.

polarization plane wave；finite element solution；Matlab

O24

A

1673-4793(2017)06-0032-04

2017-07-06

徐源(1990-)，男(汉族)，山东临沂人，中国传媒大学硕士研究生.E-mail：xuyuan@cuc.edu.cn