CN105808794A - 分析高超声速飞行目标电磁散射特性的时域积分方程方法 - Google Patents

Abstract

本发明公开了一种分析高超声速飞行目标电磁散射特性的时域积分方程方法。步骤如下：建立时域体面积分方程：根据结构散射特性，被激励目标上的总电场等于入射电场与散射电场之和，入射电场为已知激励；离散积分方程，并进行空间和时间上的测试，建立时域矩阵方程；建立等离子体中电通量密度与电场强度之间的关系；求解时域矩阵方程，分析宽频带内混合目标电磁散射特性。本发明只需离散物体模型，使用自由空间格林函数及递归卷积方法，可以分析含非均匀等离子体介质与金属混合目标，具有很高的计算精度且未知量相对较少；同时由于采用时域方法求解，可以一次计算宽频带各频点的电磁散射特性，相较于需要扫频操作的频域方法，求解消耗时间减少。

Description

分析高超声速飞行目标电磁散射特性的时域积分方程方法

一技术领域

本发明属于分析金属等离子体混合目标电磁散射特性的技术领域，特别是一种分析高超声速飞行目标电磁散射特性的时域积分方程方法。

二背景技术

在现代化战争中，武器系统的雷达散射截面积(RCS)是雷达系统对目标“可观测性”的一个重要指标。开展临近空间高超声速目标电磁特性研究，可为有效探测临近空间高超声速飞行器及武器突防提供技术支持。高超声速飞行目标由于具有很快的飞行速度(5马赫以上)以及较高的飞行高度(20Km以上)，其飞行时与空气摩擦会产生几千摄氏度的气动热，使其周围空气电离而呈离子状态存在。当电离度达到一定程度时，电离气体具有等离子体性质。此时在飞行目标表面附近的包覆流场通常被称为等离子体包覆流场、再入等离子体或等离子体壳套，从而与飞行器本体组成包含金属和等离子体的混合目标。

针对这种含色散等离子体的金属介质混合结构，在时域体面积分方法(K.Aygün,B.Shanker,andE.Michielssen,“Fasttimedomaincharacterizationoffinitesizemicrostripstructures,”Int.J.Num.ModElect.Net.Dev.&Fields,vol.15,no.6,pp.439–457,2002.)内部引入递归卷积的方法进行处理。相较于其他时域方法，时域体面积分方程方法可以精确有效的分析其瞬态电磁散射问题，而传统的时域方法，如时域有限差分方法(FDTD)(J.A.Pereda,L.A.Vielva,A.Vegas,andA.Prieto,“State-spaceapproachforthefdtdformulationfordispersivemedia,”IEEETrans.Magn.,1995,31(3):1602–1605)和时域有限元方法时域有限元方法(D.JiaoandJ.-M.Jin,“Ageneralapproachforthestabilityanalysisofthetime-domainfinite-elementmethodforelectromagneticsimulations,”IEEETrans.AntennasPropag.,vol.50,no.11,pp.1624–1632,Nov.2002.)，需要边界条件的设置，离散的未知量也相对较多，精度受到很大的限制，不能满足日益增加的对高精度的需要。而频域的方法每个频点的阻抗矩阵不同，计算一段频带的电磁特性必须扫频，而频域进行扫频非常费时，脱离实际的需求。

三发明内容

本发明的目的在于提供一种分析高超声速飞行目标电磁散射特性的时域积分方程方法，从而快速得到宽频带内电磁散射特性参数。

实现本发明目的的技术解决方案为：

一种分析高超声速飞行目标电磁散射特性的时域积分方程方法，步骤如下：

步骤1，建立时域体面积分方程：根据结构散射特性，被激励目标上的总电场等于入射电场与散射电场之和，入射电场为已知激励；

步骤2，离散积分方程，并进行空间和时间上的测试，建立时域矩阵方程；

步骤3，建立等离子体中电通量密度与电场强度之间的关系；

步骤4，求解时域矩阵方程，分析宽频带内混合目标电磁散射特性。

本发明与现有技术相比，其显著优点为：(1)可以运用时域积分方法分析含等离子体的金属介质混合目标的瞬态电磁散射特性；(2)计算精度高且未知量少，由于只需要对介质体内部及金属表面进行网格离散，不需要额外的空间未知量，则相对于时域有限差分方法和时域有限元方法可以有较少的未知量并有较高的计算精度。

四附图说明

图1为本发明的算例模型示意图。

图2为本发明算例分析含等离子体介质的混合目标在不同频率点处的双站雷达散射截面(RCS)结果图。

五具体实施方式

下面结合附图对本发明作进一步详细描述。

步骤1，建立时域体面积分方程：根据结构散射特性，被激励目标上的总电场等于入射电场与散射电场之和，入射电场为已知激励。

设待求目标为金属与等离子体混合材料，用调制高斯脉冲激励待求解目标，在混合目标的介质部分产生极化体电流J_V(r,t)，金属表面产生感应电流J_s(r,t)。由于总电场等于入射电场与散射电场之和，得到分析金属与等离子体混合目标的时域电场积分方程形式如下：

介质体上：E(r,t)＝E^inc(r,t)+E^sca(r,t)(1)

金属表面：[E^inc(r,t)+E^sca(r,t)]_tan＝0(2)

其中，E^inc(r,t)表示入射电场，E^sca(r,t)表示散射电场，E(r,t)表示总电场。

散射场具体表达式为：

$E_v^{\mathrm{sca}}(r,t)=-{\mathrm{\μ}}_0{\∫}_V\stackrel{\‾}{G}(r,r^{\′})\·{\∂}_t{J}_V(r^{\′},t),(r\∈V)---\left(3\right)$

$E_s^{\mathrm{sca}}(r,t)=-{\mathrm{\μ}}_0{\∫}_s\stackrel{\‾}{G}(r,r^{\′})\·{\∂}_t{J}_s(r^{\′},t),(r\∈S)---\left(4\right)$

其中，下标V表示体单元，S表示面单元；D(r,t)表示电通量密度、为并矢格林函数，ε₀表示自由空间的介电常数，μ₀表示自由空间的磁导率。

可得电场积分方程(EFIE):

${\∂}_t{E}^{\mathrm{inc}}(r,t)={\∂}_{t}E(r,t)+L\left\{{\∂}_{t}D(r,t)\right\}-{\mathrm{\ϵ}}_{0}L\left\{{\∂}_{t}E(r,t)\right\}+L\left\{{\∂}_t{J}_S(r,t)\right\},(r\∈V)---\left(5\right)$

$\hat{n}\×\hat{n}\×{\∂}_{t}E(r,t)=0,(r\∈S)---\left(6\right)$

其中，L{X}算子表示为：

$L\left\{X\right\}(r,t)\stackrel{\·}{=}\frac{{\mathrm{\μ}}_0}{4\mathrm{\π}}{\∫}_{\mathrm{\ξ}}\mathrm{dv}({\∂}_t^2\stackrel{\‾}{I}-c^2\▿\▿)\frac{X(r^{\′},t-\frac{R}{c})}{R}---\left(7\right)$

这里,ε₀表示自由空间的介电常数，μ₀表示自由空间的磁导率，为单位并矢，r为场点，r′为源点，R＝r-r′为场源点之间位置矢量_，D(r,t)表示电通量密度，R＝|R|＝|r-r′|，表示光速，t表示时间，ξ为积分区域。

步骤2，离散积分方程，并进行空间和时间上的测试，建立时域矩阵方程。

现对介质内部的电通量密度D(r,t)和电场强度E(r,t)在空间上用SWG基函数进行展开，金属表面的电流J_s(r,t)用RWG基函数展开:

$D(r,t)=\underset{n=1}{\overset{N_V}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_t}{\mathrm{\Σ}}}{D}_{j,n}{f}_n\left(r\right)T_j\left(t\right)=\underset{n=1}{\overset{N_V}{\mathrm{\Σ}}}{\tilde{D}}_n\left(t\right)f_n\left(r\right)---\left(8\right)$

$\begin{array}{c}E(r,t)=\underset{n=1}{\overset{N_V}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_t}{\mathrm{\Σ}}}(E_{j,n}^{+}{f}_n^{+}\left(r\right)+E_{j,n}^{-}{f}_n^{-}\left(r\right))T_j\left(t\right)\\ =\underset{n=1}{\overset{N_V}{\mathrm{\Σ}}}({\tilde{E}}_n^{+}\left(t\right)f_n^{+}\left(r\right)+{\tilde{E}}_n^{-}\left(t\right)f_n^{-}\left(r\right))\end{array}---\left(9\right)$

$J_s(r,t)=\underset{n=1}{\overset{N_S}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_t}{\mathrm{\Σ}}}{I}_{j,n}{f}_n^s\left(r\right)T_j\left(t\right)---\left(10\right)$

其中，n表示空间未知量，j表示时间步数，N_s表示空间上总的未知量，N_t表示时间上总的步数，D_j,n表示第n个电通量密度在第j时刻的系数，I_j,n表示第n个面电流在第j时刻的系数，f_n(r)表示SWG基函数，表示RWG基函数，表示上半SWG基函数，表示下半SWG基函数，T_j(t)表示时间基函数，表示在固定时间步上的第n个电通量密度未知量系数，表示第n个上半SWG基函数表示的电场未知量在第j时刻的系数，表示第n个下半SWG基函数表示的电场未知量在第j时刻的系数，表示在固定时间步上的第n个上半SWG基函数表示的电场未知量，表示在固定时间步上的第n个下半SWG基函数表示的电场未知量。

时间上用四阶拉格朗日时间基函数进行展开:

$T\left(t\right)=\left\{\begin{array}{cc}1+\frac{25}{12}t+\frac{35}{24}{t}^2+\frac{5}{12}{t}^3+\frac{1}{24}{t}^4& -1\≤t\≤0\\ 1+\frac{5}{6}t-\frac{5}{6}{t}^2-\frac{5}{6}{t}^3-\frac{1}{6}{t}^4& 0\≤t\≤1\\ 1-\frac{5}{4}{t}^2+\frac{1}{4}{t}^4& 1\≤t\≤2\\ 1-\frac{5}{6}t-\frac{5}{6}{t}^2+\frac{5}{6}{t}^3-\frac{1}{6}{t}^4& 2\≤t\≤3\\ 1-\frac{25}{12}t+\frac{35}{24}{t}^2-\frac{5}{12}{t}^3+\frac{1}{24}{t}^4& 3\≤t\≤4\\ 0& \mathrm{others}\end{array}\right.---\left(11\right)$

对(5)(6)式电场积分方程离散后在空间域上进行伽辽金测试，在时间域上进行点匹配得到时间步进的矩阵方程：

${\tilde{Z}}_0{I}_j^s=F_j^s-\underset{i=1}{\overset{j-1}{\mathrm{\Σ}}}{\tilde{Z}}_i{I}_{j-i}^s,(r^{\′}\∈S)---\left(14\right)$

其中，

$I_j^{\mathrm{vd}}=[D_{j,1},D_{j,2}...I_{j,N_V}]I_j^{\mathrm{ve}\±}=[E_{j,1}^{\±},E_{j,2}^{\±}...E_{j,N_V}^{\±}]I_j^s=[I_{j,1},I_{j,2}...I_{j,N_S}]---(18-20)$

$F_{j,m}^v=<f_m^v\left(r\right),{\&PartialD;}_t{E}^i(r,t)>,F_{j,m}^s=<f_m^s\left(r\right),{\&PartialD;}_t{E}^i(r,t)>---(21-22)$

${\tilde{Z}}_{i,\mathrm{vv}}^{m,n}=\left\{\begin{array}{cc}{Z}_{i,\mathrm{vv}}^{\mathrm{mn}}+(Z_{0,\mathrm{vv}}^{1e+,\mathrm{mn}}-{\mathrm{\&epsiv;}}_0{Z}_{0,\mathrm{vv}}^{e+,\mathrm{mn}}){\mathrm{\&beta;}}_{i,n}^{d+}+(Z_{0,\mathrm{vv}}^{1e-,\mathrm{mn}}-{\mathrm{\&epsiv;}}_0{Z}_{0,\mathrm{vv}}^{e-,\mathrm{mn}}){\mathrm{\&beta;}}_{i,n}^{d-}& 0\&le;i\&le;K\\ Z_{i,\mathrm{vv}}^{\mathrm{mn}}& i>K\end{array}\right.---\left(23\right)$

${\tilde{Z}}_{i,\mathrm{sv}}^{m,n}=\left\{\begin{array}{cc}{Z}_{i,\mathrm{sv}}^{\mathrm{mn}}+{\mathrm{\&epsiv;}}_0{Z}_{0,\mathrm{sv}}^{e+,\mathrm{mn}}{\mathrm{\&beta;}}_{i,n}^{d+}-{\mathrm{\&epsiv;}}_0{Z}_{0,\mathrm{sv}}^{e-,\mathrm{mn}}{\mathrm{\&beta;}}_{i,n}^{d-}& 0\&le;i\&le;K\\ Z_{i,\mathrm{sv}}^{\mathrm{mn}}& i>K\end{array}\right.---\left(24\right)$

$\begin{array}{cc}{Z}_{i,\mathrm{vv}}^{1e\&PlusMinus;,\mathrm{mn}}=<f_m^v\left(r\right),\left\{{\&PartialD;}_t{f}_n^{v\&PlusMinus;}\left(r\right)T_{-i}\left(t\right)\right\}>{|}_{t=0}& {\tilde{Z}}_{i,\mathrm{vs}}^{\mathrm{mn}}\end{array},=<f_m^v\left(r\right),L\left\{{\&PartialD;}_t{f}_n^s\left(r\right)T_{-i}\left(t\right)\right\}>{|}_{t=0}---(25-26)$

$\begin{array}{c}{Z}_{i,\mathrm{sv}}^{\mathrm{mn}}=<f_m^s\left(r\right),L\left\{{\&PartialD;}_t{f}_n^v\left(r\right)T_{-i}\left(t\right)\right\}>{|}_{t=0}\end{array}{Z}_{i,\mathrm{vv}}^{e\&PlusMinus;,\mathrm{mn}}<f_m^v\left(r\right),L\left\{{\&PartialD;}_t{f}_n^{v\&PlusMinus;}\left(r\right)T_{-i}\left(t\right)\right\}>{|}_{t=0}---(31-32)$

$Z_{i,\mathrm{sv}}^{e\&PlusMinus;,\mathrm{mn}}=<f_m^s\left(r\right),{\mathrm{\&epsiv;}}_{0}L\left\{{\&PartialD;}_t{f}_n^{v\&PlusMinus;}\left(r\right)T_{-i}\left(t\right)\right\}>{|}_{t=0}---\left(33\right)$

这里，V表示介质积分区域，S表示金属积分区域，v代表介质，s代表金属。Z，为建立场源之间联系的时域阻抗矩阵。F_j表示第j个时间步的激励，N_t表示时间步的个数，是第j个时间步的电通量密度系数、是第j个时间步的上下半个SWG基函数表示的电场强度的系数，是第j个时间步的面电流的系数。

此处，选用调制高斯平面波作为时域电磁散射计算时的激励：

$E^i(r,t)={\hat{u}}^i\exp [-0.5{(t-t_p-r\&CenterDot;k^i/c)}^2/{\mathrm{\&sigma;}}^2]\cos \left[2\mathrm{\&pi;}{f}_0(t-r\&CenterDot;k^i/c)\right]---\left(34\right)$

其中，kⁱ表示入射波的传播方向，为电场入射矢量，t_p表示脉冲的时延，σ表示时间标准差，f₀为计算中心频率，c为真空中的光速。

步骤3，建立等离子体中电通量密度与电场强度之间的关系。

对于等离子体介质，满足D(r,t)＝ε(r,t)*E(r,t)，则：

γ(r,t)*D(r,t)＝Eⁱ(r,t)-L{D}(r,t)+ε₀L{γ*D}(r,t)(35)

其中，

$\begin{array}{c}\mathrm{\&gamma;}(r,t)*D(r,t)=E(r,t)\\ ={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;}\left(r\right)\mathrm{\&delta;}\left(t\right)*D(r,t)+{\mathrm{\&gamma;}}_0\stackrel{\&OverBar;}{{\mathrm{\&gamma;}}_r}(r,t)*D(r,t)\\ ={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;}\left(r\right)D(r,t)+{\mathrm{\&gamma;}}_0{\&Integral;}_0^t\mathrm{d\&tau;D}(r,t-\mathrm{\&tau;})\stackrel{\&OverBar;}{{\mathrm{\&gamma;}}_r}(r,t)\end{array}---\left(36\right)$

$\mathrm{\&gamma;}(r,t)=F^{-1}\{\mathrm{\&gamma;}(r,\mathrm{\&omega;})=1/\mathrm{\&epsiv;}(r,\mathrm{\&omega;})\}(r,t)={\mathrm{\&gamma;}}_{r,\&infin;}\left(r\right)\mathrm{\&delta;}\left(t\right)+\stackrel{\&OverBar;}{{\mathrm{\&gamma;}}_r}(r,\mathrm{\&tau;})---\left(37\right)$

γ_r,∞＝1/ε_r,∞(38)

γ₀＝1/ε₀(39)

其中，F^-1表示逆傅里叶变换，ε(r,ω)为在r处ω频率下的介电常数，ε_r,∞为在r处无限大频率下的介电常数。γ(r,t)为在r处t时刻的下的磁化率，δ(t)为冲激函数。

$\mathrm{\&epsiv;}(r,\mathrm{\&omega;})={\mathrm{\&epsiv;}}_{\&infin;}\left(r\right){\mathrm{\&epsiv;}}_0+\frac{{\mathrm{\&epsiv;}}_0{\mathrm{\&omega;}}_p^2\left(r\right)}{\mathrm{j\&omega;}{\mathrm{\&upsi;}}_e\left(r\right)-{\mathrm{\&omega;}}^2}---\left(40\right)$

可得：

$\mathrm{\&gamma;}(r,t)={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_r(r,t)={\mathrm{\&gamma;}}_0({\mathrm{\&gamma;}}_{r,\&infin;}\left(r\right)\mathrm{\&delta;}\left(t\right)+\stackrel{\&OverBar;}{{\mathrm{\&gamma;}}_r}(r,t))---\left(41\right)$

现令 ${\mathrm{\&gamma;}}_{r,\&infin;}\left(r\right)=1/{\mathrm{\&epsiv;}}_{r,\&infin;}\left(r\right)=\mathrm{li}{m}_{\mathrm{\&omega;}\&RightArrow;\&infin;}{\widetilde{\mathrm{\&gamma;}}}_r(r,\mathrm{\&omega;})=1$

$\stackrel{\&OverBar;}{{\mathrm{\&gamma;}}_r}(r,t)=\mathrm{Re}\{j=\frac{{\mathrm{\&omega;}}_p^2\left(r\right)}{\sqrt{{\mathrm{\&omega;}}_p^2\left(r\right)-v_e^2\left(r\right)}}\&CenterDot;e^{[-v_e\left(r\right)+j\sqrt{{\mathrm{\&omega;}}_p^2\left(r\right)-v_e^2\left(r\right)}]t}\}U\left(t\right)---\left(42\right)$

其中，υ_e(r)为等离子体碰撞频率，ω_p(r)为等离子体频率。

$U\left(t\right)=\left\{\begin{array}{cc}0& t<0\\ 1& t\&GreaterEqual;0\end{array}\right.$ (单位阶跃函数)

通过对(36)式离散后，可得：

$\begin{array}{c}{\tilde{E}}_n^{\&PlusMinus;}\left(t\right)={\mathrm{\&gamma;}}_n^{\&PlusMinus;}\left(t\right)*{\tilde{D}}_n\left(t\right)\\ ={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;,n}^{\&PlusMinus;}{\tilde{D}}_n\left(t\right)+{\mathrm{\&gamma;}}_0{\&Integral;}_0^t\mathrm{d\&tau;}{\tilde{D}}_n(t-\mathrm{\&tau;}){\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_{r,n}^{\&PlusMinus;}\left(\mathrm{\&tau;}\right)\end{array}---\left(43\right)$

通过设定时间t＝t_j，由上式得到

$E_{j,n}^{\&PlusMinus;}={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;,n}^{\&PlusMinus;}{D}_{j,n}+{\mathrm{\&gamma;}}_0{\&Integral;}_0^{\mathrm{j\&Delta;t}}\mathrm{d\&tau;}{\tilde{D}}_n(\mathrm{j\&Delta;t}-\mathrm{\&tau;}){\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_{r,n}^{\&PlusMinus;}\left(\mathrm{\&tau;}\right)---\left(44\right)$

基于对在时间间隔[(j-1)Δt,jΔt]上的时间基函数的定义，可以表示为：

${\tilde{D}}_n\left(t\right)=\underset{k=0}{\overset{K}{\mathrm{\&Sigma;}}}\left({\mathrm{\&xi;}}_k(j-\frac{t}{\mathrm{\&Delta;t}})D_{j-k,n}\right)---\left(45\right)$

其中 ${\mathrm{\&xi;}}_k\left(\stackrel{\&OverBar;}{t}\right)=\underset{i=0}{\overset{K}{\mathrm{\&Sigma;}}}\left({\mathrm{\&alpha;}}_{k,i}{\stackrel{\&OverBar;}{t}}^i\right),$ 代入上式得，

${\tilde{D}}_n\left(t\right)=\underset{k=0}{\overset{K}{\mathrm{\&Sigma;}}}\left(\left(\underset{i=0}{\overset{K}{\mathrm{\&Sigma;}}}({\mathrm{\&alpha;}}_{k,i}\&CenterDot;{(j-\frac{t}{\mathrm{\&Delta;t}})}^i)\right)D_{j-k,m}\right)---\left(46\right)$

对于K＝4时，

$T\left(t\right)=\left\{\begin{array}{cc}1+\frac{25}{12}t+\frac{35}{24}{t}^2+\frac{5}{12}{t}^3+\frac{1}{24}{t}^4& -1\&le;t\&le;0\\ 1+\frac{5}{6}t-\frac{5}{6}{t}^2-\frac{5}{6}{t}^3-\frac{1}{6}{t}^4& 0\&le;t\&le;1\\ 1-\frac{5}{4}{t}^2+\frac{1}{4}{t}^4& 1\&le;t\&le;2\\ 1-\frac{5}{6}t-\frac{5}{6}{t}^2+\frac{5}{6}{t}^3-\frac{1}{6}{t}^4& 2\&le;t\&le;3\\ 1-\frac{25}{12}t+\frac{35}{24}{t}^2-\frac{5}{12}{t}^3+\frac{1}{24}{t}^4& 3\&le;t\&le;4\\ 0& \mathrm{others}\end{array}\right.---\left(47\right)$

对应的α_5×5如下：

$a_{5\&times;5}=\left[\begin{array}{ccccc}1& -\frac{25}{12}& \frac{35}{24}& -\frac{5}{12}& \frac{1}{24}\\ 0& 4& -\frac{13}{3}& \frac{3}{2}& -\frac{1}{6}\\ 0& -3& \frac{19}{4}& -2& \frac{1}{4}\\ 0& \frac{4}{3}& -\frac{7}{3}& \frac{7}{6}& -\frac{1}{6}\\ 0& -\frac{1}{4}& \frac{11}{24}& -\frac{1}{4}& \frac{1}{24}\end{array}\right]---\left(48\right)$

将式(46)带入到式(44)可得：

$E_{j,n}={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;,n}^{\&PlusMinus;}{D}_{j,n}+{\mathrm{\&gamma;}}_0\mathrm{\&Delta;t}\underset{v=0}{\overset{j-1}{\mathrm{\&Sigma;}}}\underset{k=0}{\overset{K}{\mathrm{\&Sigma;}}}{D}_{j-k-v,n}\underset{i=0}{\overset{K}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_{k,i}{\&Integral;}_0^1{\mathrm{\&tau;}}^i{\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_{r,n}^{\&PlusMinus;}\left((v+\mathrm{\&tau;})\mathrm{\&Delta;t}\right)\mathrm{d\&tau;}---\left(49\right)$

其中，K为拉格朗日时间基函数的阶数，α_k,i为建立时间基函数与色散介电参数之间关系的矩阵，Δt为时间步长。

对于上式中

${\&Integral;}_0^1{\mathrm{\&tau;}}^i{\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_{r,n}^{\&PlusMinus;}\left((v+\mathrm{\&tau;})\mathrm{\&Delta;t}\right)\mathrm{d\&tau;}---\left(50\right)$

由于

${\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_r\left(t\right)=\mathrm{Re}\left({\widehat{\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}}_r\left(t\right)\right)=\mathrm{Re}\{j\frac{{\mathrm{\&omega;}}_p^2}{\sqrt{{\mathrm{\&omega;}}_p^2-v_e^2}}\&CenterDot;e^{[-v_e+j\sqrt{{\mathrm{\&omega;}}_p^2-v_e^2}]t}\}U\left(t\right)---\left(51\right)$

可令， ${\widetilde{\mathrm{\&omega;}}}_c=\sqrt{{\mathrm{\&omega;}}_p^2-v_e^2},{\widetilde{\mathrm{\&delta;}}}_p=v_e-j{\widetilde{\mathrm{\&omega;}}}_c$

则

${\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_r\left(t\right)=\mathrm{Re}\left({\widehat{\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}}_r\left(t\right)\right)=\mathrm{Re}\{j\frac{{\mathrm{\&omega;}}_p^2}{{\widetilde{\mathrm{\&omega;}}}_c}\&CenterDot;e^{[-{\widetilde{\mathrm{\&delta;}}}_p]t}\}U\left(t\right)---\left(52\right)$

(50)式可化为

$\frac{k_p}{\mathrm{\&Delta;t}{\mathrm{\&epsiv;}}_{r,\&infin;}}{e}^{-\mathrm{\&zeta;v}}\frac{i!}{{\mathrm{\&zeta;}}^i}(1-e^{-\mathrm{\&zeta;}}\underset{l=0}{\overset{i}{\mathrm{\&Sigma;}}}\frac{{\mathrm{\&zeta;}}^l}{l!})---\left(53\right)$

其中, $k_p=(j{\mathrm{\&omega;}}_p^2/\left({\widetilde{\mathrm{\&omega;}}}_c{\widetilde{\mathrm{\&delta;}}}_p\right)),\mathrm{\&zeta;}=\mathrm{\&Delta;t}{\widetilde{\mathrm{\&delta;}}}_p.$ 定义：

将式(53)代入式(49)中，得

其中,

$d_k=\underset{i=0}{\overset{K}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_{k,i}\frac{i!}{{\mathrm{\&zeta;}}^i}(1-e^{-\mathrm{\&zeta;}}\underset{l=0}{\overset{i}{\mathrm{\&Sigma;}}}\frac{{\mathrm{\&zeta;}}^l}{l!})---\left(55\right)$

${\hat{S}}_j=k_p\underset{k=1}{\overset{K}{\mathrm{\&Sigma;}}}{D}_{j-k,n}{d}_k---\left(56\right)$

${\hat{C}}_j=k_p\underset{v=1}{\overset{j-1}{\mathrm{\&Sigma;}}}{e}^{-\mathrm{\&zeta;v}}\underset{k=0}{\overset{K}{\mathrm{\&Sigma;}}}{D}_{j-k-v,n}{d}_k---\left(57\right)$

则递推关系为:

${\mathrm{\&epsiv;}}_0{\mathrm{\&epsiv;}}_{r,\&infin;}{\hat{E}}_{j,n}^{\&PlusMinus;}=(1+k_p{d}_0)D_{j,n}+{\hat{S}}_j+e^{-\mathrm{\&zeta;}}\left({\mathrm{\&epsiv;}}_0{\mathrm{\&epsiv;}}_{r,\&infin;}{\hat{E}}_{j-1,n}^{\&PlusMinus;}{D}_{j-1,n}\right)---\left(58\right)$

由(此时的与为复数)，将其与式(58)(ε_r,∞＝1)进行对比，可得与推导过程如下：

$\left\{\begin{array}{c}{\hat{E}}_{j,n}^{\&PlusMinus;}=\frac{D_{j,n}-e^{-\mathrm{\&zeta;}}{D}_{j-1,\mathrm{nn}}+k_p\underset{k=0}{\overset{K}{\mathrm{\&Sigma;}}}{D}_{j-k,n}\&CenterDot;d_k}{{\mathrm{\&epsiv;}}_0}+e^{-\mathrm{\&zeta;}}{\hat{E}}_{j-1,n}^{\&PlusMinus;}\\ {\hat{E}}_{j,n}^{\&PlusMinus;}=\underset{k=0}{\overset{K}{\mathrm{\&Sigma;}}}{\mathrm{\&beta;}}_{k,m}^{d\&PlusMinus;}{D}_{j-k,n}+{\mathrm{\&beta;}}_n^{e\&PlusMinus;}{\hat{E}}_{j-1,n}^{\&PlusMinus;}\end{array}\right.---\left(59\right)$

展开后

那么，等离子体介质时域积分方程中电场强度系数E^±与电通量密度D关系如下

此时，即建立了电场强度与电通量密度之间的关系。

步骤4，求解时域矩阵方程，分析宽频带内混合目标电磁散射特性。

利用步骤3得到的电场强度与电通量密度之间的关系，带入到步骤2的(12)-(14)式中即可求解出电场强度与电通量密度的系数，进而求得电场强度与电通量密度，得到等离子内的体电流和金属表面的面电流，最终可分析宽频带内等离子体涂覆金属目标的电磁散射特性。

为了验证本发明方法的准确性与有效性，下面给出了等离子体介质涂敷金属目标的瞬态电磁特性的分析，如图1，涂敷立方体模型，金属边长0.1m，涂覆介质厚度为0.05m。涂覆的等离子体介质参数：特征频率为1e8Hz，碰撞频率为3e8Hz。该模型的双站RCS的计算结果与频域体面积分方法进行了比较，吻合得较好，如图2所示。

本算例中，入射电场采用调制高斯平面波，其表达式如下：

$E^{\mathrm{inc}}(r,t)={\hat{p}}^{\mathrm{inc}}\exp [-{\left(\frac{\mathrm{\&tau;}-t_c}{\sqrt{2}\mathrm{\&sigma;}}\right)}^2]\cos \left(2\mathrm{\&pi;}{f}_c\mathrm{\&tau;}\right)---\left(69\right)$

其中，极化方向入射方向Theta＝45deg,Phi＝0deg，σ＝6/(2πf_bw)，t_c＝3.5σ,E^inc(r,t)的频谱的中心频率为f₀＝200MHz，最高频率为300MHz，f_bw为频带宽度。时间步长Δt＝0.1lm，总时间步N_t＝300，lm是光米(lightmeter)，即光在自由空间中传播1m距离所花的时间。

综上所述，本发明方法运用时域体面积分方程分析金属等离子体介质混合目标的电磁散射特性，只需离散物体模型，使用自由空间格林函数及递归卷积方法，可以分析含非均匀等离子体介质与金属混合目标，相较于需要全空间离散的时域有限差分方法和时域有限元方法，该方法具有很高的计算精度且未知量相对较少；同时由于采用时域方法求解，可以一次计算宽频带各频点的电磁散射特性，相较于需要扫频操作的频域方法，求解消耗时间减少。

Claims (4)

1.一种分析高超声速飞行目标电磁散射特性的时域积分方程方法，其特征在于，步骤如下：

步骤1，建立时域体面积分方程：根据结构散射特性，被激励目标上的总电场等于入射电场与散射电场之和，入射电场为已知激励；

步骤2，离散积分方程，并进行空间和时间上的测试，建立时域矩阵方程；

步骤3，建立等离子体中电通量密度与电场强度之间的关系；

步骤4，求解时域矩阵方程，分析宽频带内混合目标电磁散射特性。

2.根据权利要求1所述的分析高超声速飞行目标电磁散射特性的时域积分方程方法，其特征在于，步骤1中所述建立时域体面积分方程，具体如下：

设待求解目标为金属等离子体介质混合目标，用调制高斯脉冲激励待求解目标，在混合目标的介质体内产生极化体电流J_V(r,t)，金属表面产生感应电流J_s(r,t)，由于总电场等于入射电场与散射电场之和，得到时域时间步进电场积分方程形式如下：

介质体上：E(r,t)＝E^inc(r,t)+E^sca(r,t)(1)

金属表面：[E^inc(r,t)+E^sca(r,t)]_tan＝0(2)

其中，E^inc(r,t)表示入射电场，E^sca(r,t)表示散射电场，E(r,t)表示总电场；

散射场具体表达式为：

$E_v^{\mathrm{sca}}(r,t)=-{\mathrm{\&mu;}}_0{\&Integral;}_V\stackrel{\&OverBar;}{G}(r,r^{\&prime;})\&CenterDot;{\&PartialD;}_t{J}_V(r^{\&prime;},t),(r\&Element;V)---\left(3\right)$

$E_s^{\mathrm{sca}}(r,t)=-{\mathrm{\&mu;}}_0{\&Integral;}_s\stackrel{\&OverBar;}{G}(r,r^{\&prime;})\&CenterDot;{\&PartialD;}_t{J}_s(r^{\&prime;},t),(r\&Element;S)---\left(4\right)$

其中，下标V表示体单元，S表示面单元；D(r,t)表示电通量密度、为并矢格林函数，ε₀表示自由空间的介电常数，μ₀表示自由空间的磁导率；

得电场积分方程EFIE:

${\&PartialD;}_t{E}^{\mathrm{inc}}(r,t)={\&PartialD;}_{t}E(r,t)+L\left\{{\&PartialD;}_{t}D(r,t)\right\}-{\mathrm{\&epsiv;}}_{0}L\left\{{\&PartialD;}_{t}E(r,t)\right\}+L\left\{{\&PartialD;}_t{J}_S(r,t)\right\},(r\&Element;V)---\left(5\right)$

$\hat{n}\&times;\hat{n}{\&PartialD;}_{t}E(r,t)=0,(r\&Element;S)---\left(6\right)$

其中，L{X}算子表示为：

$L\left\{X\right\}(r,t)\stackrel{.}{=}\frac{{\mathrm{\&mu;}}_0}{4\mathrm{\&pi;}}{\&Integral;}_{\mathrm{\&xi;}}\mathrm{dv}({\&PartialD;}_t^2\stackrel{\&OverBar;}{I}-c^2\&dtri;\&dtri;)\frac{X(r^{\&prime;},t-\frac{R}{c})}{R}---\left(7\right)$

这里,为单位并矢，为单位法向量，r为场点，r′为源点，R＝r-r′为场源点之间位置矢量_，D(r,t)表示电通量密度，R＝|R|＝|r-r′|，表示光速，t表示时间，ξ为积分区域。

3.根据权利要求1所述的分析高超声速飞行目标电磁散射特性的时域积分方程方法，其特征在于，步骤2中所述离散积分方程，并进行空间和时间上的测试，建立时域矩阵方程，具体如下：

现对介质内部的电通量密度D(r,t)和电场强度E(r,t)在空间上用SWG基函数进行展开，金属表面的电流Js(r,t)用RWG基函数展开:

$D(r,t)=\underset{n=1}{\overset{N_s}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{N_t}{\mathrm{\&Sigma;}}}{D}_{j,n}{f}_n\left(r\right)T_j\left(t\right)=\underset{n=1}{\overset{N_s}{\mathrm{\&Sigma;}}}{\tilde{D}}_n\left(t\right)f_n\left(r\right)---\left(8\right)$

$\begin{array}{c}E(r,t)=\underset{n=1}{\overset{N_s}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{N_t}{\mathrm{\&Sigma;}}}(E_{j,n}^{+}{f}_n^{+}\left(r\right)+E_{j,n}^{-}{f}_n^{-}\left(r\right))T_j\left(t\right)\\ =\underset{n=1}{\overset{N_s}{\mathrm{\&Sigma;}}}({\tilde{E}}_n^{+}\left(t\right)f_n^{+}\left(r\right)+{\tilde{E}}_n^{-}\left(t\right)f_n^{-}\left(r\right))\end{array}---\left(9\right)$

$J_s(r,t)=\underset{n=1}{\overset{N_s}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{N_t}{\mathrm{\&Sigma;}}}{I}_{j,n}{f}_n^s\left(r\right)T_j\left(t\right)---\left(10\right)$

其中，n表示空间未知量，j表示时间步数，N_s表示空间上总的未知量，N_t表示时间上总的步数，D_j,n表示第n个电通量密度在第j时刻的系数，I_j,n表示第n个面电流在第j时刻的系数，f_n(r)表示SWG基函数，表示RWG基函数，表示上半SWG基函数，表示下半SWG基函数，T_j(t)表示时间基函数，表示在固定时间步上的第n个电通量密度未知量系数，表示第n个上半SWG基函数表示的电场未知量在第j时刻的系数，表示第n个下半SWG基函数表示的电场未知量在第j时刻的系数，表示在固定时间步上的第n个上半SWG基函数表示的电场未知量，表示在固定时间步上的第n个下半SWG基函数表示的电场未知量；

时间上用四阶拉格朗日时间基函数进行展开:

$T\left(t\right)=\left\{\begin{array}{cc}1+\frac{25}{12}t+\frac{35}{24}{t}^2+\frac{5}{12}{t}^3+\frac{1}{24}{t}^4& -1\&le;t\&le;0\\ 1+\frac{5}{6}t-\frac{5}{6}{t}^2-\frac{5}{6}{t}^3-\frac{1}{6}{t}^4& 0\&le;t\&le;1\\ 1-\frac{5}{4}{t}^2+\frac{1}{4}{t}^4& 1\&le;t\&le;2\\ 1-\frac{5}{6}t-\frac{5}{6}{t}^2+\frac{5}{6}{t}^3-\frac{1}{6}{t}^4& 2\&le;t\&le;3\\ 1-\frac{25}{12}t+\frac{35}{24}{t}^2-\frac{5}{12}{t}^3+\frac{1}{24}{t}^4& 3\&le;t\&le;4\\ 0& \mathrm{others}\end{array}\right.---\left(11\right)$

对(5)(6)式电场积分方程离散后在空间域上进行伽辽金测试，在时间域上进行点匹配得到时间步进矩阵方程：

${\tilde{Z}}_0{I}_j^s=F_j^s-\underset{i=1}{\overset{j-1}{\mathrm{\&Sigma;}}}{\tilde{Z}}_i{I}_{j-i}^s,(r^{\&prime;}\&Element;S)---\left(14\right)$

其中， ${\tilde{Z}}_i=\left[\begin{array}{cc}{\tilde{Z}}_{i,\mathrm{vv}}& {\tilde{Z}}_{i,\mathrm{vs}}\\ {\tilde{Z}}_{i,\mathrm{sv}}& {\tilde{Z}}_{i,\mathrm{ss}}\end{array}\right],$ $Z_i^{e\&PlusMinus;}=\left[\begin{array}{c}{Z}_{i,\mathrm{vv}}^{e\&PlusMinus;}\\ Z_{i,\mathrm{sv}}^{e\&PlusMinus;}\end{array}\right],$

其中，V表示介质积分区域，S表示金属积分区域，v代表介质，s代表金属，Z、为建立场源之间联系的时域阻抗矩阵，F_j表示第j个时间步的激励，Nt表示时间步的个数，是第j个时间步的电通量密度系数、是第j个时间步的上下半个SWG基函数表示的电场强度的系数，是第j个时间步的面电流的系数。

4.根据权利要求1所述的分析高超声速飞行目标电磁散射特性的时域积分方程方法，其特征在于，步骤3中所述建立等离子体中电通量密度与电场强度之间的关系，具体如下：

对于等离子体介质，满足D(r,t)＝ε(r,t)*E(r,t)，则：

γ(r,t)*D(r,t)＝Eⁱ(r,t)-L{D}(r,t)+ε₀L{γ*D}(r,t)

(15)

其中，

$\begin{array}{c}\mathrm{\&gamma;}(r,t)*D(r,t)=E(r,t)\\ ={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;}\left(r\right)\mathrm{\&delta;}\left(t\right)*D(r,t)+{\mathrm{\&gamma;}}_0{\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_r\\ ={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;}\left(r\right)D(r,t)+{\mathrm{\&gamma;}}_0{\&Integral;}_0^t\mathrm{d\&tau;D}(r,t-\mathrm{\&tau;}){\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_r(r,\mathrm{\&tau;})\end{array}(r,t)*D(r,t)---\left(16\right)$

$\mathrm{\&gamma;}(r,t)=F^{-1}\{\mathrm{\&gamma;}(r,\mathrm{\&omega;})=1/\mathrm{\&epsiv;}(r,\mathrm{\&omega;})\}(r,t)={\mathrm{\&gamma;}}_{r,\&infin;}\left(r\right)\mathrm{\&delta;}\left(t\right)+{\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_r(r,\mathrm{\&tau;})---\left(17\right)$

γ_r,∞＝1/ε_r,∞(18)

γ₀＝1/ε₀(19)

其中，F^-1表示逆傅里叶变换，ε(r,ω)为在r处ω频率下的介电常数，ε_r,∞为在r处无限大频率下的介电常数，γ(r,t)为在r处t时刻的下的磁化率，δ(t)为冲激函数；

利用空间和时间基函数展开后，由(16)式得到：

$\begin{array}{c}{\tilde{E}}_n^{\&PlusMinus;}\left(t\right)={\mathrm{\&gamma;}}_n^{\&PlusMinus;}\left(t\right)*{\tilde{D}}_n\left(t\right)\\ ={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;,n}^{\&PlusMinus;}{\tilde{D}}_n\left(t\right)+{\mathrm{\&gamma;}}_0{\&Integral;}_0^t\mathrm{d\&tau;}{\tilde{D}}_n(t-\mathrm{\&tau;}){\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_{r,n}^{\&PlusMinus;}\left(\mathrm{\&tau;}\right)\end{array}---\left(20\right)$

设定时间t＝t_j，由上式得到

$E_{j,n}^{\&PlusMinus;}={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;,n}^{\&PlusMinus;}{D}_{j,n}+{\mathrm{\&gamma;}}_0{\&Integral;}_0^{\mathrm{j\&Delta;t}}\mathrm{d\&tau;}{\tilde{D}}_n(\mathrm{j\&Delta;t}-\mathrm{\&tau;}){\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_{r,n}^{\&PlusMinus;}\left(\mathrm{\&tau;}\right)---\left(21\right)$

进一步表示为：

$E_{j,n}^{\&PlusMinus;}={\mathrm{\&gamma;}}_0{\mathrm{\&gamma;}}_{r,\&infin;,n}^{\&PlusMinus;}{D}_{j,n}+{\mathrm{\&gamma;}}_0\mathrm{\&Delta;t}\underset{v=0}{\overset{j-1}{\mathrm{\&Sigma;}}}\underset{k=0}{\overset{K}{\mathrm{\&Sigma;}}}{D}_{j-k-v,n}\underset{i=0}{\overset{K}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_{k,i}{\&Integral;}_0^1{\mathrm{\&tau;}}^i{\stackrel{\&OverBar;}{\mathrm{\&gamma;}}}_{r,n}^{\&PlusMinus;}\left((v+\mathrm{\&tau;})\mathrm{\&Delta;t}\right)\mathrm{d\&tau;}---\left(22\right)$

其中，K为拉格朗日时间基函数的阶数，α_k,i为建立时间基函数与色散介电参数之间关系的矩阵，Δt为时间步长。