矩形单元声场波函数构造及其应用.pdf

第 43 卷 第 1 期Vol.43,No.12024 年 1 月Journal of Applied AcousticsJanuary,2024 研究报告 矩形单元声场波函数构造及其应用贺佐潦霜1,2向 宇1,2石梓玉1,2陈岩豪1,2陆 静1,2(1 广西科技大学广西汽车零部件与整车技术重点实验室柳州545006)(2 广西科技大学机械与汽车工程学院柳州545006)摘要：在波叠加法中，结构外部声场是在离散边界上对Green函数进行积分并叠加得到，但数值积分的计算效率较低。而等效源法虽然提高了计算效率，但其面源简化为点源的过程中存在较大的积分近似误差。针对上述两种方法的缺陷，构造了一种波函数以替代离散单元关于Green 函数积分的声场。首先，利用球坐标系下Helmholtz方程的解，推导了替代矩形单元积分的一般形式波函数及效率更高的内推波函数。其次，当离散单元为正方形时，将其近似成圆形域，进一步简化了内推波函数的表达式。最后，将所构造的波函数应用于声场计算。数值结果表明，在计算单个矩形单元外部辐射声场时，构造的波函数不仅保证了计算精度，而且相比于直接积分大幅度提高了计算效率。其中，矩形域一般形式和内推形式的波函数计算效率是直接积分的5 6倍，圆形域内推波函数计算效率达到了直接积分的12 13倍。在简支板声源和立方箱体辐射声源数值算例中，圆形域内推波函数在整个计算频段的声场计算精度均高于等效源法。关键词：等效源法；波叠加法；波函数；辐射声场中图法分类号:TB52文献标识码:A文章编号:1000-310X(2024)01-0100-10DOI:10.11684/j.issn.1000-310X.2024.01.013Construction and application of wave function of sound field ofrectangular elementHE Zuoliaoshuang1,2XIANG Yu1,2SHI Ziyu1,2CHEN Yanhao1,2LU Jing1,2(1 Key Laboratory of Automobile Component and Vehicle Technology in Guangxi,Guangxi University of Science andTechnology,Liuzhou 545006,China)(2 School of Mechanical and Traffic Engineering,Guangxi University of Science and Technology,Liuzhou 545006,China)Abstract:In the wave superposition method,the sound field outside the structure is obtained by integratingand superimposing the Green function on the discrete boundary,but the computational efficiency of numericalintegration is low.Although the equivalent source method improves the calculation efficiency,there is a largeintegral approximation error in the process of simplifying the area source to the point source.In view of thedefects of the above two methods,a wave function is constructed to replace the sound field of the discreteelement integration about Green function.First of all,using the solution of Helmholtz equation in sphericalcoordinates,the general wave function which replaces the rectangular element integral and the more efficientinternal wave function are derived.Secondly,when the discrete element is square,it is approximated to a circu-lar domain,which further simplifies the expression of the extrapolated wave function.Finally,the constructedwave function is applied to sound field calculation.The numerical results show that the constructed wave func-tion not only guarantees the calculation accuracy,but also greatly improves the calculation efficiency compared2022-08-12收稿;2022-11-28定稿国家自然科学基金项目(51775121)作者简介:贺佐潦霜(1999),男,江西吉安人,硕士研究生,研究方向:振动与噪声。通信作者 E-mail:第43卷 第1期贺佐潦霜等：矩形单元声场波函数构造及其应用101with direct integration when calculating the external radiation sound field of a single rectangular element.Among them,the wave function calculation efficiency of the general form and the extrapolation form in therectangular domain is 56 times higher than that of the direct integration,and the wave function calculationefficiency in the circular domain is 1213 times higher than that of the direct integration.In the numericalexamples of simply supported plate sound source and cubic box radiation sound source,the accuracy of thewave extrapolation function in the circular domain is higher than that of the equivalent source method in thewhole calculation frequency band.Keywords:Equivalent source method;Wave superposition method;Wave function;Radiated sound field0 引言在振动结构的辐射声场计算中，边界元法(Boundary element method,BEM)是一种被广泛使用的计算方法12。通常使用的边界元法包括直接边界元法(Direct-BEM)3和间接边界元法(Indirect-BEM)4，其在计算声源外部辐射声场时，由于需要数值积分的计算，因此计算效率较低，而且当积分面与振动体表面重合时，会产生奇异性问题5，需要使用额外的数值技术解决67。若将积分面设置在振动体内部(用于求解外问题)或者外部(用于求解内问题)的一个虚拟表面上，即可避免奇异积分的处理，这种方法被称为波叠加法(Wavesuperposition method,WSM)89或修正边界元法(Modified BEM,MBEM)1011。该方法避免了奇异积分的处理，但在求解振动体辐射声场时仍需对离散单元进行数值积分计算，当离散单元数量较多、求解规模较大时，单元的积分计算会耗费大量时间。WSM中的数值积分来自于对面源辐射声场的计算，若将其用单极子点源代替，则可以得到一种无需积分的高效率计算方法，即等效源法(Equivalentsource method,ESM)1213。但ESM中的单极子点源是对面源辐射声场的过度简化，其计算精度和稳定性较大程度依赖于等效源位置和数目的选择。针对以上问题，相关学者提出了优化等效源配置1415、采用具有指向性的射线波函数来代替单极子点源1617等方法，并取得了一定的研究成果。但在面源简化为点源的过程中，始终存在较大积分近似误差11，一定程度上影响了计算精度。针对上述缺陷，本文根据外问题中声源产生的空间声场，总可以展开成一系列不同阶次的球Hankel函数与球谐波函数乘积的加权和18，构造了一种波函数替代WSM中Green 函数在单元区域的积分，避免了在求解振动体外部辐射声场中复杂的积分计算。以矩形常数单元为例，推导了替代矩形单元积分一般形式和内推形式的波函数，以及当离散单元为正方形时，可将波函数进一步简化为圆形域内推波函数。文中对比了几种波函数与直接积分的计算精度和计算效率，并通过简支板声源和立方箱体辐射声源的数值仿真算例，验证了圆形域内推波函数与ESM在声场计算中的效果。1波叠加积分方程及其离散形式根据WSM理论，振动体外部r处的辐射声压可以利用式(1)计算：p(r)=Sq(rE)G(r,rE)dS(rE),(1)其中，S为振动体内部的一个封闭虚拟曲面，q(rE)为虚拟面S上rE处的源强，G(r,rE)为自由空间Green函数：G(r,rE)=eik|rrE|4|r rE|,(2)其中，|r rE|为场点r与虚拟面上rE处之间的距离，k为波数，i=1。若将虚拟面离散成N 个单元，式(1)可写为p(r)=Ni=1Siqi(rE)G(r,rE)dSi(rE).(3)当离散的单元足够小时，可以将每一个单元的源强qi(rE)均视为常数qEi，式(1)变为p(r)=Ni=1qEiSiG(r,rEi)dSi(rEi),(4)其中，qEi为i个单元的源强。式(4)就是常用的WSM，该方法虽然计算精度较高，但其在计算振动体外部不同场点r处的声压时，需要对所有单元计算Green函数的数值积分，导致整体的计算效率较低。一种更为简化的方法是用振动体内部若干个不同源强的等效源产生的声场，代替振动体辐射的声场，即将单元积分SiG(r,rEi)dSi(rEi)(面源)1022024 年 1 月直接简化为点源G(r,rEi)，得到在声学计算中广泛使用的ESM：p(r)=Ni=1q(rEi)G(r,rEi),(5)其中，rEi表示第i个等效源的位置，q(rEi)为第i个等效源的源强。ESM由于形式简单且无需积分，因此计算效率相较于WSM得到了很大提高。但该方法对物理模型过度简化，只有当等效源点位置和数目合适时才可达到较高的计算精度，并且从面源简化为点源的近似过程始终存在较大的积分误差，影响了计算精度。针对上述两种方法缺陷，为了尽可能减小积分近似误差，提高振动体外部辐射声场计算效率，本文采用波函数Ki(r,rEi)=SiG(r,rEi)dSi(rEi)来代替离散单元的数值积分。2 矩形常数单元的波函数构造2.1矩形常数单元的波函数构造在工程实际中，对于平板或圆柱等结构的外部辐射声场求解，可采用WSM设置规则的虚拟平面或圆柱面替代实际声场，这种虚拟曲面通常可以被划分为矩形单元；在近场声全息中，如果使用常规的平面矩形或正方形阵列，为了共形，也可将虚拟面划分为矩形单元。因此本文主要考虑离散单元为矩形时的波函数构造。在单元局部球坐标系下考虑一个矩形单元的积分，将该单元记为S，如图1所示。其中，单元的局部坐标系的原点位于单元质心，坐标系的x轴和y轴分别平行于矩形两边，将这两边的长度分别记为Lx和Ly。在单元区域内对Green函数的积分可用一个函数K(r,)代替：K(r,)=+Lx/2Lx/2+Ly/2Ly/2G(r,r)dS(r),(6)其中，(r,)表示场点在坐标系的位置，r为矩形单元中心到场点的距离，为对应俯仰角，为对应方位角，r为单元内部的场点。yLyxzrLxSrr图1矩形单元S 及其外域球面示意图Fig.1 Schematic diagram of rectangular elementS and its outer sphere由于式(6)是一个空间中的辐射声场，因此，可以利用Helmholtz方程在单元局部球坐标系下的解将函数K(r,)近似为如下的波函数形式：K(r,)=n=0+nm=nCnmh(2)n(kr)Ymn(,),n=0,1,N,m=n,n+1,n,(7)其中，Cnm为球面波展开系数，h(2)n(kr)为第二类n阶球Hankel 函数，k为波数，Ymn(,)为n阶m次的球谐函数，它决定了声压在不同角度的辐射属性，其表达式为Ymn(,)=(2n+1)(n m)!4(n+m)!Pmn(cos)eim,(8)其中，Pmn(cos)为连带勒让德函数。一旦确定了展开系数Cnm，矩形单元的辐射声场也就随之确定了。又由于单元积分和波函数K(r,)都满足Helmholtz方程和Sommerfeld辐射条件，因此根据微分方程的定解理论，只要两者在某一边界上等价，那么它们在整个空间中的声场分布都是等价的。为计算展开系数Cnm，在单元的外域任选一个半径为R的人工边界球面。则该人工边界上的声场分布为p(R,R,R)=+Lx/2Lx/2+Ly/2Ly/2G(R,r)dS(r),(9)第43卷 第1期贺佐潦霜等：矩形单元声场波函数构造及其应用103其中，R是人工边界球面上的点。由球谐函数正交归一性，将式(9)与式(7)联立，即可求出展开系数：Cnm=1h(2)n(kR)+Lx/2Lx/2+Ly/2Ly/2200G(R,r)Ymn(R,R)dRdxdy,(10)其 中，dR=sinRdRdR，Ymn(R,R)是Ymn(R,R)的复共轭。为了便于后文分析，令Cnm=Cnmh(2)n(kR)，再将式(10)代入式(7)，即可得到K(r,)=n=0+nm=nCnmh(2)n(kr)h(2)n(kR)Ymn(,).(11)式(11)可进一步表示为全场坐标变量形式：K(r,ri)=n=0+nm=nCnmh(2)n(k|r ri|)h(2)n(kR)Ymn(r ri)zi|r ri|),(12)其中，|r ri|为场点r与单元中心点ri之间的距离，zi为第i个单元在局域坐标系中z轴方向的单位向量。由于式(11)中人工边界球面是任意选取的，因此本文把该式构造的波函数K(r,)称为矩形域一般形式波函数。但从式(10)可以看出，在计算展开系数Cnm时需要求解四重积分，因此，当虚拟面离散的单元不一致且数目较多时，将花费大量时间计算展开系数。为缩短计算展开系数Cnm的时间，本文考虑用矩形单元的远场辐射声压解析解代替单元数值积分。2.2矩形域内推波函数的构造若将人工边界位于远场，即人工边界球面的半径RF远大于单元尺寸。此时，矩形单元的格林函数有近似解析表达式：G(rRF,r)eikRF4RFe(ikxx+ikyy),(13)其中，(kx,ky)=(ksinRFcosRF,ksinRFsinRF)。将式(13)代入矩形单元积分可得+Lx/2Lx/2+Ly/2Ly/2G(rRF,r)dS(r)=eikRF4RFLxLysinc(kxLx2)sinc(kyLy2),(14)其中，sinc(x)=sin(x)/x。将式(14)代入式(10)可得Cnm=1h(2)n(kRF)200LxLyeikRF4RFsinc(kxLx2)sinc(kyLy2)Ymn(RF,RF)dRF,(15)其中，dRF=sinRFdRFdRF，同理，将式(15)代入式(7)可得KF(r,)=n=0+nm=nCnmh(2)n(kr)h(2)n(kRF)Ymn(,).(16)式(16)可进一步表示为全场坐标变量形式：KF(r,ri)=n=0+nm=nCnmh(2)n(k|r ri|)h(2)n(kRF)Ymn(r ri)zi|r ri|).(17)由于式(16)选用的人工边界位于远场，因此本文把该式构造的波函数KF(r,)称为矩形域内推波函数。由式(15)可见，当选用远场人工边界时，仅需二重积分即可求出展开系数Cnm，而一般形式波函数在求解展开系数Cnm需要计算四重积分。当离散的单元数目较多时，可大幅度缩短计算展开系数的时间，从而提高声场计算效率。当离散的单元为正方形时，本文利用圆形单元替代正方形单元，进一步简化波函数。2.3圆形域内推波函数的构造当离散的单元为正方形单元S 时，可用一个中心点相同、面积相等的圆形域S 近似代替，则有SG(r,r)dS(r)SG(r,r)dS(r),如图2(a)所示。可以求得，边长为L的正方形单元的等效圆域半径为a=L2/。由于声场关于z轴旋转对称，因此仅在xz平面上用Legendre正交函数逼近实际指向性函数即可。如图2(b)所示，半径为R(R a)的球形人工边界上的实际指向性函数分布为(场点rp选在xz平面人工边界上)：p(R,)=20a0G(rp,r)dd,(18)其中，rp=(R,0)为人工边界xz平面上的场点，r=(,0,)为单元内部场点。1042024 年 1 月OSSxyrrpaLOxzR(a)xy?S?(b)xz平?S?图2正方形单元S 的等效圆域SFig.2 Equivalent circle regionS of square element S由于对称性，式(18)可进一步写为20a0G(rp,r)dd=20a0G(rp,r)dd.(19)人工边界xz平面上的声场可用Legendre正交多项式函数逼近展开：p(R,)=n=0CnPn(cos),(20)其中，Pn(cos)为Legendre函数。利用Legendre函数的正交性可得展开系数：Cn=(2n+1)00a0G(rp,r)Pn(cosR)sinRdRdd.(21)与式(13)类似，当选用一个远场人工边界RF时，此时圆形单元的Green函数有近似解析表达式：G(rRF,r)eikRF4RFe(ikxcos+ikysin).(22)将式(22)代入式(21)可得Cn=(2n+1)2a2eikRF0J1(kasinRF)2RFkasinRF Pn(cosRF)sinRFdRF,(23)其中，J1(kasinRF)为一阶贝塞尔函数。与前述同理，可得到圆形域外部辐射声场的等效波函数为KF(r,)=n=0Cnh(2)n(kr)h(2)n(kRF)Pn(cos).(24)式(24)可进一步表示为全场坐标变量形式：KF(r,ri)=n=0Cnh(2)n(k|r ri|)h(2)n(kRF)Pn(r ri)zi|r ri|).(25)由于式(24)的选用的人工边界位于远场，因此本文把该式构造的波函数KF(r,)称为圆形域内推波函数。可以看出，KF(r,)是在KF(r,)的基础上做了进一步的简化，构造了一个与m项无关的波函数。并且，从式(23)可以看出，展开系数Cnm由二重积分简化为一重积分Cn，进一步缩短计算声场的时间。2.4三种波函数的计算精度与计算效率对比由于本文构造的波函数与WSM类似，离散的单元位于真实边界回缩的虚拟边界上，其声压表达式不仅严格满足Helmholtz方程，且在真实边界上关于离散结点的声压分布还是严格满足形函数定义的高阶形函数19，因而相较于BEM所需的单元数更少，对于结构和声场较为简单的情况，有时每个波长仅需2个单元即可。例如，当划分的单元为正方形时，单元的边长L必须小于或等于半波长，即L 6/2，也就是kL 6。因此，本文比较3种波函数在kL=/4、/2、的情况下与单元直接积分的拟合效果。为了对比构造的3种波函数的计算精度和计算效率，算例选用位于xy平面上的正方形单元和圆形单元声源。正方形单元边长为L，圆形单元的半径为a=L2/，近场人工球面半径为L，远场人工球面半径为105L，单元中心点与人工球面中心点都位于坐标原点。计算半径为2L的球面声场，计算的场点数目为105个。算例中的积分均采用Gauss-Legendre 积分，其中，直接积分计算单元声场时，每个积分变量使用4个高斯点，计算波函数展开系数时，每个积分变量使用20个高斯点。值得指出的是，波函数的最高阶数n选择合适与否对计算结果的影响很大。图3为在kL=/4、/2、的情况下，不同阶数的矩形域内推波函数对应的重建声压相对误差，其中，相对误差由式(26)计算：Error=pn p2p2 100%,(26)其中，pn为第n阶波函数重建声压向量，p为直接积分声压向量。由图3可以看出，在kL分别为/4、/2、的情况下，当波函数的阶数n=1时，计算声压相对误差较大，但随着波函数阶数的增加，相对误差呈下降趋势；当n=4时，在不同情况下相对误差都小于0.5%，达到了足够的计算精度。进一步对比3种波函数在不同情况下的计算精度，表1为=/2时，3种波函数在kL分别为/4、第43卷 第1期贺佐潦霜等：矩形单元声场波函数构造及其应用105(a)kL=p/4(b)kL=p/2(c)kL=p1234567n00.51.01.52.0Error/%Error/%Error/%1234567n0123451234567n05101520图3矩形域内推波函数不同阶数对应的重建声压误差Fig.3Reconstructed sound pressure errors corresponding to different orders of pushwave function inrectangular domain表1kL分别为/4、/2、时，3种波函数与直接积分计算精度对比Table 1 Comparison of calculation accuracy between three wave functions and direct integrationat kL is/4、/2、030609012015018021024027030033000306090120150180210240270300330003060901201501802102402703003300030609012015018021024027030033000306090120150180210240270300330003060901201501802102402703003300030609012015018021024027030033001T100.5T100306090120150180210240270300330003060901201501802102402703003300?K?KF?KFkL=p4?kL=p2?kL=p?Error=0.0024%Error=0.0024%Error=0.1178%Error=0.0035%Error=0.0035%Error=0.3183%Error=0.0226%Error=0.0226%Error=1.0550%1T100.5T101T102T101T102T103T101T102T103T101T102T103T104T106T102T104T106T102T104T106T100.5T10(1062024 年 1 月/2、的情况下与单元直接积分的拟合图像与相对误差。当kL分别为/4、/2时，3种波函数在计算声场声压时，与直接积分的相对误差在0.5%以下，达到了足够的计算精度。当kL=时，相当于一个波长内只有两个单元，此时单元内的声场更为复杂。即便如此，矩形域一般形式波函数与内推波函数的计算精度达到了99.5%以上，即便是相对简化的圆形域内推波函数计算精度也能达到98.9%。表2记录了表1中kL=时波函数与直接积分消耗的计算时间。可以看出，在计算场点声压时，3种波函数的计算时间均远低于直接积分。但是，矩形域一般形式波函数在求展开系数消耗时间较多，导致总时间还要超过直接积分。因此，这种波函数只有在离散的矩形单元完全相同时，才具有较高的声场计算效率，否则其效率非常低。但采用远场内推之后，矩形域内推波函数求展开系数的时间大幅缩短，其计算速度约为直接积分的4 5倍，尤其是圆形域内推波函数计算声场的速度是直接积分的12 13倍。因此，即便离散单元的形状和大小不一致，内推波函数仍然具有较高的计算效率。表2当kL为时3种波函数与直接积分计算声压的CPU耗时Table 2 Three kinds of wave functions andCPU time consumption of direct integrationto calculate sound pressure at kL is 直接积分矩形域一般形式波函数K矩形域内推波函数KF圆形域内推波函数KF求展开系数时长/s4.68560.00140.0002求场点声压时长/s0.03690.00640.00640.0028总时长/s0.03694.69200.00780.0030为了进一步对比在不同数量场点下矩形域内推波函数、圆形域内推波函数与直接积分的计算效率，设置场点数量从106个增加至9106个，对比直接积分和波函数所消耗的时间，结果如图4所示。由图4可见，直接积分所消耗的计算时间约为矩形域内推波函数的5 6倍，为圆形域内推波函数的12 13倍。说明波函数的计算效率高于直接积分，且计算的场点数目越多，波函数的优势越明显。123456789?01234?/s?图4计算不同数量场点声压时直接积分和波函数的CPU耗时Fig.4 CPU time consumption of direct integra-tion and wave function when calculating soundpressure at different number of field points3圆形域内推波函数在声场计算中的应用3.1基于圆形域内推波函数的声场计算公式对于传统的ESM，空间任意r处的声压为p(r)=Ni=1qEiG(r,rEi).(27)空间任意r处质点速度为v(r)=1i0Ni=1qEiG(r,rEi),(28)其中，0为介质的密度，为角频率。在本算例中，将ESM中的G(r,rEi)采用3种波函数中精度相对最低、计算效率最高的圆形域内推波函数替代，考察波函数方法在近场声全息中的计算效果。与ESM类似，空间任意r处的声压和质点速度分别为p(r)=Ni=1WiKF(r,ri),(29)v(r)=1i0Ni=1WiKF(r,ri),(30)其中，Wi为第i个单元等效波函数的源强。3.2圆形域内推波函数在近场声全息中的应用在实际工程中，板、壳等连续分布的结构振动声源较为常见。因此，本文以四边无限大障板的简支正方形板为研究对象，对比分析圆形域内推波函数与ESM求解的重建声场。将简支板左下角的顶点作为坐标原点建立坐标系，仿真参数如表3所示。全息面位于简支板正上方0.03 m处，均匀分布第43卷 第1期贺佐潦霜等：矩形单元声场波函数构造及其应用1072020个测点。重建面在简支板正上方0.01 m处，尺寸与简支板大小相同，均匀分布4040个重建测点。虚源面尺寸与简支板大小相同，均匀划分为2020个正方形单元，单元长度为0.025 m，并用面积相同的圆形域代替正方形单元，在每一个单元中心布置一个等效源点，模型示意如图5所示。将虚源面分别置于重建面下方0.5倍、1倍、2倍单元长度位置，并在仿真中对全息面测量声压添加信噪比为30 dB的高斯白噪声，然后对比圆形域内推波函数与ESM重建声场的相对误差。表3仿真参数Table 3 Simulation parameters仿真参数数值简支板尺寸0.5 m 0.5 m 2 mm材料密度7.8 103kg/m3弹性模量2 1011Pa泊松比0.28空气密度1.29 kg/m3空气声速343 m/s简谐力幅值1 N激励点位置简支板中心-0.02-0.010.500.01z/m0.020.030.5y/mx/m00?图5仿真模型示意图Fig.5 Schematic diagram of simulation model图6为不同情况下ESM与圆形域内推波函数重建声压和振速的相对误差，其中，相对误差的定义为Error=E E2E2 100%,(31)其中，E 为重建面的解析声压或解析振速向量，E为算法重建结果向量。0100020003000?/Hz?/Hz?/Hz?/Hz?/Hz?/Hz0246810?/%?/%?/%?/%?/%?/%0100020003000051015202530350100020003000024681001000200030000246801000200030000510150100020003000051015?30 dB?30 dB?(a)0.5?(b)0.5?(c)1?(d)1?(e)2?(f)2?图6200 3000 Hz频率下ESM与圆形域内推波函数重建误差对比Fig.6 Comparison of reconstruction error between equivalent source method and circular domain extrap-olation function at 2003000 Hz frequency1082024 年 1 月由图6的计算结果可以看出，将虚源面置于重建面下方不同距离时，ESM与圆形域内推波函数的声压和振速重建误差均随着重建频率的增加而增大。但在整个计算频段，即便对全息面测量声压添加了信噪比为30 dB的高斯白噪声，波函数法的声压和振速重建误差也低于ESM，且在虚源面与重建面距离较小时优势更为明显。3.3圆形域内推波函数在声辐射计算中的应用对于无解析解的立方箱体结构，在立方箱体内部放置若干单位源强的单极子源，以这些单极子源在箱体表面单元中心点产生的质点速度作为该箱体的速度边界条件，这些单极子源在单元中心点处的辐射声压则作为标准声压，进一步验证圆形域内推波函数的有效性。以该立方箱体的中心点为坐标原点，边长为1 m，每面均匀划分为100个四边形单元，总共600个表面单元。计算时，选取每个表面单元中心点坐标乘以缩放系数进行同形缩放来获得等效源点和圆形单元波函数中心点所在位置坐标(计算时，缩放系数取为0.7)，圆形单元面积为表面四边形单元面积乘以缩放系数。立方箱体表面单元中心点与等效源点、圆形单元中心点如图7所示。y/mx/m-0.50.50.50z/m000.5-0.5-0.5?图7立方箱体表面单元中心点与等效源点、圆形单元中心点示意图Fig.7Schematic diagram of unit center point,equivalent source point and circular unit centerpoint on the surface of cubic box算例选用3个单位源强的单极子点源，分别布置在(0.3 m,0,0)、(0,0.3 m,0)、(0,0,0.3 m)处，声速及介质密度等参数与3.2节相同。本文对比了200 1200 Hz频率下ESM与圆形域内推波函数计算立方箱体表面声压相对误差。如图8所示，当分析频率等于或接近该模型对应的特征频率时，ESM与圆形域内推波函数的计算误差都较大；除特征频率外，两种方法的计算误差都很低，且随着频率的增加而增大。但是，在不同的频率下，圆形域内推波函数的计算误差始终低于ESM。2007001200f/Hz10-310-210-1100101102103?/%?图8200 1200 Hz频率下ESM与圆形域内推波函数计算声压相对误差Fig.8 Relative error of sound pressure calculatedby equivalent source method and push wave func-tion in circular domain at 2001200 Hz frequency4结论为了避免WSM在求解振动体外部辐射声场中复杂的积分计算，本文构造了一种波函数替代单元区域关于Green函数的积分。并以矩形常数单元为例，推导了替代单元积分的一般形式和内推形式波函数，以及当单元为正方形时的简化内推波函数，并对比了3种波函数与直接积分的计算精度和计算效率。最后，通过简支板声源和立方箱体辐射声源的数值仿真算例，验证了圆形域内推波函数与ESM在声场计算中的效果，主要结论如下：(1)在计算单个单元外部声场时，本文构造的3种波函数均与单元直接积分高度拟合，且矩形域一般形式和内推形式波函数的计算效率达到了直接积分的5 6倍，圆形域内推波函数的计算效率达到了直接积分的12 13倍。(2)在虚源面划分单元数目、回退距离相同的情况下，即便是3种波函数中精度相对较低的圆形域内推波函数，在无噪声和加入30 dB高斯白噪声时，重建简支板外部声压与振速误差仍低于ESM。第43卷 第1期贺佐潦霜等：矩形单元声场波函数构造及其应用109(3)对于内部放置若干单极子点源的立方箱体，在求解箱体表面声压时，当分析频率接近模型的特征频率，ESM与圆形域内推波函数的计算误差较大。除特征频率外，两种方法的计算精度较高，且在不同频率下，圆形域内推波函数的计算精度显著高于ESM。参考文献1 程昊,高煜,张永斌,等.振动体声学灵敏度分析的边界元法J.机械工程学报,2008,44(7):4551.Cheng Hao,Gao Yu,Zhang Yongbin,et al.Boundary ele-ment method for acoustic sensitivity analysis of vibratingstructuresJ.Journal of Mechanical Engineering,2008,44(7):4551.2 Liu Y.On the BEM for acoustic wave problemsJ.En-gineering Analysis with Boundary Elements,2019,107:5362.3 Tong Z,Zhang Y,Zhang Z,et al.Dynamic behavior andsound transmission analysis of a fluidstructure coupledsystem using the direct-BEM/FEMJ.Journal of Sound&Vibration,2007,299(3):645655.4 Lee J W,Chen J T,Nien C F.Indirect boundary ele-ment method combining extra fundamental solutions forsolving exterior acoustic problems with fictitious frequen-ciesJ.The Journal of the Acoustical Society of America,2019,145(5):31163132.5 Tadeu A,Antonio J.3D acoustic wave simulation usingBEM formulations:closed form integration of singularand hypersingular integralsJ.E