開式自然循環(huán)閃蒸不穩(wěn)定的線性均相流模型

侯曉凡，孫中寧

（哈爾濱工程大學(xué)核安全與仿真技術(shù)國防重點(diǎn)學(xué)科實(shí)驗(yàn)室，黑龍江哈爾濱150001）

開式自然循環(huán)閃蒸不穩(wěn)定的線性均相流模型

侯曉凡，孫中寧

（哈爾濱工程大學(xué)核安全與仿真技術(shù)國防重點(diǎn)學(xué)科實(shí)驗(yàn)室，黑龍江哈爾濱150001）

針對開式自然循環(huán)閃蒸不穩(wěn)定邊界的判定問題，本文運(yùn)用狀態(tài)空間方法對其進(jìn)行研究，提出了描述該過程的線性均相流模型，導(dǎo)出了判斷系統(tǒng)穩(wěn)定性的狀態(tài)空間表達(dá)式，并求解得到了閃蒸不穩(wěn)定的邊界。模型計算結(jié)果和實(shí)驗(yàn)結(jié)果符合較好，僅在高加熱功率工況下存在一定偏差，而形成偏差的原因是加熱段內(nèi)沸騰噴發(fā)過程對閃蒸過程產(chǎn)生強(qiáng)烈的干擾，使其無法形成穩(wěn)定的閃蒸流動。同時實(shí)驗(yàn)中通過注氣誘導(dǎo)法排除了對閃蒸過程的影響，得到新的閃蒸不穩(wěn)定邊界，與模型計算得到的下邊界符合很好，證明了本模型結(jié)果具有較高的準(zhǔn)確性。

開式自然循環(huán)；閃蒸不穩(wěn)定；線性均相流模型；狀態(tài)空間方程；注氣誘導(dǎo)法；穩(wěn)定閃蒸

流動不穩(wěn)定是制約自然循環(huán)系統(tǒng)安全穩(wěn)定運(yùn)行的重要因素［1］，特別是開式自然循環(huán)系統(tǒng)，由于其系統(tǒng)壓力較低，易于發(fā)生沸騰、流動閃蒸等劇烈的相變，同時缺乏能夠容納系統(tǒng)壓力波動，抑制空泡份額變化及流量振蕩的封閉氣空間，因此開式自然循環(huán)系統(tǒng)更容易發(fā)生流動不穩(wěn)定現(xiàn)象［2］。

在開式自然循環(huán)流動不穩(wěn)定中，閃蒸不穩(wěn)定是重要的不穩(wěn)定類型之一。閃蒸不穩(wěn)定發(fā)生時，上升段內(nèi)出現(xiàn)間歇性的閃蒸流動或者閃蒸起始點(diǎn)沿上升通道周期性地上下移動，從而導(dǎo)致自然循環(huán)流動出現(xiàn)振蕩。目前，針對開式自然循環(huán)閃蒸不穩(wěn)定的研究多為實(shí)驗(yàn)研究［3-6］，理論分析研究也集中于采用數(shù)值解析方法進(jìn)行求解［7-8］，而采用近似分析方法對閃蒸不穩(wěn)定進(jìn)行研究的文獻(xiàn)很少［9］。本文提出了描述開式自然循環(huán)系統(tǒng)閃蒸流動不穩(wěn)定的線性均相流模型，并基于系統(tǒng)控制原理中的狀態(tài)空間分析方法對閃蒸不穩(wěn)定邊界進(jìn)行分析。

1　數(shù)學(xué)模型建立

1.1開式自然循環(huán)閃蒸流動物理過程及基本假設(shè)

圖1給出了開式自然循環(huán)閃蒸流動時，加熱段和上升段內(nèi)的溫度壓力的變化情況。加熱段出口流體處于過冷狀態(tài)，隨著流體向上流動，其承受的重位壓頭和相應(yīng)的飽和溫度均不斷降低，導(dǎo)致流體運(yùn)動至上升段某一位置時達(dá)到當(dāng)?shù)仫柡蜏囟榷焖倨纬善簝上嗔鲃?，該位置稱為閃蒸起始點(diǎn)。

系統(tǒng)出現(xiàn)穩(wěn)定閃蒸流動時，閃蒸起始點(diǎn)固定在上升段的某一位置上，而系統(tǒng)呈現(xiàn)閃蒸不穩(wěn)定特性時，閃蒸起始點(diǎn)間歇性地消失或沿上升段周期性上下移動，當(dāng)閃蒸起始點(diǎn)消失（或下移）時，上升段內(nèi)空泡份額降低，自然循環(huán)驅(qū)動降低，使得自然循環(huán)流量降低，上升段內(nèi)流體溫度增加，進(jìn)而導(dǎo)致閃蒸起始點(diǎn)出現(xiàn)（或上移）。隨后上升段內(nèi)空泡份額升高，循環(huán)流量上升，上升段內(nèi)流體溫度降低，又導(dǎo)致閃蒸起始點(diǎn)再次消失（或下移）。如此反復(fù)形成了自持性的閃蒸流動不穩(wěn)定。

圖1　上升段內(nèi)閃蒸流動示意圖Fig.1 Sketch of flashing flow in upward channel

為了建立自然循環(huán)系統(tǒng)閃蒸不穩(wěn)定的線性均相流模型，本文針對上述物理過程提出如下假設(shè):

1）流體為一維流動，流型的影響可忽略不計；

2）上升段內(nèi)流動過程（包括閃蒸過程）為熱力學(xué)平衡過程，同時忽略加熱段內(nèi)過冷沸騰的影響；

3）考慮到系統(tǒng)內(nèi)溫度和壓力的變化范圍較小，故忽略溫度和壓力變化對于液相密度ρl和氣相密度ρg的影響。

1.2基本控制方程及邊界條件

針對上升段內(nèi)的閃蒸流動過程，列出均相流的質(zhì)量守恒方程和動量守恒方程:

式中:t為時間，ρ為密度，α為空泡份額，G為質(zhì)量流速，z為距加熱段入口距離，P為壓力，D為當(dāng)量直徑，?f為分液相折算系數(shù)，ξ為摩擦阻力系數(shù)，下標(biāo)g表示氣相，下標(biāo)l表示液相。同時列出能量守恒方程時，忽略了外界力做功及動能變化引起的能量變化:

式中:q為線功率密度，A為流道橫截面積，hl和hg分別為液相焓值和氣相焓值，x為干度。

另外，根據(jù)漂移流模型得到上升段兩相區(qū)內(nèi)流體的質(zhì)量含氣率x和空泡份額α之間的關(guān)系，列出輔助方程:

式中:C0為分布參數(shù)；υgj為氣相漂移速度，其具體數(shù)值的選取參照Zuber&Findlay公式［10］。

最后列出了開式自然循環(huán)上升段內(nèi)閃蒸流動的邊界方程，即

加熱段入口位置:

上升段入口位置:

上升段出口位置:

式中:下標(biāo)a為上升段出口位置，b為上升段入口位置，c表示加熱段區(qū)域，e為外界環(huán)境壓力，in和ex分別為入口和出口位置，f表示飽和狀態(tài)的參數(shù)值，xc，in為相對于上升段出口處飽和狀態(tài)的加熱段入口相對干度，hgl為汽化潛熱。

1.3方程的無量綱化

對上述方程中參數(shù)進(jìn)行無量綱化:

式中:L為長度，cp為液相定壓比熱容，ΔTin為加熱段入口過冷度。將上述無量綱參數(shù)代入方程（1）～（3），得到如下的無量綱控制方程:

為了對上述的無量綱控制方程進(jìn)行簡化，對其進(jìn)行量綱分析。結(jié)果發(fā)現(xiàn)由于Rgl數(shù)值接近0.001，導(dǎo)致αRgl等項(xiàng)與其余各項(xiàng)相比要小2～3個數(shù)量級，可以認(rèn)為是小量而忽略。此時控制方程簡化為

同理對式（4）所示的輔助方程和式（5）～（8）所示的邊界方程無量綱化和量級分析，得到無量綱的輔助方程和邊界方程，具體過程不再冗述。

2　系統(tǒng)狀態(tài)空間表達(dá)

2.1線性增量方程的提出

假定系數(shù)具有線性特性，通過外加線性小擾動的方法，對系統(tǒng)方程進(jìn)行線性化，得到描述系統(tǒng)動態(tài)特性的線性增量方程。系統(tǒng)變量采用下式表示:

將式（16）代入式（13）～（15），得到線性增量的控制方程，其中單相區(qū)域的表達(dá)形式如下

而在兩相區(qū)域內(nèi)，線性擾動方程的表達(dá)形式為

同時，將方程（16）代入無量綱的輔助方程，得到線性增量的輔助方程:

另外，將方程（16）代入無量綱的邊界方程中，得到線性增量的邊界方程:

2.2狀態(tài)空間方程組的建立

同時為了使變量和方程的數(shù)量相同，需要將線性增量的輔助方程式（24）代入控制方程，用以消去Δα這個變量，最終得到一組3×3的線性常微分方程組，即描述自然循環(huán)閃蒸不穩(wěn)定動態(tài)特性的狀態(tài)空間表達(dá)式:

其中，矩陣M和N的參數(shù)如下所示，為了方便表述，首先定義如下的參數(shù):

此時矩陣M的參數(shù)可以表示為

而矩陣N的參數(shù)則可以表示為

2.3系統(tǒng)空間狀態(tài)方程穩(wěn)定域的判斷

對于方程（28）所示的狀態(tài)空間方程組，可以通過判斷系數(shù)矩陣N（=MTK）的特征根實(shí)部的符號來判斷系統(tǒng)的穩(wěn)定性。即求解如下的特征方程:

其特征根為

當(dāng)γi均小于零時，系統(tǒng)可以維持穩(wěn)定，即表示開式自然循環(huán)系統(tǒng)能夠形成穩(wěn)定的閃蒸；而 γi存在大于零的情況時，系統(tǒng)無法維持穩(wěn)定，即此時的閃蒸流動為閃蒸不穩(wěn)定工況。

3　模型計算結(jié)果及與實(shí)驗(yàn)結(jié)果比較

由于上述計算過程繁復(fù)，較難推導(dǎo)出閃蒸不穩(wěn)定區(qū)域所對應(yīng)的數(shù)學(xué)表達(dá)形式，因此本文利用MATLAB 7.0程序?qū)Ψ匠烫卣髦颠M(jìn)行求解和對系統(tǒng)穩(wěn)定性進(jìn)行判斷，其計算過程如圖2所示。

圖3給出了模型的計算結(jié)果，圖中實(shí)線表示計算得到的閃蒸不穩(wěn)定邊界，數(shù)據(jù)點(diǎn)表示開式自然循環(huán)實(shí)驗(yàn)臺上得到的閃蒸不穩(wěn)定區(qū)域（實(shí)驗(yàn)裝置見圖4，實(shí)驗(yàn)方法見文獻(xiàn)［11］）。通過對比模型計算結(jié)果與實(shí)驗(yàn)結(jié)果發(fā)現(xiàn)，兩者之間趨勢上吻合較好。

其中在高加熱功率工況下，計算得到的不穩(wěn)定上邊界與實(shí)驗(yàn)結(jié)果存在一定偏差。此時實(shí)驗(yàn)觀察到加熱管出口過冷沸騰現(xiàn)象強(qiáng)烈而形成間歇噴泉流動，而本模型中未考慮過冷沸騰現(xiàn)象對閃蒸流動的影響，因此偏差產(chǎn)生的原因極可能是由于過冷沸騰過程的影響。

圖2　閃蒸不穩(wěn)定邊界的計算流程Fig.2 Calculation flow of flashing instability boundary

圖3　模型計算結(jié)果與實(shí)驗(yàn)結(jié)果的對比Fig.3 Comparison between the model calculation results and the experiment results

同樣在高加熱功率工況下，計算得到的不穩(wěn)定下邊界與實(shí)驗(yàn)結(jié)果也存在偏差。此時隨著加熱段入口過冷度降低，使得加熱段內(nèi)沸騰噴發(fā)過程增強(qiáng)，并與上升段內(nèi)閃蒸流動過程共同影響自然循環(huán)的流動特性［11］。而沸騰噴發(fā)過程會導(dǎo)致流量、上升段內(nèi)的溫度分布等出現(xiàn)劇烈地振蕩，進(jìn)而導(dǎo)致持續(xù)的穩(wěn)定閃蒸流動無法形成，即可以認(rèn)為加熱段內(nèi)沸騰噴發(fā)過程對于穩(wěn)定閃蒸過程起到強(qiáng)烈的干擾作用。

因此本實(shí)驗(yàn)采用注氣誘導(dǎo)法排除沸騰噴發(fā)對于穩(wěn)定閃蒸過程的干擾作用。即在閃蒸不穩(wěn)定工況時，向上升段內(nèi)持續(xù)注入氣體，使得系統(tǒng)形成穩(wěn)定的空氣-水兩相流動，此時自然循環(huán)流量快速增加，加熱段內(nèi)沸騰噴發(fā)過程被抑制。隨后逐漸減少注氣量直至停止，這一過程中閃蒸流動過程會逐漸主導(dǎo)系統(tǒng)的流動特性，而由于循環(huán)流量較高導(dǎo)致該過程中加熱段內(nèi)沸騰噴發(fā)不會出現(xiàn)。通過這種方法來排除沸騰噴發(fā)對閃蒸過程的影響，再通過觀察停止注氣后自然循環(huán)閃蒸流動的發(fā)展演化過程來判斷閃蒸流動是否穩(wěn)定。圖5給出了典型的注氣誘導(dǎo)法過程及其結(jié)果。

圖4 實(shí)驗(yàn)回路示意圖Fig.4 Sketch map of experimental loop

圖5（a）所示的實(shí)驗(yàn)工況位于模型計算的閃蒸不穩(wěn)定區(qū)域內(nèi)（圖3的不穩(wěn)定區(qū)域），注氣前系統(tǒng)處于閃蒸和沸騰共同誘發(fā)的流動不穩(wěn)定狀態(tài)。當(dāng)注氣逐漸停止后，加熱管內(nèi)沸騰始終被抑制，而上升段內(nèi)閃蒸流動出現(xiàn)振蕩，振幅不斷增大直至最終恢復(fù)為原來的流動不穩(wěn)定狀態(tài)，即說明此工況下閃蒸流動過程為發(fā)散過程。而圖5（b）所示實(shí)驗(yàn)工況則處于模型計算的穩(wěn)定閃蒸流動區(qū)域（圖3的下部穩(wěn)定區(qū)域），注氣前系統(tǒng)與工況1相似，但注氣逐漸停止后，上升段內(nèi)閃蒸流動首先出現(xiàn)較弱振蕩，隨后振蕩逐漸減弱并形成穩(wěn)定閃蒸流動，此時穩(wěn)定閃蒸驅(qū)動的流量很大，使得加熱段出口為單相流動，即工況2時自然循環(huán)閃蒸流動過程為自發(fā)收斂過程，也證實(shí)了圖3實(shí)驗(yàn)中未出現(xiàn)穩(wěn)定閃蒸流動主要是由于沸騰噴發(fā)過程的干擾作用。

圖6給出了通過重復(fù)實(shí)驗(yàn)得到的注氣誘導(dǎo)法形成的穩(wěn)定閃蒸流動邊界，本文認(rèn)為是準(zhǔn)確的閃蒸不穩(wěn)定的區(qū)域下邊界，由圖可知該邊界與計算得到的不穩(wěn)定下邊界符合很好，即證明了本文模型結(jié)果具有很好的準(zhǔn)確性。另外，可以認(rèn)為本模型計算得到下穩(wěn)定區(qū)域內(nèi)的工況，都可以通過注氣誘導(dǎo)法將閃蒸不穩(wěn)定流動轉(zhuǎn)化為持續(xù)穩(wěn)定閃蒸流動。

圖5　注氣誘導(dǎo)法誘發(fā)穩(wěn)定閃蒸流動的過程Fig.5 Flow rate variations of the method of gas injection to induce the steady flashing flow

圖6　注氣誘導(dǎo)法穩(wěn)定閃蒸邊界與模型計算結(jié)果的對比Fig.6 Comparison between model calculation results and experimental instability boundary by method of gas injection

4　結(jié)論

1）本文運(yùn)用狀態(tài)空間法對開式自然循環(huán)閃蒸不穩(wěn)定進(jìn)行研究，提出了描述該過程的線性均相流模型，并導(dǎo)出了描述系統(tǒng)穩(wěn)定性的狀態(tài)空間方程。2）通過求解狀態(tài)空間方程得到了自然循環(huán)閃蒸不穩(wěn)定的邊界，將計算結(jié)果與實(shí)驗(yàn)結(jié)果對比發(fā)現(xiàn)兩者趨勢上吻合較好。其中上邊界的偏差是由于模型未考慮過冷沸騰的影響而造成的；下邊界的偏差是由加熱段內(nèi)沸騰噴發(fā)對穩(wěn)定閃蒸的干擾作用而產(chǎn)生的。3）通過注氣誘導(dǎo)法排除了沸騰噴發(fā)干擾作用，得到了準(zhǔn)確的閃蒸不穩(wěn)定下邊界，該邊界與計算得到的不穩(wěn)定下邊界吻合很好。4）通過實(shí)驗(yàn)發(fā)現(xiàn)，計算得到的下穩(wěn)定區(qū)域內(nèi)的工況點(diǎn)，均可通過注氣誘導(dǎo)法將閃蒸不穩(wěn)定工況轉(zhuǎn)化為穩(wěn)定閃蒸流動工況。

［1］譚思超，龐鳳閣，高璞珍.自然循環(huán)過冷沸騰流動不穩(wěn)定

性實(shí)驗(yàn)研究［J］.核動力工程，2006，27（1）:18-21，93. TAN Sichao，PENG Fengge，GAO Puzhen.Experimental research on subcooled boiling flow instability of natural circulation［J］.Nuclear power engineering，2006，27（1）:18-21，93.

［2］KYUNG I S，LEE S Y.Periodic flow excursion in an open two-phase natural circulation loop［J］.Nuclear engineering and design，1996，162（2/3）:233-244.

［3］YANG S K.Stability of flashing-driven natural circulation in a passive moderator cooling system for Canadian SCWR［J］. Nuclear engineering and design，2014，276:259-276.

［4］吳莘馨，吳少融，姜勝耀，等.閃蒸現(xiàn)象及其對自然循環(huán)穩(wěn)定性的影響［R］.中國核科技報告，1991. WU Xinxin，WU Shaorong，JIANG Shengyao，et al.Flashing in riser and its effect on flow stability in natural circulation system［R］.China nuclear science and technology report，1991.

［5］FURUYA F，INADA F，VAN DER HAGEN T H J J.Flashing-induced density wave oscillations in a natural circulation BWR-mechanism of instability and stability map［J］.Nuclear engineering and design，2005，235（15）:1557-1569.

［6］SAJITH P P，VEDULA R P，MITRA S K.An Investigation of flashing-driven natural circulation［J］.International journal of green energy，2006，3（4）:369-379.

［7］GUO Xueqing，SUN Zhongning，WANG Jianjun，et al.Numerical simulation of the transient behaviors in an open natural circulation system with a large scale［J］.Annals of nuclear energy，2015，77:83-93.

［8］GUO Xueqing，SUN Zhongning，WANG Jianjun，et al.Simulating investigations on the start-up of the open natural circulation system［J］.Nuclear engineering and design，2015，289:35-48.

［9］INADA F，F(xiàn)URUYA M，TASUO A.Thermo-hydraulic instability of boiling natural circulation loop induced by flashing （analytical consideration）［J］.Nuclear engineering and design，2000，200（1/2）:187-199.

［10］ZUBER N，F(xiàn)INDLAY J A.Average volumetric concentration in two-phase flow systems［J］.Journal of heat transfer，1965，87（4）:453-468.

［11］HOU Xiaofan，SUN Zhongning，F(xiàn)AN Guangming，et al. Flow characteristics in an open two-phase natural circulation loop［C］//Proceedings of the 22nd International Conference on Nuclear Engineering（ICONE 22）.Prague，Czech Republic，2014.

本文引用格式：

侯曉凡，孫中寧.開式自然循環(huán)閃蒸不穩(wěn)定的線性均相流模型［J］.哈爾濱工程大學(xué)學(xué)報，2016，37（7）:930-935.

HOU Xiaofan，SUN Zhongning.Linear homogeneous-flow model of flashing instability in open natural circulation system［J］.Journal of Harbin Engineering University，2016，37（7）:930-935.

Linear homogeneous-flow model of flashing instability in open natural circulation system

HOU Xiaofan，SUN Zhongning
（Fundamental Science on Nuclear Safety and Simulation Technology Laboratory，Harbin Engineering University，Harbin 150001，China）

With respect to judging the flashing instability boundary in an open natural circulation system，we use a state-space method and develop a linear homogeneous-flow model to describe its physical process.From the statespace analysis results，we deduce a state space expression for estimating the stability of the system.We then determine the flashing instability boundary by solving these equations.The calculation result shows good agreement with the experimental result，with only a slight deviation under the condition of high heating rate.This deviation is mainly due to the interference effect of a boiling eruption in the steady flashing flow process in high heating power conditions，and this interference makes it difficult to establish a stable flashing flow.Next，we use a gas-injection method on the up-leg to avoid this interference effect in the experiments，and obtain an accurate lower boundary for the flashing instability.This boundary agrees well with the calculated lower boundary，which verifies the accuracy of this model.

open natural circulation；flashing instability；linear homogeneous flow model；state-space equations；gas injection method；steady flashing

10.11990/jheu.201506093

TL334

A

1006-7043（2016）07-930-06

2015-06-30.網(wǎng)絡(luò)出版日期：2016-05-06.

國家高科技研究發(fā)展計劃項(xiàng)目（2012AA050906）.

侯曉凡（1988-），男，博士研究生；孫中寧（1963-），男，教授，博士生導(dǎo)師.

孫中寧，E-mail:sunzhongning@hrbeu.edu.cn.

網(wǎng)絡(luò)出版地址：http：//www.cnki.net/kcms/detail/23.1390.u.20160506.1548.002.html