目录 目录 TOC \o "1-3" \h \z \u HYPERLINK \l "_Toc353517920" 第1章 格林函数 1 HYPERLINK \l "_Toc353517921" 1.1 位函数 1 HYPERLINK \l "_Toc353517922" 1.2 自由空间格林函数 2 HYPERLINK \l "_Toc353517923" 1.2.1 二维情况 2 HYPERLINK \l "_Toc353517924" 1.2.2 三维情况 3 HYPERLINK \l "_Toc353517925" 1.3 半空间格林函数 4 HYPERLINK \l "_Toc353517926" 第2章 无限长柱散射-CFIE 1 HYPERLINK \l "_Toc353517927" 2.1 理想导体柱散射-TMz极化 1 HYPERLINK \l "_Toc353517928" 2.2 理想导体柱散射-TEz极化 3 HYPERLINK \l "_Toc353517929" 2.3 均匀介质柱散射-TMz极化 5 HYPERLINK \l "_Toc353517930" 2.4 均匀介质柱散射-TEz极化 8 HYPERLINK \l "_Toc353517931" 第3章 理想导体目标散射-CFIE 12 HYPERLINK \l "_Toc353517932" 3.1
表
关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf
面积分方程 12 HYPERLINK \l "_Toc353517933" 3.2 线性方程组 12 HYPERLINK \l "_Toc353517934" 3.2.1 RWG基函数和试函数 12 HYPERLINK \l "_Toc353517935" 3.2.2 矩阵元素计算 13 HYPERLINK \l "_Toc353517936" 3.2.3 奇异性处理 16 HYPERLINK \l "_Toc353517937" 3.3 计算散射场 17 HYPERLINK \l "_Toc353517938" 3.4 振荡核在平面多边形上的积分 19 HYPERLINK \l "_Toc353517939" 3.5 数值算例 21 HYPERLINK \l "_Toc353517940" 第4章 均匀介质目标散射-PMCHWT 23 HYPERLINK \l "_Toc353517941" 4.1 表面积分方程 23 HYPERLINK \l "_Toc353517942" 4.2 线性方程组 24 HYPERLINK \l "_Toc353517943" 4.2.1 矩阵元素计算 24 HYPERLINK \l "_Toc353517944" 4.2.2 奇异性处理 25 HYPERLINK \l "_Toc353517945" 4.3 计算散射场 25 HYPERLINK \l "_Toc353517946" 4.4 数值算例 26 HYPERLINK \l "_Toc353517947" 第5章 非均匀介质目标散射-VIE 27 HYPERLINK \l "_Toc353517948" 5.1 体积分方程 27 HYPERLINK \l "_Toc353517949" 第6章 快速多极子 28 HYPERLINK \l "_Toc353517950" 6.1 引言 28 HYPERLINK \l "_Toc353517951" 6.2 快速多极子 29 HYPERLINK \l "_Toc353517952" 6.3 树型结构 31 HYPERLINK \l "_Toc353517953" 6.4 球谐函数 32 HYPERLINK \l "_Toc353517954" 6.5 多层快速多极子 33 格林函数 位函数 对于时谐场,将麦克斯韦方程分为只有电型源和只有磁型源的两部分 ▽× Ee = − jωμHe ▽× Em = − Jm − jωμHm ▽× He = Je + jωεEe ▽× Hm = jωεEm ▽∙ De = qe ▽∙ Dm = 0 ▽∙ Be = 0 ▽∙ Bm = qm 其中E表示电场,H表示磁场,D表示电通量,B表示磁通量,Je表示电流,Jm表示磁流。 由于▽∙ Be = 0,无散场可以用矢量的旋度表示,引入矢量位A,使得Be = ▽× A,可以得到▽× Ee = − jωμHe = − jωBe = − jω▽× A,因此▽× (Ee + jωA) = 0。无旋场可以用标量的梯度表示,引入标量位φe,则Ee + jωA = − ▽φe,因此 Ee = − jωA − ▽φe (1.1.1) 于是▽× He = Je + jωεEe = Je + jωε(− jωA − ▽φe),将▽× He = μ-1▽× (▽× A) = μ-1[▽(▽∙ A) − ▽2A]代入,整理得 ▽2A + k2A = − μJe + ▽(▽∙ A + jωμεφe) (1.1.2) 不妨令▽∙ A + jωμεφe = 0 (Lorentz Gauge Condition),于是有 ▽2A + k2A = − μJe Ee = − jωA − j(ωμε)-1▽(▽∙ A) He = μ-1▽× A 上述方程中,Ee和He可以理解为由矢量位A辐射的电场和磁场。 由于▽∙ Dm = 0,无散场可以用矢量的旋度表示,引入矢量位F,使得Dm = − ▽× F,可以得到▽× Hm = jωεEm = jωDm = − jω▽× F,因此▽× (Hm + jωF) = 0。无旋场可以用标量的梯度表示,引入标量位φm,则Hm + jωF = − ▽φm,因此 Hm = − jωF − ▽φm (1.1.3) 于是▽× Em = − Jm − jωμHm = − Jm − jωμ(− jωF − ▽φm),将▽× Em = − ε-1▽× (▽× F) = − ε-1[▽(▽∙ F) − ▽2F]代入,整理得 ▽2F + k2F = − εJm + ▽(▽∙ F + jωμεφm) (1.1.4) 令▽∙ F + jωμεφm = 0 (Gauge Condition),于是有 ▽2F + k2F = − εJm Hm = − jωF − j(ωμε)-1▽(▽∙ F) Em = − ε-1▽× F 上述方程中,Em和Hm可以理解为由矢量位F辐射的电场和磁场。 由于时谐场是线性的,因此对于既有电型源也有磁型源的时谐场,可以将上述的结论叠加得到。于是有 ▽2A + k2A = − μJe (1.1.5) ▽2F + k2F = − εJm (1.1.6) E = − jωA − j(ωμε)-1▽(▽∙ A) − ε-1▽× F (1.1.7) H = − jωF − j(ωμε)-1▽(▽∙ F) + μ-1▽× A (1.1.8) 上述方程中,E和H可以理解为由矢量位A和F共同辐射的电场和磁场。可见,只要能够表示矢量位,就能够求出空间任意一点的场量。 自由空间格林函数 在自由空间中,矢量位A满足方程▽2A + k2A = − μJe,矢量位F满足方程▽2F + k2F = − εJm。因此,可以将A和F分别理解为源Je和源Jm产生的场。不妨设自由空间中点源产生的场为g(r),r为观察点到点源的距离,则g(r)满足方程 (▽2 + k2)g(r) = − δ(r) (1.2.1) 于是矢量位A和F可以用g(r)表示为 (1.2.2) (1.2.3) 在这里,G(r)被称为格林函数,即点源产生的场。由于r = |r − r'|,r为观察点矢量,r'为源点矢量,因此,G(r)也常常写为G(r, r')。 上述结论也可以用数学公式描述。对于矢量波动方程▽2A + k2A = − μJe,可以写成标量形式 ▽2Ax + k2Ax = − μJex (1.2.4) 令 ,其中G(r, r')为格林函数。将Ax(r)代入方程 (1.2.5) 从而整理得到 ▽2G(r, r') + k2G(r, r') = − δ(r − r') (1.2.6) 对于点源,格林函数与方向无关,则G(r, r') = G(r)。 二维情况 对于二维形式,在柱坐标系中将▽2G(r)展开,得 (1.2.7) 由于格林函数与角度无关,因此 ,代入原方程得 (1.2.8) 它有两个通解,分别为CH0(1)(kr'')和CH0(2)(kr''),根据物理意义,得到格林函数的解为G = CH0(2)(kr)。为了确定系数C,我们对一个半径无限小的柱作积分 (1.2.9) 对于左边第一项,有 (1.2.10) 对于左边第二项,有 (1.2.11) 因此, ,二维标量格林函数为 (1.2.12) 二维形式的格林函数表示位于r'处点源产生的柱面波。 对于二维拉普拉斯方程,格林函数为 (1.2.13) 具体推导这里就不详细介绍。 三维情况 对于三维形式,在球坐标系中将▽2G(r)展开,得 (1.2.14) 由于格林函数与角度无关,因此 ,代入原方程得 (1.2.15) 经过整理后,有 ,它有两个通解,分别为Cejkr和Ce-jkr,根据物理意义,得到格林函数的解为G = Ce-jkr。为了确定系数C,我们对一个半径无限小的球作积分 (1.2.16) 对于左边第一项,有 (1.2.17) 对于左边第二项,有 (1.2.18) 因此, ,三维标量格林函数为 (1.2.19) 三维形式的格林函数表示位于r'处点源产生的球面波。 对于三维拉普拉斯方程,格林函数为 (1.2.20) 具体推导这里就不详细介绍。 半空间格林函数 无限长柱散射-CFIE 理想导体柱散射-TMz极化 根据等效原理 (2.1.1) 于是 (2.1.2) TMz极化波入射时,电场Einc(x, y)只有z分量,磁矢量位A(x, y)只有z分量,且 ,则 ,于是 (2.1.3) 由于无限长圆柱是二维问
题
快递公司问题件快递公司问题件货款处理关于圆的周长面积重点题型关于解方程组的题及答案关于南海问题
,任何变量均不是z方向上的函数,引入二维自由空间格林函数 ,可使体积分转变为对横截面的积分 (2.1.4) 由于PEC目标只在外表面存在电流,则对横截面的积分可变为对横截面外围的线积分,方程如下 (2.1.5) 其中 ,为场点和源点之间的距离。 用基函数将电流展开 (2.1.6) 代入积分方程得到 (2.1.7) 其中Bn表示基函数。用测试函数对方程进行测试 (2.1.8) 其中Tn表示测试函数。基函数选择脉冲函数,测试函数选择Delta函数,因此,方程可以表示为 (2.1.9) 进一步写成方程组的形式为 (2.1.10) 其中 (2.1.11) (2.1.12) 对于积分 ,当计算对角线阻抗元素时,积分存在奇异点,此时需要进行奇异性处理。当自变量很小时,零阶Hankel函数可以近似表示为 (2.1.13) 取主值部分 ,带入自作用积分项,得 (2.1.14) 其中,wm为第m条边的边长,γ = 1.781072418…。 计算半径为1 m的理想导体圆柱,入射角度为0°,TMz极化,未知量为360,用矩量法计算的双站RCS与Mie级数解比较如下: (a) (b) 图2‑1 理想导体圆柱对TMz极化波的散射:(a) f = 100 MHz;(b) f = 1 GHz 理想导体柱散射-TEz极化 根据等效原理 (2.2.1) 于是 (2.2.2) TMz极化波入射时,电场Hinc(x, y)只有z分量,J(x, y)只有t分量,于是 (2.2.3) 将 代入方程,得 (2.2.4) 由于无限长圆柱是二维问题,任何变量均不是z方向上的函数,引入二维自由空间格林函数 ,可使体积分转变为对横截面的积分,得 (2.2.5) 根据矢量恒等式A×(B×C) = (A·C)B – (A·B)C,得 (2.2.6) 由于PEC目标只在外表面存在电流,则对横截面的积分可变为对横界面外围的线积分,方程如下 (2.2.7) 其中 ,为场点和源点之间的距离。 由于 ,因此,积分方程可表示为主值积分 (2.2.8) 用基函数将电流展开 (2.2.9) 代入积分方程得到 (2.2.10) 其中Bn表示基函数。用测试函数对方程进行测试 (2.2.11) 其中Tn表示测试函数。由于基函数选择脉冲函数,测试函数选择Delta函数,因此,方程可以表示为 (2.2.12) 进一步写成方程组的形式为 (2.2.13) 其中 (2.2.14) (2.2.15) 计算半径为1 m的理想导体圆柱,入射角度为0°,TEz极化,未知量为360,用矩量法计算的双站RCS与Mie级数解比较如下: (a) (b) 图2‑2 理想导体圆柱对TEz极化波的散射:(a) f = 100 MHz;(b) f = 1 GHz 均匀介质柱散射-TMz极化 根据等效原理 (2.3.1) 于是 (2.3.2) TMz极化波入射时,电场Einc(x, y)只有z分量,磁矢量位A(x, y)只有z分量,且 ,于是 (2.3.3) 由于无限长圆柱是二维问题,任何变量均不是z方向上的函数,引入二维自由空间格林函数 ,可使体积分转变为对横截面的积分 (2.3.4) 根据矢量恒等式A×(B×C) = (A·C)B – (A·B)C,得 (2.3.5) 假定,在物体外表面产生等效电流J和磁流M,物体内表面产生等效电流– J和磁流– M。对于外表面的等效电流和磁流,理解为将整个物体对自由空间的作用等效;对于内表面的等效电流和磁流,理解为将整个自由空间及入射平面波的作用等效。因此,根据边界条件,得到方程 (2.3.6) (2.3.7) 在推导内表面方程的时候注意法线要改变方向。由于采用的是面等效电流和磁流,电流和磁流只存在于外表面,则对横截面的积分可变为对横界面外围的线积分,方程如下 (2.3.8) (2.3.9) 其中 ,为场点和源点之间的距离。 计算半径为1 m的均匀介质圆柱,相对介电常数为4,相对磁导率为1,入射角度为0°,TMz极化,未知量为720,用矩量法计算的双站RCS与Mie级数解比较如下: (a) (b) 图2‑3 均匀介质圆柱对TMz极化波的散射:(a) f = 100 MHz;(b) f = 1 GHz 计算半径为1 m的均匀介质圆柱,相对介电常数为4 − j,相对磁导率为1,入射角度为0°,TMz极化,未知量为720,用矩量法计算的双站RCS与Mie级数解比较如下: (a) (b) 图2‑4 均匀介质圆柱对TMz极化波的散射:(a) f = 100 MHz;(b) f = 1 GHz 均匀介质柱散射-TEz极化 根据等效原理 (2.4.1) 于是 (2.4.2) TEz极化波入射时,磁场Hinc(x, y)只有z分量,电矢量位F(x, y)只有z分量,且 ,于是 (2.4.3) 由于无限长圆柱是二维问题,任何变量均不是z方向上的函数,引入二维自由空间格林函数 ,可使体积分转变为对横截面的积分 (2.4.4) 根据矢量恒等式A×(B×C) = (A·C)B – (A·B)C,得 (2.4.5) 假定,在物体外表面产生等效电流J和磁流M,物体内表面产生等效电流– J和磁流– M。对于外表面的等效电流和磁流,理解为将整个物体对自由空间的作用等效;对于内表面的等效电流和磁流,理解为将整个自由空间及入射平面波的作用等效。因此,根据边界条件,得到方程 (2.4.6) (2.4.7) 在推导内表面方程的时候注意法线要改变方向。由于采用的是面等效电流和磁流,电流和磁流只存在于外表面,则对横截面的积分可变为对横界面外围的线积分,方程如下 (2.4.8) (2.4.9) 其中 ,为场点和源点之间的距离。 计算半径为1 m的均匀介质圆柱,相对介电常数为4,相对磁导率为1,入射角度为0°,TEz极化,未知量为720,用矩量法计算的双站RCS与Mie级数解比较如下: (a) (b) 图2‑5 均匀介质圆柱对TEz极化波的散射:(a) f = 100 MHz;(b) f = 1 GHz 计算半径为1 m的均匀介质圆柱,相对介电常数为4 − j,相对磁导率为1,入射角度为0°,TEz极化,未知量为720,用矩量法计算的双站RCS与Mie级数解比较如下: (a) (b) 图2‑6 均匀介质圆柱对TEz极化波的散射:(a) f = 100 MHz;(b) f = 1 GHz 理想导体目标散射-CFIE 表面积分方程 根据面等效原理,理想导体目标的表面存在等效电流,不存在等效磁流,其散射场的表达式可以表示为 (3.1.1) (3.1.2) 在导体表面存在关系 (3.1.3) (3.1.4) 则可以推出电场积分方程(EFIE)和磁场积分方程(MFIE)为 (3.1.5) (3.1.6) 其中Einc和Hinc为入射电场和入射磁场,η为波阻抗,k为波数,G(r, r’)为自由空间Green格林函数,s0表示自作用积分区域。 EFIE能够精确分析任意目标的散射问题,MFIE只能分析闭合目标的散射问题。而且,EFIE和MFIE都会遇到内谐振现象,即在某些频率点,EFIE和MFIE形成的阻抗矩阵奇异或条件数非常大,从而导致电流存在伪解。通常情况下,EFIE和MFIE的谐振频率不同,因此,为了避免谐振现象,可以引入混合场积分方程(CFIE) (3.1.7) CFIE可以真正地消除内谐振现象。一般情况下,由于MFIE只能分析闭合目标,因此,CFIE也只能分析闭合目标。 线性方程组 RWG基函数和试函数 当采用矩量法分析理想导体目标的表面积分方程时,对于任意形状物体,其表面均可以采用平面三角形贴片来模拟,表面的电流则采用平面RWG基函数,其数学表达式为 (3.2.1) 其中,An±表示相应三角形的面积,ln为边长,各符号定义如下图。 图3‑1 平面RWG基函数示意图 RWG基函数具有两个重要特性。第一个特性是棱边法向分量的连续性,它保证了电流横跨公共边时的连续;第二个特性是两个三角形的基函数散度大小相等、符号相反,该特性保证了与基函数对应的电荷的总和为零。 (3.2.2) 矩量法中测试方法的选取是多样的,较为常用的检验方法包括点匹配法和Galerkin法。点匹配就是选择δ函数作为测试函数,此时积分方程中的积分形式非常简单,并且计算过程最为简单,但计算效果具有一定的局限性。Galerkin测试就是选择基函数本身作为试函数进行检验计算,测试过程虽然繁琐,却最为稳定。当采用RWG函数作为基函数时,Galerkin测试过程常常被用于理想导体目标电磁散射的求解。 当采用Galerkin测试时,测试函数有两种选取方式,一种是电场同向离散化方程,一种是电场异向离散化方程。同向指的是测试函数与基函数同向,即选取Λm作为测试函数;异向指的是测试函数与基函数异向,即选取 作为测试函数。 矩阵元素计算 目标被平面三角形贴片离散之后,表面电流用RWG基函数模拟,并对CFIE实施Galerkin测试,从而形成线性方程组为: Z∙I = V (3.2.3) (3.2.4) 其中,Z为阻抗矩阵,I为RWG基函数的系数,V为目标表面的入射场经测试后形成的右边向量。对于电场同向离散化方程,具体表达式如下: (3.2.5) (3.2.6) , (3.2.7) 电场同向离散化方程中的EFIE和MFIE分别称为TE和TH。对于电场异向离散化方程,具体表达式如下: (3.2.8) (3.2.9) , (3.2.10) 电场异向离散化方程中的EFIE和MFIE分别称为NE和NH。参数α为CFIE系数,用来控制EFIE和MFIE的权重,一般选取范围为[0, 1]。通过求解该线性方程组,得到RWG基函数的系数,便可以计算任意理想导体目标的散射问题。 1. 电场同向离散化方程 EFIE阻抗矩阵元素的表达式可以写成 (3.2.11) 根据矢量恒等式▽∙ (gf) = g▽∙f + f∙▽g,可以得到 (3.2.12) 由于 ,于是 (3.2.13) 同理,根据矢量恒等式▽∙ (gf) = g▽∙f + f∙▽g,可以得到 , (3.2.14) 由于 ,于是 (3.2.15) MFIE阻抗矩阵元素的表达式可以写成 (3.2.16) 右边第二项根据矢量恒等式a ∙ (b × c) = b ∙ (c × a) = c ∙ (a × b),可以继续化简为 (3.2.17) 这样做的好处是可以使得格林函数的梯度仅仅积分一次,降低了填充阻抗矩阵的计算量。 2. 电场异向离散化方程 EFIE阻抗矩阵元素的表达式可以写成 (3.2.18) 根据矢量恒等式▽∙ (gf) = g▽∙f + f∙▽g,可以得到 (3.2.19) 由于 ,于是 (3.2.20) 同理,根据矢量恒等式▽∙ (gf) = g▽∙f + f∙▽g,可以得到 , (3.2.21) 由于 ,于是 (3.2.22) MFIE阻抗矩阵元素的表达式可以写成 (3.2.23) 右边第二项根据矢量恒等式a ∙ (b × c) = b ∙ (c × a) = c ∙ (a × b),可以继续化简为 (3.2.24) 这样做的好处是可以使得格林函数的梯度仅仅积分一次,降低了填充阻抗矩阵的计算量。 奇异性处理 当测试函数和基函数属于同一个三角形单元时,电场积分方程阻抗矩阵元素的填充就会出现奇异,原因是格林函数的分母出现0或者非常接近0的数值,使得数值方法无法正确计算积分值,此时,需要推导积分的解析公式。 由于 ,且 ,因此,奇异积分主要表现为以下几类 (3.2.25) (3.2.26) 代入恒等式 ,并根据散度定理得 (3.2.27) (3.2.28) 上述积分能够推导解析的
计算公式
六西格玛计算公式下载结构力学静力计算公式下载重复性计算公式下载六西格玛计算公式下载年假计算公式
,公式如下 (3.2.29) (m = n) (3.2.30) (m ≠ n) (3.2.31) 其中,a、b、c为三角形的三条边的边长,p为三角形的半周长,A为三角形的面积。 计算散射场 由等效面电流源激发的散射场的表达式为 (3.3.1) 下面推导一下自由空间并矢Green函数 的展开形式。 首先计算▽G(r, r'),令r = |r – r'|,则 (3.3.2) 然后计算▽▽G(r, r') (3.3.3) 最后计算 ,得到 (3.3.4) 对于远区散射场,r → ∞,格林函数可以近似的表示为 (3.3.5) 因此远区散射场的表达式可以简化,从而得到远区散射场的表达式如下 (3.3.6) 这是一个对目标表面电流进行积分的公式。其中,s表示散射体表面积分区域,电流J和 分别表示复振幅和相位。在一般情况下,这两项变化相对平缓,但随着波数k的引入,当k很大时,使得相位变化异常的剧烈,应用一般的数值积分方法(如高斯积分)不容易得到准确的结果。当目标表面用平面三角形单元离散且电流用RWG基函数展开后,可以推导出对每个平面三角形单元积分的解析公式。 振荡核在平面多边形上的积分 如图3‑2所示,假设S是一个具有N条边的多边形,其外形函数如下: (3.4.1) 图3‑2 任意平面多边形几何关系 设其顶点按逆时针方向分别标注为:1、2、3、…、N,其方向与我们积分时所需要的方向一致(其中顶点1同时作为最末点N + 1,同样的,顶点0 = 顶点N)。第n个顶点表示为: ,第n条边位于点n和n + 1之间(与顶点相同,边1 = 边N + 1,边0 = 边N),其向量形式表示为: In = rn+1 - rn (3.4.2) 假设第n条边的参数为: r(t) = [(1 - t)rn + (1 + t)rn+1] / 2 (3.4.3) 其中 ,则dr = [(rn+1 - rn) / 2]dt = Indt / 2。同时,我们引入第n条边的中点,表示为rnc = (rn + rn+1) / 2。 为了求解积分 ,首先将k分解为垂直分量和水平分量,与电磁波中的垂直极化和水平极化定义不同,垂直分量为垂直于多边形平面的分量,水平分量为平行于多边形平面的分量。由于k的垂直分量与多边形平面中的任意点作点乘为一常数,则积分公式可以化简如下: (3.4.4) 下面主要分析积分式 的求解。根据斯托克斯定理可知: (3.4.5) 其中 ,Σ为s的边界。下面分两种情况来讨论积分式的求解: 当k = 0时,令 ,那么,根据斯托克斯定理可知 (3.4.6) 上式便是一个多边形面积公式。 当 时,令 ,则 (3.4.7) 则积分式可以表示为: (3.4.8) 其中j0为第一类0阶球Bessel函数。 对于积分 ,其中r为三角形中任意一点,rv为三角形中基函数对应的顶点,于是 (3.4.9) 其中 (3.4.10) 其中 , (3.4.11) 当k = 0时,S(0, rv) = A(rc – rv),A为多边形面积,rc为多边形中心。 数值算例 计算半径为1 m的理想导体球的双站RCS,未知量为3072,入射电磁波为平面波,频率为300 MHz,角度为(180°, 180°),极化角为0°,矩量法计算结果与Mie级数解比较如下: (a) (b) 图3‑3 理想导体球对平面波的散射,f = 300 MHz:(a) θθ极化;(b) φφ极化。 均匀介质目标散射-PMCHWT 表面积分方程 根据面等效原理,均匀介质目标的外表面存在等效电流和等效磁流,用Je和Jm表示,其散射场的表达式可以表示为 (4.1.1) (4.1.2) 在介质外表面存在关系 , (4.1.3) 则可以推出外表面的积分方程为 (4.1.4) (4.1.5) 其中Einc和Hinc为入射电场和入射磁场,η0为介质外部空间的波阻抗,k0为介质外部空间的波数,G0(r, r’)为自由空间格林函数,s0表示自作用积分区域。 根据面等效原理,均匀介质目标的内表面存在等效电流和等效磁流,且电流和磁流的大小与外表面的相等,方向相反,用−Je和−Jm表示,其散射场的表达式可以表示为 (4.1.6) (4.1.7) 在介质内表面存在关系 , (4.1.8) 在这里注意内表面的法方向与外表面相反,则可以推出内表面的积分方程为 (4.1.9) (4.1.10) 其中η为介质内部空间的波阻抗,k为介质内部空间的波数,G(r, r’)为均匀介质内部空间格林函数,s0表示自作用积分区域。 将描述外表面的电场方程与描述内表面的电场方程相加,并将 去掉,可得 (4.1.11) 将描述外表面的磁场方程与描述内表面的磁场方程相加,并将 去掉,可得 (4.1.12) 由于Je和Jm的数量级存在差别,在实际计算均匀介质的散射时,会造成阻抗矩阵元素值差别很大,从而造成性态很差。因此,这里将得到的电场和磁场方程中的η0作简单处理,即令Je' = η0Je得 (4.1.13) (4.1.14) 这就是求解均匀介质目标电磁散射的PMCHWT方程。 线性方程组 矩阵元素计算 目标被平面三角形贴片离散之后,表面电流和磁流均用RWG基函数模拟,并对PMCHWT实施Galerkin测试,从而形成线性方程组为: (4.2.1) 其中,Z为阻抗矩阵,I为RWG基函数的系数,V为目标表面的入射场经测试后形成的右边向量。对于同向离散,具体表达式如下: (4.2.2) (4.2.3) (4.2.4) (4.2.5) (4.2.6) (4.2.7) 对于异向离散,具体表达式如下: (4.2.8) (4.2.9) (4.2.10) (4.2.11) (4.2.12) (4.2.13) 通过求解线性方程组,分别得到电流和磁流的RWG基函数的系数,便可以计算任意均匀介质目标的散射问题。 奇异性处理 PMCHWT方程的奇异性处理与CFIE中奇异性处理一致。 计算散射场 由等效面电流源激发的散射场的表达式为 (4.3.1) 对于远区散射场,格林函数可以近似的表示为 (4.3.2) 因此远区散射场的表达式可以简化,从而得到远区散射场的表达式如下 (4.3.3) 这是一个对目标表面电流和磁流进行积分的公式。其中,s表示散射体表面积分区域,电流Je、磁流Jm和 分别表示复振幅和相位。在一般情况下,这两项变化相对平缓,但随着波数k0的引入,当k0很大时,使得相位变化异常的剧烈,应用一般的数值积分方法(如高斯积分)不容易得到准确的结果。当目标表面用平面三角形单元离散且电流用RWG基函数展开后,可以推导出对每个平面三角形单元积分的解析公式。 数值算例 非均匀介质目标散射-VIE 体积分方程 快速多极子 引言 快速多极子(FMM: Fast Multi-pole Method)和多层快速多极子方法(MLFMA: MultiLevel Fast Multi-pole Algorithm)是当今最令人瞩目的积分方程数值算法,具有精度可控和高效率的优点,被广泛应用于各种复杂目标的电磁散射分析,并且被美国计算物理学会评为20世纪十大算法之一。 快速多极子方法的数学基础是矢量加法定理,即利用加法定理对积分方程中的格林函数进行处理,通过在角谱空间中展开,利用平面波进行算子对角化,最终将稠密阵与矢量的相乘计算转化为几个稀疏阵与该矢量的相乘计算。其基本原理为:将散射体表面上离散得到的子散射体分组,任意两个子散射体间的互耦根据它们所在组的位置关系而采用不同的方法计算。当场点和源点处于近邻组时,耦合作用采用传统的矩量法直接计算;而当它们处于非近邻组时,则采用聚合、转移和配置三步进行计算,实现场点和源点的分离。对于一个给定的场点组,首先将它的各个非相邻组内所有子散射体产生的的贡献“聚合”到各自的组中心表达;再将这些组的贡献由这些组的组中心“转移”至给定场点组的组中心表达;最后将得到的所有非相邻组的贡献由该组中心“配置”到该组内各子散射体。对于散射体表面上的N个子散射体,直接计算它们的互耦时,每个子散射体都是一个散射中心即为一个单极子,共需数值计算量为O(N2);而应用这种快速多极子方法,任意两个子散射体的互耦由它们所在组的组中心联系。各个组中心就是一个多极子,其数值计算量只为O(N1.5)。对于源点组来说,该组中心代表了组内所有子散射体在其非相邻组产生的贡献;对场点组来说,该组中心代表了来自该组的所有非相邻组的贡献,从而减少了散射中心的数目。 对于电大尺寸目标的散射,其未知量数目N >> 1,此时应用多层快速多极子方法将获得比快速多极子方法更高的效率。多层快速多极子方法是快速多极子方法在多层级结构中的推广。对于N体互耦,多层快速多极子方法采用多层分组计算。即对于附近区强耦合量直接计算;对于非附近区耦合量则用多层快速多极子方法实现。多层快速多极子方法基于树形结构计算,其特点是:逐层聚合、逐层转移、逐层配置、嵌套递推。对于二维情况,它将求解区域用一正方形包围,然后再细分为4个子正方形,该层记为第一层。将每个子正方形再细分为4个更小的子正方形,则得到第二层,此时共有42个正方形。依次类推得到更高层。对于三维情况,则用一正方体包围,第一层得到8个子正方体。随着层数增加,每个子正方体再细分为8个更小的子正方体。显然,对于二维,三维情况,第 层子正方形和子正方体的数目分别为4i,8i。对于散射问题,最高层的每个正方形或正方体的边长为半个波长左右,由此可以确定求解一个给定尺寸的目标散射时多层快速多极子方法所需的层数。由于每层数值计算量均为O(N)量级,共有logN层,所以多层快速多极子方法计算矩阵与矢量相乘的工作量为O(NlogN)量级。内存需求也为O(NlogN)量级。 图6‑1多层快速多极子方法中的分层结构 快速多极子 根据加法定理,标量格林函数可以展开为 (6.2.1) 其中,jl(kd)为第一类l阶球Bessel函数,hl(2)(kD)为第二类l阶球Hankel函数, 为l阶Legendre多项式,且d < D。将等式 (6.2.2) 代入格林函数表达式,得 (6.2.3) 其中 (6.2.4) (6.2.5) L为无穷求和的截断项数。 矩量法求解电磁散射时,阻抗矩阵描述的是源点对场点的作用,并假定场点为rm,源点为rn。描述源指向场的矢量rmn可以分解为 rmn = rqn + rpq + rmp (6.2.6) 其中,rp为场点对应的组中心,rq为源点对应的组中心。令d = rqn + rmp,D = rpq,于是 (6.2.7) 从公式可以看出,场点和源点从绝对值符号中分离。 图6‑2源点对场点的直接作用分解成三部分:聚合、转移、配置 对于电场同向离散的EFIE,当场点和源点处于远场时,其矩阵元素计算公式为 (6.2.8) 将展开后的格林函数代入得 (6.2.9) (6.2.10) (6.2.11) (6.2.12) 对于电场同向离散的MFIE,当场点和源点处于远场时,其矩阵元素计算公式为 (6.2.13) 将展开后的格林函数代入得 (6.2.14) (6.2.15) (6.2.16) (6.2.17) 在上述公式中,α称为转移因子,Fqn和Rmp分别称为聚合因子和配置因子,又称为辐射方向图和接收方向图。由于远场阻抗矩阵元素可以分解为三部分,转移因子α、聚合因子Fqn和配置Rmp,因此,对应的远场矩阵矢量乘操作可以转化为 (6.2.18) 其中,Np表示组p的近场。 树型结构 如图6‑1所示,在快速多极子的树型结构中,第一层有4dim个组,第二层有8dim个组,第i层有(2 × 2i) dim个组,其中dim为问题的维数。比如,对于三维散射问题,第一层有64个组,第二层有512个组,第i层有(2 × 2i) 3个组。为了方便程序访问,对于每个组,都有固定的编号,这个编号称为Morton码。Morton码的原理如图6‑3所示,它是将x方向、y方向及z方向的组的号用二进制表示,从最高位到最低位排列。然后扩展成一个3倍长度的二进制数,依次填充x、y、z从低位向高位的各位二进制数。反之,从Morton编号的组可以快速地查组所在的x、y、z方向的组号。以x、y、z方向第3,4,5组为例,假设它们是在最大组号为7的情况下,则有: 由Morton编号的组查询它的父层组的算法是将组号的二进制数右移三位(二维情况则是二位)。例如 110001101(十进制数397) >> 3 => 110001(十进制数49) 表示第i层的397组,是第i − 1层49组的子组。由Morton编号的组查询它的子层组的算法是将组号的二进制数左移三位(二维情况则是二位)后将最后三位填值为从000(0)到111(7)的数字,则可得全部子组Morton编号。 图6‑3二维情况8×8组的Morton序编号方式示意图,红色数字表示x方向的组号,粉红色数字表示y方向的组号蓝色数字表示组的编号,黑色数字则表示编号的二进制编号,蓝色则表示相应的十进制数。 构造快速多极子的树型结构,需要保存如下信息:第i层组的大小、第i层非空组的个数、非空组的Morton码、非空组的位置、非空组的父层组号、非空组的近场和远场组号、第1层非空组包含的边号和三角形号等。 球谐函数 从辐射和接收方向图的表达式中,可以看出辐射和接收方向图的内存需求为N × Nk × 2,其中N为未知量个数,乘以2表示要同时储存两个正交方向的数值,Nk = 2L2,L为多极子模式数。一般情况下,多层快速多极子的最细层组的尺寸控制在0.2~0.4λ之间,此时,模式数L一般在5~8之间,那么Nk取值在50~128之间。当未知量N非常大的时候,辐射和接收方向图需要耗费相当大的内存。 从辐射、接收方向图的公式中可以得出,方向图是波矢量的 函数,因此,可以考虑用球形调和函数来逼近方向图,达到减少内存的目的。根据最细层组的尺寸在0.2~0.4λ之间,采用球谐函数来逼近辐射和接收方向图,其最高阶数选为2阶,表达式为 (6.4.1) (6.4.2) (6.4.3) (6.4.4) (6.4.5) (6.4.6) (6.4.7) (6.4.8) (6.4.9) 于是,方向图用球谐函数表示为 (6.4.10) (6.4.11) 其中,f和r为球谐函数的系数。因此,存储方向图就转化为存储球谐函数的系数。根据球谐函数的对称性,球谐函数的系数仅仅只需要存储6个,因此,方向图的存储量就降为N × 6 × 3。在这里,最后的2之所以变成3,是因为之前存储两个正交方向的数值,此时需要存储其直角坐标系的三个分量。与不用球谐函数的存储量相比,当未知量N非常大的时候,可以减少大量的内存需求。 当采用球谐函数时,还需要额外增加一些计算量,主要在两个方面。一方面,是在计算方向图时,需要计算球谐函数的系数;另一方面,是在矩阵矢量乘时需要作球谐函数展开的计算。尽管如此,增加的计算量相对于多层快速多极子的整体计算量来说,可以忽略不计。 多层快速多极子