空气阀水力瞬变数学模型

杨开林

(流域水循环模拟与调控国家重点实验室 中国水利水电科学研究院，北京 100038)

1 研究背景

空气阀不仅可用于管道充水和放空过程中自动进排气防止管道水击破坏，而且可用于泵站事故断电水力瞬变过程中自动进排气防止液体汽化和液柱弥合水击破坏。为了兼顾两者的要求，一般需要通过水力瞬变计算分析确定空气阀的布置位置、类型和孔径。目前一般采用国际著名瞬变流专家Wylie和Streeter[1]的数学模型和算法求解空气阀水力瞬变，该数学模型成立的基本假设是：

(1)空气等熵地流进流出空气阀，即忽略气体在空气阀内的能量损失与水体和阀体的热交换；

(2)管内气体的变化遵守等温规律，认为管内空气质量通常很小，气体的温度等于液体的温度；

(3)进入管道的空气留在它可以排出的阀附近，可以采用常规特征线方法确定空气阀流量与输水管水压和流量的关系；

(4)液体表面的高度基本不变，而空气的体积和管段里的液体体积相比很小，即忽略空气阀及其配套检修阀和连接管高度的影响。

Lee等[2]假定管内气体的变化遵循多变规律，将空气阀进排气数学模型描述为多变指数k的函数，然后通过数模研究了空气阀进气流量系数Cin和排气流量系数Cout对水击压力的影响，结果表明采用较大的Cin值有利于减小输水管负压幅值，而采用较小的Cout值有利于减小排气末了液柱弥合正水击幅值。实际上，气体等熵流动就表示气体的变化遵循多变规律。鉴于Wylie和Streeter空气阀水力瞬变数学模型求解算法程序的复杂性，杨开林[3]提出用差分代替微分来解决空气阀质量流量对气压的导数不连续性的难题，然后利用牛顿-雷伏生方法求解空气阀的水力瞬变。熊水应等[4]总结了多年来液柱(水柱)分离与液柱弥合水锤综合防护的设计经验，结果表明：因空气阀进气产生的液柱分离的空腔长度常常比因液体汽化产生的大得多，由此产生的液柱弥合水锤升压有时前者比后者高。Lingireddy等[5]研究了空气阀的口径对管道压力防护的影响，认为增大进气口径和减小出气口径可以降低水柱弥合时所产生的二次水锤压力。Fontana等[6]分别在管道末端和中间高点位置布置了排气阀，然后测试了排气阀在充水过程中的排气特性，发现对于排气阀位于管道末端的情况，当充水水柱撞击排气孔口时会引发剧烈的压力波动；对于排气阀位于管道中部高点时，水流撞击排气孔口时引发的压力波动较微弱。刘梅清等[7]研究了空气阀作为事故停泵水力瞬变防护的可行性，计算结果表明，在特定条件下空气阀也能较好地抑制液柱分离的发生并有效地进行水锤防护。杨开林等[8-9]的研究表明采用较大的空气阀进排气口径对防水击危害常常是不利的，在进气过程中空气阀底部的真空度随阀口径减小而增加，但是在未设空气阀的位置可能发生水压随空气阀孔径的增加而下降到汽化压力的风险，另外在高速排气结束时的液柱弥合冲击水压随阀口径的减小而下降。胡建永等[10]从基本的空气动力学方程入手，分析推导了空气阀气体等熵流动的进排气数学模型，其最终的数学模型与Wylie和Streeter相同，然后通过算例研究了空气阀的进排气流量系数对水锤防护效果的影响。赵秀红等[11]分析了空气阀排气性能实测资料与公式计算结果的差异，提出应根据实测资料选择空气阀。张健等[12]研究了空气阀设置位置、间距、数量与管道布置的关系，提出了一种配置空气阀的方法。杨开林[13]研究了利用空气阀形成空气阀调压室，以提高空气阀控制输水管道水锤防护的效果。

需要说明的是，虽然Wylie和Streeter[1]的空气阀水力瞬变数学模型在工程计算中广泛采用，包括现有商用软件，但是该模型基本假设(1)和(2)是相互矛盾的。实际上，输水管道工程中水击过程属于瞬态过程，持续时间很短，气温T不可能等于水温，而是随气压p的变化而变化。此外，基本假设(4)不考虑空气阀及其配套检修阀和连接管高度，可能导致计算结果偏离实际，产生严重的运行安全问题。

本文首先在取消Wylie和Streeter空气阀模型基本假设(2)和(4)的条件下研究水力瞬变新数学模型，包括数值求解算法；然后通过典型算例比较现有模型和新模型水力瞬变的差异。

2 数学模型

空气阀的边界条件是相当复杂的，涉及气体动力学和水力瞬变。下面将首先建立流入和流出空气阀的质量流量与阀内压强的函数关系，然后建立气体压强、体积与输水管道水压的函数关系，最后研究空气阀水力瞬变的数值计算方法。

2.1 空气阀进排气基本方程为使问题简化，首先假定：

(1)流动为准稳态，气体为理想(完全)气体且等熵的流进流出空气阀，即忽略空气阀内沿程能量损失与流道壁和水体热交换的影响；

(2)进入管道的空气留在它可以排出的阀附近，且高程差较小，气体重力的影响可忽略不计。

在上述条件下，空气阀到下部输水管之间局部能量损失的影响采用阀孔口等效阻力系数的方法处理。

理想气体的状态方程是

p=ρRT

(1)

即

p/pa=ρT/(ρaTa)

(2)

式中：p为阀内气体的绝对压强，Pa；ρ为阀内气体的密度，kg/m3；R为气体常数，一般取287 J/(kg·K)；T为阀内气体的绝对温度，K；pa为大气的绝对压强，Pa；ρa为大气的密度，kg/m3；Ta为大气的绝对温度，K。

气体等熵流动条件是

pρ-k=const

(3)

即

(ρ/ρa)k=p/pa

(4)

式中k为多变指数，对绝热过程，k=1.4。

根据式(2)(4)可得理想气体等熵流动密度、绝对温度和压强之间的函数关系

(5)

式中pr=p/pa为压比，表示空气阀内气压与大气压之比。上式表明，等熵条件下气温T随压比pr的增加而增加，随pr的减小而减小。

根据气体动力学，任意两断面之间等熵流动的能量是守恒的，势能的影响可忽略，内能和动能之和为常数，能量守恒方程是

(6)

式中：第一项为单位质量气体的内能；第二项为单位质量气体的动能，u为任意断面的平均流速，m/s。

如图1所示，当气体通过空气阀进入输水管时，由于输水管过流截面积常常比空气阀大两个数量级，故输水管中气流速度u可忽略不计。同样，阀外为无限大空间，流速u也可忽略不计。当气体仅存在于空气阀内，可将动能项归算到阀进口孔口能量损失中。

图1 空气阀进排气示意图

对于空气阀进气，列写空气阀外大气和阀内任一断面的能量方程可得

整理得

(7)

式中：ζin为空气阀进气的总局部阻力系数，包括进出口、断面变化等局部阻力；uin为阀孔口进气流速。

prc，in=(2/(k+1))k/(k-1)

(8)

式中prc，in为阀以临界声速进气的气体压比，简称进气临界压比，典型值在0.53～0.57，例如，对于绝热流动，k=1.4，则prc，in=0.53；当k=1.2，则prc，in=0.56。空气动力学研究表明，声音是无限小的压力波，因此，只要空气阀孔口(喉部)流动达到声速，压力扰动就不可能穿过孔口，即使压比进一步减小，即pr

类似的，对于空气阀排气，可得

(9)

式中：ζout为空气阀排气的总局部阻力系数；uout为空气阀孔口排气流速。

prc，out=1/prc，in=(2/(k+1))-k/(k-1)

(10)

式中prc，out为阀以临界声速排气的压比，简称排气临界压比，典型值在1.75～1.89。例如，对于绝热流动，k=1.4，则prc，out=1.89。

(11)

用式(8)中prc，in代替式(11)中pr，可得临界声速进气的质量流量

(12)

(13)

用式(10)中prc，out代替式(13)中pr可得空气阀临界声速排气的基本方程

(14)

对式(1)两边同乘气体体积，可得用气体质量表示的理想气体状态方程

pV=MaRT

(15)

式中：V为气体体积，m3；Ma为气体质量，kg。

当假设空气阀等熵流动为绝热流动，即取多变指数k=1.4，则式(11)—(15)就是Wylie和Streeter采用的空气阀进排气基本方程。当令式(11)(13)中项k/(k-1)=1，则(11)—(15)就是Lee等采用的空气阀数学模型。由于等熵过程气体介于等温过程k=1和绝热过程k=1.4之间，所以Lee等模型缺失项k/(k-1)可能是印刷错误，指出这一点是因为一些论文采用了这一模型。

需要说明的是，上述空气阀进排气基本方程是以火箭发动机喷管气体等熵流动理论为基础的。对于火箭发动机，其喷管内气体来自燃烧室，气温T受燃烧室温度控制，可达数千K的高温，即T是火箭发动机的主要控制变量，在火箭稳定飞行状态可视为已知常数。与火箭发动机不同，输水管道工程中空气阀只在水击瞬态过程中工作，气体来自大气而不是水体，且进排气持续时间很短，而T会随气压p的瞬变而瞬变。现有空气阀排气数学模型完全照搬火箭发动机稳态模型气温T为常数的假定，且取T等于水温显然是不合理的。如果T为常数，则意味着空气阀内气体是等温流动，即多变指数k≡1，这与等熵流动的假设矛盾。如果将k=1代入式(11)(13)，则公式右边为0/0型。实际上，随着输水工程管直径的增加，空气阀规格尺寸也成正比的增加，水力瞬变过程中的进气量相当大。以南水北调中线北京段惠南庄泵站为例，输水管直径4 m，单管正常输水流量30 m3/s，一旦泵站事故断电发生水力瞬变，就会有数百，甚至上千立方米气体迅速通过空气阀进入输水管，显然，要让气温等于水温且保持不变是不可能的。下面将根据理想气体状态方程和等熵条件，把气温T变换为压比pr和大气温度Ta的函数。

(16)

(17)

及等熵条件下的气体状态方程

(18)

式(11)(12)(16)(17)(18)就是本文导出的等熵条件下的新空气阀进排气基本方程，其中参数pa、Ta、k、R、Ain、Aout、Cin、Cout为原始数据(已知)，而未知量只有压比pr=p/pa、气体体积V和气体质量Ma，它们需要联立求解输水管道的水力瞬变才能确定。美国水行业协会标准《微量进排气阀、快速进排气阀及组合式快速进排气阀》(AWWA M51)建议取Cin=Cout=0.7。

图2 空气阀关系曲线

2.2 空气阀水力瞬变的数模目前涉及到空气阀的水力瞬变计算数学模型均没有考虑空气阀进气口与输水管的高程差。随着输水管直径的增加，空气阀尺寸也成正比的增加，考虑到正常输水条件下空气阀检修的需要，在空气阀与输水管道之间常常设置检修阀和连接管，如图3所示。目前我国输水工程中空气阀标称孔径Da≤0.4 m，空气阀、检修蝶阀和连接管高度之和最大值超过2 m。

图3 空气阀安装示意图

在一般情况下，连接管直径、蝶阀孔径和空气阀孔径Da相同，而空气阀内大部分断面净空面积与进出孔口相差不大，因此，空气阀、检修蝶阀和连接管的作用可用一个具有自动进排气功能的调压室等效，如图4所示，其中：Zat为等效调压室顶高程，即空气阀顶盖高程；Z为输水管顶高程；Hs为调压室水位；HP为输水管顶测压管水头。

图4 等效空气阀调压室工作原理图

为了使下面建立的数学模型具有普遍性，假设在调压室和输水管之间有一个阻抗孔，其中等效调压室的直径可近似取空气阀孔径Da，也可根据需要将调压室截面积描述为水位Hs的函数。

空气阀等效调压室的工作原理是：在输水管道发生水力瞬变条件下，随着测压管水头HP降低到调压室顶高程Zat以下，即HP

为计算方便，在建立气体压比pr、气体体积V和质量Ma与输水管水力瞬变关系之前，先假设进入输水管的气体体积和管段里的液体体积相比很小，可采用常规特征线方法求解管道的水力瞬变。

调压室底部节点P的连续性方程为

Qs=QT-Q

(19)

式中：Qs为流入或者流出调压室的流量，m3/s；QT为流入节点P的输水管流量，m3/s；Q为流出节点P的输水管流量，m3/s。

根据特征线相容性方程[3]可得

C+：QT=CP/BP-HP/BP

(20)

C-：Q=-CM/BM+HP/BM

(21)

由式(20)和式(21)代入式(19)得

Qs=C1-C2HP

(22)

式中

(23)

在时刻t1系数C1和C2是已知量。

在正常输水工况下，空气阀完全关闭，Qs≡0。在管道系统发生水力瞬变时，如果输水管内没有气体V=0且测压管水头HP=C1/C2>Zat，则空气阀不工作，Qs=0；否则，需要采用下述方法建立调压室水位Hs、气体体积V和流量Qs的函数关系

(24)

(25)

式中：t为时间，s；As为调压室的截面积，m2。

对式(24)(25)积分并取二阶近似得

Hs=C3+C4Qs

(26)

V=V0-0.5Δt(Qs+Qs0)

(27)

式中：下标“0”表示时刻t0；Δt=t-t0。

(28)

当不考虑调压室水体惯性力、沿程水头损失和室壁的弹性，气体压力p(绝对压力)与输水管顶测压管水头HP、水位Hs和流量Qs的关系是

p/γ=HP+Ha-Hs-C5|Qs|Qs=HP+Ha-Hs-2C5|Qs0|Qs+C5|Qs0|Qs0

(29)

式中：γ为水的比重，取9807 kg/m3；Ha=pa/γ为当地大气压压头。

(30)

式中：ω为调压室底部阻抗孔面积，m2；ζ为调压室底部阻抗孔局部阻力系数。当调压室只是空气阀、检修阀和连接管的简单组合，则ω=As，这时调压室底部阻抗孔可视为突扩或突缩和90°转弯的组合，初步计算可取ζ=1.0。

把式(22)和式(26)代入式(29)，得气体流量Qs与压比pr的关系式，

papr/γ=(C1/C2-Qs/C2)+Ha-(C3+C4Qs)-(2C5|Qs0|Qs-C5|Qs0|Qs0)

整理得

Qs=C6-C7pr

(31)

式中

(32)

把式(31)代入式(27)得气体体积V与压比pr的关系式

V=C8+C9pr

(33)

式中

C8=V0-0.5Δt(C6+Qs0)，C9=0.5ΔtC7

(34)

把式(33)气体体积V代入等熵条件的气体状态方程式(18)，并对气体质量Ma取二阶近似得

(35)

由于多变指数k≠1，函数F不是压比pr的抛物线方程，不能采用Wylie和Streeter的数值算法求解。

由于式(35)函数F中只有压比pr是未知量，采用牛顿-雷伏生数值计算方法得

F+FprΔpr=0

整理得下述迭代计算求解pr的联立方程组

Δp=-F/Fp

(36a)

式中：

(36b)

(36c)

(36d)

(36e)

式中δ>0为pr的微小增量，可取δ=10-7。

当不考虑空气阀高度影响时，只需令Zat=Z，则上述空气阀水力瞬变计算数学模型同样适用。

计算是由t0=0开始，初始条件是

V0=0，Ma0=0，Qs0=0，Hs0=Zat

(37)

当空气阀排尽气体时，

Qs=0，Hs=Zat，V=0，Ma=0

(38)

当V=0、Qs=0、Hs=Zat时，空气阀接头处的边界条件就是HP、QT、Q的一般的内截面解，

QP=QT=Q=(CP-CM)/(BP+BM)，HP=CP-BpQT

(39)

采用上述计算程序，只要迭代次数足够多，计算总是收敛的。在一般情况下，经过几次迭代就可以满足计算精度要求。但在个别情况下，当空气阀排气孔径设计过小时，为了保证计算收敛，需调整收敛因子σ的大小，一般取较小的σ，例如σ=0.05甚至更小，同时增加迭代次数的控制值。

3 不同数学模型的比较

以Wylie和Streeter[1]专著中的一个典型输水管中空气阀的水力瞬变为例，比较Wylie和Streeter空气阀水力瞬变数学模型和本文新模型的差异。输水管长1219.2 m，分成两段，Δx=609.6 m；管径D=0.6096 m，对应截面积A=0.292 m2；水击波速a= 1219.2 m/s；达西-威斯巴哈沿程阻力系数f=0.02；空气阀进排气通流面积Ain=Aout=0.001858 m2，对应空气阀孔径Da=0.049 m；空气阀安装在输水管中间位置，管顶高程Z=10.36 m；大气绝对压头Ha=pa/γ=10.36 m；输水管下游是以一个水库作为边界条件，水位Z0=9.75 m；输水管进口处的强迫流量QP(1)=Q0-ΔQsin(OM×t)，Q0=0.3398 m3/s，ΔQ=0.1133 m3/s，OM=0.3 rad。计算中取空气阀流量系数Cin=Cout=0.7。

算例1：当不考虑空气阀高度，包括检修阀和连接管高度，即Zat=Hs=Z，且取大气温度为10 ℃和水温为4 ℃，即Ta=273+10=283 K和T=273+4=277 K，则在气体为绝热流动k=1.4条件下得如图5所示空气阀水力瞬变的计算结果，图中：实线和虚线分别是新模型与Wylie和Streeter空气阀进排气模型时水力瞬变的计算结果；H=HP-Z为空气阀底部压头。

图5 空气阀水力瞬变曲线(不考虑空气阀高度)

从图5可见，在空气阀排气未了(V=0)可能产生剧烈的液柱弥合水击压力，水压突然升高，然后又突然下降。采用新模型与Wylie和Streeter模型计算的水力瞬变两者最大水压Hmax和最小水压Hmin差别明显，新模型Hmax1=39.26 m和Hmin1=-2.69 m，而Wylie和Streeter模型Hmax2=36.04 m和Hmin2=-1.58 m，两者偏差相对值(Hmax1-Hmax2)/Hmax2=9%和(Hmin1-Hmin2)/Hmin2=70%。这一结果表明，新模型计算的最大水压更大，而最小负压更低。这意味着按照新模型计算结果作工程设计更为合理，是偏于安全的。在图5(c)中也示出了空气阀内气温随时间的变化：在空气阀进气V增大时，气温T随着压比pr的减小而减小，然后随着pr的增大而升高，其中最大绝对气温Tmax=293.2 K、最低绝对气温Tmin=259.7 K，或者最大气温为20.2 ℃、最低气温为-13.3 ℃。

作为进一步比较，也应用新模型计算了多变指数k=1.2时的管道水力瞬变，最大水压Hmax=37.97 m，而最小水压Hmin=-2.60 m。与绝热流动k=1.4时的结果相比，可得重要结论：输水管最大水压Hmax随多变指数k的增加而增加，但是最小水压Hmin随k的增加而减小，因此假设空气阀气体为绝热流动的计算结果作为设计依据是合理的，是偏于安全的。

算例2：当考虑空气阀高度时，以一种典型空气阀为例(1)资料来源美国GA工业公司空气阀使用说明。，孔径Da=0.049 m的空气阀高度ha≈0.327 m，当采用相同标称直径蝶阀作为检修阀时，蝶阀高度L=0.108 m，这时若取等效调压室高度h=Zat-Z=0.44 m时，则Hs0=Zat=Z+0.44 m。在此条件下，如果其他条件与算例1相同，则采用新模型得如图6所示的计算结果，图中：实线和虚线分别表示不考虑空气阀高度(h=0.0 m)和考虑空气阀高度。

图6 空气阀水力瞬变曲线(空气阀高度的影响)

观察图6可见，考虑空气阀高度影响时，Hmax3=37.92 m和Hmin3=-2.80 m；而不考虑空气阀高度影响时，Hmax1=39.26 m和Hmin1=-2.69 m；偏差相对值(Hmax1-Hmax3)/Hmax3=3.5%和(Hmin1-Hmin3)/Hmin3=-4%。这一结果表明，即使本例空气阀高度较小，仅为0.44 m，但对输水管最大和最小水压产生了明显影响。考虑空气阀高度后，计算的最大水压减小，这是有利的；但是最小水压也减小，即导致输水管负压更加严重，这是不利的。发生这一现象的主要原因是，空气阀等效调压室内水体对进排气有阻止延缓作用。

需要说明的是，在实际工程中空气阀连接管长度差别较大，不同的工程、空气阀不同的布设位置，对连接管的型式、体积大小要求也不一样，甚至没有。目前关于空气阀的要求与应用限制越来越严格，如果合理设计连接管的高度和直径将空气限制在连接管内而不进入输水管，则能够显著提高输水的安全性。例如，利用空气阀和连接管形成空气阀调压室[13]可以大大提高输水工程的水锤防护能力。

4 结论

典型算例的计算研究表明：

(1)空气阀进气V增大时，气温T随着压比pr的减小而减小，然后随着pr的增大而升高，其中最低气温可能降低到0 ℃以下。

(2)采用新空气阀水力瞬变模型和 Wylie和Streeter模型的计算结果差异明显，最大水压偏差相对值可达9%以上，而最小水压偏差相对值可达70%以上。换句话说，不考虑空气阀排气过程气温变化的影响会导致计算结果存在安全风险。

(3)空气阀及其配套检修阀和连接管高度h会对输水管最大和最小水压产生明显影响，即使h较小。

(4)合理设计连接管的高度和直径将空气限制在连接管内而不进入输水管，可以显著提高输水的安全性，例如利用空气阀和连接管形成空气阀调压室，则可大大提高输水工程的水锤防护能力。