基于动网格的三维合成射流器流场数值模拟

张 英，霍鹏飞，施坤林，李 曼，杨小会

(西安机电信息技术研究所，陕西 西安 710065)

基于动网格的三维合成射流器流场数值模拟

张 英，霍鹏飞，施坤林，李 曼，杨小会

(西安机电信息技术研究所，陕西西安710065)

针对激励出口流动模型不能够实现合成射流全场模拟这一问题，提出了基于动网格的三维合成射流器流场数值模拟方法。该方法引入了动网格的局部重构算法，通过在通常的流场数值模拟方法中插入编译的定义网格运动方式的动网格用户自定义函数程序，并设置与网格运动相关的动网格区域，实现了三维合成射流器在振动膜作周期性正弦运动下的流场求解。对比试验结果表明，采用该方法获取了射流器喷口特有的流动现象以及相应的喷口速度变化规律，说明该方法可行。

数值模拟；合成射流器；动网格；喷口速度

0 引言

合成射流器作为一种新型的流动控制器件，近年来在流动分离控制、传热传质、燃烧预混及矢量控制方面得到了广泛的应用。喷口射流速度作为合成射流器工作性能的主要指标之一，是合成射流器选型的重要依据，其大小对流体流动实施控制的效果会产生重要影响。因此，研究给定结构尺寸下的喷口射流速度与振动膜振动频率的对应关系对合成射流器结构的设计及优化有着重要意义。

近年来，有关合成射流器研究的工作主要集中在两方面：一是单独合成射流器的实验室研究，通过对合成射流器结构参数、电压参数、振动频率等进行分类，分析这些因素对合成射流器性能的影响规律并进行合成射流器的优化设计[1-4]；二是合成射流器在工程应用方面的研究，主要研究射流器对流动控制的影响效果，如机翼、叶栅分离流动控制[5-7]。

对合成射流流场的研究主要采用仿真方法，通过建立计算模型，获取给定结构尺寸下合成射流器流场的分布及射流速度。其中文献[8—9]采用激励出口流动模型，把按正弦变化规律的周期喷口平面速度分布函数作为来流的边界条件模拟射流器出口流场，这种流动模型虽然可以模拟激励器出口外流场的流动，但存在一些不足：一是不能完全真实地反映出口处流动，因为出口处流场是和时间空间变化的流场；二是出口处流动速度随时间变化规律不可能用统一的时间函数完全拟合；三是出口流动模型不能反映射流器腔体及喷口内的流动，不能够实现合成射流全场模拟，不利于研究合成射流的形成机理以及合成射流与来流的相互作用机理。文献[5,10]采用动网格方法，将振动膜每一时刻所处得位置传递给振动膜的网格节点实现实时的运动模拟，但是只是对二维合成射流器进行了数值模拟。文献[11]虽然实现了三维的动网格模拟，但是在模拟过程中对振动膜的振动作了简化，将正弦运动简化为往复运动。本文针对上述问题，提出了基于动网格的三维合成射流器流场数值模拟方法。

1 合成射流器数值模拟方法

合成射流器数值模拟方法将合成射流器振动膜、喷口内部通道及喷口处外部流场空间作为一个整体进行研究，不仅可以获取射流器结构参数及振动频率对射流器喷口速度分布的影响，还可以研究在射流过程中涡的产生及运动等流场变化，并且采用这种方法避免了实验室进行参数测量及流场观测存在的一些不利因素。

1.1 合成射流器结构

典型压电式合成射流器由腔体、喷口、压电陶瓷振动膜组成，如图1所示[12]。压电陶瓷振动膜在周期性电信号的作用下产生振动，改变腔体容积，使腔体内气体经由喷口被周期性的吸入和喷出。在振动的一个周期内，吸入和喷出的气体质量相等，净质量为零，因此是零质量射流。吸入和喷出的过程只是在向周围环境流体注入能量，并在喷口附近形成漩涡，这对漩涡随时间和空间发生变化并与周围流体发生耦合作用，从而改变喷口处流场。对飞行器而言，通过调整合成射流器安装的位置、开口方向和喷口速度大小，就可以控制飞行器表面气流分离区的位置及长度，使飞行器表面局部气流流动特性发生变化，从而影响表面的压力分布产生相应的气动力。

1.2 流场数值模拟方法

流场数值模拟步骤如下[13]：

1)建立模型，对流场进行网格划分并导入求解器；

2)设置边界条件；

3)设置湍流模型、求解器参数及收敛标准；

4)流场初始化并求解至收敛。

1.3 动网格算法

动网格算法用来计算内部网格节点的调节，动网格有三种算法，分别是铺层(Laying)、弹性光顺(Spring Smoothing)和局部重构(Local Remeshing)[14]，这三种算法可以单独使用也可以同时使用，同一个问题不同的网格区域可以使用不同的动网格算法。

铺层算法适合于四边形、六面体或三棱柱网格，铺层过程包含了边界上网格的生成和消失，这种算法会根据计算区域的扩张和收缩来相应地生成网格或消除网格，通常适用于边界有移动的情况。弹性光顺算法适用于三角形、四面体网格，网格之间的节点类似弹簧的连接，可以被压缩或拉伸，但不存在节点的生成和消除，适用于边界变形较小和运动幅度不大的情况，如果变形较大或运动幅度太大，会导致网格扭曲变形过大。局部重构算法仅适用于三角形和四面体网格，同弹性光顺法不同的是，它适合于大变形或大位移的情况，当网格的扭曲度很大时，局部网格节点和体网格就会增加或消除，即节点的数量和连接属性发生变化，通常局部重构算法和弹性光顺法联合使用。

局部重构算法的主要特点是：

1)当局部网格(体网格和面网格)的扭曲率或尺寸超过规定的范围时，对网格进行重构。

2)网格之间的连接属性发生改变，即节点的数量和连接关系都在改变。

3)适用于三角形和四面体网格。

4)适合于大变形或大位移情况。

2 基于动网格的合成射流器流场数值模拟

为了实现合成射流器全场模拟，本文引入动网格的局部重构算法，在通常流场数值模拟方法中插入了编写定义网格运动方式的动网格用户自定义函数程序(User-Defined-Function，UDF)，以及对该程序进行编译，并设置与网格运动边界相关的动网格区域两个步骤。

基于动网格技术的合成射流器流场数值模拟步骤如下：

1)建立模型，对流场进行网格划分并导入求解器；

2)编写动网格UDF程序定义网格运动方式；

3)对UDF程序进行编译并设置动网格区域；

4)设置边界条件；

5)设置湍流模型、求解器参数及收敛标准；

6)流场初始化并求解至收敛。

2.1 编写动网格UDF定义网格运动方式

本文采用局部重构算法来更新合成射流器内部网格节点随时间改变的位置。采用网格重构，需要编写用户自定义函数UDF(User Defined Function)程序，计算每个时间步长内振动膜上网格节点及腔体内网格节点在三个坐标方向上的位移，下一时刻用移动后的网格节点坐标对流场网格进行更新。

当合成射流器振动膜随着计算时间的推进按照指定的规律上下振动时，会引起振动膜附近计算区域内部网格单元发生变形，算法会标记出扭曲率超过规定扭曲率的所有网格并删除，求解器使用重构算法在原网格单元的位置上重新生成新的网格，用于计算振动膜在新的位置处时合成射流器整个流场。

2.2 对UDF进行编译并设置动网格区域

将编写好的动网格程序UDF进行编译并加载进求解器，在动网格参数设置中指定网格的最大长度尺度、最小长度尺度以及最大歪扭度等参数，并根据需要选择是否需要进行面网格重构。设置动网格区域相关参数时，需选取区域名称及运动类型，输入相邻网格的高度等参数。

3 算例及分析

3.1 建立计算模型，对流场进行网格划分并导入求解器

计算模型由合成射流器及其外部流场空间区域两部分组成，如图2所示。合成射流器包含合成射流器振动膜、腔体以及射流器喷口。合成射流器振动膜为一个边长L=8 mm的正方形薄片，喷口是一个0.1 mm×0.1 mm×0.18 mm长方体空间，喷口外部流场空间区域是一个28 mm×8 mm×20 mm的长方体。图3为合成射流器模型流场边界。

在射流器喷口喉道内部，由于不涉及网格重构，采用六面体结构化网格；在射流器腔体，振动膜的振动会引起上下两部分容积发生变化，采用与网格重构适应的四面体非结构网格。建模过程中对网格的疏密程度和网格节点之间的距离按照以下方法进行了设置。

1)采用尺寸函数法，在振动膜附近网格较密，远离振动膜的地方网格尺寸逐渐变大。

2)振动膜附近网格尺寸与沿振动膜长度每个节点的最大位移及时间步长有关。

网格划分后合成射流器内部网格如图4所示，图5是喷口内部及其附件壁面网格放大图。

3.2 动网格定义方法

模拟振动膜的周期性振动是通过移动振动膜上每一个网格节点来实现。通过编写用户自定义函数UDF程序来定义振动膜上网格节点的移动方式。本文定义振动膜按照正弦方式振动，振动函数为f=sin (ωπt)，其中ω分别取2，20，200，对应周期T为0.1 s，0.01 s，0.001 s。

已知振动膜上所有网格节点在某一时刻的坐标，在经历一个时间步长后，这些节点按照指定的正弦运动方式各自产生一个位移量，由于振动膜的运动是三维的，需要分别指定位移量在xyz三个方向上的分量。下一时刻，振动膜的新坐标通过运动之前振动膜的坐标分量和三个方向的位移量之和来进行更新，如此反复直到完成需要计算的多个周期直至流场收敛。

由于振动膜假设无厚度，并且处于合成射流器腔体中，振动膜上的网格节点被紧邻振动膜两侧腔体中的网格共用，因此当振动膜上网格节点运动时必然会引起两侧的网格变形，如果在一个时间步长后，节点的移动引起两侧网格过大的压缩和拉伸，产生歪扭度很大甚至是负体积的网格，会导致计算无法进行。因此对振动膜两侧的网格通过网格重构法进行了网格的重新划分以保证计算的顺利进行。

通过振动膜的振动和其两侧腔体内部网格的重构，就可以随时间更新振动膜节点坐标来模拟振动膜在腔体中的不同位置，实现了振动膜运动带来的腔体容积以及喷口射流速度随时间的变化。改变运动函数中的相关参数，就可以模拟出振动膜在不同参数下的喷口射流速度及合成射流器流场分布。

3.3 设置动网格区域

动网格区域可以是边界区域也可以是流场中的网格单元区域，定义了网格节点的运动方式后，求解

器自动根据运动区域的运动来调整内部节点以再生体网格。本文将振动膜设置为运动区域，运动规律按照用户自定义函数UDF定义的方式运动。

3.4 设置边界条件

本文将喷口外部空间边界、合成射流器腔体壁面及振动膜设置为壁面边界条件，对射流器喷口壁面不进行设置。无来流边界条件，流体温度设置为288 K，标准大气压P=101 325 Pa。

3.5 设置湍流模型、求解器参数及收敛标准

选用k-ε两方程湍流模型，采用耦合隐式求解器，速度与压力耦合采用同位网格上的SIMPLEC算法，对控制方程各量及湍动能均采用二阶差分格式离散。求解过程中应用多重网格技术加快收敛速度，以节省计算时间。设置收敛残差标准为两次迭代计算的绝对误差0.001。

3.6 流场初始化并求解至收敛

初始化流场选用X、Y、Z向的速度为0进行初始化，求解了振动膜振动周期T分别为0.1 s，0.01 s及0.001 s的合成射流器流场，时间步长Δt随周期不同而变化，分别为Δt=0.005，0.000 5，0.000 05。在合成射流器出口中心线上布置了四个观测点p1、p2、p3、p4。这四个点的位置如下：

p1:在喷口处中心位置，坐标(0，0.18，0)；

p2:在喷口中心线下游，坐标(0，0.23，0)；

p3:在喷口中心线下游，坐标(0，0.28，0)；

p4:在喷口中心线下游，坐标(0，0.33，0)；

3.7 计算结果

通过提取合成射流器喷口处流场速度矢量图，见图6，可以看出，与图7中的合成射流器实验纹影图[12]观察到的现象一致，在喷口两侧形成一对漩涡，这对涡关于射流喷口中心线对称，随着振动膜的周期性振动，涡对向射流器喷口下游方向运动，并在运动的过程中同喷口外部空间流场相互作用，影响喷口附件及外部空间流场参数分布。

图8、图9分别为振动膜向上和向下振动时合成射流器内部及喷口处流场矢量示意图，通过对比可以看出，两种情况下喷口内部气流矢量正好相反；进一步监测喷口出口截面上的质量流量，发现两种

情况下的气流质量流量大小相同，说明一个周期内，吸入和喷出的气流净质量为零，是零质量射流。

图10、图11、图12为射流器喷口中心线上四个监测点p1、p2、p3、p4的速度大小随振动膜振动周期的变化规律。

从图中可看出，1)对同一监测点，随着振动膜振动频率的增加，气流射流速度明显增大，当振动周期为T=0.001 s时，喷口中心处的最大射流速度可以达到60 m/s；振动频率为T=0.01 s和T=0.1 s时，射流速度分别为5 m/s和0.6 m/s，基本上射流速度和频率成正比。

2)对同一振动频率，不同监测点的射流速度不同，随着监测点距喷口处距离的增加，该监测点处射流方向的速度降低。

3)对比喷口中心处的射流速度，喷口中心线下游距喷口处0.5 mm监测点的射流速度快速下降，说明气流离开喷口后，在很短的空间内与周围流体发生能量交换，自身耗散的很快，降低了其动能。但是在距喷口距离0.5 mm后，随着距离的增加，气流的速度没有明显的降低。

4 结论

本文提出了基于动网格的三维合成射流器流场数值模拟方法。该方法引入动网格的局部重构算法，在通常流场数值模拟方法中插入了编写定义网格运动方式的动网格用户自定义函数程序，并对该程序进行编译；同时在求解器中设置了与网格运动边界相关的动网格区域。对比试验结果表明，采用这种方法模拟出了实验室测试中所观察到的射流器喷口漩涡涡对这一流动现象，同时验证了喷口射流速度随振动膜振动频率的变化规律与预测的射流速度正弦运动规律一致。证明该方法可以弥补出口流动模型无法获取腔体内部流场随振动膜振动实时变化的不足，能够有效实现合成射流全流场的数值模拟。

[1]罗振兵，夏智勋，王志吉.激励器结构参数对合成射流影响的实验研究[J].国防科技大学学报，2004，26(6)：9-14.

[2]李金波，何立东，伍伟.压电式射流激励器工作特性的实验研究[J].压电与声光，2009，31(6)：843-845.

[3]沈丹东，马炳和，邓进军，等.压电式微型合成射流器结构参数优化设计[J].航空学报，2011，32(9)：1755-1760.

[4]徐惊雷，沙江，林春峰，等.不同狭缝厚度的零质量射流PIV实验研究[J].推进技术，2007，28(2)：163-166.

[5]陈建辉，谢永慧，申仲旸，等.采用动网格技术的合成射流控制叶栅流动分离数值研究[J].中国电机工程学报，2013， 33(23)：126-133.

[6]郝礼书，乔志德.合成射流用于翼型分离流动控制的研究[J].西北工业大学学报，2006，8(4)：528-531.

[7]刘艳明，伍耐明，董金钟，等.合成射流作用机理及其对共轴射流掺混得影响[J].北京航空航天大学，2007，33(1)：5-9.

[8]何高让，汪亮. 微射流作动器出口流场数值分析[J]. 空气动力学学报, 2000,18(4):395-400.

[9]罗小兵，李志信，过增元. 不可压缩合成喷流场的数值模拟[J]. 工程热物理学报, 2001, 22(S): 56-58.

[10]唐上钦，黄长强.合成射流数值模拟方法对比[J].空军工程大学学报，2013，14(2)：1-4.

[11]马晓永，张征宇，李聪健.面向二维零质量射流的动网格数值模拟方法[J].动力学与控制学报，2007，5(3)：278-281.

[12]Smith B L, Glezer A. The formation and evolution of synthetic jets[J]. Phys. Fluids,1998,10(9):2281-2297.

[13]王福军.计算流体动力学分析[M].北京：清华大学出版社，2004.

[14]隋洪涛，李鹏飞，马世武，等.精通CFD动网格工程仿真与案例实践[M].北京：人民邮电出版社，2013.

FlowFieldNumericalSimulationof3DSyntheticJetBasedonMovingMesh

ZHANG Ying, HUO Pengfei, SHI Kunlin, LI Man, YANG Xiaohui

(Xi’an Institute of Electromechinal Information Technology, Xi’an 710065, China)

Aimed at the problem that orifice flow model cannot be used to simulate synthetic jet whole flow field,3d synthetic jet flow field numerical simulation method based on dynamic mesh was proposed. The method used the local reconstruction algorithm of dynamic grid to solve the 3d synthetic jet flow field under the condition of diaphragm vibrating at periodic sine frequency.To realize this, the compiled user defined function(UDF) describing grid movement way awas inserted into the usual flow field and dynamic zone associated with the grid movement was set up.By compared with the experiment result, the method could be used to obtain special flow phenomena of orifice jet and the corresponding velocity variation law,which demonstrates that this method is feasible.

numerical simulation; synthetic jet; dynamic mesh; orifice velocity

2017-04-13

机电动态控制重点实验室基金资助(9140C360204150C36165)

张英(1980—)，男，陕西富县人，硕士，工程师，研究方向：气动设计。E-mail:583564597@qq.com。

V211

A

1008-1194(2017)05-0052-06