离心泵叶轮内非定常流动的动态模态分解分析

陈学炳， 张人会,2， 蒋利杰， 郭广强

(1. 兰州理工大学 能源与动力工程学院，兰州 730050； 2. 兰州理工大学 甘肃省流体机械及系统重点实验室，兰州 730050)

离心泵作为一种通用机械，被广泛应用于各个领域，对国民经济具有重要作用。复杂的非定常流动是造成离心泵内压力脉动、振动、噪声等不稳定流动现象的主要原因，这些不稳定流动往往会造成泵内较大的水力损失，甚至危及泵的安全稳定运行[1-4]。因此，对离心泵内非定常流场特性进行研究具有重要意义。目前，对于离心泵内复杂的非稳态流动，国内外研究人员已开展了大量的工作[5-10]。高波等[11]基于数值模拟方法对离心泵叶轮出口沿轴向变化的圆周面流域内各节点的压力分布特性进行分析，发现在不同工况下叶轮出口圆周方向的压力系数均在隔舌附近产生极值。Pedersen等[12]和Byskov等[13]分别采用试验与数值模拟的方法对离心泵叶轮内流动进行研究，并将数值模拟结果与PIV试验测试结果进行对比分析。Zhang等[14]采用延迟分离涡模型(delay detached eddy simulation, DDES)对离心泵内流场进行模拟，详细分析了叶片扫掠蜗壳隔舌时的非定常尾流及其在不同时刻的演变规律。郭荣等[15]对射流式离心泵过流部件诱发的流动噪声和流激噪声特性进行分析，发现叶轮和导叶的动静干涉是诱发射流式离心泵内场噪声的重要因素之一。离心泵内部复杂流场具有多尺度的流动特性，目前对其复杂瞬态流动的研究大多局限于对其压力脉动和泵内部整体流场的测量及分析，但对于离心泵叶轮内非定常流动尤其是相对流动中不同频率所对应的流场结构特征及其相应的物理机理方面的研究较少，若可以将其复杂流场分解为若干不同尺度的流动结构，将对叶轮内流动机理研究提供一定的参考。

动态模态分解(dynamic modal decomposition, DMD)是一种从非定常试验测量或数值模拟流动中提取动力学信息的数据驱动算法[16-20]。Schmid[21-22]首次提出DMD方法，并将其用于对高维复杂流场进行低维近似，同时识别出脉动主频及其对应的空间结构。Rowley等[23]采用DMD方法对横向射流中的非线性流动进行分析，进而提取出不同频率的流场模态特征。洪树立等[24]采用DMD方法将包含复杂时空信息的叶栅分离流流场进行解耦，获得了流场的主要频率及与之对应的空间结构。张人会等[25]对螺旋轴流式油气混输泵导叶内流场进行模态分析，系统研究了复杂流动动态模态的流场结构及其对应的物理机理。吴亚东等[26]采用DMD方法对轴流式压气机的旋转不稳定性进行分析，将其流场按照频率分解为成对出现的不同模态，且各阶模态从频率和空间周期上与不同的空间周向模态相对应。

目前的研究表明，DMD方法在流场分解方面取得了一定的进展，但基于DMD方法分析叶轮内相对流动的非定常特性仍需进一步研究。本文基于离心泵内流场的大涡模拟(large eddy simulation, LES)数值计算结果，采用DMD方法对其叶轮内部复杂非定常流场进行动态模态分解，分析其瞬态相对速度场在不同频率下的流场结构及其随时间演化的动态特性，以研究叶轮内复杂流场中包含的流动特征及其对应的物理机制。

1 数值计算

1.1 数值计算模型

本研究以一台中比转速单级单吸离心泵为研究对象，其基本设计参数为：流量Q=200 m3/h，扬程H=20 m，转速n=1 450 r/min，效率η=81%，比转速ns=132，叶轮进口直径D1=150 mm，叶轮外径D2=260 mm，叶轮出口宽度b=33 mm，叶片数Z=6，蜗壳基圆直径D3=262 mm，蜗壳出口直径D4=125 mm。数值计算域包括进口段、叶轮、蜗壳及出口段。采用六面体结构化网格对整个数值计算域进行离散，对近壁面等流动参数变化较大的区域进行局部加密，第一层近壁面网格距离壁面0.012 mm，叶片表面平均y+值约为1.5，其他区域平均y+值约为5。采用定常模拟方法对设计工况下不同网格数量的4种方案进行数值计算，如表1所示，以验证网格无关性。随着网格数量的增加泵扬程逐渐趋于稳定，综合考虑计算资源及数值计算精度，最后确定网格总数为650万，数值计算网格如图1所示。

表1 网格无关性测试Tab.1 Grid independence test

图1 数值计算网格Fig.1 Grid of numerical calculation

1.2 模型边界条件设置

采用ANSYS Fluent 16.0对离心泵内部复杂流动进行数值模拟。定常数值计算采用SSTk-ω湍流模型，速度进口及自由出流边界条件，壁面满足无滑移固壁条件。为更加精确地描述离心泵内部复杂流动，非定常计算采用大涡模拟进行数值求解，亚格子模型选择Smagorinsky-Lilly模型，动静交界面采用滑移网格处理。采用SIMPLEC算法耦合压力及速度场，收敛精度设置为10-5，动量方程离散格式为边界二阶中心差分格式。将定常计算的结果作为非定常数值计算的初始条件，计算时间步长设置为1.15×10-4s，为保证计算结果收敛，计算20个叶轮旋转周期。

1.3 数值模拟

采用LES对离心泵各个工况进行非定常数值计算，其外特性曲线如图2所示。离心泵试验测试系统如图3所示，其包括测试泵、进出口管路、动力控制部分及各测试仪器等；动力控制部分包括电动机及变频控制柜；进口管路上安装进口蝶阀及真空表，出口管路上安装有压力表、电磁流量计及出口蝶阀等，通过调节出口阀门的开度来调整泵的运行工况；根据进出口压差计算泵扬程，同时通过变频控制柜测试其功率，结合流量及扬程计算泵的效率。通过对比数值计算与试验结果可以看出，扬程与效率的大涡模拟计算结果略高于试验测量结果。在设计工况下，扬程计算误差为3.7%，且在各个工况中最大误差不超过5%。由于小流量工况下流动复杂，其数值计算效率的误差较大，随着流量的增加，计算误差越来越小，在设计工况下效率误差为1.8%。总体上数值计算结果与试验结果近似一致，扬程及效率曲线与试验吻合良好。

图2 数值模拟与试验的外特性对比Fig.2 Comparison of pump performance curves between numerical simulation and experiment

图3 离心泵测试系统Fig.3 Test system of the centrifugal pump

2 动态模态分解

DMD方法通过提取流动中的模态，从而准确描述流动结构，是对Koopman模态分解方法的近似表达，对于线性与非线性系统均可进行分析。将n个非定常流场快照表示为一个矢量集，即

(1)

式中，第i个时刻的流场快照表示为列向量vi，且任意2个快照之间的时间间隔均为Δt。假定流场vi+1可以通过流场vi的线性映射表示

vi+1=Avi

(2)

式中，A为高维流场的系统矩阵，其将时空流场沿时间维度平移Δt，因此矩阵A能够反映动态系统的时间演化特性。

随着采集快照数量的增多，基于Arnoldi算法[27]，最后一个流场快照可以表示为之前所有快照的线性组合

(3)

式中：aT=[a1,a2,…,an-1]；r为残差向量。根据式(4)可以得到

(4)

或其矩阵形式

(5)

(6)

通过对伴随矩阵C的特征分解，得到其特征值，它是系统矩阵A特征值的近似。

(7)

(8)

因此，矩阵A可近似为

(9)

Φj=Uyj

(10)

gj=Re{lg(λj)}/Δt

(11)

fj=Im{lg(λj)}/Δt

(12)

定义模态振幅α为

α=Y-1UHv1

(13)

式中：Y为特征向量；αi为第i个模态的振幅，代表了该模态对初始流场v1的贡献。对于标准DMD算法，其模态按照该振幅进行排序，本研究中根据模态范数表示的能量对各阶模态进行排序。

根据上述公式，任意时刻的流场可以被重构为

(14)

3 结果与分析

3.1 离心泵叶轮内非定常特性分析

对离心泵全流场进行LES数值计算，并对其叶轮内非定常特性进行分析。叶轮旋转轴为z轴，从进口方向看沿逆时针方向旋转。沿z轴正向分别取三个不同截面进行流场对比分析，如图4所示。S代表与叶轮后盖板的距离比，其中S=0.5表示叶轮中间截面。

图4 离心泵叶轮内各截面位置Fig.4 Position of each section in centrifugal pump impeller

图5为设计工况下叶轮不同轴垂面上相对速度分布，可以看出叶轮内不同轴垂面上速度分布存在较大差别。在S=0.1截面上叶轮内相对速度分布基本均匀，叶片工作面存在明显的低速区，叶片背面流速由进口到出口逐渐增大，且在叶轮出口存在较大的射流区域。随着截面逐渐靠近前盖板，叶片背面的速度梯度越来越大，其流动逐渐复杂，开始发生不同程度的流动分离；叶片工作面的低速区域逐渐减小，叶片背面进口附近流速越来越大，且叶轮出口的尾迹区域逐渐增大。在靠近叶轮前盖板的截面，即S=0.9，叶轮内流动较为复杂，速度梯度较大，叶片背面出现明显的流动分离结构。

图5 设计工况下叶轮内不同截面上相对速度分布Fig.5 Relative velocity distribution on different sections in impeller under the designed flow rate condition

图6为叶轮内S=0.9截面上不同时刻相对速度的分布情况。从图6中可以看出，设计工况下叶轮内大部分区域流动基本稳定，流线分布光滑且均匀。但在叶片背面存在尺度较小的周期性分离涡结构。在小流量工况下，即0.6Qd(Qd表示设计工况流量)，叶轮流道内流动复杂，存在周期性演变的大尺度旋涡。针对CH1通道进行分析，T/6～3T/6，通道内两个尺度较小的旋涡融合成一个大尺度旋涡，其充满整个CH1通道的进口，对其造成一定的堵塞，4T/6时大尺度旋涡逐渐开始分离，且随着时间的推移逐渐耗散。

图6 S=0.9截面上不同时刻相对速度分布Fig.6 Relative velocity distribution at different moments on S=0.9 section

对叶轮叶片背面随叶轮旋转的相对运动压力脉动变化进行监测分析，如图7(a)所示，在叶片背面设置2个监测点B1及B2，对其压力数据进行快速傅里叶变换(fast Fourier transform, FFT),得到其频域变化，如图7(b)、图7(c)所示。可以看出叶轮内主要频率与叶轮旋转频率(24.17 Hz)近似一致，其余峰值基本为叶轮转频的倍数，因此叶轮、蜗壳的动静干涉是造成离心泵叶轮内不稳定流动最主要的因素。

图7 叶片背面压力频谱图Fig.7 Pressure spectrum on the suction side of blade

尽管对离心泵叶轮内部各监测点数据的频域分析能够得到叶轮内主要频率及其激励源，但是其对于不同频率所对应的流场相干结构无法解耦识别，且对于不同频率及流场结构下的流动物理机制难以进行详细分析。因此有必要对其内瞬态流场特征进行分解，分析不同频率流场结构及机理。

3.2 设计工况下速度场DMD分析

根据3.1节分析，设计工况下模型泵叶轮内流动在靠近前盖板时较为复杂，因此选取S=0.9截面上流场进行DMD分解，进一步对其非定常流动特性进行深入分析。由于旋转机械内流场具有周期性，因此选取叶轮内单个流道CH1进行分析，快照数据则选取大涡模拟计算稳定后叶轮旋转一周内180个流场(每2个时间步长选择一个流场数据)的相对速度分布数据。

图8为对设计工况下速度场进行DMD分解得到的特征值分布，从图中可以看出除少数点外基本上所有点均处于单位圆附近。若某个特征值位于单位圆上，则说明该特征值对应的流场结构处于稳定状态，若特征值位于单位圆内，表示其对应的流场收敛，若特征值位于单位圆外，则其对应的流场逐渐发散，因此叶轮内速度场DMD模态基本为稳定模态或周期性模态。由于特征值为复数，且共轭成对出现，因此图8中各模态的特征值关于Im(λ)=0上下对称分布。图9为各阶模态的能量分布，选取前4阶能量较高的模态，其对应频域分别为0，22.20 Hz，73.65 Hz，127.42 Hz。

图8 1.0Qd速度场DMD特征值分布Fig.8 DMD eigenvalue distribution of velocity field at 1.0Qd

图9 1.0Qd速度场模态的能量分布Fig.9 Energy distribution of velocity modes at 1.0Qd

图10为采用DMD模态重构的叶轮旋转周期内速度场动态演化过程。叶轮内速度场由前4阶DMD模态叠加重构，与原始速度场对比可知前4阶模态基本可以完整表达出整个流场的流动结构，对局部较大尺度流场结构能够精确刻画，最后一个时刻的流场重构平均相对误差为1.37%。因此，DMD方法能够较好地实现流场的低维近似，较为精确地反映流场的主要流动特征，通过DMD模态分析方法能够清楚地分析离心泵叶轮内复杂流场的非定常特性。

图10 DMD重构流场与原始流场对比Fig.10 Comparison between reconstructed flow field and original flow field

设计工况下的速度场前4阶模态中，1阶模态为基本模态，其频率为0，代表了非定常流场的基本稳态结构，在各阶模态中所占能量最高，表征了由于叶轮叶片型线及流动冲角引起的流场特征。基本模态的空间分布如图11所示，其能够反映叶轮出口的高速射流区与低速尾迹区及叶片背面距离进口1/3位置的流动分离区域等特征，且随着时间的变化该模态特征基本稳定。

图11 1.0Qd速度场基本模态Fig.11 Basic mode of velocity field at 1.0Qd

设计工况下速度场2阶模态的频率与叶轮旋转频率一致，反映了由于叶轮旋转，蜗壳与叶轮的相互作用对叶轮内流场结构的影响。当叶片扫掠蜗壳隔舌时叶轮出口背压发生变化，使得叶轮流道内流动发生脉动，诱导叶片背面发生流动分离。2阶模态演化过程如图12所示，可以看出叶片背面不断产生正负交替的高能流体团，在其向下游发展的过程中逐渐拉伸并扩散，且在叶轮出口形成尾迹。在叶轮旋转周期内，叶片背面共产生一对正负交替且向下游发展的流体团，因此设计工况下2阶模态表征了叶轮内流动受蜗壳干扰在叶片背面产生流动分离的流场特征。

图12 1.0Qd速度场2阶模态Fig.12 Second order mode of velocity field at 1.0Qd

图13为速度场3阶及4阶模态的演化过程，其频率分别近似为叶轮旋转频率的3倍与5倍。由图可知，3阶与4阶模态的尺度依次减小且均小于2阶模态，同时流道内高能流体团数量相比2阶模态也是依次增加。叶片背面不断交替产生逐渐拉伸、脱落的高、低速流体团，且随着时间的推移，在向下游发展的过程中逐渐耗散。因此3阶～4阶模态为动静干扰模态的高次谐波行为，表征了由于叶轮旋转及蜗壳干扰诱导产生的叶片背面不稳定涡结构分离、脱落的流场特征。

图13 1.0Qd速度场3阶及4阶模态Fig.13 Third and fourth order modes of velocity field at 1.0Qd

为进一步分析叶轮内不同频率流动结构的稳定性，将前4阶模态对应系数(第i时刻的模态系数为(λj)i-1αj)随时间的变化规律进行对比分析，如图14所示。由于1阶模态为基本模态，其不随时间发生变化，因此图中未给出。由图可知2阶模态系数振幅较高且随着时间的推移基本稳定，进一步说明离心泵叶轮内流体脉动的最主要原因是蜗壳对叶轮的干扰。3阶与4阶模态尽管初始能量较高于2阶模态，但随着时间的演化，系数的振幅逐渐减小，因此其反映了流道内脱落的涡结构在逐渐的破碎并耗散。

图14 1.0Qd速度场DMD模态系数Fig.14 DMD mode coefficients for velocity field at 1.0Qd

3.3 小流量工况下速度场DMD分析

对小流量120 m3/h(0.6Qd)、S=0.9截面上相对速度场进行DMD分析，图15为0.6Qd速度场DMD特征值分布，可以看出特征值均位于单位圆附近或单位圆内，说明流动为周期性流动。根据各阶模态的能量大小同样选取前4阶能量较高的模态进行分析，如图16所示。

图15 0.6Qd速度场DMD特征值分布Fig.15 DMD eigenvalue distribution of velocity field at 0.6Qd

图16 0.6Qd速度场模态的能量分布Fig.16 Energy distribution of velocity modes at 0.6Qd

前4阶模态流场分布如图17所示，其中1阶模态频率为0，为基本模态，同设计工况DMD基本模态一样，表征了由于叶轮流道几何形状引起的稳态流场特征。2阶模态频率为19.71 Hz，与叶轮转频近似一致，反映了由于叶轮旋转及蜗壳干扰引起叶片工作面进口流动不稳定，且随着时间推移向背面扩散，同时在叶片背面发生流动分离，与进口不稳定流体融合引起流道堵塞的失速流场特征。3阶与4阶模态频率分别为2阶模态的2倍与3倍，流体团尺度依次减小，数量依次增加，反映了流道内不稳定流体团破碎、耗散特征。图18为0.6Qd速度场DMD模态的系数变化曲线，可以看出各阶模态脉动幅值依次减小，且对比图14中1.0Qd速度场DMD模态系数可知，小流量工况下叶轮内流体脉动幅值更大。

图17 0.6Qd速度场前4阶模态Fig.17 First four modes of velocity field at 0.6Qd

图18 0.6Qd速度场DMD模态系数Fig.18 DMD mode coefficients for velocity field at 0.6Qd

4 结 论

(1) 离心泵叶轮在靠近前盖板时内部流动复杂，叶片背面易发生流动分离。对叶轮内相对速度场进行动态模态分解后其特征值基本分布在单位圆上或圆内，流动具有周期性，且其流场结构可以分解为不同能量及频率的流动特征，包括基本模态特征、动静干扰模态及其高次谐波行为流场特征。

(2) DMD方法能够剥离出反映主要流场特征的主要模态，其基本模态的频率为零，代表了非定常流场的基本稳态结构，反映出由叶轮流道几何形状引起的稳态流场特征。2阶～4阶动态模态表征了由于叶轮旋转及动静干扰在叶轮流道内产生的流动分离、旋涡脱落及耗散流场特征。

(3) 采用主要模态叠加对叶轮旋转周期内速度场演化过程进行重构，前4阶模态基本可以精确地进行流场预测，重构速度场能够保留原始速度场中重要的动态信息，实现其低维近似，且最后一个时刻的平均相对预测误差为1.37%。因此，通过DMD方法，可以提取出若干能够反映叶轮流道内流动特征的主要流场结构，为离心泵内流场稳定性的提升提供一定的理论支撑。