DOI: 10.11817/j.issn.1672-7207.2016.09.020
管道内旋转声源频域声场的数值仿真
付建，王永生，刘强，苏永生
(海军工程大学 动力工程学院，湖北 武汉，430033)
摘要：为准确计算管道内的旋转声源声场，研究基于改进的点源模型和边界元法的声场仿真方法。首先从时域上将旋转点声源离散为均匀分布在旋转轨迹上的有限的静止点声源，通过静止点声源的依次发声模拟旋转点声源的发声过程；将时域声源变换到频域，然后结合边界元方法即可计算得到管道内所有离散的静止点声源声场，该声场即为旋转点声源声场。验证基于改进点源模型计算旋转声源自由声场以及利用边界元方法计算管道内静止点声源声场的可信性，以泵喷推进器为对象，将泵喷叶片分块划分并等效为偶极源，计算分析导管内的叶轮声场。研究结果表明：以改进的点源模型和边界元法为基础的管道内旋转声源声场仿真方法是有效易行的；对导管内叶轮声场而言，入射声场与散射声场存在较大差异，导管对径向测点处叶轮声场影响较大，对轴向测点处叶轮声场的影响可以忽略。
关键词：旋转声源；管道；声场；泵喷推进器；改进的点源模型；边界元法
中图分类号：O42 文献标志码：A 文章编号：1672-7207(2016)09-3058-08
Sound filed simulation in frequency domain of rotating source in duct
FU Jian, WANG Yongsheng, LIU Qiang, SU Yongsheng
(College of Marine Power Engineering, Naval University of Engineering, Wuhan 430033, China)
Abstract: In order to simulate the sound field of rotating source in duct accurately, the method based on modified point source model and boundary element method was explored. Firstly, the rotating point source was replaced by finite stationary point sources distributing along the rotating trace in time domain, and the phonation of rotating source was simulated by the phonation of stationary sources in turn. The stationary sources were changed into frequency domain, then the sound field of stationary sources in duct was calculated by boundary element method, which is equal to the sound field of rotating source. When the creditability of modified point source model used for free sound field simulation of rotating source and BEM used for sound filed simulation of stationary source in duct was validated, the pumpjet was taken as a analysis object, the blade was divide into segments and each segment was equivalent to a dipole, then the sound field of rotating impeller in duct was calculated. The results show that the sound field simulation method of rotating source in duct based on modified point source model and BEM is effective and exercisable. There are differences between the incident sound field and scattered sound field of impeller; the influence of duct on the sound field of impeller is remarkable in radial, but the influence in axial can be neglected.
Key words: rotating source; duct; sound field; pumpjet; modified point source model; boundary element method
由于管道对噪声有明显的抑制作用，常被用作降低风扇、推进器等旋转机械噪声的有效手段。声波在管道中的传播问题相当复杂，管道声学的研究方法经历了由解析法向数值方法的转变。解析法仅适用于模型和边界条件较为简单的管道声场计算，其适用范围受到限制[1-3]；在数值方法中，由于边界元法只需在管道边界面上划分网格，不受管道形状、边界条件的限制，在管道声场的数值计算中应用最为广泛[4-8]。旋转声源如推进器、风扇的声场主要是依据声类比方程或点源模型求解。声类比方程适用于旋转机械的自由声场计算，不能考虑推进器导管或风扇管道对叶轮声场的影响[9-11]；点源模型即把旋转叶轮等效为若干个旋转点声源，所有旋转点声源声场总和即为旋转叶轮声场[6-8, 12-14]。前人对管道内旋转声源声场仿真所做研究主要侧重于管道的声传播和散射特性，没有对声源特性做出详细的分析，且多以点声源为研究对象，而本文所研究的声源类型既可以是点声源，亦可以是分布式声源。付建等[15]进一步完善了点源模型在旋转机械自由声场计算中的应用，并准确预报了旋转机械的偶极噪声。本文作者在结合改进的点源模型和边界元方法的基础上，研究了管道内旋转声源声场的仿真方法，并以泵喷推进器为对象，计算了考虑导管散射时的叶轮声场，并分析了导管对叶轮声场的影响。
1 旋转声源的频域离散方法
任何噪声源都可以看作由多个具有适当相位、幅值和位置的点声源组成。采用点源模型需要满足的条件是：声源最高频率的波长应该远大于声源的物理尺寸或者场点与声源的最近距离要远大于声源的物理尺寸[13]。将旋转叶轮表面按一定规则离散为多个网格单元，当单元尺寸满足上述要求时，可将每个表面单元看作一个小尺寸的孤立声源，叶轮总噪声即为所有孤立声源噪声的总和。
根据点源模型理论将每个旋转点声源沿运动轨迹离散为旋转圆周上均匀分布的一系列具有相位差β的静止点声源(见图1)，每个离散声源乘以图2所示的时域矩形函数(其中：T 为旋转周期，τ 代表2个相邻分布声源之间的时间差)，使得不同时刻静止点声源按照沿旋转点声源运动轨迹的排列顺序依次发声，从而模拟点声源在旋转过程中的发声过程。将所得结果进行傅里叶分解后即可得到频域声源，傅里叶变换后，点声源之间时域的延时在频域表现为相位差[15]。将离散后的点声源导入声学计算软件Virtual Lab，借助声学软件平台利用边界元方法可以完成任意边界条件下旋转声源声场仿真。上述将旋转点声源离散为静止点声源的方法称之为改进的点源模型。
图1 旋转声源离散示意图
Fig. 1 Sketch map of discrete rotating point source
图2 声源离散时域矩形函数
Fig. 2 Rectangle function for sound source discretization in time domain
2 自由场中旋转声源声场仿真与校验
在前期研究中，付建等[15]将改进的点源模型应用到自由空间旋转声源声场计算中，并进行了校核。本文仅列出自由空间旋转单级子、旋转偶极子的声场验证(见图3和图4)以及改进的点源模型在螺旋桨自由场负载噪声(亦称偶极噪声)计算中的应用结果(见图5)。螺旋桨负载噪声计算时将桨叶表面的每个单元等效为一个偶极子，文献值为基于FW-H方程计算得到的数值结果[16]。
图3 旋转单级子声压指向性
Fig. 3 Acoustic directivity of rotating monopole in free field
图4 旋转偶极子量纲一化声压指向性
Fig. 4 Nondimensionalized acoustic directivity of rotating dipole in free field
图5 单个桨叶负载噪声
Fig. 5 Loading noise induced by single blade
由旋转单级子、旋转偶极子和螺旋桨负载噪声数值仿真结果可知，基于改进的点源模型计算旋转声源自由声场的方法是准确可信的，这为下一步管道内旋转声源声场的准确计算奠定了基础。
3 管道内静止点声源声场的数值仿真
由于管道的存在使得声源对应声场不仅包括入射声场还包括散射声场，分别计算了管道内简单声源的入射声场和散射声场，并同文献值进行了对比。
3.1 边界元方法介绍
边界元法分为直接边界元法和间接边界元法。直接边界元法采用具有明确物理意义的量(如声压和振速)作为变量求解，适用于封闭空间的内声场或外声场求解；间接边界元法以振动位移和声学势函数(声压的变化量及声压梯度)作为变量求解，适用于封闭空间或非封闭空间的声场求解。间接边界元方法可以从直接边界元方法推导出来，以两端开口管道为例，通过将整个管道声场分为管道内声场和外声场2部分，用进出口边界声压和声质点速度连续的条件，将内外声场联系起来，将直接边界元的积分方程分别应用于边界表面所构成的内外域，然后2个方程相减，即可得到对应的间接边界元方程。
管道内点声源声场Pt为入射声场Pi和散射声场Ps的叠加，即：
Pt=Pi+Ps (1)
借鉴SEOL等[17]计算导管桨噪声的方法，基于Helmholtz方程预报简单声源的非自由声场，并应用边界元法求解Helmholtz方程。Helmholtz方程解的积分形式可以表示为：
(2)
(3)
其中：k =ω/c为波数，ω为声源圆频率，c为流体中声速；r为结构表面配置点y与空间场点x的距离，G(x, y)=eikr/(4πr)为自由空间格林函数的基本解，S表示结构散射边界(假设边界光滑)，E为结构外部区域，ny为结构表面点y处的外法线方向。管道作为刚性壁面，在入射声波确定后，由式(2)可以首先计算出结构表面声压，然后由式(3)可计算出对应的散射声场和总声场。
3.2 简单声源声场仿真与校验
3.2.1 管道内静止单极子声场仿真
美国NASA-Langley中心提供了验证刚性壁面导管声场预报模型的标准算例[4]，如图6所示。已知条件如下：导管的长度L和直径D均取1 m，声源S位于导管中心，坐标系xyz原点位于导管中心处，z轴与管轴重合，源项为频率已知的单极子谐波声源，并且满足ρ0c|Q|=1 N(ρ0为流体密度，Q为声源强度)。图7所示为在yoz平面上距声源2.5 m半径处750 Hz所对应的散射声场和总声场仿真值与文献值的对比(声压级
，
，P*为声压P的复数共轭，Pref = 0.2 μPa。
图6 管道内静止单极子参数定义
Fig. 6 Dimensions of stationary monopole in duct
图7 单极子的散射声场和总声场
Fig. 7 Scattering and incident sound field of monopole
3.2.2 管道内静止偶极子声场仿真
Martinez结合升力面理论求解了两端开口的直管道内点力源对应的远场噪声[3]，计算对象如图8所示。基于边界元方法计算了管道内偶极子不同频率对应的远场声压指向性(
，Pi为远场轴线处的入射声场声压)，如图9所示。
由图7和图9可知：基于边界元法计算得到的管道内简单声源声场结果与文献值吻合较好，验证了边界元法在管道内声源声场仿真中的适用性。
图8 管道内静止偶极子参数定义
Fig. 8 Dimensions of stationary dipole in duct
图9 偶极子的远场声压指向性
Fig. 9 Far field acoustic directivity of dipole
4 管道内旋转声源声场仿真
以前置导叶式泵喷推进器为对象(如图10所示，泵喷推进器叶片数为9，导叶数为11)，计算了导管内的旋转叶轮声场，并分析了导管对推进器叶轮声场的影响。叶轮是泵喷推进器的主要做功部件，且在水下航行器尾部的非均匀流场中工作，叶轮偶极噪声是声场的主要成分，因此本文仅计算导管内叶轮的偶极源声场。
图10 泵喷推进器结构示意图
Fig. 10 Schematic map of pumpjet structure
4.1 导管内叶轮声场仿真方法
推进器瞬态流场CFD计算时的叶轮壁面单元网格尺度远超过声学边界元计算对网格尺寸的要求，如果叶轮壁面的每个流体网格单元等效为一个旋转点力源，声场计算时的离散偶极子数目巨大，对计算资源的要求很高，尤其是非自由声场的仿真。因此结合点源模型理论，当叶片尺寸同分析频率对应波长相比足够小时，将单个叶片作为紧致声源处理；如果单个叶片不能作为紧致声源处理时，可以将叶片沿径向和轴向划分为若干部分，每部分作为紧致声源处理，划分方法见图11。
对叶轮偶极噪声而言，每个紧致声源对应一个点力源，其位置坐标Cjk和幅值Fijk为：
(4)
(5)
图11 叶片分块划分示意图
Fig. 11 Schematic map of blade segmentation
式中：xj为流体网格(或声网格)单元中心的不同坐标分量，j是张量编号分别对应为X，Y，Z方向；l代表不同部分编号，l=1, 2, …, m×n；NE为第l部分所包含的单元数，pi和ne分别为单元压力与法向量，i对应叶轮旋转过程的时间步编号。由式(4)和(5)可得到每个点力源的初始旋转位置以及点力源不同方向载荷分量随时间的变化曲线。叶片与点力源之间的分块等效方法在旋转叶轮声场预报中的适用性和准确性已得到了验证[15]。
在对泵喷叶轮声场进行仿真时，首先将泵喷叶轮壁面流体网格上的压力映射到叶轮声场网格，然后将叶片分块离散，每个分块等效为一个偶极子，根据瞬态流场计算时的时间步长将旋转偶极子沿运动轨迹离散为有限个有固定相位差的偶极子，最后考虑导管的声反射、散射作用，结合边界元方法即可计算得到导管内叶轮声场。叶轮声场仿真过程如图12所示。
4.2 导管内叶轮声场的计算与分析
首先以单个泵喷叶片为对象，计算分析叶片与点力源之间沿轴向和径向的不同等效方案对声场结果的影响，以确定合理的等效方案。声场计算对应的场点布置如图13所示，其中指向性场点在水平方向360°范围内，场点之间间隔3°。图14所示为不同点力源等效方案对应的单个叶片自由声场在轴向测点2处的声压频谱。
由图14可知：不同离散方案对应的单个叶片噪声在低频基本一致，随着频率升高，略有差别，但不显著，因此叶轮噪声计算时选用2×2方案。此外，单个叶片的轴向测点声压频谱在11 fA (fA为轴频)及其谐频处有线谱出现，这是因为叶轮工作于导叶后的非均匀流场中，叶轮和导叶的相互作用使得单个叶片噪声在11 fA及其谐频处出现线谱。
图12 叶轮声场仿真过程示意图
Fig. 12 Calculation progress of impeller noise in duct
图13 声场仿真时场点布置
Fig. 13 Collocation of field point
图14 不同等效方案对应的单个叶片轴向测点噪声
Fig. 14 Single blade noise in axis with different segmentation
在确定了叶片与点力源的等效方案后，计算导管内的叶轮声场。图15所示为考虑导管声散射效应前后测点1和2处的叶轮线谱噪声；图16所示为考虑导管声散射效应前后部分频率点处叶轮的远场声压指向性。
由图15可知：考虑导管散射后，径向测点1对应的叶轮噪声在叶频(fB)处变化最大，降低了9.6 dB，其余频率点多有不同量值降低；导管对轴向测点2处的噪声影响很小，变化均在1 dB之内。
由图16可知：考虑导管后，叶频处的叶轮声场在90°和276°处变化较大，其中90°处声压降低了9.6 dB；3 fB声场在100°和280°处变化较大；5 fB声场在轴向测点基本不变，其余测点声压多为降低，且在120°和345°附近变化最大。
为进一步深入分析导管声散射效应对叶轮声场的影响，将叶轮对应的每个离散偶极子按照力的矢量方向分解为FX，FY和FZ 3个分量，并计算所有偶极子的不同方向分量所对应的入射声场和散射声场声压指向性。叶频是叶轮线谱噪声的主要来源，因此仅以叶频处的远场声压指向性为对象进行分析，如图17所示。
图15 考虑导管声散射效应前后测点1和2处的叶轮线谱噪声
Fig. 15 Tonal noise of rotor at point 1 and 2 considering scattering effect of duct
图16 fB，3 fB和5 fB处的叶轮远场声压指向性
Fig. 16 Far field directivity of rotor at fB, 3 fB and 5 fB
图17 叶频处离散偶极子不同方向分量对应的声压指向性
Fig. 17 Acoustic directivity of FX, FY and FZ at blade-passing frequency
由图17可知：由于指向性场点所在平面与X轴垂直，使得FX在YZ平面的入射声场较弱，且指向性基本呈圆形；FY和FZ对应的入射声场声压指向性分别是关于Y轴和Z轴对称的“8”字形，这与自由场条件下偶极源的声场指向性一致，但对不同位置的径向测点而言，FY在90°处的声压要高于270°处。考虑导管声散射作用后，FX对应的散射声场要强于入射声场，尤其是在60°和240°附近；FY对应的散射声场在90°最强，且与入射声场相当，由于散射声场与入射声场存在相位差使得二者合成后的总声场较入射声场弱；FZ的散射声场在轴向低于入射声场，而在径向高于入射声场，使得在轴向测点散射声场对总声场的影响很小；由于径向测点1处的声场主要是由FX和FY贡献，轴向测点2处的声场主要由FZ贡献，因此考虑导管散射后，轴向测点2的声压基本不变，而径向测点1声压变化较大。
由上述计算分析结果可知：导管内叶轮的入射声场和散射声场存在较大差异，并且叶轮的入射声场在轴向测点主要是轴向力贡献，在径向测点主要是径向力贡献；导管的存在会对叶轮声场产生较大影响，尤其是在径向测点，可称其为“侧壁屏蔽效应”，由于导管长度较短其对轴向测点声场的影响可以忽略。
5 结论
1) 基于改进的点源模型将旋转点声源离散为多个静止点声源，然后利用边界元方法计算得到管道内所有离散的静止点声源声场总和，即为管道内的旋转点声源声场。所提出的管道内旋转声源声场仿真方法有效易行。
2) 以泵喷推进器为对象，将叶片分块划分等效为偶极源，计算分析了导管内的叶轮声场。对导管内叶轮声场而言，叶轮的入射声场和散射声场存在较大差异，导管对径向测点处叶轮声场影响较大，对轴向测点处叶轮声场的影响可以忽略。
3) 如何利用该方法确定叶轮在导管中的相对位置，以使得辐射能量最低，同时利用该方法指导导管的声学优化设计将是下一步工作。
参考文献：
[1] 刘秋洪, 黄曌宇, 祁大同, 等. 圆管内旋转点声源声学频域分析[J]. 西安交通大学学报, 2004, 38(5): 541-544.
LIU Qiuhong, HUANG Zhaoyu, QI Datong, et al. Frequency domain solution of a rotating point source in an infinite duct[J]. Journal of Xi’an Jiaotong University, 2004, 38(5): 541-544.
[2] 黄曌宇, 刘秋洪, 祁大同, 等. 封闭圆柱腔内旋转点声源声压的近场频域解[J]. 应用力学学报, 2004, 21(2): 37-40.
HUANG Zhaoyu, LIU Qiuhong, QI Datong, et al. Near field frequency domain solution of rotating point source in closed cylindrical cavity[J]. Chinese Journal of Applied Mechanics, 2004, 21(2): 37-40.
[3] MARTINEZ R. Diffracting open-ended pipe treated as a lifting surface[J]. AIAA Journal, 1988, 26(4): 396-404.
[4] LAKOTA B A. The radiated field generated by a monopole source in a short, rigid, rectangular duct[D]. Washington D C: George Washington University. School of Engineering and Applied Science, 1995: 36-42.
[5] 杨迎春, 周其斗, 潘雨村. 有限长导管声场预报的一种边界元方法[J]. 声学学报, 2011, 36(1): 60-65.
YANG Yingcun, ZHOU Qidou, PAN Yucun. Boundary element method for sound field prediction of finite duct[J]. Acta Acoustic, 2011, 36(1): 60-65.
[6] CHOI H L, LEE D J. Development of the numerical method for calculating sound radiation from a rotating dipole source in an opened thin duct[J]. Journal of Sound and Vibration, 2006, 295(3/4/5): 739-752.
[7] DUNN M H, TWEED J, FARASSAT F. The application of a boundary integral equation method to the prediction of ducted fan engine noise[J]. Journal of Sound and Vibration, 1999, 227(5): 1019-1048.
[8] 马晓晶. 管道内旋转声源辐射声场的数值计算方法研究[D]. 北京: 中国科学院, 2009: 19-34.
MA Xiaojing. The numerical method of the sound field radiated from rotating source in duct[D]. Beijing: Chinese Academy of Sciences, 2009: 19-34.
[9] YAMADE Y, WANG H, GUO Y, et al. Numerical prediction of sound generated from flows with a low mach number[J]. Computers & Fluids, 2007, 36(1): 53-68.
[10] KRISHNA S R, KRISHNA A R, RAMJI K. Reduction of motor fan noise using CFD and CAA simulations[J]. Applied acoustics, 2011, 72(12): 982-992.
[11] 王嘉冰, 刘敏, 吴克启. 开式轴流风扇气动噪声预测[J]. 工程热物理学报, 2007, 28(5): 778-780.
WANG Jiabing, LIU Min, WU Keqi. Aerodynamic noise prediction of the open axial flow fan[J]. Journal of Engineering Thermophysics, 2007, 28(5): 778-780.
[12] TYLER J M, SOFRIN T G. Axial flow compressor noise studies[J]. Society of Automotive Engineers, 1962, 70: 309-332.
[13] KUCUKCOSKUN K, CHRISTOPHE J, SCHRAM C, et al. Free and scattered acoustic field predictions of the broadband noise generated by a low-speed axial fan[J]. Noise Control Engineering Journal, 2013, 61(2): 180-189.
[14] 蔡建程, 祈大同, 卢傅安, 等. 离心风机基频气动偶极子噪声的数值研究[J]. 西安交通大学学报, 2009, 43(11): 66-71.
CAI Jiancheng, QI Datong, LU Fuan, et al. Numerical study on blade-passing frequency noise of aerodynamic dipole in a centrifugal fan[J]. Journal of Xi’an Jiaotong University, 2009, 43(11): 66-71.
[15] 付建, 王永生. 点源模型在旋转声源声场计算中的推广应用[J]. 哈尔滨工程大学学报, 2014, 35(6): 312-317.
FU Jian, WANG Yongsheng. Extended application of the point source model for the acoustic field solution of the rotating source[J]. Journal of Harbin Engineering University, 2014, 35(6): 312-317.
[16] TESTA C, SALVATORE F. Theoretical and numerical issues for hydroacoustic applications of the Ffowcs Williams-Hawkings equation[C]// 26th AIAA Aeroacoustics Conference. Monterey, California: American Institute of Aeronautics and Astronautics, 2005: 1-13.
[17] SEOL H, JUNG B, SUH J C, et al. Prediction of non-cavitating underwater propeller noise[J]. Journal of Sound and Vibration, 2002, 257(1): 131-156.
(编辑 赵俊)
收稿日期：2015-05-31；修回日期：2015-09-09
基金项目(Foundation item)：国家自然科学基金资助项目(51306205) (Project(51306205) supported by the National Natural Science Foundation of China)
通信作者：付建，博士，从事旋转机械水动力噪声数值预报方法研究；E-mail: fujian_qdqy@163.com
