1、管道非恒定流動數值模擬研究余祖耀 張輝意 李寶仁(華中科技大學FESTO氣動中心 武漢 430074摘要:從動力學方程出發,建立管道非恒定流動數學模型,并探討相應的數值計算方法。以簡單管道非恒定流動系統為例進行數值模擬計算,分析計算結果并得出相應結論。為管道非恒定流動數值模擬研究進行了有效探索。關鍵詞:管道 非恒定流 數值模擬中圖分類號:TB71+20 前言流體管道是流體傳輸、傳動和控制工程中用以輸送流體介質、傳遞流體動力和信息的不可或缺的元件。在實際管道系統中,由于組成系統的某一元件工作狀態的變更(如閥的開度變化、泵的脈動等或受外界干擾(如負載的變化,將不可避免地在管道內產生流量和壓力的沖擊

2、或脈動,形成非恒定流動,引起管內流動的動態過程。流體管道的動態特性對系統的穩定性和可靠性,以及系統中其它元件的正常工作有著很大影響,因此分析流體管道動態過程,對管路系統的合理設計、優化、故障診斷、消振濾波及保證管路的正常工作有重要指導意義。流體管道動態特性數值模擬具有研究周期短、經費投入少,不受模型尺寸、外界擾動、測量精度限制等優點,其作用相當于在計算機上進行復雜流體試驗3。數值模擬包含以下幾個步驟:建立反映問題(工程問題、物理問題等本質的數學模型。具體說就是要建立反映問題各量之間的微分方程及相應的定解條件。牛頓流體流動的數學模型就是納維斯托克斯方程及其相應的定解條件。尋求高效率、高準確度的計

3、算方法。計算方法不僅包括數學模型的離散化及求解方法,還包括計算網格的建立、邊界條件的處理。編制程序和進行計算。通過圖像形象地顯示計算結果。數值模擬所得出的結果可對實際流體管道系統的設計、優化等進行指導,具有重要的工程實際意義。1 數學模型1.1 管道非恒定流動動力學方程1對有壓管道中的流動進行分析,出發點仍然是運動微分方程和連續性方程。運動微分方程:21(f p vt+= (1 連續性方程:yx zvv vt x y z+=(2 式中,f質量力;哈密頓算子;p壓力;流體密度;v流體運動粘性系數;速度矢量;t時間。1.2 管道非恒定流動數學模型為簡化動力學方程,作如下假定:系統是彈性的,按分布參

4、數處理,即用流場來描述管道系統;流體是粘性的,切應力按恒定流動時的規律分布,即管道的粘性損失按恒定流動規律計算;對于運動方程,·認為流體可壓縮性對流速影響很小,可以不計;假定流動是一維的。由(1可得:(0p Qf Qx A t+=(3 式中,Q管道流量;A管道通流面積。(f Q為與管壁切應力有關的項,按恒定流動計算則取(22f Q Qf Q DS =由于在一般管道中,壓力波傳遞速度a v ,由此可將(2式簡化得:20p a Q t A x+= (4 聯立方程(3和(4,可得管道非恒定流動數學模型。2 數值模擬計算方法2.1 數學模型的離散常用的數值模擬離散方法主要有有限差分法、有限元

5、法和有限體積法。有限差分法發展較早,比較成熟,在計算流體力學中廣泛使用;有限元法因求解速度較有限差分法和有限體積法慢,在商用流體數值計算軟件中應用并不普遍;有限體積法是一種特殊的有限差分格式。考慮到方程的特殊性,使用特征線法離散數學模型。 2.2 計算方法經特征線解法變換,考慮到壓力波傳遞速度a v ,可將方程(3和(4變換為四個常微分方程(特征線方程:C +21102fQ Qdp dQ a dt S dt DS dx a dt +=+ (5C 21102fQ Qdp dQ a dt S dt DS dx a dt+=(6采用有限差分法,將管道分為N 段,每段長度/x L N =。取時間步長/

6、t x a =。于是x t 平面即分成為矩形網格,如圖1所示。 長方框初始點 黑點邊界點 圓圈界內點圖1 差分計算網格通常,0t =時刻,管道上的參數是已知的或可用恒定流動條件求得,t t =時刻各節點上的參數可通過0t =時刻各節點上的參數值而求得,2t t =時刻各節點上的參數又可通過t 時刻各節點上的參數求得,如此類推。為便于計算機運算,把x t 平面上網格的節點按位置步長序號i 和時間步長j 加以編號。則可得出如下結論(推導過程略去:沿C +11ij i j ij p C BQ +=(7 沿C11ij i j ij p C BQ +=+(8 1111(/2ij i j i j p C

7、C +=+(9 1111(/2ij i j i j Q C C B +=(101111111111i j i j i j i j i j C p BQ RQ Q +=+ (11 1111111111i j i j i j i j i j C p BQ RQ Q +=+ (12式中:/B a S =,2/2R f x DS =,S 管道的截面積。對每一已知時刻(1j t ,可按(11和(12式計算出每一節點i 的兩邊節點上的11i j C +和11i j C +值。再按(9和(10式,即可獲得j t 時刻各節點上的ij p 和ij Q 值。但對兩端的兩個邊界點而言,上游邊界點只能利用前一時刻的C

8、 特性線,下游邊界點則只能利用C +特性線,對每一邊界點都缺少(7或(8二者中的一個,故不能像內節點那樣計算,應另尋找補充方程,即邊界條件2。下面僅舉兩例邊界條件以供參考:邊界條件1:終端閥門關閉時,*0D Q =;*D p 可由(7式求得。邊界條件2:上游端恒壓時,*U p const =,*U Q 可由(10式求得。根據上面的結論,可得簡單管道瞬態響應的特征線法計算程序框圖(圖2。圖2 管道瞬態響應特征線法計算程序框圖 3 計算結果與分析按照上述方法,只要設定合理的邊界條件,可得出簡單管道在各種邊界條件下非恒定流動的數值解。考慮到恒定流動是非恒定流動充分發展的結果,將各種簡單管道及其邊界條

9、件加以組合,設定足夠長的仿真時間,可以得出串聯管道、分支管道、甚至復雜管道網絡情況下的數值解。為驗證上述計算方法的效果,下面針對圖3所示簡單管道非恒定流動系統進行數值模擬。 圖3 簡單管道非恒定流動系統工況一:恒定流動時,溢流閥設定為3MPa,試驗管道管徑為500mm,管道長為400m,入口壓力為3MPa,流量為0.2m3/s。將電磁換向閥2突然換向(流量變為0,將產生水擊現象。此時邊界條件為: 上游端:壓力恒定為3MPa。終端:流量為0。圖4和圖5分別展示試驗系統在不同參數時管道中部壓力和流量變化曲線。從圖4和圖5中可以看出:水擊過程中,產生的最大壓力隨管壁粗糙度的增大而減小,隨管徑的增大而

10、減小。但是,產生的最大流量卻隨管徑的增大而增大。圖4 不同參數時管道中部壓力變化曲線圖5 不同參數時管道中部流量變化曲線工況二:恒定流動時,溢流閥設定為3.6MPa,試驗管道管徑為500mm,管道長為400m,入口壓力為3.6MPa,流量為0.2m3/s。將節流閥2開口調大,使其穩定后流量增大25%。此時邊界條件為: 上游端:壓力恒定為3.6MPa。終端:穩定后流量為0.25 m3/s圖6和圖7分別展示試驗系統在同一參數時管道中部與管道末端壓力和流量變化曲線。 圖6 管道中部和末端壓力變化曲線 圖7 管道中部和末端流量變化曲線從圖6和圖7可以看出:在節流閥開口突然增大的瞬間,管道末端壓力突然減

11、小,但是管道中部的壓力有一段延時;管道末端的流量瞬時增大,與此同時,管道中部流量由于流體慣性,瞬時減小。分析比較發現:簡單管道非恒定流動系統兩種工況數值模擬計算結果與文獻4中相關數據和結論是吻合的。由此可見:數值模擬試驗結果是可信的,對工程實際具有較強的指導作用,且投入少、成本低、周期短。但是,數值試驗不能代替實物試驗,仍需得到實物試驗的驗證。4 結論應用數值分析方法,從流體動力學方程出發,建立起管道非恒定流動數學模型。探討如何運用特征線法將該數學模型離散化,得出相應數值計算方法,并通過實例驗證數學模型的準確性和計算方法的可行性。該數學模型和計算方法可用于編寫大型管道網絡數值模擬軟件。參考文獻

12、1 蔡亦鋼,流體傳輸管道動力學。浙江:浙江大學出版社,19902 蘇爾皇,管道動態分析及液流數值計算方法。哈爾濱:哈爾濱工業大學出版社,19853 王福軍,計算流體動力學分析。北京:清華大學出版社,20044 候國祥等,工程流體力學。武漢:華中科技大學出版社,20055 李紅,數值分析。武漢:華中科技大學出版社,20036 周明,胡斌,計算機仿真原理及其應用。武漢:華中科技大學出版社,20057 劉曉平等,管道水擊數值模擬方法分析。黃河水利職業技術學院學報,2007.7,7118 X.A.Kong, D.P.Chen, Object-oriented design of FEMprograms

13、, Computers and Structures,57(1,p.157-166, 1995. 9 G.R.Mukunda, E.D.Sotelino,S.H.Hsieh, Distributed finiteelement computations using object-oriented techniques, Engineering with Computers, 14(1,p.59-72, 1998RESEARCH ON NUMERICAL SIMULATION OF PIPELINE UNSTEADYFLOWYu Zuyao Zhang Huiyi Li Baoren(FESTO

14、 Pneumatic Centre,Huazhong University of Science and Technology, Wuhan 430074 Abstract:Mathematical model of pipeline unsteady flow is established from the fluid dynamic equation, then corresponding numerical calculation method is discussed, and the numerical simulation is done basing on a simple pipeline system, in which fluid flows unsteadily. Results are analyzed and corresponding con