CN117935837A - 一种时域多声源定位及噪声处理方法 - Google Patents

一种时域多声源定位及噪声处理方法 Download PDF

Info

Publication number
CN117935837A
CN117935837A CN202410338625.2A CN202410338625A CN117935837A CN 117935837 A CN117935837 A CN 117935837A CN 202410338625 A CN202410338625 A CN 202410338625A CN 117935837 A CN117935837 A CN 117935837A
Authority
CN
China
Prior art keywords
time
grid
boundary
derivative
point
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.)
Granted
Application number
CN202410338625.2A
Other languages
English (en)
Other versions
CN117935837B (zh
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.)
Computational Aerodynamics Institute of China Aerodynamics Research and Development Center
Original Assignee
Computational Aerodynamics Institute of China Aerodynamics Research and Development Center
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 Computational Aerodynamics Institute of China Aerodynamics Research and Development Center filed Critical Computational Aerodynamics Institute of China Aerodynamics Research and Development Center
Priority to CN202410338625.2A priority Critical patent/CN117935837B/zh
Publication of CN117935837A publication Critical patent/CN117935837A/zh
Application granted granted Critical
Publication of CN117935837B publication Critical patent/CN117935837B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G10—MUSICAL INSTRUMENTS; ACOUSTICS
    • G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
    • G10L21/00—Speech or voice signal processing techniques to produce another audible or non-audible signal, e.g. visual or tactile, in order to modify its quality or its intelligibility
    • G10L21/02—Speech enhancement, e.g. noise reduction or echo cancellation
    • G10L21/0208—Noise filtering
    • G10L21/0216—Noise filtering characterised by the method used for estimating noise
    • G10L21/0224—Processing in the time domain
    • G—PHYSICS
    • G10—MUSICAL INSTRUMENTS; ACOUSTICS
    • G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
    • G10L21/00—Speech or voice signal processing techniques to produce another audible or non-audible signal, e.g. visual or tactile, in order to modify its quality or its intelligibility
    • G10L21/02—Speech enhancement, e.g. noise reduction or echo cancellation
    • G10L21/0208—Noise filtering
    • G10L21/0216—Noise filtering characterised by the method used for estimating noise
    • G10L21/0232—Processing in the frequency domain
    • G—PHYSICS
    • G10—MUSICAL INSTRUMENTS; ACOUSTICS
    • G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
    • G10L25/00—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00
    • G10L25/27—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 characterised by the analysis technique
    • G—PHYSICS
    • G10—MUSICAL INSTRUMENTS; ACOUSTICS
    • G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
    • G10L25/00—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00
    • G10L25/48—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 specially adapted for particular use
    • G10L25/51—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 specially adapted for particular use for comparison or discrimination
    • 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)
  • Computational Linguistics (AREA)
  • Signal Processing (AREA)
  • Health & Medical Sciences (AREA)
  • Audiology, Speech & Language Pathology (AREA)
  • Human Computer Interaction (AREA)
  • Physics & Mathematics (AREA)
  • Acoustics & Sound (AREA)
  • Multimedia (AREA)
  • Quality & Reliability (AREA)
  • Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)

Abstract

本发明公开了一种时域多声源定位及噪声处理方法，涉及声源定位领域，本发明基于守恒型Euler方程，采用基于复变量微分的敏感度分析方法与时间反演相结合的实施中，能有效建立远场声信号关键频率信息与近场多声源之间的联系，分析远场声信号频率信息与近场复杂流动中的多声源的对应关系，梳理需要降噪的声信号中关键频率和关键声源，针对性的开展降噪策略研究。

Description

一种时域多声源定位及噪声处理方法
技术领域
本发明涉及声源定位领域，具体地，涉及一种时域多声源定位方法及噪声处理方法。
背景技术
随着民用、军用飞行器等航空运输的快速发展，气动噪声问题日益引起世界的广泛关注，因此亟需发展高效降噪方法。声源定位是解决噪声问题的第一步，只有精准定位声源位置，才能针对性开展噪声降噪方案研究。但是实际工程问题中大多存在多个声源且声源强度不一，这些声源可能来自不同的流动结构或物理过程，其可能辐射的多种频率的声信号，还会相互干扰，难以明确远场声信号频率信息与多个声源之间的关系，因此就无法确定哪一个是主要声源，哪一个是次要声源，导致降噪方案无法达到预期效果。
近年来出现了一种极具应用前景的声源定位方法—基于时间反演的声源定位方法（A. Deneuve, P. Druault, R. Marchiano, P. Sagaut, A coupled time-reversalcomplex differentiation method for aeroacoustic sensitivity analysis: towardsa source detection procedure, Journal of Fluid Mechanics, 2010, 642: 181–212.），通过采用伪特征形式的Euler方程（Sesterhenn J. A characteristic-typeformulation of the Navier–Stokes equations for high order upwind schemes[J].Computers&Fluids, 2001, 30(1):37-67.）进行数值求解反问题，获得声源位置信息和声源强度信息，随后采用敏感度分析将远场声信号频率信息与近场声源相关联。
然而，伪特征形式方程是一种特殊方程，针对特定问题需要进行特殊处理，难以直接应用于现有的计算流体力学体系，不具备普适性，难以推广。现有的时域声源定位方法为采用反拉克斯-温德洛夫方法与时间反演相结合，成功将基于时间反演声源定位方法从特殊方程推广到一般形式的守恒型Euler方程并提高了声源定位的精度，可以达到精准定位近场复杂流动中的噪声源，但是仍然无法解决远场程声信号频率信息与近场多声源之间的联系。此时如果盲目对近场复杂流动中的所有噪声源进行降噪处理，第一是性价比不高，实际操作难度大；第二是声源之间会相互影响，有可能某一个声源噪声确实降低，但引起其它声源噪声增大，得不偿失；第三是频率信号通常具有主频和对应得谐频信号，一般来说是由同一个声源发出，但是很多时候是多个声源发声信号之间得耦合效应导致。应此需要将远场声信号频率与近场复杂流动中多声源相关联，再进一步进行原因分析，才能找到合适的降噪方案，达到事半功倍的效果。
发明内容
本发明提供了一种时域多声源定位及噪声处理方法，本发明在时域声源定位方法的基础上，将基于复变量微分的瞬时敏感度分析方法与守恒型Euler方程、时间反演方法、反拉克斯-温德洛夫方法进一步结合，建立了一种时域多声源定位及噪声处理方法，本方法将基于复变量的敏感度分析方法应用于守恒型Euler方程；采用复变量微分方法进行敏感度分析，建立远场声信号中的频率信息与近场复杂流场中的多个声源之间的对应关系，该噪声处理方法能够从众多的噪声源中准确的找到特定频率噪声对应的一个或多个声源，进而能够准确的对某些噪声源进行准确的降噪处理。
该噪声处理方法能够从众多的噪声源中准确的找到特定频率噪声的声源，进而能够准确的对某个噪声源进行准确的降噪处理。
为实现上述发明目的，本发明提供了一种时域多声源定位及噪声处理方法，所述方法包括：
步骤1：进行数值模拟获取数据，包括：分析气动噪声问题并建立气动噪声问题对应的物理模型；在网格生成软件中建立物理模型对应的几何模型；将几何模型外表面和整个计算区域进行空间网格划分，生成计算所需要的网格，设定边界条件并导出网格文件；将网格文件导入流体仿真软件，设定初始条件进行数值模拟计算，当计算达到稳定后的时刻记为时刻t0；计算n个周期后的时刻记为时刻t1，n为整数；在时刻t0至时刻t1之间，存储网格边界处每个边界点与时刻t相关的物理量数据和时刻t1整个网格每个点的物理量数据；
步骤2：对压力信号进行处理获得远场声信号频率分布和频率强度，包括：通过傅里叶变换将时域空间内网格边界上每个网格点处的压力信号p序列转换为频域空间压力信号，获得频域空间压力信号频率分布，将频率按照幅值从大到小排序标记为fr，r=1,2,3,…，m，m为频率个数；
步骤3：对频率为f1的频域空间压力信号对应的噪声源进行定位，包括：
步骤3.1：基于步骤1中存储的数据进行变换获得时间反演所需的数据，包括：基于时刻t和时刻t1获得时间反演的时刻ttr，t大于或等于t0且小于或等于t1，基于时刻t对应的网格边界处每个边界点物理量数据，获得时刻ttr对应的网格边界处每个边界点物理量数据；
步骤3.2：对边界上的物理量数据进行处理包括：
步骤3.2.1：获得边界点的时间导数，包括：基于所有时刻的所有网格边界处每个边界点物理量数据计算获得网格边界处每个网格点对应守恒变量U在时刻ttr的时间导数；
步骤3.2.2：获得边界点的空间切向导数，包括：基于时刻ttr对应的网格边界处每个边界点物理量数据，分别计算获得时刻ttr网格边界处每个网格点在x方向无粘通量的空间切向导数和在y方向无粘通量/>的空间切向导数；
步骤3.2.3：获得边界点的空间法向一阶导数，包括：基于网格边界处每个网格点对应守恒变量U在时刻ttr的时间导数、时刻ttr网格边界处每个网格点在x方向无粘通量的空间切向导数和时刻ttr网格边界处每个网格点在y方向无粘通量的空间切向导数，计算获得网格边界处每个网格点对应守恒变量U在时刻ttr的空间法向一阶导数；
步骤3.2.4：构造多项式，包括：基于守恒变量U在网格预设位置对应的数值和守恒变量U在网格预设位置的空间法向一阶导数构造多项式；
步骤3.2.5：令虚拟点的取值为所述多项式在相应点处的取值，对频率为f1对应的压力信号成分进行傅里叶逆变换获得时域空间的单频压力信号，将单频压力信号/>以虚部的形式添加到压力信号p上获得复数形式的新压力信号/>，/>，Im为虚元；基于虚拟点的数值、网格内部点的数值和新压力信号/>进行计算获得下一时刻tk网格内部点的物理量数据；
步骤3.3：将时刻t替换为时刻tk，判断时刻tk是否小于或等于t1，若是则返回执行步骤3.1，获得时刻t0至时刻t1对应的所有网格点的物理量数据，若否则执行步骤3.4；
步骤3.4：获得频率为f1的频域空间压力信号对应的一个或多个声源位置，包括：当tk=t1 时，基于所有网格点的复数形式的新压力信号的虚部构建压强云图，基于压强云图中压力信号最大处的位置获得频率为f1的频域空间压力信号对应的声源位置信息；
步骤4：对频率为f2、f3、…、fm的频域空间压力信号依次进行步骤3所述的噪声源定位处理，获得频率为f2、f3、…、fm的频域空间压力信号对应的声源位置信息；
步骤5：基于频率为f1、f2、f3、…、fm的频域空间压力信号对应的声源位置信息分析声源产生原因并进行分类，针对不同类型的声源采用相应的降噪方案进行降噪处理。
其中，时间反演的初始条件需要网格边界处的时间导数或时序数据，但对边界处理非常敏感，容易从外部引入伪波导致声源定位分辨率不高。反Lax-Wendroff方法即反拉克斯-温德洛夫方法是一种高精度、稳定的边界处理方法，需要给定边界条件，多适用于固壁边界、动边界情况。本发明将时间反演和反Lax-Wendroff方法进行结合，实现稳定高精度高分辨率的声源定位。
其中，有多个声源时，本方法能够基于频率从压强云图中获得对应的声源位置信息。每个声源对目标的影响效果不同，本方法可以确定频率对应的声源位置信息；然后基于确定的频率对应的声源位置信息，对该位置处的声源进行降噪处理，即能够从多个声源中确定相应频率的声源，然后对不同频率的噪声进行相应的降噪处理。
其中，在一些实施例中，网格边界处每个边界点与时刻t相关的物理量数据包括：时刻t对应的密度、时刻t对应的x方向速度/>、时刻t对应的y方向速度/>和时刻t对应的压强/>，x表示横向空间坐标，y表示纵向空间坐标；
时刻ttr对应的物理量数据包括：时刻ttr对应的密度、时刻ttr对应的x方向速度/>、时刻ttr对应的y方向速度/>和时刻ttr对应的压强，其中：
。
其中，在一些实施例中，网格边界处每个网格点对应守恒变量U的时间导数包括：、/>、/>和/>，其中，网格边界处每个网格点对应守恒变量U的时间导数的计算方式为：
；
其中，U1是守恒变量U的第一分量，U2是守恒变量U的第二分量，U3是守恒变量U的第三分量，U4是守恒变量U的第四分量，是比热比。
其中，在一些实施例中，网格边界处每个网格点在x方向无粘通量的空间切向导数包括：/>，/>，/>和/>；
网格边界处每个网格点在y方向无粘通量的空间切向导数包括：/>，/>，和/>；
其中：
；
其中，是y方向无粘通量/>的第一分量，/>是y方向无粘通量/>的第二分量，是y方向无粘通量/>的第三分量，/>是y方向无粘通量/>的第四分量；
；
其中，是x方向无粘通量/>的第一分量，/>是x方向无粘通量/>的第二分量，是x方向无粘通量/>的第三分量，/>是x方向无粘通量/>的第四分量。
其中，在一些实施例中，所述步骤3.2.3具体包括：基于网格边界处每个网格点对应守恒变量U的时间导数、网格边界处每个网格点在x方向无粘通量的空间切向导数和网格边界处每个网格点在y方向无粘通量/>的空间切向导数，采用反拉克斯-温德洛夫方法计算获得网格边界处每个网格点对应守恒变量U的空间法向一阶导数。
其中，在一些实施例中，所述步骤3.2.3具体包括：
求解网格边界处每个网格点上的特征值，获得四个特征方向的特征值；
构建左特征值矩阵；
基于四个特征方向的特征值对应的特征量获得出口边界信息和入流边界信息；
基于左特征值矩阵、出口边界信息和入流边界信息构建空间法向一阶导数求解矩阵；
求解所述空间法向一阶导数求解矩阵获得网格边界处每个网格点对应守恒变量U的空间法向一阶导数。
其中，在一些实施例中，四个特征方向的特征值为：，i=1或2或3或4，其中，，/>，/>，/>，/>表示该网格点处x方向速度，/>表示该网格点的声速；
左特征值矩阵为：
；
其中，是左特征值矩阵/>的第一分量，/>是左特征值矩阵/>的第二分量，/>是左特征值矩阵/>的第三分量，/>是左特征值矩阵/>的第四分量；为左特征值矩阵/>在坐标/>处的值；
基于四个特征方向的特征值的数值大小关系判断获得出口边界对应的特征量以及入流边界对应的特征向量；
出口边界对应的特征量采用四阶拉格朗日外插获得；
入流边界对应的特征向量采用反拉克斯-温德洛夫方法获得，基于入流边界对应的特征向量获得网格边界处每个网格点对应守恒变量U的空间法向一阶导数。
其中，在一些实施例中，网格预设位置的空间法向一阶导数构造多项式为：
；
其中，是多项式/>在x=x1处的值，/>是多项式/>在x=x2处的值，是多项式/>在x=x3处的值，/>是多项式/>在x=x4处的值，/>是多项式/>在x=x1处的一阶导数值，/>是守恒变量U在x=x1处的值，/>是守恒变量U在x=x2处的值，/>是守恒变量U在x=x3处的值，/>是守恒变量U在x=x4处的值，是守恒变量U在x=x1处的空间法向一阶导数。
其中，在一些实施例中，虚拟点的取值的计算方式为：
；
其中，是多项式/>在x=xj处的值。
其中，在一些实施例中，所述方法还包括：步骤6：获得降噪结果，根据降噪结果判断是否调整当前降噪方案。
本发明提供的一个或多个技术方案，至少具有如下技术效果或优点：
由于本方法采用复变量微分方法进行瞬时敏感度分析，建立远场声信号中的频率信息与近场复杂流场中的多个声源之间的对应关系，可以明确远场声信号中的关键频率对应与近场复杂流动中的一个或者多个噪声源，然后针对该噪声源发展降噪策略研究。
附图说明
此处所说明的附图用来提供对本发明实施例的进一步理解，构成本发明的一部分，并不构成对本发明实施例的限定；
图1是本发明中时域多声源定位及噪声处理方法的整体流程示意图。
具体实施方式
为了能够更清楚地理解本发明的上述目的、特征和优点，下面结合附图和具体实施方式对本发明进行进一步的详细描述。需要说明的是，在相互不冲突的情况下，本发明的实施例及实施例中的特征可以相互组合。
在下面的描述中阐述了很多具体细节以便于充分理解本发明，但是，本发明还可以采用其他不同于在此描述范围内的其他方式来实施，因此，本发明的保护范围并不受下面公开的具体实施例的限制。
实施例一；
请参考图1，图1为本发明中时域多声源定位及噪声处理方法的整体流程示意图，所述方法包括：
步骤1：进行数值模拟获取数据，包括：分析气动噪声问题并建立气动噪声问题对应的物理模型；在网格生成软件中建立物理模型对应的几何模型；将几何模型外表面和整个计算区域进行空间网格划分，生成计算所需要的网格，设定边界条件并导出网格文件；将网格文件导入流体仿真软件，设定初始条件进行数值模拟计算，当计算达到稳定后的时刻记为时刻t0；计算n个周期后的时刻记为时刻t1，n为整数；在时刻t0至时刻t1之间，存储网格边界处每个边界点与时刻t相关的物理量数据和时刻t1整个网格每个点的物理量数据；其中，网格生成软件可以为NNW-GridStar或Pointwise或Gridgen；
步骤2：对压力信号进行处理获得远场声信号频率分布和频率强度，包括：通过傅里叶变换将时域空间内网格边界上每个网格点处的压力信号p序列转换为频域空间压力信号，获得频域空间压力信号频率分布，将频率按照幅值从大到小排序标记为fr，r=1,2,3,…，m，m为频率个数；
步骤3：对频率为f1的频域空间压力信号对应的噪声源进行定位，包括：
步骤3.1：基于步骤1中存储的数据进行变换获得时间反演所需的数据，包括：基于时刻t和时刻t1获得时间反演的时刻ttr，t大于或等于t0且小于或等于t1，基于时刻t对应的网格边界处每个边界点物理量数据，获得时刻ttr对应的网格边界处每个边界点物理量数据；
步骤3.2：对边界上的物理量数据进行处理包括：
步骤3.2.1：获得边界点的时间导数，包括：基于所有时刻的所有网格边界处每个边界点物理量数据计算获得网格边界处每个网格点对应守恒变量U在时刻ttr的时间导数；
步骤3.2.2：获得边界点的空间切向导数，包括：基于时刻ttr对应的网格边界处每个边界点物理量数据，分别计算获得时刻ttr网格边界处每个网格点在x方向无粘通量的空间切向导数和在y方向无粘通量/>的空间切向导数；
步骤3.2.3：获得边界点的空间法向一阶导数，包括：基于网格边界处每个网格点对应守恒变量U在时刻ttr的时间导数、时刻ttr网格边界处每个网格点在x方向无粘通量的空间切向导数和时刻ttr网格边界处每个网格点在y方向无粘通量的空间切向导数，计算获得网格边界处每个网格点对应守恒变量U在时刻ttr的空间法向一阶导数；
步骤3.2.4：构造多项式，包括：基于守恒变量U在网格预设位置对应的数值和守恒变量U在网格预设位置的空间法向一阶导数构造多项式；
步骤3.2.5：令虚拟点的取值为所述多项式在相应点处的取值，对频率为f1对应的压力信号成分进行傅里叶逆变换获得时域空间的单频压力信号，将单频压力信号/>以虚部的形式添加到压力信号p上获得复数形式的新压力信号/>，/>，Im为虚元；基于虚拟点的数值、网格内部点的数值和新压力信号/>进行计算获得下一时刻tk网格内部点的物理量数据；
步骤3.3：将时刻t替换为时刻tk，判断时刻tk是否小于或等于t1，若是则返回执行步骤3.1，获得时刻t0至时刻t1对应的所有网格点的物理量数据，若否则执行步骤3.4；
步骤3.4：获得频率为f1的频域空间压力信号对应的一个或多个声源位置，包括：当tk=t1 时，基于所有网格点的复数形式的新压力信号的虚部构建压强云图，基于压强云图中压力信号最大处的位置获得频率为f1的频域空间压力信号对应的声源位置信息；
步骤4：对频率为f2、f3、…、fm的频域空间压力信号依次进行步骤3所述的噪声源定位处理，获得频率为f2、f3、…、fm的频域空间压力信号对应的声源位置信息；
步骤5：基于频率为f1、f2、f3、…、fm的频域空间压力信号对应的声源位置信息分析声源产生原因并进行分类，针对不同类型的声源采用相应的降噪方案进行降噪处理。
其中，时间反演的初始条件需要网格边界处的时间导数或时序数据，但对边界处理非常敏感，容易从外部引入伪波导致声源定位分辨率不高。反Lax-Wendroff方法即反拉克斯-温德洛夫方法是一种高精度、稳定的边界处理方法，需要给定边界条件，多适用于固壁边界、动边界情况。本发明将时间反演和反Lax-Wendroff方法进行结合，实现稳定高精度高分辨率的声源定位。
其中，有多个声源时，本方法能够基于频率从压强云图中获得对应的声源位置信息。每个声源对目标的影响效果不同，本方法可以确定频率对应的声源位置信息；然后基于确定的频率对应的声源位置信息，对该位置处的声源进行降噪处理，即能够从多个声源中确定相应频率的声源，然后对不同频率的噪声进行相应的降噪处理。
其中，在一些实施例中，网格边界处每个边界点与时刻t相关的物理量数据包括：时刻t对应的密度、时刻t对应的x方向速度/>、时刻t对应的y方向速度/>和时刻t对应的压强/>，x表示横向空间坐标，y表示纵向空间坐标；
时刻ttr对应的物理量数据包括：时刻ttr对应的密度、时刻ttr对应的x方向速度/>、时刻ttr对应的y方向速度/>和时刻ttr对应的压强，其中：
。
其中，在一些实施例中，网格边界处每个网格点对应守恒变量U的时间导数包括：、/>、/>和/>，其中，网格边界处每个网格点对应守恒变量U的时间导数的计算方式为：
；
其中，U1是守恒变量U的第一分量，U2是守恒变量U的第二分量，U3是守恒变量U的第三分量，U4是守恒变量U的第四分量，是比热比。
其中，在一些实施例中，网格边界处每个网格点在x方向无粘通量的空间切向导数包括：/>，/>，/>和/>；
网格边界处每个网格点在y方向无粘通量的空间切向导数包括：/>，/>，和/>；
其中：
；
其中，是y方向无粘通量/>的第一分量，/>是y方向无粘通量/>的第二分量，是y方向无粘通量/>的第三分量，/>是y方向无粘通量/>的第四分量；
；
其中，是x方向无粘通量/>的第一分量，/>是x方向无粘通量/>的第二分量，是x方向无粘通量/>的第三分量，/>是x方向无粘通量/>的第四分量。
其中，在一些实施例中，所述步骤3.2.3具体包括：基于网格边界处每个网格点对应守恒变量U的时间导数、网格边界处每个网格点在x方向无粘通量的空间切向导数和网格边界处每个网格点在y方向无粘通量/>的空间切向导数，采用反拉克斯-温德洛夫方法计算获得网格边界处每个网格点对应守恒变量U的空间法向一阶导数。
其中，在一些实施例中，所述步骤3.2.3具体包括：
求解网格边界处每个网格点上的特征值，获得四个特征方向的特征值；
构建左特征值矩阵；
基于四个特征方向的特征值对应的特征量获得出口边界信息和入流边界信息；
基于左特征值矩阵、出口边界信息和入流边界信息构建空间法向一阶导数求解矩阵；
求解所述空间法向一阶导数求解矩阵获得网格边界处每个网格点对应守恒变量U的空间法向一阶导数。
其中，在一些实施例中，四个特征方向的特征值为：，i=1或2或3或4，其中，，/>，/>，/>，/>表示该网格点处x方向速度，/>表示该网格点的声速；
左特征值矩阵为：
；
其中，是左特征值矩阵/>的第一分量，/>是左特征值矩阵/>的第二分量，/>是左特征值矩阵/>的第三分量，/>是左特征值矩阵/>的第四分量；为左特征值矩阵/>在坐标/>处的值；
基于四个特征方向的特征值的数值大小关系判断获得出口边界对应的特征量以及入流边界对应的特征向量；
出口边界对应的特征量采用四阶拉格朗日外插获得；
入流边界对应的特征向量采用反拉克斯-温德洛夫方法获得，基于入流边界对应的特征向量获得网格边界处每个网格点对应守恒变量U的空间法向一阶导数。
其中，在一些实施例中，网格预设位置的空间法向一阶导数构造多项式为：
；
其中，是多项式/>在x=x1处的值，/>是多项式/>在x=x2处的值，是多项式/>在x=x3处的值，/>是多项式/>在x=x4处的值，/>是多项式/>在x=x1处的一阶导数值，/>是守恒变量U在x=x1处的值，/>是守恒变量U在x=x2处的值，/>是守恒变量U在x=x3处的值，/>是守恒变量U在x=x4处的值，是守恒变量U在x=x1处的空间法向一阶导数。
其中，在一些实施例中，虚拟点的取值的计算方式为：
；
其中，是多项式/>在x=xj处的值。
其中，在一些实施例中，所述方法还包括：步骤6：获得降噪结果，根据降噪结果判断是否调整当前降噪方案。
其中，本方法还包括：步骤6：验证降噪处理后的效果。依据步骤5的降噪方案调整步骤1，并进行数值仿真计算，随后采用步骤2分析频率强度是否降低。如果降低，判断能否优化步骤5中降噪方案进一步降噪，不能则结束。如果不降反增，回到步骤5进一步分析声源之间的关系和调整降噪方案，通过步骤5和步骤6的多次迭代获得最优降噪处理方案。
其中，采用守恒型Euler方程对物理问题进行数值模拟，物理问题一般指气动声学相关问题，计算时间为t0时刻到t1时刻，采用均匀计算时间步长。记录t0时刻到t1时刻网格边界处每个边界点的物理量（密度/>、x方向速度/>、y方向速度和压强/>）随时间变化的数据，记录t1时刻全流场物理量数据（密度、x方向速度/>、y方向速度/>、压强/>）。其中，x表示横向空间坐标，y表示纵向空间坐标，t表示时间。
其中，步骤1的调整方式可以为：第一种是设计层面，由于降噪方案可能是对物理模型外观上的修改，所以需要对应修改物理模型；第二种是条件上的修改，采用流动控制方式会增加一些控制条件，需要修改初始的边界条件或初始条件。具体调整方式可以根据实际需要进行设定或调整，本发明实施例不进行限定。
尽管已描述了本发明的优选实施例，但本领域内的技术人员一旦得知了基本创造性概念，则可对这些实施例作出另外的变更和修改。所以，所附权利要求意欲解释为包括优选实施例以及落入本发明范围的所有变更和修改。
显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。

Claims (10)

1.一种时域多声源定位及噪声处理方法，其特征在于，所述方法包括：
步骤1：进行数值模拟获取数据，包括：分析气动噪声问题并建立气动噪声问题对应的物理模型；在网格生成软件中建立物理模型对应的几何模型；将几何模型外表面和整个计算区域进行空间网格划分，生成计算所需要的网格，设定边界条件并导出网格文件；将网格文件导入流体仿真软件，设定初始条件进行数值模拟计算，当计算达到稳定后的时刻记为时刻t0；计算n个周期后的时刻记为时刻t1，n为整数；在时刻t0至时刻t1之间，存储网格边界处每个边界点与时刻t相关的物理量数据和时刻t1整个网格每个点的物理量数据；
步骤2：对压力信号进行处理获得远场声信号频率分布和频率强度，包括：通过傅里叶变换将时域空间内网格边界上每个网格点处的压力信号p序列转换为频域空间压力信号，获得频域空间压力信号频率分布，将频率按照幅值从大到小排序标记为fr，r=1,2,3,…，m，m为频率个数；
步骤3：对频率为f1的频域空间压力信号对应的噪声源进行定位，包括：
步骤3.1：基于步骤1中存储的数据进行变换获得时间反演所需的数据，包括：基于时刻t和时刻t1获得时间反演的时刻ttr，t大于或等于t0且小于或等于t1，基于时刻t对应的网格边界处每个边界点物理量数据，获得时刻ttr对应的网格边界处每个边界点物理量数据；
步骤3.2：对边界上的物理量数据进行处理包括：
步骤3.2.1：获得边界点的时间导数，包括：基于所有时刻的所有网格边界处每个边界点物理量数据计算获得网格边界处每个网格点对应守恒变量U在时刻ttr的时间导数；
步骤3.2.2：获得边界点的空间切向导数，包括：基于时刻ttr对应的网格边界处每个边界点物理量数据，分别计算获得时刻ttr网格边界处每个网格点在x方向无粘通量的空间切向导数和在y方向无粘通量/>的空间切向导数；
步骤3.2.3：获得边界点的空间法向一阶导数，包括：基于网格边界处每个网格点对应守恒变量U在时刻ttr的时间导数、时刻ttr网格边界处每个网格点在x方向无粘通量的空间切向导数和时刻ttr网格边界处每个网格点在y方向无粘通量的空间切向导数，计算获得网格边界处每个网格点对应守恒变量U在时刻ttr的空间法向一阶导数；
步骤3.2.4：构造多项式，包括：基于守恒变量U在网格预设位置对应的数值和守恒变量U在网格预设位置的空间法向一阶导数构造多项式；
步骤3.2.5：令虚拟点的取值为所述多项式在相应点处的取值，对频率为f1对应的压力信号成分进行傅里叶逆变换获得时域空间的单频压力信号，将单频压力信号/>以虚部的形式添加到压力信号p上获得复数形式的新压力信号/>，/>，Im为虚元；基于虚拟点的数值、网格内部点的数值和新压力信号/>进行计算获得下一时刻tk网格内部点的物理量数据；
步骤3.3：将时刻t替换为时刻tk，判断时刻tk是否小于或等于t1，若是则返回执行步骤3.1，获得时刻t0至时刻t1对应的所有网格点的物理量数据，若否则执行步骤3.4；
步骤3.4：获得频率为f1的频域空间压力信号对应的一个或多个声源位置，包括：当tk=t1 时，基于所有网格点的复数形式的新压力信号的虚部构建压强云图，基于压强云图中压力信号最大处的位置获得频率为f1的频域空间压力信号对应的声源位置信息；
步骤4：对频率为f2、f3、…、fm的频域空间压力信号依次进行步骤3所述的噪声源定位处理，获得频率为f2、f3、…、fm的频域空间压力信号对应的声源位置信息；
步骤5：基于频率为f1、f2、f3、…、fm的频域空间压力信号对应的声源位置信息分析声源产生原因并进行分类，针对不同类型的声源采用相应的降噪方案进行降噪处理。
2.根据权利要求1所述的一种时域多声源定位及噪声处理方法，其特征在于，网格边界处每个边界点与时刻t相关的物理量数据包括：时刻t对应的密度、时刻t对应的x方向速度/>、时刻t对应的y方向速度/>和时刻t对应的压强，x表示横向空间坐标，y表示纵向空间坐标；
时刻ttr对应的物理量数据包括：时刻ttr对应的密度、时刻ttr对应的x方向速度/>、时刻ttr对应的y方向速度/>和时刻ttr对应的压强，其中：
。
3.根据权利要求2所述的一种时域多声源定位及噪声处理方法，其特征在于，网格边界处每个网格点对应守恒变量U的时间导数包括：、/>、/>和/>，其中，网格边界处每个网格点对应守恒变量U的时间导数的计算方式为：
；
其中，U1是守恒变量U的第一分量，U2是守恒变量U的第二分量，U3是守恒变量U的第三分量，U4是守恒变量U的第四分量，是比热比。
4.根据权利要求3所述的一种时域多声源定位及噪声处理方法，其特征在于，网格边界处每个网格点在x方向无粘通量的空间切向导数包括：/>，/>，/>和/>；
网格边界处每个网格点在y方向无粘通量的空间切向导数包括：/>，/>，/>和；
其中：
；
其中，是y方向无粘通量/>的第一分量，/>是y方向无粘通量/>的第二分量，/>是y方向无粘通量/>的第三分量，/>是y方向无粘通量/>的第四分量；
；
其中，是x方向无粘通量/>的第一分量，/>是x方向无粘通量/>的第二分量，/>是x方向无粘通量/>的第三分量，/>是x方向无粘通量/>的第四分量。
5.根据权利要求1所述的一种时域多声源定位及噪声处理方法，其特征在于，所述步骤3.2.3具体包括：基于网格边界处每个网格点对应守恒变量U的时间导数、网格边界处每个网格点在x方向无粘通量的空间切向导数和网格边界处每个网格点在y方向无粘通量/>的空间切向导数，采用反拉克斯-温德洛夫方法计算获得网格边界处每个网格点对应守恒变量U的空间法向一阶导数。
6.根据权利要求4所述的一种时域多声源定位及噪声处理方法，其特征在于，所述步骤3.2.3具体包括：
求解网格边界处每个网格点上的特征值，获得四个特征方向的特征值；
构建左特征值矩阵；
基于四个特征方向的特征值对应的特征量获得出口边界信息和入流边界信息；
基于左特征值矩阵、出口边界信息和入流边界信息构建空间法向一阶导数求解矩阵；
求解所述空间法向一阶导数求解矩阵获得网格边界处每个网格点对应守恒变量U的空间法向一阶导数。
7.根据权利要求6所述的一种时域多声源定位及噪声处理方法，其特征在于，四个特征方向的特征值为：，i=1或2或3或4，其中，/>，/>，/>，/>，表示该网格点处x方向速度，/>表示该网格点的声速；
左特征值矩阵为：
；
其中，是左特征值矩阵/>的第一分量，/>是左特征值矩阵/>的第二分量，/>是左特征值矩阵/>的第三分量，/>是左特征值矩阵/>的第四分量；为左特征值矩阵/>在坐标/>处的值；
基于四个特征方向的特征值的数值大小关系判断获得出口边界对应的特征量以及入流边界对应的特征向量；
出口边界对应的特征量采用四阶拉格朗日外插获得；
入流边界对应的特征向量采用反拉克斯-温德洛夫方法获得，基于入流边界对应的特征向量获得网格边界处每个网格点对应守恒变量U的空间法向一阶导数。
8.根据权利要求1所述的一种时域多声源定位及噪声处理方法，其特征在于，网格预设位置的空间法向一阶导数构造多项式为：
；
其中，是多项式/>在x=x1处的值，/>是多项式/>在x=x2处的值，是多项式/>在x=x3处的值，/>是多项式/>在x=x4处的值，/>是多项式/>在x=x1处的一阶导数值，/>是守恒变量U在x=x1处的值，/>是守恒变量U在x=x2处的值，/>是守恒变量U在x=x3处的值，/>是守恒变量U在x=x4处的值，是守恒变量U在x=x1处的空间法向一阶导数。
9.根据权利要求1所述的一种时域多声源定位及噪声处理方法，其特征在于，虚拟点的取值的计算方式为：
；
其中，是多项式/>在x=xj处的值。
10.根据权利要求1所述的一种时域多声源定位及噪声处理方法，所述方法还包括：步骤6：获得降噪结果，根据降噪结果判断是否调整当前降噪方案。
CN202410338625.2A 2024-03-25 2024-03-25 一种时域多声源定位及噪声处理方法 Active CN117935837B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202410338625.2A CN117935837B (zh) 2024-03-25 2024-03-25 一种时域多声源定位及噪声处理方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202410338625.2A CN117935837B (zh) 2024-03-25 2024-03-25 一种时域多声源定位及噪声处理方法

Publications (2)

Publication Number Publication Date
CN117935837A true CN117935837A (zh) 2024-04-26
CN117935837B CN117935837B (zh) 2024-05-24

Family

ID=90763348

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202410338625.2A Active CN117935837B (zh) 2024-03-25 2024-03-25 一种时域多声源定位及噪声处理方法

Country Status (1)

Country Link
CN (1) CN117935837B (zh)

Citations (20)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2004075162A2 (en) * 2003-02-20 2004-09-02 Ramot At Tel Aviv University Ltd. Method apparatus and system for processing acoustic signals
DE102004005998B3 (de) * 2004-02-06 2005-05-25 Ruwisch, Dietmar, Dr. Verfahren und Vorrichtung zur Separierung von Schallsignalen
CN102197424A (zh) * 2008-10-24 2011-09-21 高通股份有限公司 用于一致性检测的系统、方法、设备和计算机可读媒体
CN102197422A (zh) * 2008-10-24 2011-09-21 高通股份有限公司 使用传感器陈列进行音频源接近度估计以用于减少噪音
EP2530484A1 (en) * 2011-06-01 2012-12-05 Dolby Laboratories Licensing Corporation Sound source localization apparatus and method
CN107005763A (zh) * 2014-11-25 2017-08-01 美商楼氏电子有限公司 用于非线性和时变回波消除的参考麦克风
GB201802850D0 (en) * 2018-02-22 2018-04-11 Sintef Tto As Positioning sound sources
CN109883532A (zh) * 2019-03-14 2019-06-14 哈尔滨工程大学 一种声源识别与声场预报方法
CN110100453A (zh) * 2016-11-09 2019-08-06 伯斯有限公司 在双边麦克风阵列中控制风噪声
CN110298063A (zh) * 2019-05-10 2019-10-01 北方民族大学 一种非紧致可渗透边界气动噪声数值积分计算方法
CN112182881A (zh) * 2020-09-27 2021-01-05 大连理工大学 一种基于平台现象的非线性能量阱最优刚度求解方法
CN112179669A (zh) * 2020-09-21 2021-01-05 西北工业大学 基于发动机整机试验的喷流噪声测试方法
CN113032721A (zh) * 2021-03-11 2021-06-25 哈尔滨工程大学 一种低计算复杂度的远场和近场混合信号源参数估计方法
CN114125677A (zh) * 2020-08-28 2022-03-01 奥迪康有限公司 适于定向的听力装置
CN115435891A (zh) * 2022-08-17 2022-12-06 北京市高速公路交通工程有限公司 基于矢量传声器的道路车辆声功率监测系统
CN116167109A (zh) * 2023-01-11 2023-05-26 浙江工业大学 湍流边界层激励下结构声辐射近场计算方法
CN116564265A (zh) * 2023-05-22 2023-08-08 哈尔滨工程大学 基于参量次级声源的水下结构噪声控制方法
CN116630106A (zh) * 2023-06-28 2023-08-22 沈阳创新设计研究院有限公司 一种智造实训交互教学管理方法及系统
CN116952356A (zh) * 2023-07-24 2023-10-27 中国人民解放军海军工程大学 基于浅海环境水下声全息技术的近场辐射噪声测量方法
CN117292060A (zh) * 2023-09-22 2023-12-26 长春理工大学 一种基于声波反散射的三维人体无触重建方法

Patent Citations (21)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2004075162A2 (en) * 2003-02-20 2004-09-02 Ramot At Tel Aviv University Ltd. Method apparatus and system for processing acoustic signals
DE102004005998B3 (de) * 2004-02-06 2005-05-25 Ruwisch, Dietmar, Dr. Verfahren und Vorrichtung zur Separierung von Schallsignalen
CN102197424A (zh) * 2008-10-24 2011-09-21 高通股份有限公司 用于一致性检测的系统、方法、设备和计算机可读媒体
CN102197422A (zh) * 2008-10-24 2011-09-21 高通股份有限公司 使用传感器陈列进行音频源接近度估计以用于减少噪音
CN103295579A (zh) * 2008-10-24 2013-09-11 高通股份有限公司 用于一致性检测的系统、方法和设备
EP2530484A1 (en) * 2011-06-01 2012-12-05 Dolby Laboratories Licensing Corporation Sound source localization apparatus and method
CN107005763A (zh) * 2014-11-25 2017-08-01 美商楼氏电子有限公司 用于非线性和时变回波消除的参考麦克风
CN110100453A (zh) * 2016-11-09 2019-08-06 伯斯有限公司 在双边麦克风阵列中控制风噪声
GB201802850D0 (en) * 2018-02-22 2018-04-11 Sintef Tto As Positioning sound sources
CN109883532A (zh) * 2019-03-14 2019-06-14 哈尔滨工程大学 一种声源识别与声场预报方法
CN110298063A (zh) * 2019-05-10 2019-10-01 北方民族大学 一种非紧致可渗透边界气动噪声数值积分计算方法
CN114125677A (zh) * 2020-08-28 2022-03-01 奥迪康有限公司 适于定向的听力装置
CN112179669A (zh) * 2020-09-21 2021-01-05 西北工业大学 基于发动机整机试验的喷流噪声测试方法
CN112182881A (zh) * 2020-09-27 2021-01-05 大连理工大学 一种基于平台现象的非线性能量阱最优刚度求解方法
CN113032721A (zh) * 2021-03-11 2021-06-25 哈尔滨工程大学 一种低计算复杂度的远场和近场混合信号源参数估计方法
CN115435891A (zh) * 2022-08-17 2022-12-06 北京市高速公路交通工程有限公司 基于矢量传声器的道路车辆声功率监测系统
CN116167109A (zh) * 2023-01-11 2023-05-26 浙江工业大学 湍流边界层激励下结构声辐射近场计算方法
CN116564265A (zh) * 2023-05-22 2023-08-08 哈尔滨工程大学 基于参量次级声源的水下结构噪声控制方法
CN116630106A (zh) * 2023-06-28 2023-08-22 沈阳创新设计研究院有限公司 一种智造实训交互教学管理方法及系统
CN116952356A (zh) * 2023-07-24 2023-10-27 中国人民解放军海军工程大学 基于浅海环境水下声全息技术的近场辐射噪声测量方法
CN117292060A (zh) * 2023-09-22 2023-12-26 长春理工大学 一种基于声波反散射的三维人体无触重建方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
PHILIPPE DRUAULT: ""Localization of aeroacoustic sound sources in viscous flows by a time reversal method"", 《JOURNAL OF SOUND AND VIBRATION》, vol. 33, no. 15, 15 June 2013 (2013-06-15) *
潘威振: ""材料变热物性参数辨识算法与验证"", 《中国优秀硕士学位论文全文数据库 工程科技辑》, 19 February 2016 (2016-02-19) *

Also Published As

Publication number Publication date
CN117935837B (zh) 2024-05-24

Similar Documents

Publication Publication Date Title
JP4891910B2 (ja) ノイズ源の遠距離音場分析
JP5496192B2 (ja) 音響場を再構成するための方法
EP2297557B1 (en) Reconstructing an acoustic field
Lee et al. Source characterization of a subsonic jet by using near-field acoustical holography
Lee et al. Patch near-field acoustical holography in cylindrical geometry
Zhang et al. An alternative hybrid time-frequency domain approach based on fast iterative shrinkage-thresholding algorithm for rotating acoustic source identification
Zamponi et al. 3D Generalized Inverse Beamforming in wind tunnel aeroacoustic testing: application to a Counter Rotating Open Rotor aircraft model
Zhongming et al. A monotonic two-step iterative shrinkage/thresholding algorithm for sound source identification based on equivalent source method
Lehmann et al. Beamforming for measurements under disturbed propagation conditions using numerically calculated Green’s functions
Ruhala et al. Planar near-field acoustical holography in a moving medium
Zhao et al. Improved technique for evaluation of wall pressure fluctuations from turbulent boundary layer
CN106017837B (zh) 一种等效声模拟源的模拟方法
CN117935837B (zh) 一种时域多声源定位及噪声处理方法
Pan et al. Multipole orthogonal beamforming combined with an inverse method for coexisting multipoles with various radiation patterns
Huijssen et al. Numerical evaluation of source–receiver transfer functions with the Fast Multipole Boundary Element Method for predicting pass-by noise levels of automotive vehicles
CN117932199B (zh) 一种时域声源定位方法
Li et al. Functional generalized inverse beamforming based on the double-layer microphone array applied to separate the sound sources
CN114996902A (zh) 一种考虑车辆运动影响的通过噪声室内测量方法
Thor et al. Source identification of a bladeless fan by using sonah in cylindrical coordinates
Xu et al. Denoising weighting beamforming method applied to sound source localization with airflow using microphone array
Arcondoulis et al. Experimental application of an acoustic beamforming array pairing method using CLEAN-SC
Yue et al. Method for estimating the radiated sound power of moving sound sources with linear array beamforming
Tsuchiya et al. Two-dimensional finite-difference time-domain simulation of moving sound source and receiver with directivity
CN115782480B (zh) 一种滚动轮胎结构噪声时域预测方法
Pignier et al. Aeroacoustic study of a submerged air inlet using an IDDES/FW-H approach and sound source modelling through direct numerical beamforming

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
GR01 Patent grant
GR01 Patent grant