双螺杆压缩机热力结构特性分析

支明宇，李志峰，2，万旭东，赵永强，2

(1.陕西理工大学 机械工程学院，陕西 汉中 723001；2.陕西省工业自动化重点实验室，陕西 汉中 723001)

引言

螺杆式压缩机在石油化工、食品处理加工、大型冷库等应用领域具有显著优势，现阶段受到越来越多的重视[1]。随着螺杆式压缩机应用领域的逐渐加大，螺杆转子在运转过程中的受热变形直接影响压缩机的安全运行，其自身的热力变形结构特性研究也被广大学者更加关注。

目前MUSTAFIN TN等[2]提出了利用数值计算法计算转子工作的温度场仿真模型，根据压缩机结构参数和工作条件仿真转子的温度场。DING H等[3]提出一种新式求压缩机共轭传热的混合时间尺度耦合方法并应用于双螺杆压缩机的建模。MARWAN C等[4]利用Modelica建立一个种新型双螺杆压缩机模型，并描述在高温下的压缩循环。SUN-SEOK B等[5]开发了基于MATLAB的转子型线程序来改变压缩机性能。CAO F等[6]提出工作腔压力分布数学模型，通过传感器实验验证准确性。HSIEH S H等[7]求解无内热瞬态热传导偏微分方程组，得到转子在周期对流边界条件下的瞬态热传导问题解。

国内在转子热应力方面的研究如魏静等[8]通过弱耦合方法实现双螺杆捏合机转子的热力耦合、流固耦合分析得到流场是导致转子应力变化的主要因素。孙瑾亭等[9]建立变螺距螺杆真空泵有限元模型，进行热分析仿真模拟。邓亦攀[10]建立压缩机热力学模型，对不同转速下的压缩机工作腔容积效率进行分析。赵永强等[11-12]通过对三螺杆泵进行热力耦合分析得到主、从动螺杆温度和热力变形分布规律；对主、从动螺杆进行静力学扭矩、流场压力两种载荷的耦合分析，得到啮合间隙是影响三螺杆泵容积效率的关键。王智博等[13]通过热力耦合仿真得到双螺杆泵转子的温度分布规律。操建平等[14]通过傅里叶级数拟合螺杆泵定子应力应变，并对常规螺杆泵及等壁厚螺杆泵进行对比。张远深等[15]通过建立数学模型对干式螺杆真空泵进行热力耦合，并对冷却系统进行设计优化。陈宗武等[16]对螺杆真空泵转子进行温度场分析，模拟推导出泵腔内冷配合间隙对转子形状的改进，优化降低了泵的返流量。何雪明等[17]探索并验证出一种替换动网格技术的新方法用于双螺杆压缩机的CFD计算仿真模拟。仝继钢等[18]建立微型高压压缩机工作腔热力学模型，得出压缩腔内容积效率变化规律。吴晓等[19]建立微小齿轮泵流场模型，对不同顶隙下的工作流量脉动率进行分析。张占东等[20]基于功率平衡原理建立流量调节策略。

以上研究对压缩机及螺杆泵等结构性能进行研究，虽得到螺杆转子温度场和应力场分布结果，但大多均以介质温度、热对流系数等为影响因素，未考虑在实际工作中由转速不同导致压缩机热力变形结构特性的影响，为准确模拟双螺杆压缩机在实际工作中的受力变形情况，本研究考虑以不同转速带来的性能差异，通过间接耦合方法将温度场与应力场联系在一起；在ANSYS软件中，通过求解静力方程，对比螺杆转子在静应力、稳态温度场和热力耦合场3种工况得到的应力、变形分布规律，为实现双螺杆压缩机阴、阳转子间隙合理分布，优化压缩机性能提供理论依据。

1 数学模型

1.1 耦合理论

对螺杆压缩机阴、阳转子采用间接耦合方法，先得到转子温度场分析结果，再通过施加约束及力载荷，其中力载荷主要为转子受到的转矩。通过数值插值技术施加到应力场，实现螺杆转子热力耦合。

1.2 温度场理论

稳态热力学平衡方程：

KT1=Q

(1)

式中，K—— 传导矩阵，包含导热系数、对流系数、辐射率等

T1—— 节点温度向量

Q—— 节点热流率向量

双螺杆压缩机工作时转子运送的介质具有进出口温差大、高流动性特点，易与转子发生强制热对流现象，即流体与固体表面之间发生热的交换过程。由牛顿冷却方程可知：

q=h(Ts-Tb)

(2)

式中，q—— 表示螺杆转子的热流密度

h—— 对流换热系数(或称给热系数)

Ts—— 固体表面温度，即螺杆转子表面温度

Tb—— 周围流体温度，即内部介质温度

1.3 力学理论

通过对经典力学理论的动力学通用方程的联合，建立双螺杆压缩机阴、阳转子力学模型的运动方程：

(3)

式中，M—— 表示螺杆转子质量矩阵

C—— 螺杆转子阻尼系数矩阵

K—— 螺杆转子刚度系数矩阵

f(t) —— 力矢量

1.4 热力学理论

根据温度场理论和应力学理论，对螺杆转子进行热力耦合分析得到螺杆转子的温度分布和应力、变形结果。考虑实际工作状态下，稳态热分析中温度不随时间变化，其中热对流换热系数与热流量也不随时间变化，流出系统的流量与进入系统的流量相同。在三维坐标系中稳态热控制分析方程：

(4)

在考虑边界条件和初值条件情况下，热力耦合数学模型：

(5)

式中，i—— 区分介质的不同项

Γ—— 研究区域的外边界

μ—— 热交换系数

K1—— 热传导矩阵

式（1）中：TCDab为两国的贸易结合度指数；Xab为a国对b国的出口额；Xa、Mb分别为a国和b国的出口总额和进口总额；Mw为世界进口总额；若TCDab＞1，说明a国和b国贸易联系紧密；若TCDab＜1，则a国和b国贸易联系松散。

n—— 区域的界法线向量

Ti—— 系统初始温度分布

Tα(t) —— 与系统相互作用的外界环境温度

2 有限元模型前处理

2.1 几何模型建立

螺杆转子端面型线由多段复杂曲线组成，需在现有数据点通过MATLAB生成局部端面型线基础上导入三维CAD软件，通过阵列命令将导入的型线以数据拟合方法得到完整端面型线，进行螺旋线扫略得到阴、阳转子三维模型如图1所示，表1给出了阴、阳转子主要参数。考虑网格精度、计算时间等综合因素影响，采用网格单元结构简单、计算时间短的四面体单元进行网格划分，在阴、阳转子啮合齿面进行网格局部加密设置。划分后的阴、阳转子网格数218 066，节点数374 686。

表1 阴、阳转子基本参数 mm

图1 阴、阳转子几何模型

考虑到双螺杆空压机的实际工况，结合齿轮、齿轮轴、蜗杆等耐磨且受力较大零件为参考依据，本研究所选用的材料为30CrMnTi合金，具有弯曲强度较大、耐磨性能高等优点。30CrMnTi合金的材料力学性能如表2所示。

表2 30CrMnTi合金材料参数

2.2 设置边界条件

螺杆转子在工作过程中，转子表面温度低于物料温度，与周围物料发生热对流，施加热对流换热系数，并设置环境温度进行稳态热分析。考虑压缩机实际工况，在转子进气端面轴向只保留绕Z轴转动和移动，在排气端面保留Z轴转动，以确保排气端具有最小轴向间隙并避免磨损，进气端具有最大轴向间隙使其扩张。为确保得到的数据更加具有收敛性，打开若弹簧命令。施加进气口温度25 ℃，排气口温度95 ℃，室温25 ℃；在螺杆转子啮合齿面施加热流密度为25 W/mm2，取对流换热系数100 W/(m2·k) 如图2所示。进行热力耦合分析时因阴、阳转子存在转矩作用，依次设置不同转速来确保转矩大小，具体流程如图3所示。

图2 稳态热分析边界条件

图3 热力耦合流程图

对阴、阳螺杆转子取3，4，5 mm不同尺寸进行网格划分并进行稳态热分析计算。由表3可见，不同尺寸网格的阴、阳转子网格数、节点数相差较大但最大形变量相差较小，验证了网格敏感性与收敛性，保证仿真结果正确性。

表3 不同网格尺寸差异对比量

3 模拟结果

依照双螺杆压缩机的功率130 kW，设定阴、阳转子啮合间隙为0.2 mm，阳转子转速依次设为2650，2950，3250，3550，3850 r/min的不同工况。通过数值模拟分析双螺杆压缩机阴、阳转子在静应力、稳态温度场和热力耦合场3种情况下的变形、应力变化规律。以阳转子转速2650 r/min为例，分析转子在稳态温度场下温度分布和变形、应力变化情况如下。

3.1 稳态温度场分析

螺杆转子稳态温度场下温度分布如图4所示，由于把室温加载在转子两端轴上，故阴、阳转子进气端的温度与轴段的温度相同，由于存在对流换热影响，转子温度从进气端向排气端逐渐升高，在排气端温度达到最大。图5为阴、阳转子稳态温度场轴向截面云图，可知变化趋势大致相似，以阳转子为例，在靠近进气口处进气温度设置与室温相同，因齿面可直接与热源接触，从热源吸收热量，而轴心热量温度均通过齿面传导至内部，吸收的热量相对于齿面来说较慢、几乎不产生热对流换热，温度散失较少，故靠近进气端转子内部温度大于螺杆表层温度；在出口处由于压力作用使大部分热量聚集于啮合齿表面，齿面从热源吸收热量的趋势大于散热的趋势，导致出气端螺杆表层温度大于转子内部温度。由图6知阴、阳转子在啮合齿面同一径向层温度不同，同一温度层大致以螺杆中部为基准向两端各呈现椭圆面形状分布。

图4 螺杆转子稳态温度场

图5 阴、阳转子稳态温度场轴向截面云图

图6 螺杆转子同温度层

阴、阳转子进、出口轴段采用一端游支，另一端固定。在螺杆的入口端保留Z轴方向的位移和绕Z轴转动，出口端仅保留绕Z轴转动，螺杆转子在稳态温度场下变形结果如图7所示。

图7 稳态温度场变形图

由于转子一端固定，一端游支导致螺杆发生热膨胀现象，最大变形量发生在螺杆进气端轴段。阳转子最大变形量为0.066422 mm，应力为744.92 MPa，螺杆工作时主要以螺旋齿面转动传送介质，其啮合面间最大变形量为0.060979 mm。阴转子最大形变量为0.04978 mm，发生在螺旋齿面上的最大形变量为0.046454 mm。由此知阴、阳转子螺旋齿面之间总变形极限0.107433 mm，小于阴、阳转子啮合间隙0.2 mm，符合安装理论要求且满足实际工况，可通过变形了解其啮合状态，确保安装间隙设置安全。

3.2 耦合场分析

如果只加约束，不加载螺杆转子受到的转矩载荷，螺杆转子耦合分析得到的只为温度场单独作用下热应力、热变形结果。通过温度模块与应力模块结合形式，在加上约束条件的同时，也加上转矩载荷，得到的将是应力-温度共同作用下的热力耦合结果，针对约束、载荷的不同，对双螺杆压缩机阴、阳转子进行有限元数值模拟分析，阴、阳转子在热力耦合下的变形如图8和图9所示。

图8 耦合下阳转子变形图

图9 耦合下阴转子变形图

在稳态温度场基础上施加转矩进行热力耦合，以主要工作面进行分析，阴、阳转子变形量在径向上由齿根向齿顶逐渐增大。由于螺杆转子依靠一组同步齿轮进行动力传递，转矩通过阳转子传递到阴转子，最大变形发生在阳转子进气轴端，形变量为0.068407 mm，比热变形增加0.001985 mm，应力749.16 MPa，阴转子最大形变量为0.065331 mm，与温度场相比增大0.015551 mm。

3.3 稳态静力学变形

设置螺杆转子在稳态下的变形情况，对比分析，阴转子变形量大于阳转子，最大变形0.015996 mm，阳转子最大变形0.0052431 mm。最大应力187.77 MPa，与稳态温度场、热力耦合场相比相差较大，可以看出温度场变化对螺杆转子之间的变形有较大影响。同时对转速在2650 r/min时3种不同工况下的阴、阳螺杆转子数值仿真结果进行对比，结果如表4所示。

表4 阴、阳转子不同场变形量对比

静应力下，阴转子变形量大于阳转子，在X，Y，Z三方向转子变形情况一致，均为在Y方向变形量最大，在Z方向变形量最小；稳态温度场下的最大应力、最大形变均为阳转子大于阴转子，且最大形变量均发生在轴向位置；热力耦合工况下，阴阳转子仿真数值略大于温度场，最大变形也均在Z轴方向，阴转子在Y方向变形小于X方向，相差0.012 μm，大小可忽略不计。由于温度的存在，导致温度场与耦合场的变形量及应力数值明显大于静力场下的应力、变形，可知温度是导致螺杆转子变形的主要因素，耦合场由于施加转矩作用得到数值比温度场更加准确可靠。

因阴、阳转子发生咬合、磨损的危险截面均发生在螺旋啮合齿面，故选取3种工况下不同转速的阴、阳转子螺旋齿面最大形变量进行对比分析，如图10和图11所示。

图10 阳转子螺旋齿面最大变形量

对比以上3种工况发现，静应力下的齿面变形较温度场及耦合场相差较大，且阳转子整体变形量大于阴转子；由于温度场的边界条件无转矩的施加，故最大变形量是一条水平直线；温度场与耦合场相比，由于耦合场存在转矩的施加，其齿面最大形变在所设转速范围内均大于温度场的单独作用，由图10和图11趋势可知，在某一转速时耦合场的最大变形量等于温度场的变形量，当双螺杆压缩机的转速大于该转速时，耦合场的最大变形量将小于温度场的最大变形量。因此，在实际工作中，随转速的增大，螺旋齿面最大变形量将逐步减小，且阴转子比阳转子受到的转速影响更加明显，所以在实际工况中，应多加注意阴转子的齿面变形对螺杆压缩机工作性能的影响。

图11 阴转子螺旋齿面最大变形量

4 结论

通过对双螺杆压缩机阴、阳转子进行有限元分析，利用ANSYS中温度场和应力场模块的相结合数值模拟，得到双螺杆压缩机的阴、阳转子在静力场、稳态温度场和热力耦合场条件下的热变形分布规律，并且在不同工况条件下进行了结构特性分析。

(1) 从稳态温度场计算结果可知，从进气端向排气端温度逐渐升高，进气端温度较低，排气端温度最高，且进气端齿面温度低于轴心温度，排气端齿面温度高于轴心温度；

(2) 在静力场结构下,阳转子变形为0.0052413 mm，阴转子变形为0.015996 mm，且主要变形发生在Y方向；耦合场下阳转子变形为0.068407 mm, 阴转子变形为0.065331 mm，静应力下的转子变形与耦合场相比差异明显，说明单独静应力分析与实际工况相差过大，应考虑多重因素的多场耦合模拟；

(3) 热变形是导致阴、阳转子变形的主要原因，阳转子的整体变形量在稳态温度场和热力耦合场均大于阴转子的整体变形量，但齿面间变形在耦合场下正好相反，在实际使用过程中应多加注意阴转子齿面变形对压缩机的安全运行影响；

(4) 选取的转速不同导致转矩差异，在稳态温度场和热力耦合情况下转子变形大致相同，最大变形量均发生在轴端位置。在实际工况下，通过选取适当的转速来调节转矩影响最大变形量，从而改善双螺杆压缩机工作性能。

因此，在进行双螺杆压缩机阴、阳转子配合间隙设计时必须充分考虑到转速对阴、阳转子变形的影响，可避免阴、阳转子之间磨损加剧甚至抱死现象，为双螺杆式压缩机阴、阳转子工作间隙合理选择提供了理论依据。