CN111339688A - 基于大数据并行算法求解火箭仿真模型时域方程的方法 - Google Patents

Abstract

本发明公开了一种基于大数据并行算法求解火箭仿真模型时域方程的方法，包括步骤：对火箭模型建模，获取WLP‑FDTD时域方程组；对火箭电磁仿真模型中WLP‑FDTD时域方程组进行多map划分；对所述map块内若干待求电场值进行“键‑值”对设置；对所述“键‑值”进行分组迭代计算，直至计算到目标拉盖尔阶数时的待求电场值。本发明通过将大数据并行技术引入到火箭模型的WLP‑FDTD时域方程的计算求解中，解决了三维柱坐标下超大尺寸细微结构的电磁问题，提高计算效率，降低硬件成本。

Description

基于大数据并行算法求解火箭仿真模型时域方程的方法

技术领域

本发明涉及数据处理技术领域，尤其涉及一种基于大数据并行算法求解火箭仿真模型时域方程的方法。

背景技术

近年来，时域有限差分(FDTD)法被广泛应用于电磁波问题的建模当中。尽管如此，其时间步长受CFL(Couran-Friedrich-Lecy)稳定性条件的限制成为困扰学者们的一个难题。同时，柱坐标系中FDTD算法的时间步长同样受CFL稳定性条件限制，为了消除此种限制，学者们做了大量的研究，大致可以分为两个方向。一是交替方向隐式FDTD(ADI-FDTD)方案。二是无条件稳定加权拉盖尔多项式FDTD(WLP-FDTD)方案。这两种方案都可以有效的消除传统柱坐标系FDTD算法的CFL稳定性条件限制，使仿真更加有效。尽管如此，当时间步长较大时，ADI方案会导致较大的数值色散误差，而无条件WLP方案在计算精度上能够较好的降低数值色散误差，但会导致大量的内存消耗，降低计算效率。尤其是一些超大尺寸细微结构模型会直接导致无法计算。

WLP-FDTD方法与传统FDTD方法不同在于，WLP按拉盖尔基函数的阶数由低到高计算，待求量为基函数的展开系数，其隐式差分需要展开一个大型的高度稀疏矩阵，其求解过程虽然精度较高，但是非常耗费内存。为避免大型稀疏矩阵的求解，在二维空间或小尺寸三维空间内，通常采用LU分解来实现对WLP的高效求解。但在一些超大尺寸细微结构中，由于矩阵的阶数较高，通常会导致无法计算。

发明内容

技术目的：针对现有技术中当前柱坐标系WLP-FDTD算法对三维问题进行计算时矩阵阶次较高导致计算时间很长甚至无法计算的缺陷，本发明公开了一种基于大数据并行技术求解火箭模型WLP-FDTD的方法，通过将大数据并行技术引入到火箭模型的WLP-FDTD时域方程的计算求解中，解决了三维柱坐标下超大尺寸细微结构的电磁问题，提高计算效率，降低硬件成本。

技术方案：针对上述技术问题，本发明提供以下技术方案。

一种基于大数据并行算法求解火箭仿真模型时域方程的方法，包括以下步骤：

步骤S1、对火箭模型建模，获取WLP-FDTD时域方程组：在雷电云层中，雷电以亚纳秒脉冲波打入火箭载体，以柱形体系结构建立火箭模型，基于ρ，

z三个坐标轴建立火箭电磁仿真模型；

步骤S2、对火箭电磁仿真模型中WLP-FDTD时域方程组进行多map划分：将所述模型导入大数据的HDFS平台，利用所述HDFS平台对时域方程组限定边界后进行网格拆分，划分为若干map块；

步骤S3、对所述map块内若干待求电场值进行“键-值”对设置：对待求电场量矩阵进行展开，将待计算的电场的行号作为“键”数据，将柱坐标系WLP-FDTD内各网格点的电场仿真数据作为“值”数据；

步骤S4、对所述“键-值”进行分组迭代计算，直至计算到目标拉盖尔阶数时的待求电场值：将若干map块并行计算所得电场值和磁场值，进行是否为边界点判断，并且将所得分组结果在的HDFS平台上构建两组数据A、B；其中将属于本map块内的电场值和磁场值归于B组，并将其直接代入柱坐标系WLP-FDTD时域方程组中，进行逐步迭代；将属于相邻块内边界点的电场值和磁场值存储于A组，在需计算B组内电场值时进行调用。

优选地，所述电磁仿真模型中时域方程组为：

其中，

为i坐标轴方向的电场矢量，

为i坐标轴方向的磁场矢量，

为沿着i坐标轴一阶中心差分算子，这里

为三个坐标轴，ε为介质介电系数，μ为磁导系数，s＞0，为时间尺度因子，q为拉盖尔阶数。

优选地，所述步骤S2中划分为若干map块具体步骤包括：

步骤S21、将所述WLP-FDTD时域方程组合并为矩阵公式，公式为：

其中，q为拉盖尔阶数，[A]为稀疏矩阵，

为待求电场量，

为已知电流源矢量，

为已知磁流源矢量，[β^q-1]为0～q-1阶电场和磁场分量累加项，所述公式按拉盖尔阶数由低到高迭代计算待求电场量；

步骤S22、矩阵降维：根据所述火箭模型为柱形体系结构，即火箭模型为绕z轴旋转的对称体，因此不考虑

方向的网格划分；

步骤S23、限定边界：对待求电场量矩阵限定其计算的左右矩阵：

其中，Δρ、Δz为ρ、z方向的单位网格尺寸，I_max，K_max为ρ、z方向的最外边界，即为Mur/PML边界；

步骤S24、划分map块：利用HDFS平台对限定边界后的待求电场量矩阵进行网格拆分，划分为若干map块。

优选地，所述步骤S3中柱坐标系WLP-FDTD内各网格点的map“键-值”对设置具体步骤包括：

步骤S31、展开所述待求电场量矩阵[E^q]：为简化稀疏矩阵[A]的排列，并将其设置为便于程序运行的分块带状稀疏矩阵，拓展并展开所述待求电场量矩阵；

其中，I_max，K_max为ρ、z方向的Mur/PML边界，q为拉盖尔阶数，

这里i，k分别代表ρ，z方向的空间位置，

为虚设为零的待求电场分量；

步骤S32、确定所述待求电场量矩阵[E^q]排列顺序：对所述电场量矩阵[E^q]的展开值，确定所述待求电场量矩阵[E^q]排列顺序，即

其中，I_max，K_max为ρ、z方向的Mur/PML边界，q为拉盖尔阶数，

表示结果取整；

步骤S33、确定“键-值”对：将待计算的电场值行号编号，即相应的ρik、

和zik作为“键”数据，将相应网格点的电场仿真数据作为“值”数据，即计算获得的

优选地，所述步骤S4中是否边界点判断方法为：首先根据FDTD吸收边界原理划分外边界线，即将

和

这四条边上的

电场值单独设为Mur/PML吸收边界，该边界上的电场值单独计算；其次将若干map块的分界线作为内边界线；最后在HDFS平台上判断所有并行计算的“键-值”对落于区域上的位置，“键-值”对计算电场值过程中所需变量值均落于内边界上时，该组数据为边界点，数据存储于A组，否则该组数据为非边界点，数据存储于B组。

优选地，步骤S4中A、B两组内电场值的具体计算过程包括：

步骤S41、A组数据为若干边界点的电场值计算，计算过程中电场值的计算由相邻块边界点的磁场量进行传递计算；

步骤S42、B组数据为若干非边界点的电场值计算，计算过程中所述电场值由自身磁场量进行独立迭代计算；

步骤S43、电场值由低阶到高阶迭代更新计算：根据WLP-FDTD时域方程组公式及所述步骤S41和S42对所有待计算电场值进行迭代更新，按照拉盖尔阶数q，由低阶到高阶迭代计算。

有益效果：本发明将大数据并行技术引入到火箭模型的WLP-FDTD时域方程的计算求解中，解决了三维柱坐标下超大尺寸细微结构的电磁问题，提高超大尺寸细微结构模型在电磁环境下的计算效率，降低硬件成本，同时使得电磁仿真结果更加贴近实际工程应用环境。

附图说明

图1为本发明的流程图；

图2为传统三维柱坐标WLP-FDTD算法划分map块示意图；

图3为本发明中火箭载体示意图；

图4为本发明中火箭载体在雷电云层中以亚纳秒脉冲波打入后示意图；

图5为本发明中柱坐标WLP-FDTD算法划分map块示意图；

图6为实施例中本发明的仿真结果图。

具体实施方式

以下结合附图对本方案做进一步的说明和解释。

本发明解决的关键问题是构建大数据多map并行技术中的“键-值”对，求解WLP的大型稀疏矩阵，从而实现WLP的由低到高阶的迭代计算，因此本发明提供了一种基于大数据并行算法求解火箭仿真模型时域方程的方法。

本发明以传统三维柱坐标WLP-FDTD算法为基础，其稀疏矩阵，如附图2所示，这里给出一个100*100的非零元素(黑色部分)分布矩阵，可见WLP的系数矩阵为一个高度稀疏的带状矩阵。以往的诸多柱坐标WLP算法中，均采用三元组法对其进行求解。但事实证明，对于超大尺寸细微结构，即使三元组高效求解方案仍然会导致计算困难的情况发生，因而在实际的工程项目中无法得到广泛应用。

应对实际工程项目中，往往存在超大尺寸，精细结构的三维柱坐标模型需要仿真计算的情况，如火箭模型，为使仿真计算的模型尺寸和网格剖分与实体一致，必须通过某种手段将现有模型拆分后并行计算，否则以现有的硬件条件，想要实现其快速计算，只能利用天河机/银河机，这会大大的增加计算的运行成本。本发明利用大数据Hadoop的多map技术来并行求解三维柱坐标系下WLP-FDTD算法问题，如附图1的处理过程。

为了更好的解释说明本发明的执行方案，这里给出本发明的总体设计流程。

火箭模型建模

附图4所示为火箭模型载体的计算环境，实则为雷电云层环境中亚纳秒脉冲在火箭载体中的传播情况研究。理论上，火箭模型的柱坐标系下WLP-FDTD计算方程组如下所示：

其中，

和

分别为电场矢量和磁场矢量，

分别是沿ρ，

Z坐标轴的是一阶中心差分算子，ε为介质介电系数，μ为磁导系数，s＞0是时间尺度因子。

方程组(1)所示的FDTD时域方程，适用于包括火箭模型在内的所有柱形体系结构的电磁仿真计算。方程组(1)在仿真计算中内存消耗方面具有较大的劣势，极有可能导致大尺寸细微结构模型仿真过程中无法计算的情况发生(即使当前比较前沿的高效处理方法也无法避免)。

显然，如果对实际工程项目中需计算如图4所示的火箭载体进行仿真计算，在该尺寸下以图5的方式进行网格剖分(每1cm一个网格)，势必会造成由于模型尺寸过大，网格剖分过细而无法建模和计算。这里需利用大数据Hadoop的多map并行技术进行处理，因此，如图1所示，本发明的核心内容共分为三个部分：(1)多map块划分；(2)map“键-值”对设置；(3)内外边界计算数据分组。下面对此三部分进行具体阐述说明。

(1)多map块划分

WLP-FDTD的计算原理中，柱坐标系电场方程组(1)合在一起的矩阵为：

[A][E^q]＝[J^q]+[M^q]+[β^q-1] (2)

其中，q为拉盖尔的阶数，[A]为稀疏矩阵(每行有9个非零元素，且与q无关)，

为待求电场量，

为已知电流源矢量，

为已知磁流源矢量，[β^q-1]为0～q-1阶电场和磁场分量累加项，显然式(2)等式右侧为已知量，即式(2)的求解过程按照拉盖尔基函数阶数q，由低到高渐进迭代计算。因此，在FDTD的电磁场计算中，相邻网格节点数据存在着较强的相互依赖性，其单位网格上的电磁场计算必须利用相邻网格点上的数据进行迭代计算。WLP-FDTD电磁场计算中虽然将传统的时域变换为拉盖尔域，但其相邻网格点上的电磁场数据的依赖性是不变的。因此WLP-FDTD算法单纯的套用传统大数据方法是不合理的，还需考虑到相邻网格点间的相互依存关系。即图5所示内容，对实际尺寸火箭载体的多map块划分必须以WLP-FDTD的计算原理为基准，否则会出现计算结果有误。事实上，本发明在火箭载体进行多map块划分时，采用的方式是对计算模型通过在大数据Hadoop的HDFS平台上进行网格拆分的方式。但为了保证仿真结果的准确，必须对后期的数据结果进行内外边界分组，以便WLP-FDTD算法中相邻网格点上的电磁场值间的迭代运算(详见(3)内外边界计算数据分组)。

(2)map“键-值”对设置

式(2)的电场方程组的矩阵形式可简写为如下形式：

[A][E^q]＝[Q] (3)

如附图5所示，令柱坐标下WLP-FDTD计算区域的左右边界为：

其中，Δρ，Δz为ρ，z方向的单位网格尺寸，I_max，K_max为图5所示的ρ，z方向的最外边界，即为Mur/PML边界。这里并未考虑

方向的网格划分，是因为如图2的Part2所示，火箭载体可以理解为绕z轴旋转的对称体，因此，可以不作为计算标准，对矩阵降维处理。

进而，电场方程组[E^q]可展开为：

其中，

这里i，k分别代表ρ，z方向的空间位置，如果

为整网格点，则e_ρ和e_z则分别在ρ(i)和z(k)方向为半网格，参考图5左侧的电磁场分布图。

为虚设的待求场分量，实际计算过程只需将其值赋零，火箭载体可以理解为绕z轴旋转的对称体，因此，可以不作为计算标准，因此

方向不做考虑。这样做的目的是为了简化[A]矩阵的排列，使得矩阵行列号和空间坐标之间的关系有规律可循，并使[A]矩阵成为分块带状稀疏矩阵，从而便于程序实现。

设[e^q]的第n个元素为

则有

式中，

表示结果取整。

一旦矩阵[E^q]的排列顺序确定下来，那么[Q]矩阵元素和[A]矩阵中的非零元素也就随之确定。

因此本发明中符合式(5)矩阵分布的map块的“键-值”对设置应为：

“键”：待计算电/磁场的行号，以此为编号；“值”：一个网格点的电/磁场仿真数据。

即，如果“键”为z11，则对应的“值”为

(3)内外边界计算数据分组

经(1)中所述分析，可知map块划分必须考虑相邻map块上的网格点上的电磁场值的相互依赖性，如附图5所示，黑色实线为计算模型的外吸收边界Mur’s/PML，虚线为map块之间相互依赖的内边界，事实上，这里在进行HDFS平台划分之前既已确定了内外边界，其内外边界划分原理可参考FDTD吸收边界原理。即附图5右侧所示黑色实线的外边界即不属于A、B的任何一种情况的判别方法是利用电磁场的空间位置进行判别的，以

为例，其所属的外边界场量为

和

这四条边上的

场必须单独设为Mur/PML吸收边界，这是FDTD处理吸收边界问题的基本法则，此时，此四条边界上个的

并不能进行多个map块间的传递计算。传统的FDTD方法对于处于Mur/PML边界上的场值也是单独计算。因此，据柱坐标系下WLP-FDTD的计算原理和吸收边界设定原理，这里必须将内外边界的计算进行分组，如此，能够更好的使得外边界上(Mur/PML吸收边界)的场量进行单独计算，而属于内边界情况的场量可以进行map块上的互相依托计算。这里Mur/PML吸收边界上的场量计算可参考论文“Improved ADI iterativealgorithm with two-step Gauss-Seidel procedure for efficient Laguerre-basedBOR-FDTD method”。

这里根据前文柱坐标系下WLP-FDTD计算方程式(1)中

和

的计算进行解释说明，参照图5，以下情况是不包含黑色的Mur/PML外边界。

(3-1)

计算规则

式(1)中的

需要利用到其中的

在考虑扣除掉黑色外边界的情况下，属于虚线内边界和其内部的场量计算规则如下：

①如果计算

所需的磁场分量

均在附图5的map2块区域，则

可以独立计算。

②如果计算

所需

在map1块的内边界上，则需map1块将

值从map1块传递至map2块，方可正确求解

(3-2)

和

计算规则

式(1)的

需要利用到其中的

和

式(1)的

需要利用到其中的

和

区别于

在于不同的电磁场量间的循环迭代运算，其规则可参考(3-1)中的规则①和②。

因此，本发明利用分组方案，将多个map块并行计算所得电/磁场值，进行是否为内外边界点判断，并且将所得分组结果在Hadoop的HDFS平台上构建两组组数据A、B。如果属于本块区域内，则属于B组，并将其直接代入柱坐标系WLP-FDTD方程中，进行逐步迭代。如果属于相邻块内边界点的场值，则属于A组，在需计算某个电/磁场时，方便随时调用。

完成以上三部分，并顺利读取A，B数据后，可利用柱坐标系WLP-FDTD迭代方程有低阶到高阶更新电磁场值。

实施例

现给出一个具体仿真案例作为补充说明。如图4所示的火箭载体在雷电云层中以亚纳秒脉冲波打入后，其传播情况的研究。如附图3所示，云层底边长为100.0m，高位20.0m，云层中介质层电导率σ_s＝1.0×10³s/m，相对介电常数ε_rs＝10.0。火箭载体的尺寸为箭头长为10.0m，箭体长为70.0m，箭头和箭体半径为4.0m，尾翼长8.0m，尾翼半径为8.0m，亚纳秒脉冲波为：E_i(t)＝E₀cos(2πf₀t)exp(-4π((t-t₀)/τ)²)，其中，E₀＝1000V/m，f₀＝10GHz，τ＝1.5/f₀，t₀＝1.75/f₀。WLP的时间影响因子s＝8×10¹¹，阶为q＝210。

表1

图6为本发明和传统柱坐标WLP-FDTD和利用MPI+CUDA快速处理技术后的WLP-FDTD在观察点处的时域波形对比。可见，三种计算方式的仿真结果吻合度较高，可以证明本发明方法的正确性。此外，表1为本发明和传统柱坐标WLP-FDTD和利用MPI+CUDA技术后的WLP-FDTD计算时间对比，表1表明，本发明的方法与传统柱坐标WLP-FDTD和当前较流行的快速并行算法对比，其计算效率有非常显著的提升。至此，可以证明，本发明算法的不仅能够使超大尺寸细微结构模型得以简便计算，节约科研成本，还可以大大提高计算效率。

以上所述仅是本发明的优选实施方式，应当指出：对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (6)

1.一种基于大数据并行算法求解火箭仿真模型时域方程的方法，其特征在于，包括以下步骤：

步骤S1、对火箭模型建模，获取WLP-FDTD时域方程组：在雷电云层中，雷电以亚纳秒脉冲波打入火箭载体，以柱形体系结构建立火箭模型，基于ρ，

z三个坐标轴建立火箭电磁仿真模型；

步骤S2、对火箭电磁仿真模型中WLP-FDTD时域方程组进行多map划分：将所述模型导入大数据的HDFS平台，利用所述HDFS平台对时域方程组限定边界后进行网格拆分，划分为若干map块；

步骤S3、对所述map块内若干待求电场值进行“键-值”对设置：对待求电场量矩阵进行展开，将待计算的电场的行号作为“键”数据，将所述柱坐标系WLP-FDTD内各网格点的电场仿真数据作为“值”数据；

步骤S4、对所述“键-值”进行分组迭代计算，直至计算到目标拉盖尔阶数时的待求电场值：将若干map块并行计算所得电场值和磁场值，进行是否为边界点判断，并且将所得分组结果在的HDFS平台上构建两组数据A、B；其中将属于本map块内的电场值和磁场值归于B组，并将其直接代入柱坐标系WLP-FDTD时域方程组中，进行逐步迭代；将属于相邻块内边界点的电场值和磁场值存储于A组，在需计算B组内电场值时进行调用。

2.根据权利要求1所述的一种基于大数据并行算法求解火箭仿真模型时域方程的方法，其特征在于，所述电磁仿真模型中时域方程组为：

其中，

为i坐标轴方向的电场矢量，

为i坐标轴方向的磁场矢量，

为沿着i坐标轴一阶中心差分算子，这里

为三个坐标轴，ε为介质介电系数，μ为磁导系数，s为大于零的时间尺度因子，q为拉盖尔阶数。

3.根据权利要求1所述的一种基于大数据并行算法求解火箭仿真模型时域方程的方法，其特征在于，所述步骤S2中若干map块划分的具体步骤包括：

步骤S21、将所述WLP-FDTD时域方程组合并为矩阵公式，公式为：

其中，q为拉盖尔阶数，[A]为稀疏矩阵，

为待求电场量，

为已知电流源矢量，

为已知磁流源矢量，[β^q-1]为0～q-1阶电场和磁场分量累加项，所述公式按拉盖尔阶数由低到高迭代计算待求电场量；

步骤S22、矩阵降维：根据所述火箭模型为柱形体系结构，即火箭模型为绕z轴旋转的对称体，因此不考虑

方向的网格划分；

步骤S23、限定边界：对待求电场量矩阵限定其计算的左右矩阵：

其中，Δρ、Δz为ρ、z方向的单位网格尺寸，I_max，K_max为ρ、z方向的最外边界，即为Mur/PML边界；

步骤S24、划分map块：利用HDFS平台对限定边界后的待求电场量矩阵进行网格拆分，划分为若干map块。

4.根据权利要求1所述的一种基于大数据并行算法求解火箭仿真模型时域方程的方法，其特征在于，所述步骤S3中柱坐标系WLP-FDTD内各网格点的map“键-值”对设置具体步骤包括：

步骤S31、展开所述待求电场量矩阵[E^q]：为简化稀疏矩阵[A]的排列，并将其设置为便于程序运行的分块带状稀疏矩阵，拓展并展开所述待求电场量矩阵；

其中，I_max，K_max为ρ、z方向的Mur/PML边界，q为拉盖尔阶数，

这里i，k分别代表ρ，z方向的空间位置，

为虚设为零的待求电场分量；

步骤S32、确定所述待求电场量矩阵[E^q]排列顺序：对所述电场量矩阵[E^q]的展开值，确定所述待求电场量矩阵[E^q]排列顺序，即

其中，I_max，K_max为ρ、z方向的Mur/PML边界，q为拉盖尔阶数，

表示结果取整；

步骤S33、确定“键-值”对：将待计算的电场值行号编号，即相应的ρik、

和zik作为“键”数据，将相应网格点的电场仿真数据作为“值”数据，即计算获得的

5.根据权利要求1所述的一种基于大数据并行算法求解火箭仿真模型时域方程的方法，其特征在于，所述步骤S4中是否边界点判断方法为：

根据FDTD吸收边界原理划分外边界线，即将

和

这四条边上的

电场值单独设为Mur/PML吸收边界，该边界上的电场值单独计算；利用所述若干map块的分界线作为内边界线；

在HDFS平台上判断所有并行计算的“键-值”对落于区域上的位置，“键-值”对计算电场值过程中所需变量值均落于内边界上时，该组数据为边界点，数据存储于A组，否则该组数据为非边界点，数据存储于B组。

6.根据权利要求1所述的一种基于大数据并行算法求解火箭仿真模型时域方程的方法，其特征在于，步骤S4中A、B两组内电场值的具体计算过程包括：

步骤S41、A组数据为若干边界点的电场值计算，计算过程中电场值的计算由相邻块边界点的磁场量进行传递计算；

步骤S42、B组数据为若干非边界点的电场值计算，计算过程中所述电场值由自身磁场量进行独立迭代计算；

步骤S43、电场值由低阶到高阶迭代更新计算：根据WLP-FDTD时域方程组公式及所述步骤S41和S42对所有待计算电场值进行迭代更新，按照拉盖尔阶数q，由低阶到高阶迭代计算。

Images