董曉梅，董俊華，余雛麟，高炳軍

(河北工業(yè)大學(xué) 化工學(xué)院，天津 300130)

全球許多核電廠報(bào)告的事件表明，管道中的熱分層會(huì)導(dǎo)致管道材料產(chǎn)生過度的熱應(yīng)力和疲勞，可能導(dǎo)致嚴(yán)重的熱疲勞開裂和管道系統(tǒng)的泄漏事故[1]。盡管如此，老式核電站的設(shè)計(jì)并沒有考慮這種熱效應(yīng)[2]。導(dǎo)致意外熱疲勞的熱分層現(xiàn)象是由冷熱流體之間的密度差異產(chǎn)生的[3，4]。一種典型的熱分層常發(fā)生在連接到反應(yīng)堆冷卻系統(tǒng)(RCS)的安全注入系統(tǒng)(RIS)中，是由于RCS向安注管線的湍流滲入以及支管止逆閥內(nèi)漏造成的[5，6]。安全注入系統(tǒng)僅在核電廠發(fā)生事故時(shí)向核反應(yīng)堆供應(yīng)應(yīng)急水冷卻時(shí)才運(yùn)行，在正常工作狀態(tài)下，由于RIS中的所有閥門都關(guān)閉，RIS的冷卻劑是隔離的。但是，如果連接硼水箱和一回路主管道支管的閥門可能會(huì)出現(xiàn)微量的滲漏，冷水將進(jìn)人高溫水側(cè)的管道[7]，并在一定范圍內(nèi)形成熱分層。泄漏流量的準(zhǔn)確確定對(duì)于安注管線的熱疲勞分析至關(guān)重要，而由于核管道安注管線流體流動(dòng)的特殊性，對(duì)于存在冷熱流體熱分層、湍流滲入、冷熱交替等復(fù)雜流體流動(dòng)與傳熱現(xiàn)象的問題，目前尚缺乏有效的泄漏監(jiān)測(cè)手段。本研究中，通過流固耦合計(jì)算分析評(píng)估在RIS管線內(nèi)漏條件下的熱分層特征，提取止逆閥閥前管道外壁溫度并定義熱分層特征溫度參數(shù)，得出該特征溫度參數(shù)與泄漏量的定量關(guān)系式。 在此，僅以熱安注管線為例進(jìn)行分析討論。

1 數(shù)值計(jì)算模型

1.1 物理模型

一回路熱段管內(nèi)徑為736.6 mm，壁厚為67 mm，長度為4 904 mm，水平放置。安注管內(nèi)徑為131.78 mm，壁厚為18.26 mm，通過管嘴與水平面呈30°和主管道熱段連接，其上安裝有止逆閥。利用SolidWorks建立含有主管道熱段和安注管的幾何模型，支管末端閥門泄漏的區(qū)域是研究的重點(diǎn)，在閥門中建立月牙形的開孔模擬閥門泄漏，如圖1所示，其中閥門開度與泄漏流速的關(guān)系由公式(1)表示[8]。

(1)

圖1 幾何模型Fig.1 Geometric model

式中：v——泄漏流體速度，m/s；

x——閥門開度。

1.2 網(wǎng)格劃分與質(zhì)量檢查

用ANSYS ICEM劃分流體域和固體域的網(wǎng)格，流體域和固體域網(wǎng)格導(dǎo)入Fluent進(jìn)行流場(chǎng)、溫度場(chǎng)計(jì)算。止逆閥和止逆閥內(nèi)的流體域使用四面體劃分非結(jié)構(gòu)網(wǎng)格，其他流體域和固體域均采用六面體網(wǎng)格，劃分結(jié)構(gòu)化網(wǎng)格。為了盡可能真實(shí)地模擬流體混合過程，本文對(duì)網(wǎng)格模型進(jìn)行了網(wǎng)格無關(guān)性驗(yàn)證，計(jì)算時(shí)，先對(duì)模型進(jìn)行初步網(wǎng)格劃分，在此基礎(chǔ)上，對(duì)網(wǎng)格進(jìn)行加密，直到網(wǎng)格的大小對(duì)計(jì)算結(jié)果影響忽略不計(jì)時(shí)，才最終的確立網(wǎng)格大小和分布。網(wǎng)格Quality質(zhì)量低于0.2的網(wǎng)格處在結(jié)構(gòu)和非結(jié)構(gòu)網(wǎng)格的交界面處，對(duì)流場(chǎng)和溫度場(chǎng)影響較小。最終網(wǎng)格如圖2所示，網(wǎng)格的數(shù)量為75萬個(gè)，已具有良好的網(wǎng)格無關(guān)性。

圖2 網(wǎng)格模型Fig.2 Meshof the model

1.3 計(jì)算方法

采用三維全尺寸穩(wěn)態(tài)流固耦合傳熱的計(jì)算流體動(dòng)力學(xué)(CFD)方法對(duì)安注管線發(fā)生熱分層現(xiàn)象時(shí)的流動(dòng)與傳熱進(jìn)行數(shù)值研究。與以往的k-ε湍流模型相比，將計(jì)算控制區(qū)域擴(kuò)展到安注管管壁所在的固體區(qū)域，采用更適合于求解具有流動(dòng)分離情況的剪切應(yīng)力傳輸模型Transition SST湍流模型[9]。選擇三維、雙精度、基于壓力隱式求解器，求解控制方程包括連續(xù)性方程，質(zhì)量、能量守恒方程及湍流輸運(yùn)方程。并將密度處理為溫度的線性函數(shù)來估算浮力項(xiàng)的影響。

1.3.1 控制方程

(1)全浮力模型

管內(nèi)流體溫差較大，且工作在高溫、高壓下，密度隨溫度的變化明顯。因此，本文采用全浮力模型進(jìn)行計(jì)算。在包括浮力的計(jì)算中，動(dòng)量方程需增加1個(gè)源項(xiàng)[10]。全浮力模型則根據(jù)求解的實(shí)際密度通過公式(2)直接對(duì)浮力進(jìn)行評(píng)估，進(jìn)而求解全場(chǎng)浮力。

SM,buoy=(ρ-ρref)g

(2)

式中：SM,buoy——浮力源項(xiàng)，N；

ρ——密度，kg/m3；

ρref——參考密度，kg/m3；

g——重力加速度，m/s2。

(2)連續(xù)性方程

(3)

式中：t——時(shí)間，s；

xj(j=1，2，3)——直角坐標(biāo)分量，m；

uj(j=1，2，3)——速度分量，m/s。

(3)動(dòng)量守恒方程

(4)

(5)

式中：μ——黏性系數(shù)，kg/(m·s)；

p——壓力，Pa；

(4)能量守恒方程

(6)

式中：cp——比定壓熱容，J/(kg·K)；

λ——導(dǎo)熱系數(shù)，W/(m·K)；

S——能量方程的源項(xiàng)。

(5)k和ω輸運(yùn)方程

(7)

(8)

式中：Gk，Gω——名義速度梯度下k的產(chǎn)生項(xiàng)，ω產(chǎn)生項(xiàng)；

Γk，Γω——k和ω的有效擴(kuò)散率；

Yk，Yω——k和ω的有效耗散率；

Sk，Sω——用戶自定義源項(xiàng)。

Transition SST數(shù)學(xué)方程模型在k-ωSST 的基礎(chǔ)上額外增加了兩個(gè)輸運(yùn)方程，分別是間歇度γ和瞬態(tài)起始標(biāo)準(zhǔn)，可表示為：

(9)

(10)

有效擴(kuò)散率Γk和Γω可表示為

(11)

式中：σk和σω——k和ω的湍流普朗特?cái)?shù)；

μt——湍流黏度。

(12)

式中：ɑ*可在低雷諾數(shù)時(shí)限制湍流粘度的大小，ɑ*的定義如下：

(13)

(14)

1.3.2 邊界條件

由于主管和熱安注管外側(cè)都有保溫層，所以主管道和熱安注管的外壁設(shè)為絕熱邊界。忽略管道端部橫截面?zhèn)鳠崃?，因此也設(shè)為絕熱邊界。流體域和固體域的交界面設(shè)為耦合壁面邊界條件。

熱安注管進(jìn)口設(shè)為速度進(jìn)口邊界；主管道的流量按功率運(yùn)行時(shí)的流量確定，主管道進(jìn)口也設(shè)為速度入口邊界條件；壓力出口，壓強(qiáng)為15.5 MPa。為了避免進(jìn)口段對(duì)混合過程模擬的影響，消除邊界條件對(duì)流場(chǎng)的干擾，同時(shí)減小計(jì)算量，通過fluent的二次開發(fā)UDF，采用DEFINE_ PROFILE編寫，將發(fā)展完全后的速度分布添加到主、支管的進(jìn)口邊界條件中。

1.3.3 物性參數(shù)

本文通過自定義材料方式實(shí)現(xiàn)材料參數(shù)的定義，熱安注管和主管材料為Z2 CN 18-10，假定流體為不可壓縮流體，物性參數(shù)采用變物性。金屬材料屬性取自RCC-M規(guī)范[11]，對(duì)于流體的物理性質(zhì)，使用基于技術(shù)數(shù)據(jù)中的物性數(shù)據(jù)將相關(guān)性作為溫度的函數(shù)進(jìn)行計(jì)算，材料屬性隨溫度變化規(guī)律如表1所示。

表1 材料參數(shù)表Table 1 Material parameters

1.3.4 模型驗(yàn)證

為了驗(yàn)證計(jì)算模型的準(zhǔn)確性，利用本文的數(shù)值計(jì)算方法來模擬Nobuo Nakamori試驗(yàn)工況[12]。得到水平支管某剖面上豎直方向高度H與溫度分布的無量綱關(guān)系的模擬數(shù)據(jù)與試驗(yàn)數(shù)據(jù)對(duì)比結(jié)果如圖 3 所示，模擬溫度與試驗(yàn)數(shù)據(jù)較好吻合，相對(duì)誤差小于 2%，表明本文所采用的數(shù)值計(jì)算方法的準(zhǔn)確性。圖3中D為支管直徑；T為剖面上某點(diǎn)的溫度；Tmin為剖面上的最低溫度；Tmax為剖面上的最高溫度。

圖3 某剖面高度與溫度分布的無量綱關(guān)系Fig.3 No-dimensional relation of temperature and Height in section plane

2 計(jì)算方案與結(jié)果分析

2.1 計(jì)算方案

核電廠主管道內(nèi)熱水參數(shù)受用電負(fù)荷影響，是不穩(wěn)定的。根據(jù)閥前主管道流體溫度Thot、流量Qhot及閥后管內(nèi)泄漏工質(zhì)的溫度Tcold范圍，主管道內(nèi)熱流體溫度Thot取564～601 K；流量Qhot范圍為90%～ 110 %Q額(Q額=2.379×107L/h)；安注管泄漏的冷流體溫度Tcold范圍為293～328 K；泄漏量QLeak取5～400 L/h。得出不同工況下的溫度場(chǎng)并獲取閥前監(jiān)測(cè)位置管外壁溫度數(shù)據(jù)，每一個(gè)工況對(duì)應(yīng)一個(gè)監(jiān)測(cè)位置的管壁溫度分布。采用單一因子變量法著重分析和探討了一回路主管道與安注管工質(zhì)參數(shù)是如何影響止逆閥前管外壁溫度的。依據(jù)數(shù)學(xué)建模方法，采用非線性回歸模型建立管外壁溫度與影響因素之間的函數(shù)關(guān)系式。

2.2 計(jì)算結(jié)果與分析

Thot=599.15 K、Qhot=2.379×107L/h、Tcold=327.15 K、QLeak=400 L/h工況下管外壁溫度分布云圖如圖4 (a)所示，可見近主管道區(qū)域受湍流滲入的影響管外壁溫度均勻，而靠近閥門較大的區(qū)域內(nèi)由于流體熱分層管壁溫度上下存在溫差。不同圓周位置上管外壁溫度沿管長的變化如圖4(b)所示，可見管外壁溫度自主管道開始沿管長都呈下降趨勢(shì)，但湍流滲入?yún)^(qū)圓周方向無明顯溫差，熱分層區(qū)有溫差，較大范圍內(nèi)溫差分布規(guī)律相似。

圖4 管外壁溫度分布Fig.4 Temperature distribution of the pipe outer surface

不同泄漏流量下閥前50 mm位置橫截面流體與管壁的溫度云圖如圖5所示。由圖可見在所討論的泄漏流量范圍內(nèi)流體均存在較為明顯的熱分層現(xiàn)象。上部溫度較高，下部溫度較低，最高溫度在最上部，最低溫度在最下部。管壁受分層流體的作用，也呈現(xiàn)相應(yīng)的溫度變化，即上部溫度高下部溫度低，最高溫度在最上部，最低溫度在最下部。

圖5 不同泄漏流量下截面溫度分布Fig.5 Temperature distribution of cross section under different leakage flow rates

圖6為截面直徑線位置溫度分布曲線，最高最低位置金屬管壁內(nèi)外溫差很小。直徑線位置流體溫度從高到低是逐漸變低的，熱分層是漸變的，這是由于流體沿垂直方向的擴(kuò)散與熱傳導(dǎo)造成的。管壁與流體最高最低點(diǎn)溫差均隨泄漏流量增加而增大，泄漏流量為400 L/h時(shí)，最大溫差達(dá)到105 ℃。

圖6 截面直徑線上溫度分布Fig.6 Vertical temperature distribution of Cross Section

不同泄漏流量時(shí)，閥前50 mm管外壁面沿圓周方向溫度分布如圖7所示，溫度分布左右基本對(duì)稱，QLeak越大，管外壁溫度越低。

圖7 不同泄漏量下管外壁溫度分布Fig.7 Pipe wall temperatures at various leakage rate

3 泄漏預(yù)測(cè)模型的建立

3.1 熱分層特征溫度參數(shù)

將溫度監(jiān)測(cè)點(diǎn)設(shè)在止逆閥前50 mm處的位置AT6B—7，在測(cè)量位置的一側(cè)沿圓周方向30°等間距安裝7個(gè)熱電偶，如圖8所示。根據(jù)7個(gè)測(cè)點(diǎn)溫度可整理如下熱分層溫度特征參數(shù)。

圖8 溫度測(cè)點(diǎn)的布置Fig.8 Arrangement of temperature measuring points

令

(15)

(16)

(17)

(18)

(19)

(20)

則熱分層致溫差為

ΔT=T0°-T180°

(21)

熱分層致當(dāng)量面平均溫度為

(22)

(23)

式中：Si——兩相鄰測(cè)溫點(diǎn)間監(jiān)測(cè)截面的面積。

3.2 熱分層特征溫度參數(shù)單變量分析

熱分層特征溫度參數(shù)隨主管道流體溫度的變化關(guān)系如圖9所示，當(dāng)量面平均溫度TA隨著主管道流體溫度的增加而增加，且呈線性增加趨勢(shì)。這是因?yàn)橹鞴艿纼?nèi)流體的溫度越高，則管壁溫度會(huì)直接增大，熱流體的傳熱速率加快，熱流密度的增加從而管壁的溫度相應(yīng)的增加。主管道溫度增大時(shí)，管外壁的溫差ΔT逐漸增大，呈非線性，運(yùn)用最小二乘法進(jìn)行擬合時(shí)，管壁溫差ΔT與主管道溫度成三次關(guān)系。

圖9 溫度特征值隨主管道溫度變化關(guān)系Fig.9 Variation relationship of characteristic temperature with the main pipe temperature

熱分層特征溫度參數(shù)隨主管道流量的變化關(guān)系如圖10所示，當(dāng)主管道內(nèi)流量逐漸增大時(shí)，當(dāng)量面平均溫度TA呈非線性降低，運(yùn)用最小二乘法進(jìn)行擬合時(shí)，當(dāng)量平均溫度TA與主管道流量成冪函數(shù)關(guān)系。主管道流量增大時(shí)，溫差ΔT非線性降低，運(yùn)用最小二乘法進(jìn)行擬合時(shí)，溫差ΔT與主管道流量成一階指數(shù)衰減函數(shù)關(guān)系。

圖10 溫度特征值隨主管道流量變化關(guān)系Fig.10 Variation relationship of characteristic temperature with the main pipe flow rate

熱分層特征溫度參數(shù)隨泄漏流體溫度的變化關(guān)系如圖11所示，當(dāng)泄漏工質(zhì)溫度逐漸增大時(shí)，當(dāng)量面平均溫度TA線性增加。管道泄漏工質(zhì)溫度增大時(shí)，管壁溫差ΔT線性降低。

圖11 溫度特征值隨泄漏工質(zhì)溫度變化關(guān)系Fig.11 Variation relationship of characteristic temperature with the leakage flow temperature

熱分層特征溫度參數(shù)隨泄漏流量的變化關(guān)系如圖12所示，當(dāng)泄漏工質(zhì)流量逐漸增大時(shí)，監(jiān)測(cè)位置的當(dāng)量面平均溫度TA有明顯的減小趨勢(shì)，且成二次函數(shù)關(guān)系。當(dāng)管道泄漏流體的流量增大時(shí)，管壁溫差ΔT非線性增大。運(yùn)用最小二乘法進(jìn)行擬合時(shí)，溫差ΔT與泄漏流量成二次關(guān)系。

圖12 溫度特征值隨泄漏量變化關(guān)系Fig.12 Variation relationship of characteristic temperature with the leakage flow rate

3.3 熱分層特征溫度參數(shù)多變量擬合

根據(jù)不同工況條件進(jìn)行的流固耦合計(jì)算結(jié)果，利用MATLAB統(tǒng)計(jì)工具箱中的regress命令，對(duì)熱分層特征溫度參數(shù)關(guān)于主管道內(nèi)流體溫度Thot、主管道內(nèi)流體流速Q(mào)hot、安注流體溫度Tcold、泄漏流量QLeak進(jìn)行多變量回歸計(jì)算可得

(24)

(25)

式中：α0……α5，β0……β7——系數(shù)。

TA和ΔT統(tǒng)計(jì)計(jì)算的R2分別為0.9967和0.9957，表明所定義的熱分層特征溫度參數(shù)與回歸變量之間有良好的相關(guān)性。

從而可建立內(nèi)漏流量QLeak的閥門內(nèi)漏量化預(yù)測(cè)模型

QLeak,TA=g(Thot,Qhot,Tcold,TA)

(26)

QLeak,ΔT=f(Thot,Qhot,Tcold,ΔT)

(27)

現(xiàn)場(chǎng)泄漏流量計(jì)算時(shí)，通過現(xiàn)場(chǎng)在閥前管外壁布置的熱電偶，得到測(cè)溫點(diǎn)處的溫度，根據(jù)公式(15)～公式(23)計(jì)算出所述測(cè)溫點(diǎn)的當(dāng)量面平均溫度TA測(cè)，將TA測(cè)代入到公式(26)中，得到QLeak,TA，根據(jù)公式(21)計(jì)算出所述測(cè)溫點(diǎn)的溫度差ΔT測(cè)，將ΔT測(cè)代入到公式(27)中，得到QLeak,ΔT，若滿足公式(28)，那么取它們兩個(gè)的算數(shù)平均值如公式(29)所示，若不滿足，則取兩個(gè)之間的較大值如公式(30)所示。

(28)

(29)

QLeak=max{QLeak,TA,QLeak,ΔT}

(30)

3.4 誤差分析

擬合關(guān)聯(lián)式結(jié)果與數(shù)值計(jì)算結(jié)果吻合較好，最大誤差為6.471 7%，可滿足工程閥門內(nèi)泄漏流量需要。

4 結(jié)論

以熱安注管線止逆閥泄漏為例，探討了基于閥前管壁溫度場(chǎng)測(cè)量的閥內(nèi)泄漏流量監(jiān)測(cè)方法。通過流固耦合計(jì)算對(duì)熱安注管線由于閥門內(nèi)泄漏造成的熱分層現(xiàn)象進(jìn)行了詳細(xì)分析，定義了閥前監(jiān)測(cè)截面管外熱分層特征溫度參數(shù)，通過多變量回歸計(jì)算獲取了熱分層特征溫度參數(shù)與主管道內(nèi)流體溫度Thot、主管道內(nèi)流體流速Q(mào)hot、泄漏流體溫度Tcold和泄漏流量QLeak之間的關(guān)系式。只要根據(jù)閥前監(jiān)測(cè)截面管外壁7點(diǎn)監(jiān)測(cè)溫度獲取熱分層特征溫度參數(shù)，即可根據(jù)電廠已知參數(shù)(主管道內(nèi)流體溫度Thot、主管道內(nèi)流體流速Q(mào)hot、泄漏流體溫度Tcold)較為準(zhǔn)確地預(yù)測(cè)止逆閥的泄漏流量。