程洪霞 吳建華

（太原理工大學(xué)水利科學(xué)與工程學(xué)院，太原 030024）

1 蝶閥關(guān)閉角間隔點(diǎn)的計(jì)算機(jī)存儲(chǔ)及直線插入法

在國(guó)內(nèi)外大量的供水系統(tǒng)水錘分析文獻(xiàn)中，給出了蝶閥關(guān)閉角每隔5°時(shí)，相對(duì)應(yīng)的阻力系數(shù)和過流面積與閥門面積比值的離散值，見表1。

表1 阻力系數(shù)值及過流面積與閥門面積的比

在水力過渡過程的計(jì)算中，將關(guān)閉角度每隔5°的阻力系數(shù)和過流面積值輸入計(jì)算機(jī)中，并求出相對(duì)應(yīng)的水頭損失系數(shù)。停泵后，任意時(shí)刻閥門關(guān)閉角對(duì)應(yīng)的水頭損失系數(shù)，用直線插入法計(jì)算。

在表1給出的數(shù)值及所采用的直線插入法的基礎(chǔ)上，運(yùn)用多元線性回歸模型，給出過流面積與閥門面積的比值關(guān)于關(guān)閉角度的函數(shù)關(guān)系式，以及阻力系數(shù)和關(guān)閉角度的函數(shù)關(guān)系式，使水力過渡過程的數(shù)值模擬更為精確、易行。

2 多元線性回歸模型

在實(shí)際問題中，影響結(jié)果y的因素往往不止一個(gè)。一般地，設(shè)有 x1，x2，…，xk，k 個(gè)因素，通?？紤]如下的線性關(guān)系式[1]：

其中，β0，β1，…，βk為 k+1 個(gè)未知數(shù)，ε 是隨機(jī)變量，一般假設(shè) E（ε）=0，D（ε）=σ2，σ2＞0。為了估計(jì)未知參數(shù) β0，β1，…，βk及 σ2，對(duì) y 與 x1，x2，…，xk同時(shí)作 n 次獨(dú)立觀察，得 n 組觀察值（yt，xt1，xt2，…，xtk），t=1，2，…，n（n＞k+1），滿足關(guān)系式：

其中，ε1，…εn互不相關(guān)且均是與ε同分布的隨機(jī)變量。為用矩陣形式表示上式，令：

于是，多元線性回歸模型可表示為：

其中，In為n階單位矩陣。

用最小二乘法，得β的最小二乘估計(jì)為：

3 數(shù)學(xué)模型的建立與求解

3.1 過流面積和閥門面積比值Y1與關(guān)閉角度x之間的函數(shù)關(guān)系

根據(jù)關(guān)閉角度x與過流面積和閥門面積比值Y1在表1中的數(shù)值，由MATLAB軟件可繪制出（xi，Y1i）（其中i=1，2，…，19）的散點(diǎn)圖，如圖1所示。

根據(jù)圖1，對(duì)過流面積和閥門面積比值Y1與關(guān)閉角度x的關(guān)系，可用二次函數(shù)進(jìn)行回歸分析。經(jīng)回歸分析的實(shí)驗(yàn)驗(yàn)證：①當(dāng)只采用一個(gè)二次函數(shù)進(jìn)行回歸分析時(shí)，誤差較大，所以采用分段函數(shù)的形式；②分段區(qū)間選擇［0，50），［50，80），［80，85），［85，90］上較為精確；③在區(qū)間段［0，50），［50，80）適合用二次函數(shù)進(jìn)行回歸分析，這時(shí)誤差很小，在區(qū)間［80，85］，［85，90］上適合用原來的直線插入法，這樣更符合實(shí)際。下面給出在區(qū)間［0，50），［50，80）上，用二次函數(shù)進(jìn)行回歸分析的計(jì)算方法。

圖1 （xi，Y1i）的散點(diǎn)圖

圖2 Y1與x的函數(shù)關(guān)系圖

在用二次函數(shù)y=β0+β1x+β2x2進(jìn)行回歸分析的計(jì)算時(shí)，首先令：xk=xk，（k=1，2），從而轉(zhuǎn)換為多元線性回歸模型，運(yùn)用上述多元線性回歸模型原理，并由MATLAB軟件進(jìn)行矩陣計(jì)算（MATLAB中矩陣計(jì)算方法見文獻(xiàn)2），可以求得過流面積和閥門面積的比值 Y1與關(guān)閉角度 x 在區(qū)間［0，50），［50，80）上的函數(shù)關(guān)系為：

由MATLAB軟件可繪制出（5）式以及在區(qū)間［80，85］，［85，90］上運(yùn)用直線插入法得到的函數(shù)圖象，如圖2?？梢钥闯?，圖2與圖1很接近。

為進(jìn)一步判斷由（5）式確定的函數(shù)關(guān)系的準(zhǔn)確性，可以計(jì)算閥門關(guān)閉角度每隔5°時(shí)，由（5）式得到的過流面積和閥門面積的比值與原來的過流面積和閥門面積比值的相對(duì)誤差，如表2。

表2 過流面積和閥門面積比的相對(duì)誤差

由表2可知，由（5）式得到的過流面積和閥門面積的比值與原來的過流面積和閥門面積比值的相對(duì)誤差是很小的，即由多元線性回歸模型得到的過流面積和閥門面積比值Y1與關(guān)閉角度 x在區(qū)間［0，50），［50，80）上的函數(shù)關(guān)系是合理的。

3.2 阻力系數(shù)Y2與關(guān)閉角度x的函數(shù)關(guān)系

根據(jù)表1中關(guān)閉角度x和相對(duì)應(yīng)的阻力系數(shù)Y2的數(shù)值，可得到關(guān)閉角度x和對(duì)應(yīng)的阻力系數(shù)的自然對(duì)數(shù)1nY2的數(shù)值，如表3。由MATLAB軟件可繪制出（xi，1nY2i）（其中i=1，2，…，19）的散點(diǎn)圖，如圖 3 所示。

表3 關(guān)閉角度所對(duì)應(yīng)阻力系數(shù)的自然對(duì)數(shù)

根據(jù)圖3，對(duì)阻力系數(shù)的自然對(duì)數(shù)1nY2與關(guān)閉角度x的關(guān)系，可用二次函數(shù)進(jìn)行回歸分析。經(jīng)回歸分析的實(shí)驗(yàn)驗(yàn)證：①當(dāng)只采用一個(gè)二次函數(shù)進(jìn)行回歸分析時(shí)，誤差較大，所以需要采用分段函數(shù)的形式；②分段區(qū)間選擇［0，25），［25，50），［50，80），［80，85］，［85，90］ 較為精確；③在區(qū)間段［0，25），［25，50），［50，80）用二次函數(shù)進(jìn)行回歸分析，這時(shí)誤差很小，在區(qū)間［80，85］，［85，90］上適合用原來的直線插入法。

對(duì)阻力系數(shù)的自然對(duì)數(shù)1nY2和關(guān)閉角度x，用二次函數(shù)1nY2=β0+β1x+β2x2進(jìn)行回歸分析，通過轉(zhuǎn)換為多元線性回歸模型，運(yùn)用多元線性回歸模型的原理，并由MATLAB軟件進(jìn)行矩陣計(jì)算，得：

Y2可由（6）式推出。

圖3 （xi，1nY2i）的散點(diǎn)圖

圖4 Y2與x函數(shù)關(guān)系圖

由MATLAB軟件可繪制出（6）式及在區(qū)間［80，85］，［85，90］上運(yùn)用直線插入法得到的函數(shù)圖象，如圖4?？梢钥闯?，圖4與圖3較為接近。

表4 阻力系數(shù)的相對(duì)誤差

為進(jìn)一步判斷由（6）式所確定的阻力系數(shù)與關(guān)閉角度的函數(shù)關(guān)系的準(zhǔn)確性，可以計(jì)算閥門關(guān)閉角度每隔5°時(shí)，由（6）式計(jì)算得到的阻力系數(shù)與原阻力系數(shù)的相對(duì)誤差，如表4。

由表4可知，由（6）式所得阻力系數(shù)與原阻力系數(shù)的相對(duì)誤差較小，所以，用多元線性回歸模型得到的阻力系數(shù)Y2與關(guān)閉角度 x 在區(qū)間［0，25），［25，50），［50，80）上的函數(shù)關(guān)系是可行的。

4 結(jié)束語

根據(jù)蝶閥關(guān)閉角每隔5°時(shí)，相對(duì)應(yīng)的阻力系數(shù)和過流面積與閥門面積比的離散數(shù)值，首先，用MATLAB軟件分別繪制出過流面積與閥門面積的比值關(guān)于關(guān)閉角度的離散圖，以及阻力系數(shù)和關(guān)閉角度的離散圖；其次，運(yùn)用多元線性回歸分析原理，給出了在區(qū)間［0，80）上過流面積與閥門面積的比值關(guān)于關(guān)閉角度的函數(shù)關(guān)系式，以及阻力系數(shù)和關(guān)閉角度的函數(shù)關(guān)系式；最后，分別對(duì)兩個(gè)函數(shù)關(guān)系式進(jìn)行了相對(duì)誤差分析，證實(shí)了兩個(gè)函數(shù)關(guān)系的合理性和有效性。通過文中得到的兩個(gè)函數(shù)關(guān)系式，可以使蝶閥關(guān)閉規(guī)律在水力過渡過程數(shù)值模擬的運(yùn)用中變得更為精確、易行。

［1］孫榮恒.應(yīng)用數(shù)理統(tǒng)計(jì)［M］.北京：科學(xué)出版社，2003：195-198.

［2］龔劍，朱亮.MATLAB入門與提高［M］.北京：清華大學(xué)出版社，2000：30-39.