基于柱面等效源区域分解的电磁散射特性仿真方法

摘要

本发明公开了一种基于柱面等效源区域分解的电磁散射特性仿真方法，步骤如下：建立模型，划分子区域并进行剖分；展开各个子等效柱面和子散射体上的电磁流，并确定每个子等效柱面的等效入射电磁流基函数系数；将非旋转对称体的子等效柱面等效入射电磁流基函数系数转换成RWG基函数系数；确定各个子等效柱面的等效入射电磁流在对应子散射体表面产生的入射电场，确定各个子区域内子散射体表面的散射电流；确定子等效柱面上的散射电磁流，将其作为其它各子区域的入射流，更新各子等效柱面上的散射电磁流系数直至达到平衡；确定远场雷达散射截面积，完成电磁散射特性仿真。本发明提供了一种稳定高效、适用于任意形状金属目标的电磁仿真方法。

说明书

技术领域

本发明涉及电磁仿真技术领域，特别是一种基于柱面等效源区域分解的电磁散射特 性仿真方法。

背景技术

旋转对称体是指绕一轴线旋转对称的三维目标，由于它的结构特性，可以由局部全 域基函数展开，使用旋转对称体矩量法快速分析计算(M.Andreasen,"Scattering from  bodies of revolution,"Antennas and Propagation,IEEE Transactions on,vol.13,pp.303-310, 1965.)，相比较传统的分域RWG基函数的矩量法(Rao,S.M.;Wilton,D.;Glisson,A.W., Electromagnetic scattering by surfaces of arbitrary shape,"Antennas and Propagation, IEEE Transactions on,vol.30,no.3,pp.409,418,May1982)大大节约了计算资源。但实际目 标中有很多目标结构不具有或者仅部分具有旋转对称体结构，不适合单独使用旋转对称 体矩量法分析。Durham,T.E.，Christodoulou，C.G使用旋转对称基函数分析旋转对称 部分，RWG基函数分析非旋转对称部分，旋转对称部分和非旋转对称部分的相互耦合 使用的时两者基函数相互直接耦合的方法仿真（Durham,T.E.;Christodoulou,C.G., "Electromagnetic radiation from structures consisting of combined body of revolution and  arbitrary surfaces,"Antennas and Propagation,IEEE Transactions on,vol.40,no.9, pp.1061,1067,Sep1992），该方法的缺点是两种基函数直接相互作用计算时间和计算内 存的消耗大严重的制约了该方法仿真问题的规模，对于完全不具有旋转对称体结构的目 标仍然不能使用旋转对称体矩量法。

区域分解方法是近年来电磁仿真技术中热门的研究方向。W.-D.Li,W.Hong,and  H.-X.Zhou提出了积分方程的重叠型区域分解方法（W.-D.Li,W.Hong,and H.-X.Zhou, “Integral equation-based overlapped domain decomposition method for the analysis of  electromagnetic scattering of3D conducting objects,”Microwave and Optical Technology  Letters,vol.49,no.2,pp.265-274,2007.），M.K.Li,and W.C.Chew提出了基于等效原 理的区域分解方法（M.K.Li,and W.C.Chew,“Wave-Field Interaction With Complex  Structures Using Equivalence Principle Algorithm,”Antennas and Propagation,IEEE  Transactions on,vol.55,no.1,pp.130-138,2007.），P.Zhen,W.Xiao-Chuan,and L.Jin-Fa 提出了非重叠型区域分解方法（P.Zhen,W.Xiao-Chuan,and L.Jin-Fa,“Integral Equation  Based Domain Decomposition Method for Solving Electromagnetic Wave Scattering From  Non-Penetrable Objects,”Antennas and Propagation,IEEE Transactions on,vol.59,no.9, pp.3328-3338,2011.），这些方法都是建立在积分方程方法的基础上的，它们将原先的 大问题分解成一个个小问题来求解，但是由于问题的维数较多，增加了计算资源的开销， 且操作步骤比较复杂。

发明内容

本发明的目的在于提供一种高效稳定的基于柱面等效源区域分解的电磁散射特性 仿真方法，能够对任意形状目标的电磁散射特性进行快速的仿真。

实现本发明的技术解决方案是：一种基于柱面等效源区域分解的电磁散射特性仿真 方法，步骤如下：

第1步，建立散射体模型及等效柱面模型，划分子区域，对每个子区域内部的子散 射体和等效柱面进行剖分；

第2步，根据第1步中子区域的剖分方式，定义相应的基函数，展开各个子等效柱 面和子散射体上的电磁流，并确定每个子等效柱面的等效入射电磁流旋转对称体基函数 系数；

第3步，对于非旋转对称体子区域的子等效柱面，将对应的等效入射电磁流旋转对 称体基函数系数转换成RWG基函数系数；

第4步，确定各个子等效柱面的等效入射电磁流在对应子散射体表面产生的入射电 场；

第5步，由第4步求得的入射电场，确定各个子区域内子散射体表面的散射电流；

第6步，由第5步得到的各个子区域内子散射体上的散射电流向外辐射，确定对应 子等效柱面上的散射电磁流；

第7步，将第6步的各个子区域内子等效柱面上的散射电磁流，作为除自身之外其 它各子区域的入射流，更新各个子区域内子等效柱面上的入射电磁流系数，并重新计算 第3～7步，更新各子区域内子等效柱面上的散射电磁流系数，直至散射电磁流系数达到 平衡，进入下一步；

第8步，由第7步所得等效柱面的表面散射电磁流，确定远场雷达散射截面积，完 成电磁散射特性仿真。

本发明与现有技术相比其显著效果是：（1）利用了目标结构特性进行区域分解，很 好的结合了传统的旋转对称体方法和多层快速多极子方法的优越性；（2）节约了计算资 源的开销，对于各种类导弹结构的散射特性的分析本发明方法尤为适合；（3）整个方法 理论可靠，易于实现，对于任意形状金属目标的电磁散射特性的快速仿真着重要的意义。

附图说明

图1为散射体柱面等效源区域分解示意图，其中(a)为任意散射体柱面等效源区域 分解示意图，(b)为含部分旋转对称结构的散射体柱面等效源区域分解示意图。

图2为等效柱面剖分示意图，其中(a)为各子散射体块及等效柱面剖分示意图，(b) 为非旋转对称部分三角形面片剖分示意图，(c)为旋转对称部分母线剖分示意图。

图3为RWG基函数示意图。

图4为实施例中导弹模型示意图，(a)为侧视图，(b)为主视图。

图5为实施例中导弹模型柱面等效源区域划分示意图。

图6为实施例中双站雷达散射截面积。

具体实施方式

下面结合附图，以一部分含有旋转对称体结构的散射体电磁散射特性的仿真为例， 对本发明作进一步详细说明。

本发明基于柱面等效源区域分解的电磁散射特性仿真方法，包含以下步骤：

第1步，建立散射体模型及其等效柱面，划分子区域，对每个子区域内部的子散射 体和子等效柱面进行剖分。

（1.1）建立一个可以完全包围待求散射体的封闭等效柱面，若散射体结构中部分 为旋转对称体，那么等效柱面的轴向和散射体旋转对称部分的轴向一致；

（1.2）划分子区域：若模型不包含旋转对称体时，将等效柱面分解成尺寸完全相 同的子等效柱面，每个子等效柱面及其包围的子散射体形成对应的子区域，如图1（a） 所示；若模型是由旋转对称体部分和非旋转对称部分组合而成的散射体时，如图1（b） 所示，则将等效圆柱面等分成两种子等效柱面：一种子等效柱面包围旋转对称体部分， 为旋转对称体子区域；另一种子等效柱面包围非旋转对称部分，为非旋转对称体子区域。

（1.3）对每个子区域进行剖分，如图2（a）所示；对于非旋转对称体子区域的子 等效柱面，分别进行三角形面片剖分和母线线段剖分，对于非旋转对称体子区域的子散 射体进行三角形面片剖分，图2(b)为三角形面片剖分示意图；对于旋转对称体子区域的 子等效柱面和子散射体的母线进行线段剖分，图2(c)为母线线段剖分示意图；两个子等 效柱面相交的面使用一致的线段剖分，保证两个子区域电磁流的连续性，记录各个子区 域的剖分信息。

第2步，根据第1步中子区域的剖分方式，定义相应的基函数，展开各个子等效柱 面和子散射体上的电磁流，并确定每个子等效柱面的等效入射电磁流旋转对称体基函数 系数。

（2.1）定义基函数：

（a）在三角形面片剖分单元上定义RWG基函数

其中n为RWG基函数编号，即为三角形剖分单元中第n条边，l_n为第n条边的长 度，表示第n条边对应的上三角形的面积，表示第n条边对应的下三角形的 面积，表示第n条边对应的上三角形中r点指向自由顶点（即除去边n两个顶点剩下 的三角形顶点）的方向向量，表示第n条边对应的下三角形中自由顶点指向r点的方 向向量，如图3(a)所示。

则散射体表面的电磁流用RWG基函数展开为：

$J\left(r\right)\approx \underset{n=1}{\overset{N^{\mathrm{rwg}}}{\Sigma }}{a}_n^{\mathrm{rwg}}{f}_n^{\mathrm{rwg}}\left(r\right)$              （2）

$M\left(r\right)\approx \underset{n=1}{\overset{N^{\mathrm{rwg}}}{\Sigma }}{b}_n^{\mathrm{rwg}}{f}_n^{\mathrm{rwg}}\left(r\right)$

其中J(r)表示散射体表面任意一点r点的电流，M(r)表示r点的磁流，表示 电流J(r)使用RWG基函数展开对应于第n个基函数的展开系数，表示磁流 M(r)使用RWG基函数展开对应于第n个基函数的展开系数，N^rwg表示等效柱 面使用RWG基函数展开后未知量的个数。

（b）在母线剖分线段上定义旋转对称体基函数，该基函数在母线方向上为三角 形基函数形式，在周向方向上为一个指数函数形式，如图3（b）所示：

$f_{\alpha n^{\prime }}^{\mathrm{BoR},t}\left(r\right)=\frac{T_{n^{\prime }}\left(t\right)}{\rho \left(r\right)}{e}^{\mathrm{j\alpha \phi }}\hat{t}\left(r\right)$，n'＝1,…,N^BoR（3）

$f_{\alpha n^{\prime }}^{\mathrm{BoR},\phi }\left(r\right)=\frac{T_{n^{\prime }}\left(t\right)}{\rho \left(r\right)}{e}^{\mathrm{j\alpha \phi }}\hat{\phi }\left(r\right)$

其中上标BoR表示旋转对称体，T_n'(t)表示三角基函数，是一维局部基函数，T_n'(t) 定义在两条相连接的剖分线段上，我们分别称这两条线段为前段和后段，其表达式为 T_n'(t)表示三角基函数，其表达式为：

其中n'为BoR三角基函数编号，t表示r点的母线方向分量；ρ(r)表示r点到旋转 对称体旋转轴的垂直距离；φ表示r点的周向角；e^jαφ表示傅里叶级数展开对应于第α个 模式的指数项；表示r点的母线方向；表示r点的周向方向；N^BoR表示旋转对称 体基函数对应的未知量的个数；表示第n个三角基函数对应的前段的起点切向分量， 表示第n个三角基函数对应的前段的终点切向分量即后段的起点切向分量，表示 第n个三角基函数对应的后段的终点切向分量；Δ_n表示前段的长度，Δ_n+1表示后段的长 度；

$J\left(r\right)=\underset{\alpha =-\infty }{\overset{\infty }{\Sigma }}\underset{n^{\prime }=1}{\overset{N^{\mathrm{BoR}}}{\Sigma }}[a_{\alpha n^{\prime }}^{\mathrm{BoR},t}{f}_{\alpha n^{\prime }}^t\left(r\right)+a_{\alpha n^{\prime }}^{\mathrm{BoR},\phi }{f}_{\alpha n^{\prime }}^{\phi }\left(r\right)]$            （5）

$M\left(r\right)=\underset{\alpha =-\infty }{\overset{\infty }{\Sigma }}\underset{n^{\prime }=1}{\overset{N^{\mathrm{BoR}}}{\Sigma }}[b_{\alpha n^{\prime }}^{\mathrm{BoR},t}{f}_{\alpha n^{\prime }}^t\left(r\right)+b_{\alpha n^{\prime }}^{\mathrm{BoR},\phi }{f}_{\alpha n^{\prime }}^{\phi }\left(r\right)]$

其中表示电流J(r)使用旋转对称体基函数展开对应于第α个模式的第n' 个基函数母线方向的展开系数；表示电流J(r)使用旋转对称体基函数展开 对应于第α个模式的第n'个基函数周向方向的展开系数；表示磁流M(r)使用旋 转对称体基函数展开对应于第α个模式的第n'个基函数母线方向的展开系数； 表示磁流M(r)使用旋转对称体基函数展开对应于第α个模式的第n'个基函 数周向方向的展开系数。

（2.2）计算每个子等效柱面的等效入射电磁流旋转对称体基函数系数：

均匀平面波照射在子等效柱面：

$E^{\mathrm{inc}}\left(r\right)=\hat{E}{e}^{-\mathrm{jk}·r}$

（6）

$H^{\mathrm{inc}}\left(r\right)=\frac{1}{\eta }\hat{k}\times E\left(r\right)$

E^inc为入射电场，H^inc为入射磁场，为平面波传播方向单位矢量，k为传播常数， 表示电场方向单位矢量，η为平面波阻抗。下面使用c标记子等效柱面，p标记子散 射体，则柱面i上等效入射电流等效入射磁流分别为：

$J_{c,i}^{\mathrm{inc}}=-\hat{n}\times H_c^{\mathrm{inc}}$           （7）

$M_{c,i}^{\mathrm{inc}}=\hat{n}\times E_c^{\mathrm{inc}}$

定义旋转对称体测试函数为基函数的共轭：

其中下标m'表示旋转对称体测试函数编号，β表示旋转对称体测试函数模式 数编号，表示对应于切向方向第β个模式的第m'个测试函数，表示 对应于周向方向的第β个模式第m'个测试函数，表示第i个子等效柱面上旋 转对称体基函数对应的未知量的个数，Mod表示需要的总模式数；

使用公式（8）对公式（7）测试，生成子等效柱面的激励向量

V_c,i＝[V_c,i,-Mod,…,,V_c,i,β…,V_c,i,Mod-1,V_c,i,Mod]^T：

$V_{c,i,\beta }=\left(\begin{array}{c}⟨w_{m^{\prime }\beta ,i}^{\mathrm{BoR},t},J_{c,i}^{\mathrm{inc}}⟩\\ ⟨w_{m^{\prime }\beta ,i}^{\mathrm{BoR},\phi },J_{c,i}^{\mathrm{inc}}⟩\\ ⟨w_{m^{\prime }\beta ,i}^{\mathrm{BoR},t},M_{c,i}^{\mathrm{inc}}⟩\\ ⟨w_{m^{\prime }\beta ,i}^{\mathrm{BoR},\phi },M_{c,i}^{\mathrm{inc}}⟩\end{array}\right)---\left(9\right)$

其中i＝1,…,N，N表示子区域总数，β＝-Mod,…,Mod，确定所有子区域的子 等效柱面的等效入射电磁流。

确定第i个子区域的子等效柱面的等效入射电磁流旋转对称体基函数系数求解矩 阵(U^BoR_i)^-1，其中U^BoR_i的元素由下式给出：

${[U^{\mathrm{BoR}}}_{i,m^{\prime }{n}^{\prime }}]=[⟨w_{\beta m^{\prime },i}^{\mathrm{BoR}},f_{\alpha n^{\prime },i}^{\mathrm{BoR}}⟩]---\left(10\right)$

其中α表示基函数模式数编号，β表示测试 函数模式数编号，化简（10）式：

$\left[{U^{\mathrm{BoR}}}_{i,m^{\prime }{n}^{\prime }}\right]=<w_{\beta m^{\prime },i},f_{{\mathrm{\alpha n}}^{\prime },i}>={\int }_{t_m}{\int }_0^{2\pi }{\mathrm{\rho w}}_{m^{\prime },i}\left(t\right)f_{n^{\prime },i}\left(t\right)e^{\mathrm{j\phi }(\alpha -\beta )}\mathrm{d\phi dt}$    （11）

$=\underset{q=1}{\overset{M_p}{\Sigma }}\frac{\pi }{4{\rho }_{i,q}}{T}_{m^{\prime }q,i}{T}_{n^{\prime }q,i}{\Delta }_{q,i}$

其中T_m'q,i、T_n'q,i为三角基函数，由公式（4）给出，T_m'q,i表示第i个子区域对应于 第q条线段的第m'个测试函数、T_n'q,i表示第i个子区域对应于第q条线段的第n'个三角 基函数，M_q表示一个基函数包含的剖分线段数，这里取2，Δ_q,i表示第i个子区域对应 的的第q条剖分线段的长度，由于三角基函数只和相邻的两个基函数有重叠，因此矩阵 是一个三带条型矩阵，可以通过追赶法求解得到电流系数求解矩阵

则子区域i的子等效柱面模式β对应的等效入射电磁流系数向量为：

$\left(\begin{array}{c}{a}_{c,i,\beta }^{t,\mathrm{inc}}\\ a_{c,i,\beta }^{\phi ,\mathrm{inc}}\\ b_{c,i,\beta }^{t,\mathrm{inc}}\\ b_{c,i,\beta }^{\phi ,\mathrm{inc}}\end{array}\right)=\left(\begin{array}{c}{\left(U_i^{\mathrm{BoR}}\right)}^{-1}⟨w_{m^{\prime }\beta ,i}^{\mathrm{BoR},t},J_{c,i}^{\mathrm{inc}}⟩\\ {\left(U_i^{\mathrm{BoR}}\right)}^{-1}⟨w_{m^{\prime }\beta ,i}^{\mathrm{BoR},\phi },J_{c,i}^{\mathrm{inc}}⟩\\ {\left(U_i^{\mathrm{BoR}}\right)}^{-1}⟨w_{m^{\prime }\beta ,i}^{\mathrm{BoR},t},M_{c,i}^{\mathrm{inc}}⟩\\ {\left(U_i^{\mathrm{BoR}}\right)}^{-1}⟨w_{m^{\prime }\beta ,i}^{\mathrm{BoR},\phi },M_{c,i}^{\mathrm{inc}}⟩\end{array}\right)---\left(12\right)$

其中表示第i个子等效柱面c上的母线方向t对应的模式β等效入射电流系数 向量；表示第i个子等效柱面c上的周向φ对应的模式β等效入射电流系数向量； 表示第i个子等效柱面c上的母线方向t对应的模式β等效入射磁流系数向量；表示第i个子等效柱面c上的周向φ对应的模式β等效入射磁流系数向量。

第3步，对于非旋转对称体子区域的子等效柱面，将对应的等效入射电磁流旋转对 称体基函数系数转换成RWG基函数系数；

（3.1）对于非旋转对称体子区域的子等效柱面，柱面等效入射电磁流使用RWG 基函数表示，通过旋转对称体基函数系数可以确定第i个子等效柱面c上的任意一点r处 的电流J_c,i(r)和磁流M_c,i(r)：

$J_{c,i}\left(r\right)=\underset{\beta =-\infty }{\overset{\infty }{\Sigma }}\underset{n^{\prime }=1}{\overset{{N_i}^{\mathrm{BoR}}}{\Sigma }}[a_{c,i,{\mathrm{\beta n}}^{\prime }}^{\mathrm{BoR},t}{f}_{{\mathrm{\beta n}}^{\prime }}^t\left(r\right)+a_{c,i,{\mathrm{\beta n}}^{\prime }}^{\mathrm{BoR},\phi }{f}_{{\mathrm{\beta n}}^{\prime }}^{\phi }\left(r\right)]$                （13）

$M_{c,i}\left(r\right)=\underset{\beta =-\infty }{\overset{\infty }{\Sigma }}\underset{n^{\prime }=1}{\overset{{N_i}^{\mathrm{BoR}}}{\Sigma }}[b_{c,i,{\mathrm{\beta n}}^{\prime }}^{\mathrm{BoR},t}{f}_{{\mathrm{\beta n}}^{\prime }}^t\left(r\right)+b_{c,i,{\mathrm{\beta n}}^{\prime }}^{\mathrm{BoR},\phi }{f}_{{\mathrm{\beta n}}^{\prime }}^{\phi }\left(r\right)]$

其中表示旋转对称体基函数母线方向方向的转对称体基函数，表示旋转 对称体基函数周向方向的转对称体基函数，表示第i个区域柱面c上的电流J_c,i(r) 使用旋转对称体基函数展开对应于第β个模式的第n'个基函数母线方向的展开 系数，表示表示第i个区域柱面c上的电流J_c,i(r)使用旋转对称体基函数展开 对应于第β个模式的第n'个基函数周向方向的展开系数，表示第i个区域柱面c上 的磁流M_c,i(r)使用旋转对称体基函数展开对应于第β个模式的第n'个基函数母线 方向的展开系数，表示表示第i个区域柱面c上的磁流M_c,i(r)使用旋转对称体基函 数展开对应于第β个模式的第n'个基函数周向方向的展开系数。

（3.2）由于无论使用哪种基函数表示，任意一点的电磁流的值是固定的，RWG基 函数展开电磁流的表达式是：

$J_{c,i}\left(r\right)\approx \underset{n=1}{\overset{N_i^{\mathrm{rwg}}}{\Sigma }}{a}_{c,i,n}^{\mathrm{rwg}}{f}_n^{\mathrm{rwg}}\left(r\right)$               （14）

$M_{c,i}\left(r\right)\approx \underset{n=1}{\overset{{N_i}^{\mathrm{rwg}}}{\Sigma }}{b}_{c,i,n}^{\mathrm{rwg}}{f}_n^{\mathrm{rwg}}\left(r\right)$

其中J_c,i(r)、M_c,i(r)与公式（13）中含义相同，表示第i个子区域的子等效柱面c上 电流RWG基函数展开系数，表示第i个子区域的子等效柱面c上磁流RWG基函数 展开系数，表示第i个子等效柱面使用RWG基函数展开后未知量的个数。

（3.3）选择RWG基函数为测试函数对公式（14）两边进行测试，得到如下方程：

$\underset{i=1}{\overset{N}{\Sigma }}{a}_{c,i,n}^{\mathrm{rwg}}⟨f_{n,i}^{\mathrm{rwg}}\left(r^{\prime }\right),f_{m,i}^{\mathrm{rwg}}\left(r\right)⟩=⟨J_{c,i}\left(r^{\prime }\right),f_{m,i}^{\mathrm{rwg}}\left(r\right)⟩$       （15）

$\underset{i=1}{\overset{N}{\Sigma }}{b}_{c,i,n}^{\mathrm{rwg}}⟨f_{n,i}^{\mathrm{rwg}}\left(r^{\prime }\right),f_{m,i}^{\mathrm{rwg}}\left(r\right)⟩=⟨M_{c,i}\left(r^{\prime }\right),f_{m,i}^{\mathrm{rwg}}\left(r\right)⟩$

将公式（1）带入方程（15）并写出矩阵方程形式如下：

${U^{\mathrm{rwg}}}_{c,i}{I}_{c,i}={\bar{V}}_{c,i}---\left(16\right)$

I_c,i为待求RWG基函数系数向量、为已知激励向量，(U^rwg_c,i)^-1为RWG基函数 系数求解矩阵，U^rwg_c,i具体元素为：

$\left[{U^{\mathrm{rwg}}}_{c,i,m,n}\right]=[⟨f_{n,i}^{\mathrm{rwg}}\left(r^{\prime }\right),f_{m,i}^{\mathrm{rwg}}\left(r\right)⟩]---\left(17\right)$

m表示RWG基函数测试函数编号，n表示RWG基函数编号，矩阵U_c,i,m,n为稀疏 阵，每行每列包含5个非零元素，但对于未知量很大的情况，计算时间长和内存消耗大。 为此我们利用RWG基函数是分域基函数，分别在每个三角形面片上建立方程（17），逐 个三角形求解RWG基函数的系数，则U_c,i^-1被划分成N_triangle个3×3的子矩阵求逆的问题， N_triangle为三角形个数。将这个求逆问题计算复杂度降低到了O(N)。则RWG基函数系 数可以由求得。

同理，我们也可以通过同样的理论由RWG基函数的系数计算出旋转对称体基函数 的系数，使用旋转对称体测试函数测试公式（13），并将公式（14）带入即能得到系数 求解方程，其系数求解矩阵在第2步中已经给出，不再赘述。

第4步，确定各个子等效柱面的等效入射电磁流在对应子散射体表面产生的入射电 场。

电流和磁流在空间中产生的电场E(r)和磁场H(r)可以由下式计算得到，其中k₀为 自由空间的传播常数：

$E\left(r\right)=-{\mathrm{jk}}_0\mathrm{\eta A}\left(J\right)+\frac{\eta }{{\mathrm{jk}}_0}▿\Psi \left(J\right)-C\left(M\right)$           （18）

$H\left(r\right)=-\frac{{\mathrm{jk}}_0}{\eta }A\left(M\right)+\frac{1}{{\mathrm{jk}}_0\eta }▿\Psi \left(M\right)+C\left(J\right)$

其中

$A\left(x\right)=\underset{\partial \Omega }{\int }\mathrm{xgd}{s}^{\prime }$

$\psi \left(x\right)=\underset{\partial \Omega }{\int }(▿·x)\mathrm{gd}{s}^{\prime }$                   （19）

$C\left(x\right)=\mathrm{pv}\underset{\partial \Omega }{\int }x\times {▿}^{\prime }\mathrm{gd}{s}^{\prime }$

其中x代表电流J或磁流M，A(x)表示矢量位，ψ(x)表示标量位，C(x)表示旋 度场，pv表示主值积分。g表示源点r'到场点r的格林函数，其表达式为：

$g(r,r^{\prime })=\frac{e^{-\mathrm{jk}|r-r^{\prime }|}}{|r-r^{\prime }|}---\left(20\right)$

定义微积分算子L、K：

$L\left(x\right)={\mathrm{jk}}_0\mathrm{\eta A}\left(x\right)-\frac{\eta }{{\mathrm{jk}}_0}▿\Psi \left(x\right)---\left(21\right)$

K(x)＝C(x)

则：

$\left(\begin{array}{c}E\\ H\end{array}\right)=\begin{array}{c}\left(\begin{array}{cc}-L& -\mathrm{\eta K}\\ K& -\frac{1}{\eta }L\end{array}\right)\end{array}\left(\begin{array}{c}J\\ \frac{1}{\eta }M\end{array}\right)---\left(22\right)$

根据公式（22），我们可以计算出子等效柱面的等效入射电磁流在其内部子散射体 上产生的电场：

$E_{P,i}=\left(\begin{array}{cc}-L_{\mathrm{pc},i}& -K_{\mathrm{pc},i}\end{array}\right)\left(\begin{array}{c}{J}_{c,i}\\ M_{c,i}\end{array}\right)---\left(23\right)$

E_P,i表示第i个子区域的子等效柱面等效入射电流J_c,i和等效入射磁流M_c,i在对应 内部子散射体P上产生的电场。

旋转对称体子区域的子散射体表面使用旋转对称体测试函数离散公式（23），得到 由旋转对称体基函数表示的入射场E_PBoR,i；非旋转对称体子区域的子散射体表面使用 RWG测试函数确定入射场E_PRWG,j。

第5步，由第4步求得的入射电场，确定各个子区域内子散射体表面的散射电流。

由电场积分公式：

L_pp,i(J)|_t＝-E_p,i|_t    （24）

L_pp,i(J)具体形式由公式（21）给出，下标pp表示积分算子L作用在第i个子区域 内部的子散射体p上。

若内部为旋转对称体，则使用旋转对称体矩量法建立阻抗矩阵方程，母线剖分段数 在500以内可以使用LU分解直接计算各个模式阻抗矩阵方程的逆，大于500则可以使 用矩阵行列抽取技术加速迭代求解方程，求得旋转对称子散射体表面的散射电流系数 散射电流系数由Mod个子模式电流组成β表示模式数；

若内部为非旋转对称散射体，则使用基于RWG基函数的矩量法求解，对每个子区 域以柱面大小建立八叉树，利用多层快速多极子技术快速迭代求解子散射体表面的电流 系数。多层快速多极子求解散射体电磁散射技术是一个较为成熟的技术，这里不再给出 详细过程，求解到的非旋转对称子散射体表面的散射电流系数

第6步，由第5步得到的各个子区域内子散射体上的散射电流向外辐射，确定对应 子等效柱面上的等效散射电磁流。

由公式（7）和公式（18）可以得到如下关系：

$\left(\begin{array}{c}{J}_{c,i}^{\mathrm{sca}}\\ \frac{1}{\eta }{M}_{c,i}^{\mathrm{sca}}\end{array}\right)=\begin{array}{c}\end{array}\left(\begin{array}{c}-\hat{n}\times H_{c,i}^{\mathrm{sca}}\\ \frac{1}{\eta }\hat{n}\times E_{c,i}^{\mathrm{sca}}\end{array}\right)=\left(\begin{array}{c}-\hat{n}\times K_{\mathrm{cp}}\\ \frac{1}{\eta }\hat{n}\times L_{\mathrm{cp}}\end{array}\right)\left[j_{p,i}^{\mathrm{sca}}\right]---\left(25\right)$

其中表示的是内部子散射体上的散射电流系数，若内部为旋转对称结构子散射体则 该电流系数对应为若内部为非旋转对称结构子散射体则该电流系数对应为

（1）确定具有旋转对称结构子散射体对应的子等效柱面上的散射电磁流：

旋转对称体子散射体上的散射电流可以由第5步求得的散射电流系数表示 $J_{\mathrm{PBoR},i}^{\mathrm{sca}}\left(r\right)=\underset{\beta =-\infty }{\overset{\infty }{\Sigma }}\underset{n^{\prime }=1}{\overset{N^{\mathrm{BoR}}}{\Sigma }}[a_{\mathrm{PBoR},i,\beta n^{\prime }}^{\mathrm{sca},t}{f}_{\beta n^{\prime }}^t\left(r\right)+a_{\mathrm{PBoR},i,\beta n^{\prime }}^{\mathrm{sca},\phi }{f}_{\beta n^{\prime }}^{\phi }\left(r\right)],$为的元素，β 为模式数，n'为基函数编号，见公式（3）。

使用旋转对称体测试函数对（25）两边测试：

$\left(\begin{array}{c}{S}_{c,i,\beta ,m^{\prime }}^{\mathrm{jt}}\\ S_{c,i,\beta ,m^{\prime }}^{\mathrm{j\phi }}\\ S_{c,i,\beta ,m^{\prime }}^{\mathrm{mt}}\\ S_{c,i,\beta ,m^{\prime }}^{\mathrm{m\phi }}\end{array}\right)=\left(\begin{array}{c}⟨w_{\beta m^{\prime }}^m\left(r\right),J_{c,i,\beta ,m^{\prime }}^{\mathrm{sca}}⟩\\ ⟨w_{\beta m^{\prime }}^{\phi }\left(r\right),J_{c,i,\beta ,m^{\prime }}^{\mathrm{sca}}⟩\\ ⟨w_{\beta m^{\prime }}^t\left(r\right),\frac{1}{\eta }{M}_{c,i,\beta ,m^{\prime }}^{\mathrm{sca}}⟩\\ ⟨w_{\beta m^{\prime }}^{\phi }\left(r\right),\frac{1}{\eta }{M}_{c,i,\beta ,m^{\prime }}^{\mathrm{sca}}⟩\end{array}\right)$                      （26）

$=\left(\begin{array}{cc}⟨w_{\beta ,m^{\prime }}^t\left(r\right),-\hat{n}\times K_{\mathrm{cp}},f_{\beta ,n^{\prime }}^t\left(r\right)⟩& ⟨w_{\beta ,m^{\prime }}^t\left(r\right),-\hat{n}\times K_{\mathrm{cp}},f_{\beta ,n^{\prime }}^{\phi }\left(r\right)⟩\\ ⟨w_{\beta ,m^{\prime }}^{\phi }\left(r\right),-\hat{n}\times K_{\mathrm{cp}},f_{\beta ,n^{\prime }}^t\left(r\right)⟩& ⟨w_{\beta ,m^{\prime }}^{\phi }\left(r\right),-\hat{n}\times K_{\mathrm{cp}},f_{\beta ,n^{\prime }}^{\phi }\left(r\right)⟩\\ ⟨w_{\beta ,m^{\prime }}^t\left(r\right),\frac{1}{\eta }\hat{n}\times L_{\mathrm{cp}},f_{\beta ,n^{\prime }}^t\left(r\right)⟩& ⟨w_{\beta ,m^{\prime }}^t\left(r\right),\frac{1}{\eta }\hat{n}\times L_{\mathrm{cp}},f_{\beta ,n^{\prime }}^{\phi }\left(r\right)⟩\\ ⟨w_{\beta ,m^{\prime }}^{\phi }\left(r\right),\frac{1}{\eta }\hat{n}\times L_{\mathrm{cp}},f_{\beta ,n^{\prime }}^t\left(r\right)⟩& ⟨w_{\beta ,m^{\prime }}^{\phi }\left(r\right),\frac{1}{\eta }\hat{n}\times L_{\mathrm{cp}},f_{\beta ,n^{\prime }}^{\phi }\left(r\right)⟩\end{array}\right)\left(\begin{array}{c}{a}_{\mathrm{PBoR},i,\beta ,n^{\prime }}^{\mathrm{sca},t}\\ a_{\mathrm{PBoR},i,\beta ,n^{\prime }}^{\mathrm{sca},\phi }\end{array}\right)$

则对应子等效柱面上的等效散射电磁流系数为：

$\left(\begin{array}{c}{a}_{c,i,\beta }^{t,\mathrm{sca}}\\ a_{c,i,\beta }^{\phi ,\mathrm{sca}}\\ b_{c,i,\beta }^{t,\mathrm{sca}}\\ b_{c,i,\beta }^{\phi ,\mathrm{inc}}\end{array}\right)=\left(\begin{array}{c}{\left(U_i^{\mathrm{BoR}}\right)}^{-1}{S}_{c,i,\beta }^{\mathrm{jt}}\\ {\left(U_i^{\mathrm{BoR}}\right)}^{-1}{S}_{c,i,\beta }^{\mathrm{j\phi }}\\ {\left(U_i^{\mathrm{BoR}}\right)}^{-1}{S}_{c,i,\beta }^{\mathrm{mt}}\\ {\left(U_i^{\mathrm{BoR}}\right)}^{-1}{S}_{c,i,\beta }^{\mathrm{m\phi }}\end{array}\right)---\left(27\right)$

其中系数求解矩阵由第二步中求得。表示第i个柱面c上的母线方向对应的模式β等效散射电流系数向量；表示第i个柱面c上的周向对应的模式β 等效散射电流系数向量；表示第i个柱面c上的母线方向对应的模式β等效散射磁 流系数向量；表示第i个柱面c上的周向对应的模式β等效散射磁流系数向量

（2）确定具有非旋转对称结构子散射体对应的子等效柱面上的散射电磁流：

非旋转对称结构子散射体上的散射电流可以由第五步求得的散射电流系数 表示$J_{\mathrm{PRWG},j}^{\mathrm{sca}}\left(r\right)=\underset{n=1}{\overset{N^{\mathrm{rwwg}}}{\Sigma }}{a}_{\mathrm{PRWG},j,n}^{\mathrm{sca}}{f}_n^{\mathrm{rwg}}\left(r\right),$为的元素，n为基函数编号， 见公式（1）。

使用RWG测试函数对（25）两边测试：

$\left(\begin{array}{c}{S}_{c,j,n^{\pm }}^j\\ S_{c,j,n^{\pm }}^m\end{array}\right)=\left(\begin{array}{c}⟨w_{\mathrm{rwg},n^{\pm }}\left(r\right),J_{c,j}^{\mathrm{sca}}⟩\\ ⟨w_{\mathrm{rwg},n^{\pm }}\left(r\right),\frac{1}{\eta }{M}_{c,i,\beta }^{\mathrm{sca}}⟩\end{array}\right)$             （28）

$=\left(\begin{array}{c}⟨w_{\mathrm{rwg},n^{\pm }}\left(r\right),-\hat{n}\times K_{\mathrm{cp}},f_{\mathrm{rwg},n}\left(r\right)⟩\\ ⟨w_{\mathrm{rwg},n^{\pm }}\left(r\right),\frac{1}{\eta }\hat{n}\times L_{\mathrm{cp}},f_{\mathrm{rwg},n}\left(r\right)⟩\end{array}\right)\left[a_{\mathrm{PRWG},j}^{\mathrm{sca}}\right]$

其中n^±表示第j个区域的第n条边对应于上、下三角形各定义一个系数，分别表示上、 下三角形面片上的电磁流系数。

按照三角形循环，使用公式（17）各个三角形求解表示电磁流系数：

$\left(\begin{array}{c}{a}_{c,j,\mathrm{ntr}}^{\mathrm{sca}}\\ b_{c,j,\mathrm{ntr}}^{\mathrm{sca}}\end{array}\right)=\left(\begin{array}{c}{\left(U_{\mathrm{ntr}}^{\mathrm{rwg}}\right)}^{-1}{S}_{c,j,\mathrm{ntr}}^j\\ {\left(U_{\mathrm{ntr}}^{\mathrm{rwg}}\right)}^{-1}{S}_{c,j,\mathrm{ntr}}^m\end{array}\right)---\left(29\right)$

其中表示三角形ntr对应的电磁流系数求解矩阵，表示子等效柱面上三角 形ntr对应的电流、磁流散射向量，由（28）计算得到，表示子等效柱面上三角形 ntr对应的等效散射电流系数，表示子等效柱面上三角形ntr对应的等效散射磁流 系数。

接下来使用电磁流不同基函数系数转换算法将（29）式求得的RWG基函数电磁流 系数转换为旋转对称体基函数系数，方法参考第3步。

第7步，将第6步的各个子区域内子等效柱面上的散射电磁流，作为除自身之外其 它各子区域的入射流，更新各个子区域内子等效柱面上的入射电磁流系数，并重新计算 第3～7步，更新各子区域内子等效柱面上的散射电磁流系数，直至散射电磁流系数达到 平衡，进入下一步；

根据公式（18），各个子区域的子等效柱面上的散射电磁流会在除自身之外的其它 子区域表面形成新的入射场：

$\left(\begin{array}{c}{J}_j^{\mathrm{inc}}\\ \frac{1}{\eta }{M}_j^{\mathrm{inc}}\end{array}\right)=\underset{j,i\ne j}{\overset{N_d}{\Sigma }}{T}_{\mathrm{ji}}^{\mathrm{cc}}·\left(\begin{array}{c}{J}_i^{\mathrm{sca}}\\ \frac{1}{\eta }{M}_i^{\mathrm{sca}}\end{array}\right)=\left(\begin{array}{cc}-\hat{n}\times K& \hat{n}\times \frac{1}{\eta }L\\ -\frac{1}{\eta }\hat{n}\times L& -\hat{n}\times K\end{array}\right)·\left(\begin{array}{c}{J}_i^{\mathrm{sca}}\\ \frac{1}{\eta }{M}_i^{\mathrm{sca}}\end{array}\right)---\left(30\right)$

其中分别表示第6步中求得的第i个子区域柱面上的等效散射电、磁流，分别表示第j个子区域更新的等效入射电、磁流，叠加上第2步得到的 外部均匀平面波产生的等效入射电磁流。重复第3～7步，更新各个子区域内子等效柱面 上的散射电磁流系数，直至各个子区域内子等效柱面表面的散射电磁流达到平衡，不再 跟新，进入下一步。当两个子区域相接时，要保证两边剖分线段一致且两边电磁流相同。

第8步，由第7步所得等效柱面的表面散射电磁流，确定远场雷达散射截面积，完 成电磁散射特性仿真。

下面结合具体的实施例对本发明做进一步详细描述。

实施例1

根据本发明所述方法对模型号导弹模型使用本发明方法建模仿真其电磁散射特性。

某型号导弹模型，外形和尺寸如图4所示。弹径2m，弹长24m，弹身部分旋转对 称体结构，方向舵部分为任意结构部分。工作频率为：300MHz，入射角度为弹头方向 入射(0⁰,0⁰) ，垂直极化。建立一个半径为4.1m，高度为25m的柱面将弹体全部包围。 以5.0m间隔横切柱面成5等分，如图5所示。每个子圆柱具有相同的结构，对其中一 个建立母线的线段剖分和三角形基函数剖分，并将信息公用给其它柱面。内部子散射体 根据内部结构选择线段剖分或三角形贴片剖分，区域1内子散射体使用RWG基函数计 算，区域2到区域5使用旋转对称体基函数计算。按照本发明的计算步骤实施，最终计 算得到双站雷达散射截面积，如图6中所示，结果很好的跟商用仿真软件FEKO吻合， 证明本发明方法的正确性。

综上可知，本发明利用了目标结构特性进行区域分解，很好的结合了传统的旋转对 称体方法和多层快速多极子方法的优越性；节约了计算资源的开销，对于各种类导弹结 构的散射特性的分析本发明方法尤为适合；整个方法理论可靠，易于实现，对于任意形 状金属目标的电磁散射特性的快速仿真着重要的意义。