基于Godunov格式的有壓管道中水柱分離的模擬方法
【技術(shù)領(lǐng)域】
[0001 ]本發(fā)明涉及一種基于Godunov格式的有壓管道中水柱分離的模擬方法，屬于水電 站(栗站)水力學(xué)數(shù)值模擬計(jì)算技術(shù)領(lǐng)域。
【背景技術(shù)】
[0002] 當(dāng)輸水管道系統(tǒng)中液體壓力降低到汽化壓力時(shí)，可能發(fā)生瞬變蒸汽空化現(xiàn)象。這 些空化泡可能發(fā)育長大到占據(jù)整個(gè)管道截面的程度，破壞液柱的連續(xù)性。Bergant和 Simpson等學(xué)者在文獻(xiàn)中指出分離的液柱再彌合可能產(chǎn)生突升的高壓脈沖，導(dǎo)致水力系統(tǒng) 破壞。因此，為避免水柱分離現(xiàn)象，需對該現(xiàn)象進(jìn)行嚴(yán)格的瞬變分析。
[0003] 幾十年來，很多學(xué)者采用各種方法來預(yù)測水柱分離過程中的瞬變壓力，應(yīng)用最廣 的是由Streeter和Wylie基于特征線法(M0C，Method of Characteristies)提出的離散蒸 汽一空穴模型(DVCM,discrete vapor-cavity model)，該模型允許汽穴在壓力達(dá)到或低于 汽化壓力的任何計(jì)算節(jié)點(diǎn)形成，并且作為一內(nèi)部邊界占據(jù)整個(gè)管道截面，因該模型易于嵌 入到標(biāo)準(zhǔn)的水錘商業(yè)包中，且再現(xiàn)了水柱分離過程中很多的物理特性，成為目前模擬水柱 分離現(xiàn)象最常用的模型。
[0004] 經(jīng)典M0C-DVCM的主要缺點(diǎn)是當(dāng)采用相對細(xì)密的網(wǎng)格時(shí)，計(jì)算結(jié)果會出現(xiàn)不真實(shí)的 壓力脈沖。Chaiko在FVM方法基礎(chǔ)上采用均質(zhì)連續(xù)模型來描述水柱分離過程中的汽穴過程， 然而，其算法相比M0C而言過于復(fù)雜，實(shí)現(xiàn)困難。因此，如果能夠合理保留各種方法的優(yōu)點(diǎn)， 引入新的方法和技巧消除缺陷，就能有助于水柱分離過程中瞬變壓力的數(shù)值模擬方法的不 斷改進(jìn)。

【發(fā)明內(nèi)容】

[0005] 發(fā)明目的：為彌補(bǔ)現(xiàn)有M0C-DVCM方法存在的虛假震蕩等不足，本發(fā)明基于FVM，提 供了一種算法簡單，易于實(shí)現(xiàn)的模擬方法，可以滿足計(jì)算網(wǎng)格細(xì)密時(shí)不會出現(xiàn)人為的脈沖 峰值。
[0006] 技術(shù)方案:一種基于Godunov格式的有壓管道中水柱分離的模擬方法，具體步驟如 下：
[0007] 步驟1:建立瞬變流基本微分方程，根據(jù)工程實(shí)例確定計(jì)算域、初始條件以及邊界 條件；
[0008] 步驟2:根據(jù)FVM劃分計(jì)算網(wǎng)格，并建立離散方程；
[0009] 步驟3:采用Godunov方法求解純對流時(shí)的離散方程，得到純對流時(shí)控制單元界面 處的數(shù)值通量，并取得二階精度；
[0010]步驟4:通過基于二階Runge-Kutta離散格式的時(shí)間算子分裂法，在純對流控制方 程的解中引入源項(xiàng)，從而得到最終解的二階顯式FVM-Godunov格式；
[0011]步驟5:根據(jù)計(jì)算得到的壓力判定是否形成水柱分離:若已形成，則進(jìn)行水柱分離 模塊計(jì)算;反之，進(jìn)入下一時(shí)步計(jì)算。
[0012] 進(jìn)一步地，步驟2中，將空間域x離散為長度為Δ χ，時(shí)間域t離散為間隔為Δ t的控 制單元，數(shù)量為N。對第i個(gè)控制單元，定義其上、下游界面編號分別為i_l/2、i+l/2。
[0013] 進(jìn)一步地，步驟2中，對控制單元i，建立的流動變量u的積分方程為：
[0014]
(I)
[0015] 其中，上標(biāo)η和n+1分別代表t和t+Δ t時(shí)步；:1? 為u在整個(gè)控制體的 平均值；u = Η是測壓管水頭，V是平均截面速率;f為單元界面處的通量： 為源項(xiàng);f為達(dá)西-威斯巴哈摩阻系數(shù);D為管徑。
[0016] 進(jìn)一步地，步驟3包含以下子步驟：
[0017] 步驟3.1:求解內(nèi)部控制單元界面處通量。
[0018] 首先，基于黎曼問題，根據(jù)Godunov格式，對任一內(nèi)部控制單元i (l〈i〈N)，界面i+1/ 2處的通量為：
[0019]
[0020]其=
?測壓管水頭;V是平均截面速率；P是V的平均 值，為一常數(shù)；砒、邱為在η時(shí)步時(shí)，u分別到界面i+1/2左、右側(cè)兩側(cè)的平均值。
[0021 ] 接著，通過引入MUSCL-Hancock格式計(jì)算內(nèi)部單元通量f i+i/2，從而取得空間和時(shí)間 上的二階精度。
[0022] 更進(jìn)一步地，步驟3.1中，在引入MUSCL-Hancock格式計(jì)算內(nèi)部單元通量fi+i/2過程 中，需選擇斜率限制器，以保證解中不出現(xiàn)虛假振蕩。
[0023] 步驟3.2:求解邊界控制單元界面處通量。為在邊界面處也取得二階精度，分別在 起始控制單元1上游側(cè)、終點(diǎn)控制單元N下游側(cè)構(gòu)建兩個(gè)虛擬控制單元Ι-^Ιο,以及I N+1、IN+2， 并假定在虛擬單元處的流動信息與邊界處是一致的。從而可求解邊界黎曼問題，且相應(yīng)的 Godunov通量f 1/2和ft+i/2也可像內(nèi)部單元那樣進(jìn)行計(jì)算。
[0024] 更進(jìn)一步地，步驟3中，純對流時(shí)離散方程中的對流項(xiàng)滿足CFL條件（(：〇1^&1^-Friedr^-T η'ιτ,τ …1·μ、）肚一最大時(shí)間步長 Δ tmax,cFL:
[0025] ?3)
[0026] 其中，Cr為柯朗數(shù)，1為矩陣S的特征值。
[0027] 進(jìn)一步地，步驟4中，引入源項(xiàng)后，瞬變流基本微分方程解的二階FVM-Godunov格式 為：
[0028] (4)
[0029] 其中，砰+1為n+1時(shí)步，控制單元i在純對流時(shí)，流動變量u的通量；0f+1為采用時(shí)間 分裂法第一次更新后的通量。
[0030] 更進(jìn)一步地，步驟4中，源項(xiàng)滿足以下穩(wěn)定性約束，并可推得適用于源項(xiàng)的最大時(shí) 間步長八tmax,s:
[0031]
(5)
[0032] 更進(jìn)一步地，包含對流項(xiàng)和源項(xiàng)的最大允許時(shí)間步長為：
[0033] Δ tmax = min( Δ tmax,CFL，A tmax,s) (6)
[0034] 進(jìn)一步地，步驟5中，水柱分離計(jì)算模塊中，假定汽穴僅在控制體的中部形成，為計(jì) 算汽穴體積，將每個(gè)控制體等分為兩個(gè)相同的控制體。當(dāng)為純單相流動時(shí)，第i個(gè)控制體內(nèi) 兩部分的壓力和流量均由FVM計(jì)算。一旦任一小控制體壓力達(dá)到或低于液體汽化壓力，則兩 部分壓力設(shè)為汽化壓力水頭，計(jì)算體內(nèi)流量則依然采用Godunov格式求得。汽穴體積由汽穴 上下游側(cè)的流量差值計(jì)算：
[0035] (7)
[0036]共T，設(shè)、Vn別73n+iH、」2FH、」，第i個(gè)控制體內(nèi)，上、下游小控制體內(nèi)的流量。
[0037] 更進(jìn)一步地，步驟5中，一旦汽穴體積小于或等于0，則汽穴潰滅，流動恢復(fù)為純單 相流動，壓力和流量再次均由二階FVM-Godunov格式計(jì)算。
[0038] 有益效果:與現(xiàn)有技術(shù)相比，本發(fā)明具有如下優(yōu)點(diǎn)：
[0039] (1)本發(fā)明提供的FVM-DVCM方法成功地克服了 FVM方法在追蹤汽穴并預(yù)測其發(fā)育、 潰滅方面的難題，并且同M0C-DVCM方法一樣簡單且易于實(shí)現(xiàn)；（2)當(dāng)采用相對細(xì)密的網(wǎng)格 時(shí)，結(jié)果中不會出現(xiàn)人為的脈沖峰值；（3)非線性對流項(xiàng)很容易加入到解中，而這些項(xiàng)在 M0C-DVCM中通常是忽略的，因此它可將該模型的應(yīng)用范圍擴(kuò)展到包含馬赫數(shù)的問題；（4)該 計(jì)算方法通過時(shí)間算子分裂方法，為模擬多維水錘瞬變提供了框架，而M0C解較難擴(kuò)展到包 含多維瞬變的流體問題。
【附圖說明】
[0040]圖1為本發(fā)明的基本流程圖；
[0041 ]圖2為實(shí)施例的水庫-管道-閥門-水庫系統(tǒng)示意圖；
[0042] 圖3為實(shí)施例的水錘區(qū)的網(wǎng)格離散系統(tǒng)圖；
[0043] 圖4為實(shí)施例的水柱分離區(qū)的網(wǎng)格離散系統(tǒng)圖；
[0044]圖5為實(shí)施例下，網(wǎng)格數(shù)為32時(shí)，閥門末端的壓力曲線圖；
[0045]圖6為實(shí)施例下，網(wǎng)格數(shù)為256時(shí)，閥門末端的壓力曲線圖；
[0046]圖2中：1_上游水庫;2-輸水管;3-球閥；4-下游水庫。
【具體實(shí)施方式】
[0047] 下面結(jié)合具體實(shí)施例，進(jìn)一步闡明本發(fā)明，應(yīng)理解這些實(shí)施例僅用于說明本發(fā)明 而不用于限制本發(fā)明的范圍，在閱讀了本發(fā)明之后，本領(lǐng)域技術(shù)人員對本發(fā)明的各種等價(jià) 形式的修改均落于本申請所附權(quán)利要求