刘益才，颜鹏，马卫武，陈思明，陈丽新，李照龙

(中南大学 能源科学与工程学院，湖南 长沙，410083)

近年来国内外学者从多个角度讨论热声振荡的产生机理[1-6]。研究表明热声热机的声学特性和频率存在一定关系，在某一频率下声学特性表现最佳[7-8]。当某个力作用于结构时，会引起结构自身的振动，这是结构声产生的一个重要原因即振动激励[9-11]。根据结构声理论，结构声是热声热机所固有的声学特性，在发生结构振动时机体是以多个模态振动。因此，当机体振动频率与热声系统中气体介质的振动发生耦合时，势必会对系统效率造成影响。黄谦等[12-13]利用Solidworks和ANSYS软件对平板和回热器的热声机理从结构动力学的角度进行了模拟。结果显示，回热器的振动中存在多种模态，其中在低频段各点应力应变基本吻合，而在高频段各点差异逐渐增大；在一阶固有频率3 300 Hz附近，回热器内各点的应力和应变达到了最大。然而将平板或回热器从热声系统中独立出来孤立分析，会产生诸多误差。如回热器的套筒在压力载荷下发生了较大应变，而实际上，当回热器在谐振管的约束下发生应变的幅度很小；且就一般热声系统而言，系统工作时的谐振频率无法达到如此之高。本文作者对热声整机进行了三维建模，利用商业软件ANSYS对其进行了结构动力学分析。

1 模型概述

结构中的声波指的是某一物理量(力、位移、速度或能量)以结构为媒体向某一方向传播，结构本身作为一个整体并不运动[9]。在一维情况下，假设结构的质量元仅在x方向分布，且波动沿x方向以速度c传播，则结构声可用下式表示：

式中：t为时间；ω为角频率；A为常数，k=ω/c为波数，表示声波相位沿空间的变化规律。结构中声波的表达式与流体中声波的表达式一致，各物理量意义也相似。结构中的声波存在纵波、横波和弯曲波，纵波是指结构中微小质量单元的运动方向和声波的传播方向一致，表达式为：

式中：B为材料刚度的相关参数，与材料弹性模量E及泊松比σ有关，

固体中，质量微元振动的方向可以和声波传播方向垂直，即横波。横波导致质量微元发生形变导致切应力，故被称为横向切变波：

动力学基本控制方程如下：

其中：M为结构质量矩阵；a为节点加速度矢量；C为结构阻尼矩阵；c为节点速度矢量；K为结构刚度矩阵；u为节点位移矢量。当进行模态分析时，F(t)为零，包含流-固耦合自由振动，有无阻尼自由振动根据具体情况而定；当进行谐响应分析时，F(t)和u(t)都为谐函数。

取初始频率f0=100 Hz时，氮气的半波长1.7 m作为谐振管的长度。固体介质为锻造不锈钢，其密度ρs=7 800 kg/m3，弹性模量E=200 GN/m2，泊松比μ=0.28。板叠回热器长度为60 mm，外径为44 mm，平板厚度和板间距均为0.6 mm，谐振管厚度为2 mm。谐振管各封闭端和连接处利用Glue命令进行焊接。由于平板是直接焊接在谐振管的内表面，因此在平板和谐振管的相交面设置了DOF面约束，并在谐振管内表面施加了内压力。

根据热声热机结构特性和谐振管内部的受力情况，单元类型采用PLANE82和SOLID95。模型总共包含31 607个单元及70 889个节点。该模型基于以下几种假设：(1)谐振管内均匀分布气体介质，不考虑周围流体介质的黏性，较之在实际流体中板叠的振动幅度衰减较慢；(2)固体材料的物性参数如弹性模量、密度等不随外界因素的变化而变化；(3)设热声热机处于恒压环境中。

2 数值研究

2.1 结构静力分析

分别选取了0.3，0.5和0.7 MPa 3种充气压力分析热声热机内部的结构响应。静力分析可以模拟充气开始时到热声热机起振前系统的结构参数，此时假定载荷和响应是固定不变的，即载荷和结构的响应随时间的变化非常缓慢。图1所示为系统中最大应力和最大应变随充气压力的变化关系。随着压力的提升，应力应变基本呈现线性增大的趋势，且峰值点一一对应。

图2和3所示为0.5 MPa压力下热声热机不同角度的应力云图。分析结果表明系统内部的应力、应变分布是基本吻合的，总体看来系统中外表面受力微小，应力和应变主要分布在内表面。从图2可以看到：由于填充同性多孔类固体介质，在板叠回热器所处位置的谐振管外部应力和应变极其微弱，而之前单独对回热器的结构动力分析中，该位置谐振管的应变不容忽视[13]；谐振管其他位置所受应力较大，且同样分布均匀；管左边封闭端的不锈钢圆板上应力和应变呈同心圆状由外及内逐渐增大，并在圆板中心位置达到峰值。谐振管中段呈现了上凸下凹的弯曲状态，且整个谐振管径向应力由里至外逐渐加强。由图3可知：从宏观上看板叠整体各个面的中部位置几乎没有受到应力作用，应力载荷集中分布在四周边缘处，各个面应力等值图呈现出由里及外逐渐增大的特点；而对单独一块平板而言，应力载荷在轴向和径向均呈现出“中间小，两侧大”的特点，且随着平板表面积的增大，应力载荷的极大值与极小值之比也逐渐增大，这与之前的研究相符[12]。

图1 3种不同压力下热声热机最大应力和最大应变的分布Fig.1 Distributions of maximum stress point and strain point in thermoacoustic engine under three different pressures

图2 0.5 MPa时热声热机的应力云图Fig.2 Stress contour of thermoacoustic engine under 0.5 MPa

图3 0.5 MPa时回热器的应力云图Fig.3 Stress contour of regenerator under 0.5 MPa

固体介质应力和应变最终影响回热器和谐振管的声学特性，进一步影响了回热器中气体微团的运动。由结构声辐射理论，在热声热机的固体介质(如平板)中弯曲波所引起的振动幅度最大，而弯曲波的相速度与频率有关，因此都存在色散现象。任意形状的弯曲波经过一段时间的传播后，该声波的各个频率分量的相对相位的变化有可能导致其形状的变化，从而进一步影响到结构和流体之间的声耦合。因而分析系统的模态和谐响应尤为重要。

2.2 模态以及谐响应分析

模态分析可以确定热声热机固体介质的固有频率和振型，同时也是谐响应分析的出发点。根据热声系统情况，选取充气压力为0.5 MPa，相位角为0°。表1列出了系统1～8阶模态的频率值。在不同振动模态时，热声系统各位置的应变不同，说明静力下热声热机的应变起主要作用，而动力下回热器的应变是静力和波动两者的叠加。可以看到：热声系统的3阶和4阶模态下固有频率接近气体介质的振荡频率，因而可以重点观察3阶和4阶模态时热声系统的结构动力参数。

表1 热声热机各阶模态频率Table1 Frequencies of each modal in thermoacoustic engine

谐响应分析是用于确定线性结构在承受随时间按正弦(简谐)规律变化的载荷时稳态响应的一种算法。分析的目的是计算一定材料和结构在一定频率范围内的周期性载荷作用下的响应情况。该技术计算结构的稳态受迫振动时，不考虑结构发在激励开始时的瞬态振动。谐响应分析是一种线性分析，非线性特性被忽略。在热声热机内部重点研究位置选取了监测点来研究结构振动的响应情况，其中：点a设定在板叠内中间平板的上端面几何中心；点b位于该板的几何中心；点c位于该板下端面的几何中心；点d位于谐振管封闭端端面中心，也就是系统最大应力位置处；点e位于谐振管上部内表面。图4～8所示为上述点处应变随频率的变化曲线，其中对点d和点e处曲线进行了拟合。

图4 a点处应变随频率的变化Fig.4 Curves of strain versus frequency at point a

图5 b点处应变随频率的变化Fig.5 Curves of strain versus frequency at point b

图6 c点处应变随频率的变化Fig.6 Curves of strain versus frequency at point c

图4～6分别代表了板叠中间板x，y和z方向上应变随频率的变化关系。由图4～6可以看到：位于平板上、中和下部的3个监测点的应变变化高度一致：x和z方向的应变均非常微小，且随着频率的增大亦无十分明显的变化，只是当频率超过500 Hz以后，上端面和下端面分别出现了应变减小和应变增大的趋势，不过变化幅度不大；反观监测点y方向的应变在300 Hz以后明显上升，直至500 Hz附近达到峰值，随后随着频率的增大应变不断减小直至变为负值。

可以看到：随着结构频率的上升，平板的应变变化幅度总体是不断增大的，且各处受振动的影响相似，在垂直于平板的y方向应变最大而其他2个方向在低频段几乎没有发生应变，这也与实际情况相吻合。在频率为500 Hz时平板的应力应变达到极大值，且随着频率的进一步提高而增大，说明回热器作用是激发和维持高阶热声振荡。

图7和8所示为谐振管内表面监测点的应变随频率变化的拟合曲线。在热声系统这一有限结构内，结构声的传播在遇到内表面时会发生折射、反射和衍射等现象，其中边界的反射现象尤为突出，这些边界反射导致结构中的声波在稳态条件下仅在某些频率处存在较大的响应。图中曲线以一种类似正弦或余弦的形式发展，因而，应变在某些频率处达到极值，这些频率就是热声系统的本征频率。点d位于谐振管封闭端，由之前的应力云图可知此处谐振管的应变达到系统峰值，尤其是z方向也就是轴向应变最为强烈。该点处各个方向的应变在5阶、6阶模态即300 Hz附近同时达到最大值点。值得注意的是位于谐振管内表面的e点，在3阶、4阶模态即100 Hz附近应变幅度最大，这与气体振荡频率一致。如前述，高密度的气体介质的振动在固体介质中传播，会引起固体介质的结构振动，结构振动会反作用于气体介质。当固体介质在某一模态中达到其本征频率，且该频率与声激励频率相近时，声波与固体介质会在结构中产生以此阶模态为基础的相互作用，从而产生共振和耦合[14-15]。

图7 d点处应变随频率变化拟合曲线Fig.7 Fitting curves of strain versus frequency at point d

图8 e点处应变随频率变化拟合曲线Fig.8 Fitting curves of strain versus frequency at point e

3 结论

(1)利用商业软件ANSYS建立了热声热机模型，并对其进行了结构静力、模态和谐响应分析。在分别施加了0.3，0.5和0.7 MPa 3种压力载荷后得到了热声系统内部位移、应力和应变等参数的随压力和频率的变化关系。随着压力的提升，应力应变基本呈现线性增大的趋势，且峰值点一一对应。

(2)在板叠回热器位置处，谐振管外部应力和应变极其微弱，而其他位置所受应力较大；谐振管径向应力由里至外逐渐加强，而封闭端在圆板中心位置达到峰值并呈同心圆状由内及外逐渐减小；板叠处各个平板前后两端面几乎没有受到应力作用，上下端面发生了较大应变。

(3)随着结构频率的上升，平板的应变变化幅度总体是不断加强的，且各处受振动的影响相似，在垂直于平板的y方向应变最大而其他2个方向在低频段几乎没有发生应变，与实际情况相吻合。在频率为500 Hz时平板的应力应变达到极大值，且随着频率的进一步提高而增大，说明回热器激发和维持高阶热声振荡。

(4)热声系统的3阶和4阶模态下固有频率接近气体介质的振荡频率，根据检测点d的结果验证了高密度的气体介质的振动会在回热器中产生振动激励和声激励并作用于谐振管，使声波和固体介质产生共振和耦合。

[1]张晓青, 郭方中.热声热机中回热器的数值计算与分析[J].华中科技大学学报, 2001, 29(3): 87-89.ZHANG Xiao-qing, GUO Fang-zhong.Numerical calculation and analysis of regenerator in thermoacoustic heat engines[J].Journal of Huazhong University of Science and Technology,2001, 29(3): 87-89.

[2]HU Zhong-jun, LI Qing, LI Qiang, et al.A high frequency cascade thermoacoustic engine[J].Cryogenics, 2006, 46(11):771-777.

[3]陈熙, 李青, 李正宇, 等.热声热机网络模型演化与机理初探[J].低温工程, 2004(1): 27-31.CHEN Xi, LI Qing, LI Zheng-yu, et al.Network model of thermoacoustic engine and study on the mechanism[J].Cryogenics, 2004(1): 27-31.

[4]罗二仓, 戴巍, Ray R.交变流动回热器的热声功能和回热功能[J].工程热物理学报, 2006, 27(1): 1-4.LUO Er-cang, DAI Wei, Ray R.Thermoacoustic and recuperative functions of a cyclic-flow regenerator[J].Journal of engineering thermophysics, 2006, 27(1): 1-4.

[5]刘益才, 张明研, 黄谦, 等.热声热机板叠式回热器结构数值计算[J].中南大学学报: 自然科学版, 2010, 41(3): 1186-1189.LIU Yi-cai, ZHANG Ming-yan, HUANG Qian, et al.Numerical calculation of the structure of flat-plate regenerator in thermoacoustic engine[J].Journal of Central South University:Science and Technology, 2010, 41(3): 1186-1189.

[6]DAI Wei, YU Guo-yao, LUO Er-cang, et al.300 Hz thermoacoustically driven pulse tube cooler for temperature below 100 K[J].Applied Physics Letters, 2007, 90(2): 104-124.

[7]谢秀娟, 李青, 禹智斌.热声系统中的回热器声学特性的理论及实验研究[J].低温工程, 2004(6): 33-37.XIE Xiu-juan, LI Qing, YU Zhi-bing.Theoretical and experimental research of acoustic characteristics of the regenerator in thermoacoustic system[J].Cryogenics, 2004(6):33-37.

[8]刘益才, 杨智辉, 刘振利, 等.热声热机谐振管截止频率选择机理[J].中南大学学报: 自然科学版, 2006, 37(4): 759-762.LIU Yi-cai, YANG Zhi-hui, LIU Zhen-li, et al.Mechanism of cutoff frequency of resonator tube in thermoacoustic engine[J].Journal of Central South University: Science and Technology,2006, 37(4): 759-762.

[9]杜功焕, 朱哲民, 龚秀芬.声学基础[M].2版.南京: 南京大学出版社, 2001: 259-267.DU Gong-huan, ZHU Zhe-min, GONG Xiu-fang.Acoustics foundation[M].2nd ed.Nanjing: Nanjing University Press, 2001:259-267.

[10]马大猷, 沈豪.声学手册[M].修订版.北京: 科学出版社,2006: 164-176.MA Da-qiu, SHEN Hao.Acoustics handbook[M].Revised Edition.Beijing: Science Press, 2006: 164-176.

[11]何祚镛.结构振动与声辐射[M ].哈尔滨: 哈尔滨工程大学出版社, 2001: 67-98.HE Zuo-yong.Structure vibration and radiation[M].Harbin:Harbin Engineering University Press, 2001: 67-98.

[12]黄谦, 刘益才, 张明研, 等.热声回热器结构频率的数值研究[J].低温与超导, 2008, 36(11): 1-6.HUANG Qian, LIU Yi-cai, ZHANG Ming-yan, et al.Numerical simulation of structure frequency in thermoacoustic regenerator[J].Cryogenics and Superconductivity, 2008, 36(11):1-6.

[13]刘益才, 辛天龙, 陈思明, 等.热声板叠式回热器结构动力学分析[J].中南大学学报: 自然科学版, 2011, 42(11):3551-3555.LIU Yi-cai, XIN Tian-long, CHEN Si-ming, et al.Analysis of structural dynamics in thermoacoustic flat-plate regenerator[J].Journal of Central South University: Science and Technology,2011, 42(11): 3551-3555.

[14]刘益才, 李青, 杨智辉, 等.微型回热器等效直径与振荡频率的理论研究[J].低温工程, 2006, 151(3): 55-59.LIU Yi-cai, LI Qing, YANG Zhi-hui, et al.Theoretical research on equivalent diameter and oscillation frequency of miniature regenerator[J].Cryogenics, 2006, 151(3): 55-59.

[15]CHEN Yan-yan, LUO Er-cang, DAI Wei.Heat transfer characteristics of oscillating flow regenerator filled with circular tubes or parallel plates[J].Cryogenics, 2007, 47: 40-48.