雷 恒，張 兵

（1.黃河水利職業(yè)技術(shù)學(xué)院 水利系，河南 開封 475003；2.小流域水利河南省高校工程技術(shù)研究中心，河南 開封 475003）

0 引言

在應(yīng)用計算機輔助設(shè)計技術(shù)進行水輪機的選型設(shè)計、經(jīng)濟運行和過渡過程數(shù)字仿真等工作時，常采用大量離散數(shù)據(jù)點表示水輪機特性曲線。而水輪機特性曲線數(shù)據(jù)是由轉(zhuǎn)輪模型試驗獲得的，這些特性曲線均是一系列大量的離散數(shù)據(jù)點，即以單位流量Q11和單位轉(zhuǎn)速n11作為坐標(biāo)系的一系列等值平面曲線數(shù)據(jù)點［1］。實質(zhì)上是用二維的圖形表述三維的信息，表達不直觀。本文針對不規(guī)則離散的水輪機特性數(shù)據(jù)，利用MATLAB強大的數(shù)據(jù)處理功能［2］，通過三次樣條插值和Delaunay三角剖分求取水輪機特性曲線。

1 水輪機特性離散數(shù)據(jù)采集

水輪機流量特性為Q11＝f（α，n11）（α為水輪機導(dǎo)葉開度），實際上是水輪機模型綜合特性在三維坐標(biāo)系（Q11，n11，α）中的三維曲面投影到（n11，Q11）平面上的一系列等高線。在MATLAB程序中，首先要為X、Y、Z三個變量進行賦值，再把在Excel表格中處理好的數(shù)據(jù)列復(fù)制、粘貼到MATLAB的變量值表中，這樣能夠把大量的列項量數(shù)值快速地復(fù)制給MATLAB中的變量。為了表示各采樣點的地理位置，需要繪制采樣點的分布圖，如圖1所示。

2 數(shù)據(jù)插值和等值線的繪制［3］

水輪機特性采樣點多，呈不規(guī)則分布，其生成的等值線或等值面（如圖2所示）需要進行數(shù)據(jù)插值和數(shù)據(jù)網(wǎng)格化。常用的插值算法有線性插值、最近鄰插值、三次多項式插值和樣條插值4種。其中線性插值和最近鄰插值結(jié)果所構(gòu)成的曲面不平滑，精度較低；三次多項式插值結(jié)果所構(gòu)成的曲面較平滑；樣條插值結(jié)果所構(gòu)成的曲面最平滑，精度高。

在任意插值計算區(qū)間［a，b］的每一個子區(qū)間［xk，xk＋1］上，三次樣條函數(shù)為：

其中：x為水輪機的單位轉(zhuǎn)速；hk為相鄰水輪機單位轉(zhuǎn)速的差值，記為hk＝xk＋1－xk；mk為三次樣條插值函數(shù)的二階導(dǎo)數(shù)；Ak、Bk為積分常數(shù)sk（xk＋1）。y為水輪機的單位流量；mk可通過下面方程組解出：

三次樣條插值的計算步驟如下：首先確定端點條件，然后構(gòu)成相應(yīng)的方程式組（2），求出mk（k＝1，2，3，…，n），最后根據(jù)式（1）求得各子區(qū)間［xk，xk＋1］上的三次樣條插值函數(shù)。

圖1 水輪機流量特性曲線離散數(shù)據(jù)樣本分布

圖2 等值線圖

3 插值三維網(wǎng)格曲面繪制和Delaunay三角網(wǎng)

通過MATLAB可以實現(xiàn)三維網(wǎng)格曲面圖形的繪制［4］。插值數(shù)據(jù)生成的三維網(wǎng)格見圖3，生成的三維表面圖見圖4。

圖3 插值數(shù)據(jù)形成的三維網(wǎng)格圖

不規(guī)則三角網(wǎng)可減少規(guī)則格網(wǎng)方法帶來的數(shù)據(jù)冗余，同時計算效率也較高。如何把一個散點集合剖分成不均勻的三角形網(wǎng)格，這就是散點集的三角剖分問題。不規(guī)則三角網(wǎng)模型是根據(jù)區(qū)域有限個點集將區(qū)域劃分為相連的三角面網(wǎng)格，區(qū)域中任意點落在三角面的頂點、邊上或三角形內(nèi)。在實際中，運用最多的三角剖分是Delaunay三角剖分［5］。

要滿足Delaunay三角剖分必須符合兩個重要的準(zhǔn)則：一是最大－最小化角特性，即在散點集可能形成的三角剖分中，Delaunay三角剖分所形成的三角形最小角最大；二是空圓特性，即Delaunay三角網(wǎng)是唯一的（任意4點不能共圓），在Delaunay三角網(wǎng)中，任一三角形的外接圓范圍內(nèi)不會有其他點存在。

圖4 插值數(shù)據(jù)形成的三維表面圖

生成的三維三角網(wǎng)格圖見圖5，生成的三維三角網(wǎng)曲面圖見圖6。

圖5 三維三角網(wǎng)格圖

圖6 三維三角網(wǎng)曲面圖

4 結(jié)果比較分析

通過計算，對三維曲面中的開度值α與原圖的等開度值α作比較，處理結(jié)果見圖7。

圖7 處理結(jié)果比較

在圖7中取誤差最大的一條等開度線（α＝20 mm）作誤差分析，其數(shù)據(jù)見表1。通過比較可知，最大絕對誤差不超過0.2 mm。由此可見，本文提出的求取水輪機特性曲線中隨機插值點的方法是可行的，能滿足工程的實際要求。

表1 水輪機特征值計算結(jié)果及與原值對比（α＝20 mm）

5 結(jié)論

本文針對不規(guī)則離散的水輪機特性，利用MATLAB強大的數(shù)據(jù)處理功能，通過三次樣條插值和Delaunay三角剖分確保了其計算精度和唯一性，對水輪機的正確選型和優(yōu)化設(shè)計是有利的。

［1］ 張昌期.水輪機－原理與數(shù)學(xué)模型［M］.武漢：華中工學(xué)院出版社，1988.

［2］ 蔡旭輝，劉衛(wèi)國，蔡立燕.MATLAB基礎(chǔ)與應(yīng)用教程［M］.北京：人民郵電出版社，2009.

［3］ 張蓉生，劉志鵬，屈波.水輪機模型特性圖計算機輔助數(shù)據(jù)采集及擬合與運轉(zhuǎn)特性曲線的生成［J］.機械工程學(xué)報，2006（4）：222－226.

［4］ 王宣懷，沈祖詒，孫涌.基于主曲線方法的水輪機特性曲線的數(shù)值擬合［J］.水力發(fā)電學(xué)報，2009，28（3）：181－186.

［5］ 劉巖，關(guān)振群，張洪武，等.面向大規(guī)模科學(xué)計算的三維Delaunay快速插點算法［J］.中國科學(xué)：物理學(xué) 力學(xué) 天文學(xué)，2012（2）：192－198.