CN114444353A - 基于有限元方法的二维声场模型的声场获取方法 - Google Patents
基于有限元方法的二维声场模型的声场获取方法 Download PDFInfo
- Publication number
- CN114444353A CN114444353A CN202210054039.6A CN202210054039A CN114444353A CN 114444353 A CN114444353 A CN 114444353A CN 202210054039 A CN202210054039 A CN 202210054039A CN 114444353 A CN114444353 A CN 114444353A
- Authority
- CN
- China
- Prior art keywords
- sound field
- unit
- node
- sound
- stiffness matrix
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 88
- 239000011159 matrix material Substances 0.000 claims abstract description 61
- 239000013598 vector Substances 0.000 claims abstract description 51
- 238000004364 calculation method Methods 0.000 claims abstract description 26
- 238000004800 variational method Methods 0.000 claims abstract description 5
- 230000006870 function Effects 0.000 claims description 45
- 230000008569 process Effects 0.000 claims description 14
- 238000012360 testing method Methods 0.000 claims description 12
- 239000013535 sea water Substances 0.000 claims description 11
- 230000014509 gene expression Effects 0.000 claims description 10
- 230000005855 radiation Effects 0.000 claims description 6
- 230000009466 transformation Effects 0.000 claims description 6
- 238000010521 absorption reaction Methods 0.000 claims description 3
- 238000013507 mapping Methods 0.000 claims description 3
- 238000012876 topography Methods 0.000 claims description 3
- 230000005540 biological transmission Effects 0.000 abstract description 2
- 239000000243 solution Substances 0.000 description 17
- 230000006872 improvement Effects 0.000 description 6
- 230000010354 integration Effects 0.000 description 6
- 230000008859 change Effects 0.000 description 3
- 230000007613 environmental effect Effects 0.000 description 3
- 238000010586 diagram Methods 0.000 description 2
- 239000012088 reference solution Substances 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 239000012086 standard solution Substances 0.000 description 2
- 230000003044 adaptive effect Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000012821 model calculation Methods 0.000 description 1
- 230000005477 standard model Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
Abstract
本发明属于水下声传播技术领域，具体地说，涉及一种基于有限元方法的二维声场模型的声场获取方法，包括：根据获取目标区域的海洋环境参数，将该目标区域离散为有限个计算单元；根据单元形状和节点数，构造每个计算单元的基函数；使用标准Galerkin变分法，推导出Helmholtz方程的弱形式；根据节点的各局部编号和各总体编号的一一对应的关系，将所有计算单元的单元刚度矩阵集成总体刚度矩阵；根据节点的各局部编号和各总体编号的一一对应的关系，将所有计算单元的单元外载荷向量集成总体外载荷向量；将上述结果组成线性代数方程组，作为二维声场模型，得到声场，即各节点的声压；再通过选择的基函数和插值方法，得到空间各节点的声压值。
Description
技术领域
本发明属于水声物理与水声工程中的水下声传播技术领域，具体地说，涉及一种基于有限元方法的二维声场模型的声场获取方法。
背景技术
声波是目前所知的唯一一种可以在海水介质中远距离传播的物理形式，是探测海洋环境、开发海洋资源、实现水下信息传输的重要载体。海水中声波的传播特性研究是海洋学研究的重要领域之一，也是研究其他海洋声学问题的基础。
目前常用的声场计算模型难以准确计算复杂环境中的声场。例如：射线法仅适用于计算高频声波的传播；抛物方程方法对于水平变化剧烈的环境中的声场计算精度较差；未经特殊处理的简正波方法和波数积分方法无法计算水平变化环境中的声场；耦合简正波方法无法精确计算非对称的地形起的声场变化。
因此，现有的方法无法精确计算二维复杂海洋环境中的声场，无法解决二维复杂环境中的水下声传播问题。
发明内容
为解决现有技术存在的上述缺陷，本发明提出了基于有限元方法的二维声场模型的声场获取方法，是一种用于计算二维水下声传播问题的有限元方法，该方法包括：
根据获取目标区域的海洋环境参数，离散该目标区域；将该目标区域离散为有限个计算单元，各个计算单元之间以节点相连；根据单元形状和节点数，构造每个计算单元的基函数；
使用标准Galerkin变分法，推导出Helmholtz方程的弱形式；
建立各计算单元的单元刚度矩阵，根据节点的各局部编号和各总体编号的一一对应的关系，将所有计算单元的单元刚度矩阵集成总体刚度矩阵；
建立各计算单元的单元外载荷向量，根据节点的各局部编号和各总体编号的一一对应的关系，将所有计算单元的单元外载荷向量集成总体外载荷向量；
将总体刚度矩阵、声压向量和总体外载荷向量组成线性代数方程组，作为二维声场模型，求解该二维声场模型得到声场，即各节点的声压；
再通过选择的基函数和插值方法，得到空间各点的声压值。
作为上述技术方案的改进之一，所述海洋环境参数包括：海水海底声速、海水海底密度、该目标区域的海底地形、声源强度和声源位置。
作为上述技术方案的改进之一，所述根据获取目标区域的海洋环境参数，将该目标区域离散为有限个计算单元，各个计算单元以节点相连；根据单元形状和节点数，构造每个计算单元的基函数；其具体过程包括：
使用三节点三角形单元作为计算单元，根据获取目标区域的海洋环境参数，对该目标区域进行剖分，将该目标区域离散为有限个三节点三角形单元，每个计算单元的基函数为Nn(x,z),n＝1,2,3；其中，n表示节点的局部编号，(x,z)为所建立的直角坐标系；z轴通过声源且垂直向下；x轴平行于海面；
假设三个节点表示为An(xn,yn),n＝1,2,3，并且定义Jacobi矩阵的行列式|J|＝(x2-x1)(y3-y1)-(x3-x1)(y2-y1)，其中，xn为第n个节点的横坐标；yn为第n个节点的纵坐标；得到计算单元的三个节点对应的三个基函数的表达式：
作为上述技术方案的改进之一，所述使用标准Galerkin变分法，推导出Helmholtz方程的弱形式；其具体过程包括：
对于无限长线源问题，建立直角坐标系(x,z)，让z轴通过声源且垂直向下，x轴平行于海面；此时不均匀介质中的亥姆霍兹方程表示为
其中，ρ表示海水密度，k(x,z)表示波数，记为
其中，c(x,z)表示声速，f为声源频率；Sω表示声源强度，zs表示声源深度；
为梯度算符；p为声压；δ为狄拉克函数；
使用标准Galerkin变分方法推导控制方程的弱形式，具体过程如下：
对上述亥姆霍兹方程的两边同时乘以试验函数q(x,z)，并在该目标区域的求解域Ω上进行积分，得到如下表达式：
对上述表达式应用高斯定理，得到Helmholtz方程的弱形式：
其中，Γ表示边界，Ω表示目标区域的求解域；
应用三节点三角形单元对应的拉格朗日插值函数对该目标区域的声压和试验函数进行离散，得到离散的声压p(x,z)和离散的试验函数q(x,z)；
q(x,z)＝Nm(x,z),m＝1,2,3
推导出离散后的不考虑边界条件的Helmholtz方程：
其中，pn是待求的第n个节点的声压值；
是pn的系数，将其作为单元刚度矩阵的元素kmn；n和m都是节点的局部编号，n＝1,2,3，m＝1,2,3；Nn和Nm分别为节点n和m对应的基函数；
和
为Nn和Nm的梯度；
对于单元外载荷向量，将离散后的试验函数q(x,z)带入Helmholtz方程的弱形式的右侧，得到
将其作为单元外载荷向量的第m个元素。
作为上述技术方案的改进之一，所述建立各计算单元的单元刚度矩阵，根据节点的各局部编号和各总体编号的一一对应，将所有计算单元的单元刚度矩阵集成总体刚度矩阵；其具体过程包括：
建立各计算单元的单元刚度矩阵，遍历所有计算单元，根据计算单元中每个节点的局部编号与总体刚度矩阵中的各元素的对应关系，将单元刚度矩阵的各个元素分别与总体刚度矩阵的各个对应元素相加，进而将单元刚度矩阵的元素集成到总体刚度矩阵中的对应位置。
作为上述技术方案的改进之一，所述将总体刚度矩阵、声压向量和总体外载荷向量组成线性代数方程组，作为二维声场模型，并求解该二维声场模型，得到声场；再通过选择的基函数和插值方法，得到空间各节点的声压值；其具体过程包括：
对于无限大区域，声场满足辐射条件：当r→∞时，声压p→0；
此时需要对无限大区域进行截断，采用人工边界条件模拟出射条件，选取完美匹配层方法模拟辐射条件；其中，完美匹配层为认为设置的吸收层，用于吸收各个角度的出射波；
使用完美匹配层模拟无限大区域，进行坐标变换，将实坐标轴映射为复坐标轴，具体过程如下：
对于深度有限的介质，假设需要计算的物理域的水平距离为R，在物理域的右侧设置一个厚度为δx的完美匹配层，完美匹配层中的声速和密度与相邻的左侧介质相同；
其中，
为坐标变换后的复坐标；i为虚数单位；ξ为积分变量；σx(ξ)为人为规定的衰减函数，
其中，σxm是决定完美匹配层最大衰减量的常数；
对于完美匹配层，将pn的系数记为kmn，kmn为完美匹配层中单元刚度矩阵的元素：
利用上述得到的针对完美匹配层的单元刚度矩阵的元素，建立总体刚度矩阵；再结合总体外载荷向量和声压向量组成线性代数方程组，作为二维声场模型，求解该二维声场模型得到声场；
Kp＝F
其中，K表示总体刚度矩阵，F表示总体外载荷向量，p表示声场；
再根据得到的声场，通过选择的基函数和插值方法，得到空间各节点的声压值。
本发明与现有技术相比的有益效果是：
1、本发明的方法适用范围广，与射线法相比，适用于高频和低频问题，与抛物方程方法相比，适用于计算水平变化剧烈的海洋环境中的声场，与简正波和波数积分方法相比，适用于计算水平变化环境中的声场；
2、本发明的方法使用有限元方法模拟水下声传播，精度高，普适性强，能够解决现有模型无法解决的二维复杂环境中的水下声传播问题；从理论上讲，有限元方法可以通过不断细分网格，获得理想精度的计算结果；
3、本发明的方法易于并行，随着计算机性能的不断提升与并行技术的不断发展，可以有效地提升有限元模型的计算效率。
附图说明
图1是本发明的一种基于有限元方法的二维声场模型的声场获取方法的方法流程图；
图2a是深度有限的波导中声源频率为100Hz时使用波数积分方法计算的双层介质的传播损失示意图；
图2b是深度有限的波导中声源频率为100Hz时使用有限元模型计算的双层介质的传播损失示意图；
图3是声源频率100Hz时接收深度60m处的传播损失示意图；
图4a是声源频率为30Hz，声源深度20m，水平距离为50m时使用简正波模型DGMCM计算的可穿透海底斜坡波导中的声场传播损失；
图4b是声源频率为30Hz，声源深度20m，水平距离为50m时使用有限元模型计算的可穿透海底斜坡波导中的声场传播损失；
图4c是声源频率为50Hz，声源深度20m，水平距离为50m时使用简正波模型DGMCM计算的可穿透海底斜坡波导中的声场传播损失；
图4d是声源频率为50Hz，声源深度20m，水平距离为50m时使用有限元模型计算的可穿透海底斜坡波导中的声场传播损失。
具体实施方式
现结合附图和实例对本发明作进一步的描述。
如图1所示，本发明提供了一种基于有限元方法的二维声场模型的声场获取方法，是一种用于计算二维水下声传播问题的有限元方法，该方法包括：
根据获取目标区域的海洋环境参数，离散该目标区域；将该目标区域离散为有限个计算单元，各个计算单元以节点相连；根据单元形状和节点数，构造每个计算单元的基函数；
其中，所述海洋环境参数包括：海水海底声速、海水海底密度、该目标区域的海底地形、声源强度和声源位置。由于三角形单元对复杂海洋环境边界的适应性更强，尤其是海面和海底会呈现出复杂的边界情况，因此，使用三角形单元对物理域(目标区域)进行剖分。
具体地，使用三节点三角形单元作为基本计算单元，根据获取目标区域的海洋环境参数，对该目标区域进行剖分，将该目标区域离散为有限个三节点三角形单元，每个基本计算单元的基函数为Nn(x,z),n＝1,2,3；其中，n表示节点的局部编号，(x,z)为所建立的直角坐标系；z轴通过声源且垂直向下；x轴平行于海面；
假设三个节点的坐标分别表示为An(xn,yn),n＝1,2,3，并且定义Jacobi矩阵的行列式|J|＝(x2-x1)(y3-y1)-(x3-x1)(y2-y1)，其中，xn为第n个节点的横坐标；yn为第n个节点的纵坐标；得到计算单元的三个节点对应的三个基函数的表达式；
使用标准Galerkin变分法，推导出Helmholtz方程的弱形式；
具体地，对于无限长线源问题，建立直角坐标系(x,z)，让z轴通过声源且垂直向下，x轴平行于海面；此时不均匀介质中的Helmholtz方程(亥姆霍兹方程)表示为
其中，ρ表示海水密度，k(x,z)表示波数，记为
其中，c(x,z)表示声速，f为声源频率；Sω表示声源强度，zs表示声源深度；
为梯度算符；p为声压；δ为狄拉克函数；
使用标准Galerkin变分方法推导控制方程的弱形式，具体过程如下：
对上述亥姆霍兹方程的两边同时乘以试验函数q(x,z)，并在该目标区域的求解域Ω上进行积分，得到如下表达式：
对上述表达式应用高斯定理，得到下述表达式，即Helmholtz方程的弱形式(待求控制方程)：
其中，Γ表示边界，Ω表示目标区域的求解域；
应用三节点三角形单元对应的拉格朗日插值函数(基函数)对该目标区域的声压和试验函数进行离散，得到离散的声压p(x,z)和离散的试验函数q(x,z)；
q(x,z)＝Nm(x,z),m＝1,2,3
推导出离散后的不考虑边界条件的Helmholtz方程
其中，pn是待求的第n个节点的声压值；
是pn的系数，将其作为单元刚度矩阵的元素kmn；n和m都是节点的局部编号，n＝1,2,3，m＝1,2,3；Nn和Nm分别为节点n和m对应的基函数；
和
为Nn和Nm的梯度；
对于单元外载荷向量，我们将离散后的试验函数q(x,z)带入Helmholtz方程的弱形式的右侧，就可以得到
将其作为单元外载荷向量的第m个元素。
建立各计算单元的单元刚度矩阵，根据节点的各局部编号和各总体编号的一一对应的关系，将所有计算单元的单元刚度矩阵集成总体刚度矩阵；
具体地，建立各计算单元的单元刚度矩阵，遍历所有计算单元，根据计算单元中每个节点的局部编号与总体刚度矩阵中的各元素总体编号的对应关系，将单元刚度矩阵的各个元素分别与总体刚度矩阵的各个对应元素相加，进而将单元刚度矩阵的元素集成到总体刚度矩阵中的对应位置；
例如，我们假设单元刚度矩阵为k，对于三节点三角形单元，它包含九个元素：k11,k12,k13,k21,k22,k23,k31,k32,k33。若单元中各节点的局部编号1,2,3分别对应总体编号的23,34,35,则将单元刚度矩阵的九个元素分别与总体刚度矩阵的九个对应元素相加，设总体刚度矩阵为K，则这九个元素分别是K23,23,K23,34,K23,35,K34,23,K34,34,K34,35,K35,23,K35,34,K35,35。
建立各计算单元的单元外载荷向量，根据节点的各局部编号和各总体编号的一一对应的关系，将所有计算单元的单元外载荷向量集成总体外载荷向量；
具体地，将单元外载荷向量集成为总体外载荷向量的方法与刚度矩阵相同。
建立各计算单元的单元外载荷向量，遍历所有计算单元，根据计算单元中每个节点的局部编号与总体外载荷向量中的各元素的总体编号的对应关系，将单元刚度外载荷向量的各个元素分别与总体外载荷向量的各个对应元素相加，进而将单元刚度外载荷向量的元素集成到总体外载荷向量中的对应位置。
例如，我们假设单元外载荷向量为f，对于三节点三角形单元，该单元外载荷向量f包含三个元素：f1,f2,f3；若单元节点的局部编号1,2,3分别对应总体编号的23,34,35,则将单元外载荷向量的三个元素分别与总体外载荷向量的三个对应元素相加，设总体外载荷向量为F，则这三个元素分别是F23,F34,F35。
将边界条件的影响带入上述的总体刚度矩阵和总体外载荷向量，基于推导出的Helmholtz方程的弱形式，将总体刚度矩阵、声压向量和总体外载荷向量组成线性代数方程组，作为二维声场模型(二维水下声传播模型)，并求解该二维声场模型，得到声场，即各节点的声压；再通过选择的基函数和插值方法，得到空间各点的声压值。
具体地，边界条件的影响包括：第一类边界条件和第二类边界条件；其中，第一类边界条件为Dirichlet边界条件；
对于第一类边界条件，声压等于常数，使用“对角元素置一”法；
对于第二类边界条件，常见的边界条件是声压的法向导数为零，此时Helmholtz方程左侧的第一项
为零，不需要做特别的处理。
对于无限大区域，声场满足辐射条件，即熄灭条件：当r→∞时，声压p→0；
此时需要对无限大区域进行截断，采用人工边界条件模拟出射条件，我们选取完美匹配层方法模拟辐射条件；其中，完美匹配层为认为设置的吸收层，用于吸收各个角度的出射波；
使用完美匹配层模拟无限大区域，进行坐标变换，将实坐标轴映射为复坐标轴，使人工边界处的反射系数为0，在完美匹配层中吸收各个角度的出射波；
在完美匹配层中，将实数坐标变换到复数域。例如，对于深度有限的介质，我们设需要计算的物理域的水平距离为R，在物理域的右侧设置一个厚度为δx的完美匹配层，完美匹配层中的声速和密度与相邻的左侧介质相同；
其中，
为坐标变换后的复坐标；i为虚数单位；ξ为积分变量；σx(ξ)为人为规定的衰减函数，
其中，σxm是决定完美匹配层最大衰减量的常数；
对于完美匹配层，将pn的系数记为kmn，kmn为完美匹配层中单元刚度矩阵的元素：
利用上述得到的针对完美匹配层的单元刚度矩阵的元素，建立总体刚度矩阵；再结合总体外载荷向量和声压向量组成线性代数方程组，将其作为二维声场模型，并根据该二维声场模型，求解该二维声场模型得到声场，即各节点的声压；
当得到总体刚度矩阵和总体外载荷向量时，我们就可以求解未知的节点声压。可以写成如下方程
Kp＝F
其中，K表示总体刚度矩阵，F表示总体外载荷向量，p表示声场，即各节点的声压向量。对于上述的有限元方法，总体刚度矩阵K是一个对称稀疏阵。因此，我们使用常用的线性代数软件包LAPACK来存储矩阵K并求解线性代数方程组，提升内存的使用率，提高线性方程组的计算速度。
根据选取的基函数和插值方法，计算空间任一点的声场，得到空间各点的声压值，之后可以通过基函数与插值得到空间各点的声压值。根据与标准问题的解析解或现有标准模型的计算结果进行比对，进行误差分析。具体地，空间任一点的声场可以表示为，
我们经过上述操作和计算可以得到每个节点的声压pn，利用各节点声压和基函数就可以得到空间任意一点的声压
我们可以得到我们需要的任意一点的声压，对它们进行需要的处理，得到需要的数据。
对于存在解析解的标准问题，我们可以将有限元模型的计算结果与解析解进行比对，计算误差；对于不存在解析解的标准问题，我们可以使用常用的声场计算模型与有限元模型的计算结果进行比对，计算误差，从而验证有限元模型的正确性。其中，有限元模型的结果与解析解或者参考解相减，就是误差。
如图2a和2b所示，表示深度有限的波导中声源频率为100Hz时使用不同方法计算的双层介质的传播损失；其中，图2a表示波数积分方法计算的标准解；图2b表示使用有限元模型计算的有限元解；通过图2a和2b所示，对于深度有限的波导，半解析的波数积分方法计算的结果和有限元模型计算的结果吻合得非常好，有限元模型的精度非常高。
如图3所示，该图表示声源频率100Hz时接收深度60m处的传播损失，其中，点划线表示使用完美匹配层处理出射边界的有限元解，虚线表示波数积分方法计算的标准解。可以看到，两种方法计算的声场具有非常高的一致性，最大误差小于0.1dB，说明有限元模型的精度非常高。
如图4a-4d所示，声源频率分别为30Hz和50Hz，声源深度20m，水平距离为50m时使用不同方法计算的可穿透海底斜坡波导中的声场传播损失；
其中，图4a表示声源频率30Hz时简正波模型DGMCM计算的传播损失；
图4b表示声源频率30Hz时有限元模型计算的传播损失；
图4c表示声源频率50Hz时简正波模型DGMCM计算的传播损失；
图4d表示声源频率50Hz时有限元模型计算的传播损失。
这四张图表明在不同声源频率下，有限元模型计算的结果和简正波模型DGMCM计算的参考解具有非常高的一致性，精度非常高。同时说明，有限元模型不仅适用于计算不同频率的声场，也适用于计算水平变化的海洋环境中的声场。
最后所应说明的是，以上实施例仅用以说明本发明的技术方案而非限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行修改或者等同替换，都不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。
Claims (6)
1.一种基于有限元方法的二维声场模型的声场获取方法，该方法包括：
根据获取目标区域的海洋环境参数，离散该目标区域；将该目标区域离散为有限个计算单元，各个计算单元之间以节点相连；根据单元形状和节点数，构造每个计算单元的基函数；
使用标准Galerkin变分法，推导出Helmholtz方程的弱形式；
建立各计算单元的单元刚度矩阵，根据节点的各局部编号和各总体编号的一一对应的关系，将所有计算单元的单元刚度矩阵集成总体刚度矩阵；
建立各计算单元的单元外载荷向量，根据节点的各局部编号和各总体编号的一一对应的关系，将所有计算单元的单元外载荷向量集成总体外载荷向量；
将总体刚度矩阵、声压向量和总体外载荷向量组成线性代数方程组，作为二维声场模型，并求解该二维声场模型得到声场，即各节点的声压；
再通过选择的基函数和插值方法，得到空间各点的声压值。
2.根据权利要求1所述的基于有限元方法的二维声场模型的声场获取方法，其特征在于，所述海洋环境参数包括：海水海底声速、海水海底密度、该目标区域的海底地形、声源强度和声源位置。
3.根据权利要求1所述的基于有限元方法的二维声场模型的声场获取方法，其特征在于，所述根据获取目标区域的海洋环境参数，将该目标区域离散为有限个计算单元，各个计算单元以节点相连；根据单元形状和节点数，构造每个计算单元的基函数；其具体过程包括：
使用三节点三角形单元作为计算单元，根据获取目标区域的海洋环境参数，对该目标区域进行剖分，将该目标区域离散为有限个三节点三角形单元，每个计算单元的基函数为Nn(x,z),n＝1,2,3；其中，n表示节点的局部编号，(x,z)为所建立的直角坐标系；z轴通过声源且垂直向下；x轴平行于海面；
假设三个节点表示为An(xn,yn),n＝1,2,3，并且定义Jacobi矩阵的行列式|J|＝(x2-x1)(y3-y1)-(x3-x1)(y2-y1)，其中，xn为第n个节点的横坐标；yn为第n个节点的纵坐标；得到计算单元的三个节点对应的三个基函数的表达式：
4.根据权利要求1所述的基于有限元方法的二维声场模型的声场获取方法，其特征在于，所述使用标准Galerkin变分法，推导出Helmholtz方程的弱形式；其具体过程包括：
对于无限长线源问题，建立直角坐标系(x,z)，z轴通过声源且垂直向下，x轴平行于海面；此时不均匀介质中的亥姆霍兹方程表示为
其中，ρ表示海水密度，k(x,z)表示波数，记为
其中，c(x,z)表示声速，f为声源频率；Sω表示声源强度，zs表示声源深度；
为梯度算符；p为声压；δ为狄拉克函数；
使用标准Galerkin变分方法推导控制方程的弱形式，具体过程如下：
对上述亥姆霍兹方程的两边同时乘以试验函数q(x,z)，并在该目标区域的求解域Ω上进行积分，得到如下表达式：
对上述表达式应用高斯定理，得到Helmholtz方程的弱形式：
其中，Γ表示边界，Ω表示目标区域的求解域；
应用三节点三角形单元对应的拉格朗日插值函数对该目标区域的声压和试验函数进行离散，得到离散的声压p(x,z)和离散的试验函数q(x,z)；
q(x,z)＝Nm(x,z),m＝1,2,3
推导出离散后的不考虑边界条件的Helmholtz方程：
其中，pn是待求的第n个节点的声压值；
是pn的系数，将其作为单元刚度矩阵的元素kmn；n和m都是节点的局部编号，n＝1,2,3，m＝1,2,3；Nn和Nm分别为节点n和m对应的基函数；
和
为Nn和Nm的梯度；
对于单元外载荷向量，将离散后的试验函数q(x,z)带入Helmholtz方程的弱形式的右侧，得到
将其作为单元外载荷向量的第m个元素。
5.根据权利要求1所述的基于有限元方法的二维声场模型的声场获取方法，其特征在于，所述建立各计算单元的单元刚度矩阵，根据节点的各局部编号和总体编号的一一对应，将所有计算单元的单元刚度矩阵集成总体刚度矩阵；其具体过程包括：
建立各计算单元的单元刚度矩阵，遍历所有计算单元，根据计算单元中每个节点的局部编号与总体刚度矩阵中的各元素总体编号的对应关系，将单元刚度矩阵的各个元素分别与总体刚度矩阵的各个对应元素相加，进而将单元刚度矩阵的元素集成到总体刚度矩阵中的对应位置。
6.根据权利要求1所述的基于有限元方法的二维声场模型的声场获取方法，其特征在于，所述将总体刚度矩阵、声压向量和总体外载荷向量组成线性代数方程组，作为二维声场模型，并求解该二维声场模型得到声场，即各节点的声压；再通过选择的基函数和插值方法，得到空间各点的声压值；其具体过程包括：
对于无限大区域，声场满足辐射条件：当r→∞时，声压p→0；
此时需要对无限大区域进行截断，采用人工边界条件模拟出射条件，选取完美匹配层方法模拟辐射条件；其中，完美匹配层为认为设置的吸收层，用于吸收各个角度的出射波；
使用完美匹配层模拟无限大区域，进行坐标变换，将实坐标轴映射为复坐标轴，具体过程如下：
对于深度有限的介质，假设需要计算的物理域的水平距离为R，在物理域的右侧设置一个厚度为δx的完美匹配层，完美匹配层中的声速和密度与相邻的左侧介质相同；
其中，
为坐标变换后的复坐标；i为虚数单位；ξ为积分变量；σx(ξ)为人为规定的衰减函数，
其中，σxm是决定完美匹配层最大衰减量的常数；
对于完美匹配层，将pn的系数记为kmn，kmn为完美匹配层中单元刚度矩阵的元素：
利用上述得到的针对完美匹配层的单元刚度矩阵的元素，建立总体刚度矩阵；再结合总体外载荷向量和声压向量组成线性代数方程组，作为二维声场模型，并求解该二维声场模型得到声场；
Kp＝F
其中，K表示总体刚度矩阵，F表示总体外载荷向量，p表示声场；
再根据得到的声场，通过选择的基函数和插值方法，得到空间各点的声压值。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210054039.6A CN114444353A (zh) | 2022-01-18 | 2022-01-18 | 基于有限元方法的二维声场模型的声场获取方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210054039.6A CN114444353A (zh) | 2022-01-18 | 2022-01-18 | 基于有限元方法的二维声场模型的声场获取方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN114444353A true CN114444353A (zh) | 2022-05-06 |
Family
ID=81368456
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210054039.6A Pending CN114444353A (zh) | 2022-01-18 | 2022-01-18 | 基于有限元方法的二维声场模型的声场获取方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114444353A (zh) |
-
2022
- 2022-01-18 CN CN202210054039.6A patent/CN114444353A/zh active Pending
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110135052B (zh) | 浅海信道下弹性结构辐射声场的计算方法 | |
| CN107576388A (zh) | 一种浅海信道下三维结构声源辐射声场预报方法 | |
| CN113259034B (zh) | 一种并行的耦合海洋声学预报系统及运行方法 | |
| CN112257311B (zh) | Pekeris波导下结构声振计算的FE/WSM方法 | |
| CN109341845A (zh) | 一种海洋环境稳态声场空间实时仿真的方法及装置 | |
| CN114779170A (zh) | 一种浅海近场声源定位方法 | |
| CN112577592A (zh) | 基于空间傅里叶变换的有限空间平面近场声全息测量方法 | |
| CN114286279A (zh) | 基于全局矩阵耦合简正波的三维声场模型的声场获取方法 | |
| CN116879952B (zh) | 点源弹性波海底反射系数的计算方法、存储介质和设备 | |
| CN108875234B (zh) | 应用于船舶三维声弹性分析的浅海声传播计算方法 | |
| Tu et al. | Application of a spectral method to simulate quasi-three-dimensional underwater acoustic fields | |
| CN113253248A (zh) | 一种基于迁移学习的小样本垂直阵目标距离估计方法 | |
| Brick et al. | Fast direct solution of 3-D scattering problems via nonuniform grid-based matrix compression | |
| CN114444353A (zh) | 基于有限元方法的二维声场模型的声场获取方法 | |
| CN114841892A (zh) | 一种基于全连接网络的稀疏导波数据恢复方法 | |
| Wu et al. | Dispersion error control for underwater acoustic scattering problems using a coupled cell-based smoothed radial point interpolation method | |
| Hospital-Bravo et al. | Numerical modeling of undersea acoustics using a partition of unity method with plane waves enrichment | |
| Martinelli | An application of the level set method to underwater acoustic propagation | |
| Duan et al. | PINNs for Sound Propagation and Sound Speed Field Estimation Simultaneously | |
| Gao et al. | An integrated calculation model for vibroacoustic radiation and propagation from underwater structures in oceanic acoustic environments | |
| CN116522699B (zh) | 使用波叠加-有限元-pml联合计算浅海下体目标辐射噪声的方法 | |
| Meng et al. | Adaptive Interference Simulation of Environment Noise | |
| Sturm et al. | Accurate treatment of a general sloping interface in a finite-element 3D narrow-angle PE model | |
| Deo et al. | Predicting Wave Propagation for Varying Bathymetry Using Conditional Convolutional Autoencoder Network | |
| Jiawen et al. | Solving differential equations with neural networks: application to the normal-mode equation of sound field under the condition of ideal shallow water waveguide |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination |