CN116522699A - 使用波叠加-有限元-pml联合计算浅海下体目标辐射噪声的方法 - Google Patents

使用波叠加-有限元-pml联合计算浅海下体目标辐射噪声的方法 Download PDF

Info

Publication number
CN116522699A
CN116522699A CN202211106264.6A CN202211106264A CN116522699A CN 116522699 A CN116522699 A CN 116522699A CN 202211106264 A CN202211106264 A CN 202211106264A CN 116522699 A CN116522699 A CN 116522699A
Authority
CN
China
Prior art keywords
sound
source
field
sound source
equivalent
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
Application number
CN202211106264.6A
Other languages
English (en)
Inventor
唐宇航
李海超
刘宝
王雪仁
高晟耀
张�浩
彭子龙
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
People's Liberation Army 92578
Original Assignee
People's Liberation Army 92578
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by People's Liberation Army 92578 filed Critical People's Liberation Army 92578
Priority to CN202211106264.6A priority Critical patent/CN116522699A/zh
Publication of CN116522699A publication Critical patent/CN116522699A/zh
Pending legal-status Critical Current

Links

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
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/04—Constraint-based CAD
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00—Details relating to the application field
    • G06F2113/08—Fluids
    • 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/10—Noise analysis or noise optimisation
    • 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
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Computing Systems (AREA)
  • Fluid Mechanics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)

Abstract

本发明提供一种使用波叠加‑有限元‑PML联合计算浅海下体目标辐射噪声的方法，选取典型的体积声源—球形声源作为研究对象，将球形声源的辐射声场等效为其对称轴上若干等效点源辐射声场的线性叠加，采用有限元法建立球形声源、浅海环境、边界条件的联合动力学模型，计算相应源强下每个等效点源诱发的声场，将其线性叠加快速求得球形声源诱发的远场声辐射。本发明用于计算的海水中近、远场处的频域声压级相比有限元发在计算精度和效率上都具有显著优势。

Description

使用波叠加-有限元-PML联合计算浅海下体目标辐射噪声的 方法
技术领域
本发明涉及浅海环境声辐射的计算方法，具体涉及一种使用波叠加-有限元-PML联合计算浅海下体目标辐射噪声的方法，属于数据处理技术领域。
背景技术
波叠加法在研究声学中体积声源辐射声场方面有很多例子，类似球形声源，圆柱声源的辐射问题，通过波叠加法求解得到的球形声源及圆柱声源的辐射声压和解析解吻合的很好；但是将波叠加法运用于海洋声学中的研究相对较少，在浅海地震波场领域更是鲜有研究。
发明内容
发明目的：针对浅海中体积声源诱发的远场声辐射频域预报问题，文中提出一种基于波叠加-有限元联合法的快速预报方法。选取典型的体积声源—球形声源作为研究对象，将球形声源的辐射声场等效为其对称轴上若干等效点源辐射声场的线性叠加，采用有限元法建立球形声源、浅海环境、边界条件的联合动力学模型，计算相应源强下每个等效点源诱发的声场，将其线性叠加快速求得球形声源诱发的远场声辐射。
技术方案：为解决上述技术问题，本发明提供的使用波叠加-有限元-PML联合计算浅海下体目标辐射噪声的方法，包括以下步骤：
步骤1，基于波叠加原理，将球形声源的辐射声场等效为其对称轴上若干等效点源辐射声场的线性叠加；
步骤2，根据与球形声源表面上等效点源数量一致的场点振速匹配原则，确定等效点源的强度；
步骤3，采用有限元法建立球形声源、浅海环境、边界条件的联合动力学模型，计算相应源强下每个等效点源诱发的声压；
步骤4，将步骤3计算结果线性叠加，快速求得球形声源诱发的声辐射。
有益效果：从计算精度和效率两个方面，比较本发明方法和有限元法，分别计算海水中近、远场处的频域声压级，本发明方法具有更高的计算精度和效率。
附图说明
图1是球形声源辐射声场的波叠加法求解示意；
图2是海面为压力释放边界时的等效源链；
图3是PML建模原理；
图4是有限元模型的建立；
图5是有限元模型网格的划分；
图6是球形声源处网格的局部细化；
图7是两种方法求解的硬海底时脉动球和摆动球激发的海水中近、远场频域声压级对比。
具体实施方式
实施例1：
1.波叠加法求解
波叠加法是在浅海波导中，将球形声源的辐射声场等效为内部若干点源辐射声场的线性叠加，如图1所示，r轴、z轴分别表示水平距离向和垂直向，R为球形声源的半径，将个数为n的场点设置在球形声源外圆周上，以空心点表示，个数为N的等效源点设置在球形声源对称轴上，如实心点所示，为了区分场点和源点，令第l个场点和源点的标号分别为rl与Nl，θl为第l个场点矢量与r轴的夹角，球心的垂直坐标为z0。
球形声源上场点r处声场p(r)看做是其对称轴上所有等效源辐射声场的线性叠加，其表达式为：
式中：i为模为-1的虚数，ρ为球形声源所在介质中的密度，c为球形声源所在介质中的声速，k＝2πf/c为波数，f为声波频率，r为对应场点到坐标原点O的距离，rl为第l个等效源到坐标原点O的距离，s(rl)为等效声源分布在处的等效源强度，G(r,rl)为声场Green函数，求解浅海海底地震波场，选取海底为弹性介质的浅海波导中的Green函数，在计算G(r,rl)时，虚源通常只取前四阶即可，其表达式为:
式中的R01、R02、R03、R04的计算式分别为：
式(2.1.2)和式(2.1.3)中，r'为场点到源点的水平距离，m为镜像反射的次数，γ1为海面反射系数，通常考虑海面为压力释放边界，故γ1通常为-1；γ为海底反射系数，h为海深，h1，h2为接收点到海面和海底的距离，如图2所示。
式(2.1.3)中的海底反射系数γ通常由海底介质的属性决定，在海底为液态海底的条件下，海底反射系数γ为：
其中，θ为入射角，α＝ρ1/ρ0，β＝c0/c1。海水与海底介质的速度分别为c0与c1，海水和海底介质的密度分别为ρ0、ρ1。由此，可以得到弹性海底下浅海声场的Green函数G(r,rl)。由声压和振速的关系，可以得到表面振速的表达式为：
在球形声源表面振速已知的情况下，通过式(2.1.5)得到对应的等效源强度s(rl)。
2.有限元求解
在求解海底地震波场时，由于没有类似无限大声场中球形声源辐射声压的解析解，故在对其仿真得到的方法进行验证时，需要与计算精度较高以及技术较为成熟的有限元法得到的结果进行对比，下面，对有限元方法的基本原理进行介绍：
2.1轴对称有限元
在二维轴对称坐标系中，声源位于海水中时，声学单元控制方程即Helmholtz方程可以表示为如下的形式，即
式中，p表示单元上的声压。zs为等效声源到海底的距离。对声压采用与结构相同的形函数进行插值处理，则每个单元上的声压可以表示为
p＝Nfp (2.2.2)式中，Nf表示由形函数N(ξ,η)构成的矩阵，表示为
Nf＝[N1(ξ,η) N2(ξ,η) … Nm(ξ,η)]。
p为单元节点声压矢量，p＝[p1 p2 … pm]T。
类似于结构单元的推导，根据Hamilton原理可知，声学轴对称单元对应的刚度矩阵和质量矩阵可以分别表示为：
式中，
矩阵与/>均采用Guass积分进行计算。进一步将各轴对称结构单元根据各个离散节点进行组合，即可获得海水相应矩阵。最后可得声学有限元方程为：
(Kf+jωCf-ω2Mf)p＝{Fi} (2.2.2)
式中：{Fi}为声学激励向量，Ma、Ka、Ca分别为声学系统的质量矩阵、刚度矩阵以及阻尼矩阵。
2.2无限域空间的PML截断
解决有限元方法中截断无限域问题的一种常用方法是采用完全匹配层(PML)。PML通常为一个人工吸收层，用于截断波动方程中的计算域，通过附加PML到计算领域，PML可以吸收波传播的内部不反射波回到域有界介质的界面和PML，这是一个属性，区别于PML与其他吸收边界条件。虽然PML可以建立在任何坐标系中，如图3所示，但本节研究的重点是在笛卡尔坐标系中的PML的发展。在图中，红色的图域是原始的计算图域，外部的附加图层是PML。
PML的概念最初是由贝伦格通过时域有限差分引入电磁波建模的。大多数关于PML的应用都与有限差分时域(FDTD)有关，并且在频域中只有很少的使用有限元的应用方法。虽然完全匹配的层(PML)被确定适用于许多线性系统和弹性动力问题，弹性波的首次实际应用是通过使用速度应力公式。这些公式随后被扩展到各向异性介质的情况。
经典PML的控制方程，包括将位移场u分成四部分(u＝u1+u2+u3+u4)。
式中，和/>I是3×3的恒等式张量。阻尼函数d和当前位置与PML域和原始计算域界面的距离有关，阻尼函数d与当前位置到PML域和原始计算域界面的距离有关：
其中，R是理论上的反射系数，cp是P波速度，lp是PML的厚度。
3.球形声源诱发的远场声辐射建模
基于有限元法建立水下体积声源诱发远场声辐射的有限元模型，首先通过有限元软件COMSOL Multiphysics中的声学模块建立有限元模型，如图4所示。将模型划分为海水和海底两个部分。其中，海水和海底均采用用压力声学单元代替，这类单元能够求解海水中声压及声压级等物理量，能很好地模拟海洋环境；由于需要验证的为浅海远场声辐射中的声传播损失，在构建海水-海底分界面时，需对其边界条件进行进一步约束，海水-海底分界面处的位移以及应力应满足：(1)法向应力连续；(2)法向位移连续；为了模拟无限域问题，采用2.2节中的完美匹配层(PML)来对无限域进行截断，其厚度设置为对应区域最小波长的6倍。
以球形声源为例，将球的半径设置为R＝2.5m，模型水平距离设置为2000m，海水深度设置为100m，海底深度设置为200m，将球形声源设置在深度为50m的位置，如图4所示。由于球形声源的高度轴对称性，在建模时采用轴对称模型将三维的问题简化为二维的问题以减小计算量。
为了更好地验证此方法的正确性，在对球形声源施加激励时分别采用脉动和摆动的形式。其中，脉动球表面法向振速为v0，摆动球表面法向振速为vi，v0＝1m/s，vi＝v0 sinθi，其中，sinθi＝|z-z0|/R，z0为球心的z方向坐标，z为对应圆周上场点的坐标，如图1所示；海底底质参数选取表2.1中粉砂以及泥质岩。
表2.1介质类型及参数
Tab.2.1 Media types and parameters
在建立的有限元模型上，根据2.1节中的内容，通过波叠加法将等效源点加载在图1中球形声源的母线上。假设场点的个数与波叠加等效源个数一致均为Ni，故场点坐标(ri,zi)如式(2.3.1)所示：
划分网格为有限元方法前处理中一个重要环节，网格类型的选取以及网格数量的设置往往会影响计算的时间以及计算结果的精度。在划分网格的过程中，为提高计算精度，通常选取为规则网格，例如二维网格中的三角形网格，四边形网格，三维网格中的四面体网格以及六面体网格；对于网格数量较多的复杂模型，可以考虑将需要求解的区域进行局部细化，其余部分粗化的方式来提高计算效率；在对网格大小进行选取时，需满足其大小要不大于该区域中所有波动成分最小波长的1/6。
根据上述原则，在进行球形声源及波场的网格划分时，为了保证模型本身的精确性以及两者对比结果的准确性，将图1中球形声源边界上的网格数量均设置为Ni，以保证波叠加法中源点和场点数量的一致以及两者对比下网格质量的一致；海水-海底分解面上的网格数量设置为500，以确保在绘制传播损失时图像的高对比度，在选取网格时，选取二维下精度较高的三角形网格，并将球形声源处的网格局部细化，结果如图5、图6所示。
4.有效性验证
在完成建模后，对其有效性进行验证。故将分析频率设置为20-100Hz，间隔为3Hz，选取表1中的两种海底介质验证联合法建模正确性。
为了更好地验证联合法对远场声压的快速预报效果，从海水中声场、海底两个不同的方面进行验证，从计算精度和效率两个方面对结果进行分析。
在对声学有效性的验证中，通常以声压、声压级为出发点进行对比验证，为了更好地分析联合法的精度，从两个不同方面进行声学有效性的分析，分别为：海水中近场及远场长点处的频域声压级，其中，近场场地选取为距离声源中心5m处，远场则是选取在1500m处，分析频率选为20-50Hz，间隔为1Hz。为了定量对比，将声压级定义为：
式中：p1为对应场点处的声压，参考声压pref＝10-6Pa。
将等效源设置为单极子源，等效源点和场点数量为n＝N＝20，选取式(2.1.2)中的Green函数进行计算，利用式(2.1.6)计算得到不同频率时等效源对应的强度，并采用有限元法计算对应的频率及源强下每个等效源点所对应的近场与远场处的声压级；同时，采用有限元法计算添加表面振速时球形声源上对应处的声压级，将两种方法的结果进行对比分析，具体结果如图7所示。图中，实线表示有限元法求解得到的海水中近、远场场点处频域声压级，“*”代表联合法计算得到的海水中近、远场场点处频域声压级。
从图7中可以看出，波叠加-有限元-PML联合法对海水中近场及远场场点处频域声压级的拟合结果与有限元法得到的结果数值及变化趋势基本一致，
为了更精确地衡量两者之间的误差，利用统计学中的相对误差与相关系数对每条曲线的拟合程度进行评估，得到结果见表2.2。
表2.2不同仿真条件下两种方法近场频域声压级计算结果的相关系数及相对误差
Tab.2.2 Correlation coefficients and relative errors of near-fieldfrequency-domain sound pressure levels calculated by the two methods underdifferent simulation conditions
从表2.2中的数据可以看出，联合法对弹性海底下浅海波导的声压传播损失的拟合精度较高，相对误差保持在0.01％左右，相关系数均接近1，证明两者具有很高的计算精度。
实施例2
高效性验证：继续比较两者计算速度差异，在相同硬件环境下对两种方法的计算速度进行测试，测试采用实施例1有效性验证的数值模型，其余仿真条件均一致，取计算平均值进行对比，结果如表2.6所示。
表2.6两种方法平均运行时间对比结果
Tab.2.6 Compariswiths of the average running time of two methods
从表2.6中数据可以看到，对于相同仿真条件的有限元方法，本方法可在保证计算精度时极大地缩短计算时间，提高运算效率。

Claims (3)

1.一种使用波叠加-有限元-PML联合计算浅海下体目标辐射噪声的方法，其特征在于包括以下步骤：
步骤1，基于波叠加原理，将球形声源的辐射声场等效为其对称轴上若干等效点源辐射声场的线性叠加；
步骤2，根据与球形声源表面上等效点源数量一致的场点振速匹配原则，确定等效点源的强度；
步骤3，采用有限元法建立球形声源、浅海环境、边界条件的联合动力学模型，计算相应源强下每个等效点源诱发的声压；
步骤4，将步骤3计算结果线性叠加，快速求得球形声源诱发的声辐射。
2.根据权利要求1所述的使用波叠加-有限元-PML联合计算浅海下体目标辐射噪声的方法，其特征在于：所述步骤1中，球形声源上场点r处辐射声场p(r)：
式中：i为模为-1的虚数，ρ为球形声源所在介质中的密度，c为球形声源所在介质中的声速，k＝2πf/c为波数，f为声波频率，r为对应场点到坐标原点O的距离，rl为第l个等效源到坐标原点O的距离，s(rl)为等效声源分布在处的等效源强度，G(r,rl)为声场Green函数。
3.根据权利要求2所述的使用波叠加-有限元-PML联合计算浅海下体目标辐射噪声的方法，其特征在于：所述步骤2中，等效源强度s(rl)计算公式为：
CN202211106264.6A 2022-09-10 2022-09-10 使用波叠加-有限元-pml联合计算浅海下体目标辐射噪声的方法 Pending CN116522699A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211106264.6A CN116522699A (zh) 2022-09-10 2022-09-10 使用波叠加-有限元-pml联合计算浅海下体目标辐射噪声的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211106264.6A CN116522699A (zh) 2022-09-10 2022-09-10 使用波叠加-有限元-pml联合计算浅海下体目标辐射噪声的方法

Publications (1)

Publication Number Publication Date
CN116522699A true CN116522699A (zh) 2023-08-01

Family

ID=87394587

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211106264.6A Pending CN116522699A (zh) 2022-09-10 2022-09-10 使用波叠加-有限元-pml联合计算浅海下体目标辐射噪声的方法

Country Status (1)

Country Link
CN (1) CN116522699A (zh)

Similar Documents

Publication Publication Date Title
CN110135052B (zh) 浅海信道下弹性结构辐射声场的计算方法
CN107576388A (zh) 一种浅海信道下三维结构声源辐射声场预报方法
CN113259034B (zh) 一种并行的耦合海洋声学预报系统及运行方法
CN110750934A (zh) 深海弹性结构与环境耦合声辐射预报方法
CN112257311B (zh) Pekeris波导下结构声振计算的FE/WSM方法
CN110399680A (zh) 一种浅海弹性结构辐射声场计算方法
CN115630446B (zh) 用于潜航器结构低频辐射噪声实时仿真的快速推演方法
Petris et al. A numerical method for the solution of the three-dimensional acoustic wave equation in a marine environment considering complex sources
CN114779170A (zh) 一种浅海近场声源定位方法
CN113330439A (zh) 浅海多层海底地声参数反演方法、装置、计算机设备及存储介质
Huang et al. Study on the integrated calculation method of fluid–structure interaction vibration, acoustic radiation, and propagation from an elastic spherical shell in ocean acoustic environments
Jiang et al. Combining dynamic XFEM with machine learning for detection of multiple flaws
CN116879952B (zh) 点源弹性波海底反射系数的计算方法、存储介质和设备
CN114286279A (zh) 基于全局矩阵耦合简正波的三维声场模型的声场获取方法
Xu et al. An efficient ray-tracing method and its application to Gaussian beam migration in complex multilayered anisotropic media
CN113758559A (zh) 基于散射声场分离算法在浅海信道中测量方法
Sessarego et al. Acoustic scattering by an elastic spherical shell near the seabed
CN112926231A (zh) 一种基于等效源法的有限空间中近场声全息测量方法
Schroeder On the interaction of elastic waves with buried land mines: an investigation using the finite-difference time-domain method
CN116522699A (zh) 使用波叠加-有限元-pml联合计算浅海下体目标辐射噪声的方法
Wu et al. A meshfree radial point interpolation coupled with infinite acoustic wave envelope element method for computing acoustic fields
Wu et al. A novel hybrid superposition method for predicting ship seismic wave in shallow sea
Li et al. A 3D reflection ray-tracing method based on linear traveltime perturbation interpolation
Huang et al. Study on calculation methods for acoustic radiation of axisymmetric structures in finite water depth
CN114252148A (zh) 一种基于长椭球波叠加的声场重建方法

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