一种计算三维高速平动目标电磁散射场的方法与流程

本发明涉及电磁场计算领域，尤其是一种计算电磁散射场的方法。
背景技术：
：随着航空器件的飞行速度越来越快，关于高速运动物体的电磁研究也逐渐成为了热点问题。然而，随着物体的高速运动，根据相对论原理，物体在运动状态下的电场散射场也发生了改变。针对运动目标电磁散射计算的数值方法开始出现，截止目前，部分学者提出了关于运动目标电磁散射的数值方法，比如等效面电流法，运动边界法等。以往的方法因坐标系限制、边界条件或者求解结果单一或不能适应于三维目标复杂多样的形状等缺点，更不能适用于高速乃至超高速运动物体的电磁场值研究。因此，非常有必要设计一种能够解决复杂目标高速平动电磁场值的技术。技术实现要素：为了克服现有技术的不足，本发明提供一种计算三维高速平动目标的电磁散射场的技术方案。本发明以杏仁体为目标，设计建立相对于目标静止的运动坐标系。使用fortran语言处理，生成杏仁体网格模型，实现时间步长和空间步长的自由选取，然后在相对运动坐标系中引入入射波，用fdtd计算相对运动坐标系中的电磁散射场，最后将其转换为静止坐标系下的电磁散射场。本发明解决其技术问题所采用的技术方案包括以下步骤：步骤一：建立相对于高速平动目标静止的运动坐标系k'，并确定其特征；步骤1.1：建立相对运动坐标系为k'，相对运动坐标系k'相对于高速平动的目标是静止的，即相对运动坐标系k'与目标有相同的运动速度和运动方向；步骤1.2：建立相对运动坐标系k'的空间步长δx′,δy′,δz′和静止坐标系k的空间步长δx,δy,δz之间的关系，以及相对运动坐标系k'的时间步长δt′和静止坐标系的δt的关系；目标与相对运动坐标系k'有相同的速度大小v以及运动角度θ′v,其中θ′v代表运动方向与z轴正半轴的夹角，范围为0°-180°，代表运动方向投影到xoy面后与x轴正半轴的夹角，范围0°-360°，角度根据公式(1)可得到空间步长δx′,δy′,δz′和δx,δy,δz的关系，以及时间步长δt′和δt的关系：其中，c为光速，x,y,z和t是静止坐标系k中的空间分量和时间分量，x′,y′,z′和t′是相对运动坐标系k'中的空间分量和时间分量，v为目标运动速度，速度在x、y、z方向上的分量分别为vz＝vcosθv，式中步骤1.3：根据courant稳定性条件且δx′＝δy′＝δz＝δ′，确定相对运动坐标系k'中空间步长δx′,δy′,δz′以及时间步长δt′的大小，至此，相对于高速平动目标静止的运动坐标系k'建立完成；步骤二：建立杏仁体网格模型；步骤2.1：利用计算机建模，本发明的三维高速平动目标选用杏仁体，建立数学模型的控制方程如式(2)所示：其中，ψ为控制杏仁体外形在yoz平面内轨迹的角度量，取值范围-π≤ψ≤π；d为杏仁体沿x方向上的总长度；a1和a2为杏仁体在x方向上正负半轴分别占总长度的比例系数，且满足a1+a2＝1；b1和b2分别为控制杏仁体在x负半轴和x正半轴在y方向和z方向长度的系数；n为控制杏仁体在z方向上长度的缩减系数，且n为大于零的整数；c1和c2为控制杏仁体x正半轴沿x方向长度的缩减系数，满足且0＜c1＜1，从而在x正半轴形成尖端，b1和b2满足|b1-b2(1-c1)|＜10-5；根据式(2)所建立数学模型，杏仁体在x方向的尺寸范围为x∈(a1×d,a2×d)，y方向的尺寸范围为y∈(-b1×d,b1×d)，z方向的尺寸范围为步骤2.2：根据步骤2.1中杏仁体在x，y，z方向的尺寸范围和相对运动坐标系k'空间步长的大小，对网格数采用整数位取整，确定三维空间的散射边界，吸收边界，连接边界和总场边界，其中散射边界为i∈(ismin,ismax)；j∈(jsmin,jsmax)；k∈(ksmin,ksmax)；分别代表x方向上散射边界的范围为(ismin,ismax)，即x方向的最小值和最大值；y方向上散射边界的范围为(jsmin,jsmax),即y方向的最小值和最大值；z方向上散射边界的范围为(ksmin,ksmax)，即z方向的最小值和最大值；步骤2.3：将绝对尺寸换算为网格数，若空间某网格(i,j,k)的位置满足步骤2.1中的数学模型，即坐标为(i,j,k)的网格，该网格的实际坐标为(x,y,z)＝(i×δx,j×δy,k×δz)，其中x，y，z分别代表该网格的实际位置坐标值，δx,δy,δz分别代表空间步长，i,j,k分别为该网格在x,y,z方向上的网格数，将坐标值(x,y,z)代入步骤2.1中的数学模型，若x,y,z的值满足公式(2)的关系式，则此网格为杏仁体内部，即设置为金属；若x,y,z的值不满足公式(2)的关系式，则此网格为杏仁体外部，设为空气，依次类推，即可得到空间内的杏仁体模型，并已将杏仁体模型按网格剖分为长方体元胞模型，得到txt格式的平面模型数据；步骤三：在相对运动坐标系k'中引入入射平面波；步骤3.1：确定静止坐标系k下的入射波，定义幅值，频率以及入射角度θi，θi代表入射方向与z轴正半轴的夹角，范围为0°-180°，代表入射方向投影到xoy面后与x轴正半轴的夹角，范围0°-360°，入射波满足其中，是入射波的电场强度，e0为幅度初值，t为时间，ωi是角频率，等于2πf，f是入射波的频率,ki是波数，等于ωi/c；步骤3.2：应用洛伦兹变换得到相对运动坐标系k'和静止坐标系k中入射波的幅值e0和e′0，频率f和f′，以及入射角度θi,与θ′i,的转换关系，在总场边界处将入射波引入到相对于杏仁体目标静止的相对运动坐标系k'中；步骤四：运用fdtd算法计算得到运动坐标系k'中电场与磁场数据；步骤4.1：输入运行时间步数timestop；步骤4.2：输入杏仁体的运动速度v及运动角度θ′v,步骤4.3：输入入射波的入射角度θi,并将其转化为θ′i,步骤4.4：运用三维fdtd算法的差分计算，电场和磁场相差半个时间步，从而在时间上迭代求解，直至时间步数从1到timestop完成，求得各个时刻空间电场和磁场迭代公式如下：其中和分别代表网格节点(i′,j′,k′)在n+1时刻和n时刻在相对运动坐标系k'中的电场值的e′x分量，和分别代表在时刻网格节点(i′,j′,k′)和(i′,j′-1,k′)的相对运动坐标系k'中的磁场值的h′z分量，和分别代表在时刻网格节点(i′,j′,k′)和(i′,j′,k′-1)的相对运动坐标系k'中的磁场值的h′y分量，ε(m)为空间位置m处的介电系数，σ(m)为空间位置m处的电导率；步骤五：将步骤4.4得到的电场e和磁场h转换到静止坐标系k，并插值，即得到杏仁体高速平动的电磁散射场，详细步骤如下：步骤5.1：结合麦克斯韦方程和洛伦兹变换将相对运动坐标系k'下的电磁场值变换到静止坐标系k中；按照场值变换为相对运动坐标系k′中电场的平行分量，为相对运动坐标系k′中电场的垂直分量，为磁场的平行分量，为磁场的垂直分量；步骤5.2：由于速度的影响，相对运动坐标系k′和静止坐标系k的网格大小不一致，由步骤5.1得到静止坐标系k中的电磁场值是由相对运动坐标系k'中的场值变换而来，这些场值对应的网格大小为相对运动坐标系k'中的δ′，因此要重新划分静止坐标系k中的网格大小，即在静止坐标系k中，按δx＝δy＝δz＝δ′插值得到静止坐标k系下各空间网格的电磁场值，重新划分网格之后，进行三维插值，使用fortran调用matlab引擎，matlab含有三维插值的函数meshgrid和interp3，插值方式为linear，即线性插值，即可得到电磁场值；步骤5.3：对步骤5.2得到的电磁场值，采用fft算法，得到其频率分布，即可同时提取其幅度和相位。本发明针对三维高速平动目标的计算包含但不仅限于正方体、长方体、球体、圆柱体、圆锥体、棱柱体和圆台，以及正方体、长方体、球体、圆柱体、圆锥体、棱柱体和圆台的组合；还包括锥球和多面体。本发明的有益效果在于采用fortran语言进行网格建模，引入相对于杏仁体目标静止的运动坐标系，引入入射波，先计算得到相对于杏仁体目标静止的运动坐标系中的场值，然后再进行场值变换和插值，最后得到静止坐标系中的电磁散射场，包括了幅值分布，频率分布和相位分布；本发明方法解决了两个坐标系间的转换关系，利用这些转换关系，得到了杏仁体目标在高速平动状态下的电磁散射场值；本发明方法新颖，不涉及格林函数、矩阵、渐进函数和基函数等其他数值方法必须采用的数学公式，简单易操作，适于处理更多三维空间复杂目标高速平动电磁散射场的优点。附图说明图1是本发明的相对平动目标静止的运动坐标系k'示意图。图2是本发明的杏仁体建模截面展示图。图3是本发明的边界示意图。图4是本发明的计算过程图。图5是本发明的插值示意图。图6是本发明的xoy截面的电场分布图。图7是本发明的xoy截面的频率分布图。图8是本发明的xoy截面的相位分布图。具体实施方式下面结合附图和实施例对本发明进一步说明。一种计算三维高速平动目标电磁散射场的方法，包括以下步骤：步骤一：建立相对于高速平动目标静止的运动坐标系k'，并确定其特征；步骤1.1：建立相对运动坐标系为k'，相对运动坐标系k'相对于高速平动的目标是静止的，即相对运动坐标系k'与目标有相同的运动速度和运动方向；步骤1.2：建立相对运动坐标系k'的空间步长δx′,δy′,δz′和静止坐标系k的空间步长δx,δy,δz之间的关系，以及相对运动坐标系k'的时间步长δt′和静止坐标系的δt的关系；目标与相对运动坐标系k'有相同的速度大小v以及运动角度θ′v,其中θ′v代表运动方向与z轴正半轴的夹角，范围为0°-180°，代表运动方向投影到xoy面后与x轴正半轴的夹角，范围0°-360°，角度根据公式(1)可得到空间步长δx′,δy′,δz′和δx,δy,δz的关系，以及时间步长δt′和δt的关系：其中，c为光速，x,y,z和t是静止坐标系k中的空间分量和时间分量，x′,y′,z′和t′是相对运动坐标系k'中的空间分量和时间分量，v为目标运动速度，速度在x、y、z方向上的分量分别为vz＝vcosθv，式中步骤1.3：根据courant稳定性条件且δx′＝δy′＝δz＝δ′，确定相对运动坐标系k'中空间步长δx′,δy′,δz′以及时间步长δt′的大小，至此，相对于高速平动目标静止的运动坐标系k'建立完成；如图1所示，表示与高速平动目标相对静止的运动坐标系k'，目标与该坐标系有相同的速度大小v，以及运动角度θ′v,此处设v＝0.1c，其中角度根据公式(1)可得到空间步长δx′,δy′,δz′和δx,δy,δz的关系，时间步长δt′和δt的关系，根据courant稳定性条件，取空间步长δx′,δy′,δz′＝λ′/120，以及时间步长δt′＝δx′/2c。步骤二：建立杏仁体网格模型；步骤2.1：使用fortran语言建模，本发明的三维高速平动目标选用杏仁体，建立数学模型的控制方程如式(2)所示：其中，ψ为控制杏仁体外形在yoz平面内轨迹的角度量，取值范围-π≤ψ≤π；d为杏仁体沿x方向上的总长度；a1和a2为杏仁体在x方向上正负半轴分别占总长度的比例系数，且满足a1+a2＝1；b1和b2分别为控制杏仁体在x负半轴和x正半轴在y方向和z方向长度的系数；n为控制杏仁体在z方向上长度的缩减系数，且n为大于零的整数；c1和c2为控制杏仁体x正半轴沿x方向长度的缩减系数，满足且0＜c1＜1，从而在x正半轴形成尖端；为了保证模型左右两部分在x＝0平面吻合，b1和b2需要满足|b1-b2(1-c1)|＜10-5；本发明中，取a1＝-0.416667，a2＝0.58333，b1＝0.193333，b2＝4.83345，c1＝0.96，c2＝2.08335，d＝0.252374，n＝3；根据式(2)所建立数学模型，杏仁体在x方向的尺寸范围为x∈(a1×d,a2×d)＝(-0.105156,0.147217)，y方向的尺寸范围为y∈(-b1×d,b1×d)＝(-0.04879,0.04879)，z方向的尺寸范围为步骤2.2：根据步骤2.1中杏仁体在x，y，z方向的尺寸范围和相对运动坐标系k'空间步长的大小，对网格数采用整数位取整，确定三维空间的散射边界，吸收边界，连接边界和总场边界，其中散射边界为i∈(ismin,ismax)；j∈(jsmin,jsmax)；k∈(ksmin,ksmax)；分别代表x方向上散射边界的范围为(ismin,ismax)，即x方向的最小值和最大值；y方向上散射边界的范围为(jsmin,jsmax),即y方向的最小值和最大值；z方向上散射边界的范围为(ksmin,ksmax)，即z方向的最小值和最大值；如图3所示，规定了吸收边界，散射边界以及总场边界，各值如下表：imin,imaxjmin,jmaxkmin,kmax吸收边界-112,129-90,90-76,76散射边界-42,59-20,20-6,6总场边界-62,79-40,40-26,26步骤2.3：将绝对尺寸换算为网格数，若空间某网格(i,j,k)的位置满足步骤2.1中的数学模型，即坐标为(i,j,k)的网格，该网格的实际坐标为(x,y,z)＝(i×δx,j×δy,k×δz)，其中x，y，z分别代表该网格的实际位置坐标值，δx,δy,δz分别代表空间步长，i,j,k分别为该网格在x,y,z方向上的网格数，将坐标值(x,y,z)代入步骤2.1中的数学模型，若x,y,z的值满足公式(2)的关系式，则此网格为杏仁体内部，即设置为金属；若x,y,z的值不满足公式(2)的关系式，则此网格为杏仁体外部，设为空气，依次类推，即可得到空间内的杏仁体模型，并已将杏仁体模型按网格剖分为长方体元胞模型，得到txt格式的平面模型数据；步骤三：在相对运动坐标系k'中引入入射平面波；步骤3.1：确定静止坐标系k下的入射波，定义幅值，频率以及入射角度θi,θi代表入射方向与z轴正半轴的夹角，范围为0°-180°，代表入射方向投影到xoy面后与x轴正半轴的夹角，范围0°-360°，入射波满足其中，是入射波的电场强度，e0为幅度初值，t为时间，ωi是角频率，等于2πf，f是入射波的频率,ki是波数，等于ωi/c；本发明在总场边界处引入入射波，坐标系k中规定入射波幅值e0＝1，频率f＝1ghz，角度k′中入射波的参数参考公式(1)及相位不变性得到幅值e′0＝1.106，频率f′＝1.105ghz，以及入射角度步骤3.2：应用洛伦兹变换得到相对运动坐标系k'和静止坐标系k中入射波的幅值e0和e′0，频率f和f′，以及入射角度θi,与θ′i,的转换关系，在总场边界处将入射波引入到相对于杏仁体目标静止的相对运动坐标系k'中；步骤四：运用fdtd算法计算得到运动坐标系k'中电场与磁场数据；步骤4.1：输入运行时间步数timestop；步骤4.2：输入杏仁体的运动速度v及运动角度θ′v,步骤4.3：输入入射波的入射角度θi,并将其转化为θ′i,步骤4.4：运用三维fdtd算法的差分计算，电场和磁场相差半个时间步，从而在时间上迭代求解，直至时间步数从1到timestop完成，求得各个时刻空间电场和磁场迭代公式如下：其中和分别代表网格节点(i′,j′,k′)在n+1时刻和n时刻在相对运动坐标系k'中的电场值的e′x分量，和分别代表在时刻网格节点(i′,j′,k′)和(i′,j′-1,k′)的相对运动坐标系k'中的磁场值的h′z分量，和分别代表在时刻网格节点(i′,j′,k′)和(i′,j′,k′-1)的相对运动坐标系k'中的磁场值的h′y分量，ε(m)为空间位置m处的介电系数，σ(m)为空间位置m处的电导率；如图4所示，表示计算过程的流程图，先进行建模及初始化运动坐标系k′中的δx′,δy′,δz′,δt′，然后输入时间步timestop＝8692，以及各角度参数，随着时间步利用公式(3)进行迭代计算；步骤五：将步骤4.4得到的电场e和磁场h转换到静止坐标系k，并插值，即得到杏仁体高速平动的电磁散射场，详细步骤如下：步骤5.1：结合麦克斯韦方程和洛伦兹变换将相对运动坐标系k'下的电磁场值变换到静止坐标系k中；按照场值变换为相对运动坐标系k′中电场的平行分量，为相对运动坐标系k′中电场的垂直分量，为磁场的平行分量，为磁场的垂直分量；步骤5.2：由于速度的影响，相对运动坐标系k′和静止坐标系k的网格大小不一致，由步骤5.1得到静止坐标系k中的电磁场值是由相对运动坐标系k'中的场值变换而来，这些场值对应的网格大小为相对运动坐标系k'中的δ′，因此要重新划分静止坐标系k中的网格大小，即在静止坐标系k中，按δx＝δy＝δz＝δ′插值得到静止坐标k系下各空间网格的电磁场值，重新划分网格之后，进行三维插值，使用fortran调用matlab引擎，matlab含有三维插值的函数meshgrid和interp3，插值方式为linear，即线性插值，即可得到电磁场值；步骤5.3：对步骤5.2得到的电磁场值，采用fft算法，得到其频率分布，即可同时提取其幅度和相位。本发明针对三维高速平动目标的计算包含但不仅限于正方体、长方体、球体、圆柱体、圆锥体、棱柱体和圆台，以及正方体、长方体、球体、圆柱体、圆锥体、棱柱体和圆台的组合；还包括锥球和多面体。如图5所示，插值过程中，坐标系k将按照δx＝δy＝δz＝δ′的格式重新插值，其中δ′为坐标系k'中的离散间隔；如图6所示，得到txt数据文件后，打开origin画图软件，关掉book1，选择newmatric，导入数据文件，即可得到xoy截面的电场场值分布图，不仅如此，也可以得到空间任一个截面的场值分布，即杏仁体目标的电磁散射场；如图7所示，对上述场值进行fft算法，杏仁体高速平动目标xoy面的频率分布也被得到，前向散射大于后向散射，同时空间任一点的频率值都可以被得到；如图8所示，对fft后的数据提取相位，可得到杏仁体高速平动目标xoy面的相位分布图。当前第1页12