活塞是一种应用广泛的声源，其声辐射阻抗和辐射特性对于扬声器和水下换能器设计很有帮助[1]．常见的圆形活塞声源有两种特殊情况，即自由空间中的圆形活塞和无限大障板上的圆形活塞．长管管端活塞的辐射阻抗一定范围内居于这两种特殊情况之间[2]，其更加接近实际的声源设备．对于一个音箱来说，扬声器可以简化为一个振动活塞而箱体被当作有限的障板．在高频处，箱体的作用接近无限大障板；在低频处，箱体更接近有限长的管道，管道的存在使活塞只一端辐射声波．对于有限长管端圆形活塞的声辐射阻抗，很难通过瑞利积分得到其理论解．Beranek等[3]直接用半无限长开口圆管终端的辐射阻抗代替长管管端活塞的辐射阻抗．莫喜平等[4]用有限元软件ANSYS计算了无障板圆形活塞换能器的特性，研究了与无限障板活塞的区别．本研究认为这种无障板圆形活塞换能器即是一种管道较短的管端活塞声源．纪伟等[5]和Valacas等[6]分别研究了管侧圆形活塞与矩形活塞的情况．由于管道声学的发展和管状乐器的精细设计需求，因此已有研究大多集中在开口管道终端的辐射问题[7-9]．回顾开口管道终端的声辐射研究工作，Levine等[10]通过Wiener-Hopf方法从理论上得到了半无限长圆管的管端修正和反射系数，在此基础上可进一步解决带法兰圆管的声辐射问题[11]．对于形状复杂的管道，Dalmont等[12]通过边界元法和有限差分法计算了多种管道的管端修正值．Selamet等[13]用子域边界元法分析了长径比对管端修正的影响，认为长径比大于3后管端修正值与半无限长情况一致．目前对于有限长管端活塞问题研究较少，已有的文献大都直接或间接沿用Beranek等[3]的处理方式，即用半无限长开口管道辐射阻抗代替管端活塞的辐射阬；但是两个问题是否完全一致尚无文献进行分析，当管端活塞声源的管道较短时，其阻抗能否继续用半无限长管道阻抗等效须要进一步明确，有限长管道端部活塞的辐射特性和影响因素须要进行研究．1 理论公式1.1　半无限长开口圆管的终端辐射问题图1(a)为半无限长开口圆管终端辐射问题的简化模型，图中：a为管道半径；Ψ为速度势．图1(b)为管端活塞辐射模型，图中：l为管道长度；v为活塞振速．10.13245/j.hust.239536.F001图1两类问题简化模型图1(a)的管道中有一平面声波从左端向右端传播，在右端管口处存在向外辐射的声波及反射声波，当管道足够长(l→∞)，管中声压速度势表示为ψ=Aexp(ikz)+Bexp(-ikz), (1)式中：z为管道中的轴向位置，z的正向为声波沿管道传播方向；k为声波数；A和B为待定系数．在右端口的辐射阻抗体现为管内声波向外辐射的能力．在Su等[10]的研究基础上推导出了半无限长开口圆管终端的辐射阻抗理论解．该问题存在两个重要参数：反射系数R和管端修正值δ．反射系数可以描述管口反射的强弱，管端修正值可理解为辐射中的附加质量等效增加的管道长度，分别可表示为R=exp-2kaπ∫0kaarctan(-J1(x)/N1(x))x[(ka)2-x2]1/2dx; (2)      δa=1π∫0kalogπJ1(x)[(J1(x))2+(N1(x))2]1/2x[(ka)2-x2]1/2∙dx+1π∫0∞log1/(2I1(x)K1(x))x[x2+(ka)2]1/2dx, (3)式中：x为积分变量；J1(x)为第一类贝塞尔函数；N1(x)为第二类贝塞尔函数；I1(x)为第一类修正贝塞尔函数；K1(x)为第二类修正贝塞尔函数．设RM为辐射阻，XM为辐射抗，已知辐射阻抗为ZM=RM+iXM，则推导得到的归一化辐射阻ℜ0和归一化辐射抗X0分别为：R0=RMπa2ρ0c=1-|R|21+|R|2+2|R|cos(2kδ); (4)X0=XMπa2ρ0c=2|R|sin(2kδ)1+|R|2+2|R|cos(2kδ). (5)在有限长开口管道终端辐射问题数值计算中，Dalmont等[12]认为：无法从管道开口处直接得到开口管的辐射阻抗，须要在左端加振动活塞制造管中传播的平面波，并在左端施加吸收边界条件．文献[12-13]都用此方法先计算活塞面上的阻抗，然后用阻抗转移公式得到右端口处的辐射阻抗．然而管端活塞的辐射阻抗可直接用活塞表面的声压与法向振速来计算．1.2　管端活塞辐射问题如图1(b)所示，圆管右端装有一个刚性活塞沿法向振动，管左端封闭，整个管道为刚性，有v=v0exp(-iωt), (6)式中：v0为活塞振动速度幅值；ω为圆频率；t为时间．该问题与开口管道终端辐射相同的是，管壁上的刚性边界条件一致．在管道右端有所不同，对于开口管道终端的辐射问题，假定管道足够长且管中是平面声波，这表明在管道开口截面上各点的声压及法向振速相等．然而在活塞辐射问题中，虽然活塞上各点的法向振速相同，但活塞上各点的声压并不一定相等．推测只有当频率较低时，活塞辐射面上声压几乎相等，这时可近似为开口管道终端辐射问题．在数值仿真中，活塞被离散为Ne个有限单元，pi为单元声压，Si为单元的面积，则辐射阻抗为ZM1=∑i=1NepiSi/v．2 数值计算模型2.1　几何模型为了研究管端活塞声源的长径比、矩形活塞宽长比对辐射阻抗的影响，讨论常用的2种管端活塞声源，即圆管端圆形活塞声源和矩形管端活塞声源．为使所有活塞面积一致，将矩形活塞尺寸(长为L、宽为W)进行等效，满足关系a=WL/π，其中a为圆形活塞半径(矩形活塞的等效半径）．2.2　FE-FE-BE方法对于声源辐射问题，不少文献使用有限元软件进行仿真[4，14]．由于在流域中实施了人工截断，在截断边界处不满足远场辐射条件，因此阻抗和远场声压计算会存在误差．边界元法(BEM)是解决结构声辐射问题比较准确的方法，文献[12-13]也都对开口管道问题使用BEM求解．然而，为准确计算辐射阻抗，活塞上往往须要划分很细的单元，开口管道问题中管道内外都须要用边界单元．这造成边界元数目很大，计算效率低．本研究使用一种结构有限元耦合内域流体有限元、耦合外域流体边界元(FE-FE-BE)方法[15-16]对管端活塞辐射进行仿真．如图2所示，FE-FE-BE法在结构外部流体域中选取人工边界ГA，将流体域划分为内域流体Ω1和外域流体Ω，对结构和内域流体应用有限元法，对外域流体实施边界元法．传统的BEM应用在结构表面ГS上，FE-FE-BE法应用在ГA上．这种方法保持了传统BEM的优点，可以自动满足远场辐射条件．该方法不必在复杂的结构表面上划分边界单元，边界单元网格尺度也不必很细，只须能够描述声波长．因此，相比传统BEM可以减少边界元数目，提高计算效率．10.13245/j.hust.239536.F002图2FE-FE-BE法几何模型采用FE-FE-BE方法、人工截断的有限元法(FE-FE)计算了无限大障板上圆形活塞的辐射阻抗．图3为FE-FE-BE法计算模型，底部圆面代表刚性活塞，半球部分为内域流体有限元网格，人工边界为包裹内域流体的半球壳，底部流体边界施加刚性条件来模拟无限大障板．有限元法模型与图3相似，但采用了更大的流体域，外边界为人工截断并施加阻抗吸收边界．图4为无限大障板上圆形活塞归一化辐射阻抗，图中Z0为归一化辐射阻抗，有10.13245/j.hust.239536.F003图3无限大障板辐射阻抗计算模型10.13245/j.hust.239536.F004图4无限大障板上圆形活塞归一化辐射阻抗Z0=ℜ0+iX0．由图4可见：FE-FE-BE数值结果与理论值一致，最大误差小于2%，验证了FE-FE-BE方法计算活塞声源辐射问题的正确性．有限元法存在一定误差，尤其在低频和高频段误差较为明显．对比140个频率点的平均计算时间，FE-FE-BE法为32.2 s，有限元法为85.5 s．由于人工截断的吸收边界并不严格满足远场辐射条件，因此在低频处须要应用范围足够大的流体域才能得到更准确的结果，高频则须要网格尺度够小．相比之下，FE-FE-BE法在保证精度同时改善了计算效率．3 辐射阻抗影响因素3.1　开口管道和管端活塞声辐射问题的区别及长径比的影响分析共计算了圆管端活塞长径比l/a=0.2，1，5三种工况下的辐射阻抗，并且与半无限长开口管道终端的辐射阻抗理论值进行对比．图5结果表明：只有在低频下管端活塞辐射阻抗与开口管道终端辐射阻抗一致，中高频段二者出现一定差别．当ka0.5时，辐射阻二者几乎相等；当ka1时，两者的辐射抗相等．在低频段，声波长相比于活塞尺度足够长，活塞上各个位置的声压几乎相等，这与开口管道终端的平面波条件比较接近．随着频率增大，二者在管端的边界条件不能近似相等，所以在中高频段的辐射阻抗出现差别．10.13245/j.hust.239536.F005图5不同长径比圆管端活塞声源辐射阻抗图5中不同长径比下的结果对比得出：不同长径比下的归一化辐射阻抗在低频和高频段几乎一致，在中频段的曲线峰值附近表现出差别．当l/a＞1时，辐射阻抗基本不再改变．当频率很低时，声波长远大于管道长度，所以长径比变化对低频段辐射阻抗影响很小．同样地，高频声波长很小，辐射阻抗对长径比变化不敏感．总体上看：管道长径比对管端活塞辐射阻抗影响有限，当管道比较短时须要考虑这种影响．3.2　管端矩形活塞宽长比影响文献[17]研究了无限大障板上矩形活塞不同宽长比下的辐射阻抗，较少针对矩形管端活塞宽长比变化对辐射阻抗影响规律的研究．计算矩形活塞的宽长比W/L=0.10，0.25，1.00三种工况，a=WL/π，长径比l/a均为5，活塞形状与矩形管道截面形状一致．图6显示：矩形活塞宽长比对辐射阻抗大小影响显著，当频率比较低时，不同宽长比下辐射阻大小基本相同，但辐射抗明显变化．整体上看，宽长比增大使辐射阻抗曲线的峰值变小且峰值的位置向右偏移．在ka2范围内，宽长比为1即方形活塞的辐射阻抗最大，证明方形活塞比矩形活塞辐射效率高．10.13245/j.hust.239536.F006图6管端不同宽长比的矩形活塞的辐射阻抗4 长径比和矩形活塞宽长比对指向性的影响4.1　长径比对指向性影响为研究管道长径比对声源指向特性的影响，计算同一频率下3种长径比(l/a=0.2，1.0，5.0)管端活塞的远场声压．为覆盖图5显示的三个主要频段，低频段(阻抗上升段)、中频段(阻抗峰值附近)及高频段(辐射阻波动、辐射抗下降)，共设置5组计算频率(ka=0.5，1.0，2.0，3.0，5.0)．其中，在距离活塞中心10 m的圆周上均匀分布72个远声场点．图7为不同长径比的管端活塞声源指向性图谱，图中Lp为声压级．10.13245/j.hust.239536.F007图7不同长径比的管端活塞声源指向性图谱图7结果显示：当频率比较低时，长径比对指向特性影响很小．当频率增高时，活塞前方的指向性不随长径比而改变；活塞后方的指向性随长径比改变较大，且波瓣随长径比增大而增多．原因主要是不同长度管道会影响边界条件．另外，尽管长径比对前方指向性的影响很小，但声压强度在中频段因管道长径比而不同．在低频和高频处，不同长径比下活塞前方的声压强度几乎不变，而在ka=1处出现一定区别．这与图5中辐射阻抗与长径比的规律对应，即管道长径比对中间频段辐射阻抗的影响较明显．4.2　矩形活塞宽长比对指向性影响为研究管端矩形活塞宽长比对声源指向性的影响，分别计算了当W/L=0.10，0.25，0.50，1.00时的远场声压．图8为以活塞中心为原点的球面上的远场声压分布云图，图中：p为声压；Z轴垂直于活塞表面；Y轴与矩形活塞短边(W)平行；X轴与长边(L)平行；YOZ参考面为水平面；XOZ参考面为垂向面．当ka=5，W/L=0.10时，高声压区主要在水平面附近呈带状分布，声源前方区域声压强度最大．10.13245/j.hust.239536.F008图8远场球面上声压分布云图(ka=5，W/L=0.10)图9和图10分别为管端不同宽长比的矩形活塞声源及方形活塞声源的远场声压级随方位角变化的特征曲线，图中：Lp1为相对声压级；θ为角度．结果表明，同一W/L下的指向特性随着无因次频率(ka)增加而增强．当W=L时，水平面内和垂向面内的结果相同；当W与L不相等时，水平面内和垂向面内的结果有明显区别．当W/L比较小(W/L=0.10)时，水平面内的指向性很弱，曲线较为平坦；而垂向面内的指向性强，产生较多的波峰．随着W/L增加，同一ka下的水平指向特性有增强趋势，最终在W=L处水平和垂向的曲线趋于一致．矩形活塞声源设计中必须考虑活塞宽长比对指向特性的影响．只有当ka比较小并且W与L都远小于声波长时，宽长比对指向性的影响可以忽略．10.13245/j.hust.239536.F009图9矩形活塞声源的远场声压级随方位角变化特征10.13245/j.hust.239536.F010图10方形活塞声源的远场声压级随方位角变化特征(W/L=1.00)5 结语a．当ka0.5时，圆管端部活塞辐射阻抗几乎与半无限长开口圆管终端辐射阻抗完全相等，其他频段有一定差别．在高频段，将二者的辐射阻抗等效处理会带来较大误差．b．长径比变化主要影响辐射阻抗曲线峰值附近，对高频和低频的辐射阻抗影响很小．长径比大于1，管道长度对阻抗影响可近似忽略．管道长径比对辐射指向性作用主要体现在活塞后方．c．W/L增大使矩形活塞辐射阻抗曲线峰值减小，峰值位置右移．W/L对于矩形活塞声源空间指向性影响显著，水平面内和垂向面内的辐射特征差异较大，与矩形活塞短边垂直的垂向面内的指向性更强．