丁 聰

(上?？睖y設(shè)計研究院有限公司，上海 200434)

1 概述

我國沿江濱湖地區(qū)已建許多大型低揚程泵站。停泵過渡過程的研究是大型軸流泵、泵裝置、泵站重要研究內(nèi)容之一[1-2]。泵站在起動、停機等水力瞬變過程中常出現(xiàn)最不利工況，在運行安全可靠性方面出現(xiàn)問題，因此，停泵過渡過程的計算是校核泵組運行安全性的重要內(nèi)容。低揚程泵裝置停泵動態(tài)特性計算需分析確定機組轉(zhuǎn)動力矩平衡方程、流道的非恒定流方程、閘門水力損失表達式，并合理的引入水泵的全特性，最終建立停泵過渡過程的數(shù)學(xué)模型[3]。水泵的全特性較為復(fù)雜，為方便計算，可利用無因次新坐標參數(shù)化得到的Suter特性曲線。為有效減少Suter特性曲線擬合取值誤差，可使用較成熟的移動最小二乘法進行特性曲線的擬合插值[4]。

建立的停泵過渡過程的數(shù)學(xué)模型為微分方程組，因此，停泵過渡過程的仿真求解實質(zhì)為解微分方程。Matlab中的Simulink采用模塊化思路搭建數(shù)學(xué)模型，可方便進行微分方程的求解，并較容易修改模塊參數(shù)。因此，本文擬利用Simulink進行低揚程水泵裝置停泵過渡過程的仿真計算。

2 停泵過渡過程數(shù)學(xué)模型建立

2.1 機組轉(zhuǎn)動力矩平衡方程

停泵過渡過程中，低揚程泵機組的轉(zhuǎn)動力矩平衡方程為[5]：

(1)

式中GD2為回轉(zhuǎn)部分的飛輪力矩；n為機組轉(zhuǎn)速；Mz為的阻力矩，包括水泵水力矩、水泵慣性附加力矩、泵組摩擦力矩、風(fēng)阻力矩。

由于泵組摩擦力矩及風(fēng)阻力矩較小，可忽略不計，僅考慮水泵水力矩、水泵慣性附加力矩。

水泵水力矩為水泵固有特性，詳見2.3節(jié)。水泵慣性附加力矩可使用下式計算。

(2)

其中θ為泵葉片平面包角；α為泵葉輪輪轂比；βr為根部翼型安放角；βR為外緣翼型安放角。

2.2 水頭平衡方程

低揚程水泵裝置停泵時流道為滿水狀態(tài)。取進、出水池兩個斷面，依據(jù)伯努利方程可得：

(3)

水頭平衡方程可簡寫為[6]：

(4)

式中Hsta為出、進水池水位差，Hsta=Z2-Z1。

2.3 水泵水力特性

低揚程水泵裝置停泵過程中水泵必然經(jīng)歷水泵工況、泵制動工況、水輪機工況。數(shù)值模擬停泵過渡過程必然要使用水泵的全特性曲線。水泵全特性曲線是以流量為橫坐標、揚程和轉(zhuǎn)矩為縱坐標繪制不同轉(zhuǎn)速下的四象限特性曲線；或以相對速度為橫坐標、相對流量為縱坐標，以兩簇等相對揚程和等相對轉(zhuǎn)矩曲線表示的性能關(guān)系。水泵全特性曲線形式比較復(fù)雜，不便于數(shù)值計算時使用。

Suter、Marchal和Flesh[7]根據(jù)泵的相似原理，提出x=π+arctan(v/β)、WH(x)=h/(β2+v2)、WB(x)=m/(β2+v2)新體系，將水泵全特性曲線變換為2條無因次曲線。

(5)

式中h=H/HN，β=n/nN，v=Q/QN，m=M/MN，為水泵揚程、轉(zhuǎn)速、流量、轉(zhuǎn)矩?zé)o量綱量。

由于停泵過渡過程中水泵特性連續(xù)變化，因此，為保證模擬精度需使用實測的連續(xù)的水泵性能曲線。但試驗所得的特性曲線也均為離散點，若直接采用離散點插值計算，會大大降低模擬精度。移動最小二乘法(MLS)是一種基于點的近似方法，具有擬合精度高、通用性強的特點。因此，為較好地使用已知Suter特性曲線，使用移動最小二乘法(MLS)對Suter特性曲線進行擬合，其基本原理如下[8]。

設(shè)待求函數(shù)u(x)在域Ω內(nèi)個節(jié)點xI(I=1，2，…N)處數(shù)值已知，即uI=u(xI)。計算點x的域Ωx內(nèi)函數(shù)u(x)的函數(shù)uh(x)可項式表示為：

(6)

其中p(x)為式基向量，a(x)為向量，m代的項數(shù)。通常選取單項式作為基函數(shù)，此處選擇二次基函數(shù)，即m=3，pT(x)=[1，x，x2]。LS是通過令uh(x)與節(jié)點的函數(shù)值uI之加權(quán)平方和最小構(gòu)造近似函數(shù)，即使下式取最小值。

(7)

其中n為點x的求解域Ωx內(nèi)含的節(jié)點數(shù)。w(x-xI)為點xI相權(quán)函數(shù)，此處取高斯函數(shù)。令：

W(x)=diag(w1(x)，w2(x)，…，wn(x))，

wi(x)=w(x-xi)，

B(x)=pTW(x)，

由式(7)取最小值解得系數(shù)向量a(x)。

a(x)=A-1(x)B(x)u

(8)

將式(8)代入式(6)中得：

uh(x)=?(x)u

(9)

2.4 閘門關(guān)閉特性數(shù)學(xué)模型[9]

水流經(jīng)過快速閘門屬于淹沒孔口出流，則流量可表示為：

(10)

式中k為閘門相對開度，k=Ak/Az=(t0-t)/t0，t0為全關(guān)需要時間，Ak為開啟面積；Az為面積；μ為流量系數(shù)，μ=0.6-0.176e。

閘門的水頭損失可表示為：

(11)

式中ζz為局部損失系數(shù)，ζz=1.2(k-2.63-1)。

2.5 停泵過渡過程數(shù)學(xué)模型

結(jié)合式(1)(2)(5)可得機組轉(zhuǎn)動力矩平衡綜合方程：

(12)

結(jié)合式(4)(11)可得停泵水頭平衡綜合方程：

(13)

3 停泵過渡過程Simulink模型搭建

Matlab軟件的Simulink提供了各種數(shù)學(xué)模塊，僅需根據(jù)仿真需要將不同模塊組合搭建即可。但Simulink模塊搭建有其自有的規(guī)則。為方便停泵過渡過程數(shù)學(xué)模型的搭建，需將式(12)進行變形。

(14)

根據(jù)停泵過渡過程數(shù)學(xué)模型表達式(13)(14)，搭建停泵過渡過程模型原理圖(見圖1所示)。

圖1 停泵過渡過程模型原理示意

根據(jù)停泵過渡過程模型原理圖在Simulink中進行模塊搭建，搭建仿真模塊如圖2所示。

圖2 停泵過渡過程Simulink模塊示意

由于Suter特性曲線以流量、轉(zhuǎn)速的無因次量為橫坐標，而計算中需利用Suter特性曲線求取水泵揚程、轉(zhuǎn)矩。因此，可理解為給定流量、轉(zhuǎn)速，求取揚程、轉(zhuǎn)矩，如此可利用Simulink模塊中的自定義MATLAB Function進行計算，并可將移動最小二乘法(MLS)對Suter特性曲線擬合計算內(nèi)置其中。

4 停泵過渡過程動態(tài)特性

根據(jù)上述資料，計算得數(shù)學(xué)模型中的各項系數(shù)如下：

圖3 閘門未關(guān)閉停泵過渡過程示意

在正轉(zhuǎn)正流、正轉(zhuǎn)逆流階段，泵的流量、轉(zhuǎn)速下降急劇，2.63 s左右泵開始倒流，轉(zhuǎn)速降至0.423nN；85 s左右泵開始倒轉(zhuǎn)，反轉(zhuǎn)流量為0.769QN；泵組趨于穩(wěn)定時反轉(zhuǎn)流量為1.552QN，轉(zhuǎn)速為1.426 8nN，轉(zhuǎn)速已超反轉(zhuǎn)1.2倍要求。因此，泵站停泵時需配備事故閘門，閘門應(yīng)做到快速關(guān)閉。

考慮閘門關(guān)閉，假設(shè)閘門線性關(guān)閉，20 s內(nèi)全部關(guān)閉，仿真計算的結(jié)果如圖4所示。

圖4 閘門線性關(guān)閉停泵過渡過程示意

由于停泵初期閘門關(guān)閉較小，泵的流量、轉(zhuǎn)速下降急劇。2.63 s左右泵開始倒流，轉(zhuǎn)速降至0.423nN；85 s左右泵開始倒轉(zhuǎn)，反轉(zhuǎn)流量為0.761QN。但隨著閘門孔口的減小，泵組反轉(zhuǎn)流量、反轉(zhuǎn)轉(zhuǎn)速下降速度減緩，并開始減小，最終將趨于0流量和0轉(zhuǎn)速。

若計算閘門不同關(guān)閉規(guī)律下的泵站停泵過渡特性，僅需對gate close模塊進行修改。若想更換水泵葉輪，只需將MATLAB Function模塊中的水泵Suter特性曲線進行更換。若模擬其他泵站的停泵過渡特性，僅需修改相應(yīng)的水泵物理參數(shù)、流道特性參數(shù)以及泵站基本水位特性。由此可見，僅需簡單修改Simulink模型中的模塊參數(shù)便可完成不同閘門關(guān)閉規(guī)律的停泵特性計算，以及其他新工程的模擬計算。

5 結(jié)語

本文分析確定了低揚程水泵裝置停泵過渡過程的機組轉(zhuǎn)動力矩平衡方程、流道的非恒定流方程、閘門水力損失表達式，并采用移動最小二乘法，較合理地擬合了水泵Suter特性曲線。然后使用Simulink模塊化搭建了泵站停泵過渡過程模型，并模擬了閘門未動作和閘門動作時的停泵特性。通過模型搭建及模擬計算可以發(fā)現(xiàn)，Simulink模塊可方便搭建泵站停泵數(shù)學(xué)模型，并且簡單修改模塊參數(shù)，便可比較方案的合理性和完成新工程的計算。Simulink可方便地將模擬計算通用化。