反應(yīng)堆堆內(nèi)構(gòu)件流體激振力的數(shù)值模型研究

馮志鵬，黃 旋，劉 帥，沈平川，蔡逢春

中國(guó)核動(dòng)力研究設(shè)計(jì)院核反應(yīng)堆系統(tǒng)設(shè)計(jì)技術(shù)重點(diǎn)實(shí)驗(yàn)室，四川 成都 610213

引言

反應(yīng)堆堆內(nèi)構(gòu)件（簡(jiǎn)稱堆內(nèi)構(gòu)件）是核電站反應(yīng)堆的主要設(shè)備，位于反應(yīng)堆壓力容器內(nèi)，在高溫、高壓、強(qiáng)輻照下運(yùn)行，并且長(zhǎng)期承受冷卻劑的沖刷[1]。冷卻劑的循環(huán)流動(dòng)會(huì)引起堆內(nèi)構(gòu)件的長(zhǎng)期振動(dòng)[2]，且總是伴隨著反應(yīng)堆的運(yùn)行而存在，極易導(dǎo)致由流致振動(dòng)引起的疲勞磨損[3]、噪聲[4]等問(wèn)題，為確保結(jié)構(gòu)在使用壽期內(nèi)的完整性，只能通過(guò)優(yōu)化結(jié)構(gòu)設(shè)計(jì)和調(diào)整流速，將流致振動(dòng)引起的結(jié)構(gòu)振動(dòng)響應(yīng)控制在可接受的范圍之內(nèi)。因此，堆內(nèi)構(gòu)件流致振動(dòng)理論分析與計(jì)算是核電設(shè)備安全分析不可缺少的重要內(nèi)容[5]。

針對(duì)核工程中的流致振動(dòng)問(wèn)題，國(guó)內(nèi)外學(xué)術(shù)界和工程界開(kāi)展了大量的理論、實(shí)驗(yàn)和數(shù)值分析研究。但由于流致振動(dòng)問(wèn)題十分重要但又十分復(fù)雜，且堆內(nèi)構(gòu)件結(jié)構(gòu)復(fù)雜，包含了下部堆內(nèi)構(gòu)件、上部堆內(nèi)構(gòu)件、壓緊彈簧和U 形嵌入件等，同時(shí)流致振動(dòng)經(jīng)常出現(xiàn)在高雷諾數(shù)、非定常、分離流和比較柔性的結(jié)構(gòu)中[6-8]，很難從理論模型上來(lái)分析，設(shè)計(jì)主要依賴于經(jīng)驗(yàn)系數(shù)和實(shí)驗(yàn)。Au-Yang 等[9]在反應(yīng)堆1：9 水力回路上進(jìn)行了流致振動(dòng)實(shí)驗(yàn)，基于實(shí)驗(yàn)數(shù)據(jù)提出了適用于圓筒結(jié)構(gòu)的隨機(jī)湍流激勵(lì)的功率譜密度（PSD）。喻丹萍等[10]針對(duì)秦山核電二期工程的堆內(nèi)構(gòu)件，做了大量的流致振動(dòng)實(shí)驗(yàn)，其研究結(jié)果表明，吊籃流致振動(dòng)主要是由湍流引起的各態(tài)歷經(jīng)的平穩(wěn)隨機(jī)振動(dòng)，吊籃可能出現(xiàn)梁式和殼式振動(dòng)。梁志堅(jiān)等[11]應(yīng)用極值載荷分布法分析了秦山600 MW 核電站堆內(nèi)構(gòu)件流致振動(dòng)實(shí)驗(yàn)的實(shí)測(cè)載荷。張可豐等[12]以堆內(nèi)構(gòu)件防斷組件為研究對(duì)象，研究了隨機(jī)湍動(dòng)力、擬靜態(tài)湍動(dòng)力、漩渦脫落、吊籃運(yùn)動(dòng)等流體激勵(lì)機(jī)理對(duì)防斷組件流致振動(dòng)的影響。楊翊仁等[13-14]采用隨機(jī)振動(dòng)的模態(tài)分析方法研究了吊籃結(jié)構(gòu)縮比模型的流致振動(dòng)響應(yīng)。姚偉達(dá)等[15]采用實(shí)驗(yàn)與計(jì)算相結(jié)合的分析方法，對(duì)秦山核電廠反應(yīng)堆吊籃的流致振動(dòng)響應(yīng)進(jìn)行了計(jì)算。席志德等[16-17]用大渦模擬（LES）方法計(jì)算了秦山核電站二期1：5 模型實(shí)驗(yàn)的壓力容器和吊籃之間環(huán)腔內(nèi)的流場(chǎng)，獲取了流體作用在吊籃上的流體壓力，Simoneau 等[18]詳細(xì)介紹了大渦模擬在核能領(lǐng)域用于計(jì)算湍流方面的應(yīng)用，從平均流場(chǎng)、波動(dòng)的幅值和頻率等方面來(lái)看，運(yùn)用大渦模擬與實(shí)驗(yàn)結(jié)果吻合很好，但需要非常精細(xì)的網(wǎng)格，因此，極大地限制了其在實(shí)際工程中的應(yīng)用。

綜上所述，流致振動(dòng)是關(guān)系到核反應(yīng)堆結(jié)構(gòu)設(shè)計(jì)安全性和可靠性的重要問(wèn)題，目前要在理論上求解作用于堆內(nèi)構(gòu)件表面的脈動(dòng)壓力載荷十分困難，因此，需要采用實(shí)驗(yàn)與計(jì)算相結(jié)合的方法來(lái)研究堆內(nèi)構(gòu)件的流致振動(dòng)問(wèn)題。本文首先建立堆內(nèi)構(gòu)件流場(chǎng)區(qū)域的幾何模型和網(wǎng)格模型，然后建立適用于工程的三維流場(chǎng)分析模型，最后研究各關(guān)鍵部件的流體激振力特性，為堆內(nèi)構(gòu)件流致振動(dòng)分析評(píng)價(jià)提供關(guān)鍵參數(shù)和輸入載荷。

1 數(shù)值模型

1.1 控制方程

本文中的流體為水，且流速較低，當(dāng)作不可壓縮流體處理?？刂品匠淌荖-S 方程，在直角坐標(biāo)系下表示為[19]

連續(xù)性方程

下標(biāo)：i，j=1,2,3。

1.2 湍流模型

本文主要采用大渦模擬方法[20]對(duì)湍流流場(chǎng)進(jìn)行求解。對(duì)N-S 方程在物理空間經(jīng)濾波即可得到大渦的控制方程

需要使用不同亞格子尺度模型進(jìn)行建模，本文采用Boussinesq 假設(shè)[21]，由渦黏形式計(jì)算亞格子湍流應(yīng)力

下標(biāo)：k=1,2,3。

此外，在研究流場(chǎng)模型時(shí)，還用到了k-ε 及SST兩種湍流模型，可參見(jiàn)文獻(xiàn)[20]。

1.3 幾何模型

堆內(nèi)構(gòu)件結(jié)構(gòu)非常復(fù)雜，在建立流場(chǎng)分析模型時(shí)，須對(duì)實(shí)際結(jié)構(gòu)進(jìn)行適當(dāng)?shù)暮?jiǎn)化，如忽略局部結(jié)構(gòu)特征、刪除連接螺栓、刪除螺栓孔、消除過(guò)小縫隙、平滑突出結(jié)構(gòu)表面的小構(gòu)件等，在合理降低計(jì)算規(guī)模的同時(shí)，又充分保留其流場(chǎng)特征。在反應(yīng)堆壓力容器內(nèi)壁包圍的總體空間中，通過(guò)布爾運(yùn)算扣除堆內(nèi)構(gòu)件固體區(qū)域，即為堆內(nèi)構(gòu)件流場(chǎng)分析所需的流體區(qū)域，如圖1 所示，其中，堆芯段內(nèi)布置燃料組件，流場(chǎng)十分復(fù)雜，由于本文不關(guān)注堆芯段的流場(chǎng)細(xì)節(jié)，因此，通過(guò)多孔介質(zhì)的方式模擬，以減小計(jì)算規(guī)模。

圖1 流體域幾何模型Fig.1 Geometry model of fluid domain

1.4 網(wǎng)格模型

流場(chǎng)網(wǎng)格采用ICEM CFD 軟件來(lái)完成。為了使劃分的網(wǎng)格能夠盡可能捕捉幾何上的精細(xì)結(jié)構(gòu)，對(duì)關(guān)鍵的面和線創(chuàng)建部件，并對(duì)不同部件分別設(shè)置網(wǎng)格尺寸。本文將上部堆內(nèi)構(gòu)件分為17個(gè)部件、下部堆內(nèi)構(gòu)件分為46 個(gè)部件、堆芯曲面分為3 個(gè)部件。對(duì)不同的部件分別設(shè)置最大尺寸、最小尺寸和網(wǎng)格生長(zhǎng)率等參數(shù)。參數(shù)的設(shè)置應(yīng)保證該部件上的面網(wǎng)格能夠較為精確地描述部件的幾何邊界，特別是對(duì)于重點(diǎn)關(guān)注區(qū)域或?qū)α鲌?chǎng)影響較大的區(qū)域應(yīng)該設(shè)置較小的網(wǎng)格尺寸，以捕捉流場(chǎng)特征。

上部堆內(nèi)構(gòu)件、下部堆內(nèi)構(gòu)件（不包含堆芯）的流體域結(jié)構(gòu)復(fù)雜，采用非結(jié)構(gòu)化網(wǎng)格劃分，通過(guò)一體化生成方式生成網(wǎng)格，網(wǎng)格細(xì)節(jié)如圖2a、圖2b 所示。堆芯部分由于采用多孔介質(zhì)來(lái)模擬，因此，采用六面體網(wǎng)格劃分，如圖2c 所示。

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

在綜合考慮計(jì)算資源和計(jì)算效率的前提下，進(jìn)行了一系列網(wǎng)格敏感性測(cè)試，以期在計(jì)算精度和計(jì)算效率之間取得平衡，最終確定了本文流場(chǎng)網(wǎng)格的劃分方式和網(wǎng)格尺寸，網(wǎng)格總數(shù)為36 050 000。

1.5 計(jì)算參數(shù)

根據(jù)運(yùn)行壓力和進(jìn)口溫度，流體介質(zhì)的密度為742 kg/m3，動(dòng)力黏度為9.242×10-5Pa·s。

經(jīng)過(guò)大量計(jì)算，本文選取預(yù)測(cè)結(jié)果比較好的3種湍流模型進(jìn)行對(duì)比研究：大渦模擬模型、SST 模型、k-ε 模型。

入口采用速度入口邊界條件，其值18.645 m/s。當(dāng)湍流模型為SST 模型或k-ε 模型時(shí)，入口的湍流強(qiáng)度設(shè)置為1.5%，水力直徑為697.8 mm。3 個(gè)出口采用壓力出口，設(shè)置壓強(qiáng)為0。重力方向?yàn)閦軸負(fù)方向，大小為9.81 m/s2。

對(duì)堆芯模擬體的壓降采用多孔介質(zhì)模型模擬，由堆芯壓降和堆芯模擬題長(zhǎng)度可得到多孔介質(zhì)模型的慣性阻力系數(shù)為16.62 m-1。

壓力-速度耦合采用SIMPLEC 算法，動(dòng)量離散為邊界中心差分，瞬態(tài)求解格式是二階隱式，采用雙精度求解器。時(shí)間步長(zhǎng)取為0.001 s。

2 三維流場(chǎng)模型分析

2.1 流場(chǎng)特征

以LES 的計(jì)算結(jié)果為例來(lái)分析堆內(nèi)構(gòu)件的流場(chǎng)特征。整體反應(yīng)堆流場(chǎng)的流線分布如圖3 所示，可以看出，在上部堆內(nèi)構(gòu)件和下部堆內(nèi)構(gòu)件的二次支承組件區(qū)域流線最為紊亂。

圖3 反應(yīng)堆結(jié)構(gòu)流線圖Fig.3 Streamline of reactor

過(guò)出、入口管的中心沿反應(yīng)堆軸向作切面，并在該面上列出速度矢量分布圖，分別見(jiàn)圖4 和圖5，可以看出，整個(gè)二次支承區(qū)域的流場(chǎng)最為紊亂，并有較大的漩渦產(chǎn)生，而對(duì)于上部堆內(nèi)構(gòu)件區(qū)域，出、入口處的流場(chǎng)較為紊亂，下降段環(huán)腔和堆芯部分的流場(chǎng)較為規(guī)整，主要是由于堆芯部分采用了簡(jiǎn)化的多孔介質(zhì)來(lái)模擬，僅關(guān)注其壓降，而忽略了堆芯的結(jié)構(gòu)細(xì)節(jié)。從圖4 和圖5 還可以看出，吊籃與壓力容器環(huán)腔的流體在兩列入口交匯的區(qū)域內(nèi)流體流速較大，作用在吊籃外表面的流體主要沿吊籃的軸向分布，吊籃外表面主要受到軸向流的作用。入口正對(duì)的吊籃壁面區(qū)域?yàn)轳v點(diǎn)，壓力最大。

圖4 出口軸向切面上速度矢量圖Fig.4 Velocity vectors in axial section of outlet

圖5 入口軸向切面上速度矢量圖Fig.5 Velocity vectors in axial section of inlet

2.2 壓力脈動(dòng)

為了對(duì)比研究不同湍流模型對(duì)壓力脈動(dòng)的預(yù)測(cè)能力，取位于吊籃、流量分配板、二次支承柱連接板上的3 個(gè)測(cè)點(diǎn)的脈動(dòng)壓力進(jìn)行對(duì)比分析，其中，監(jiān)測(cè)點(diǎn)如圖6 所示。測(cè)點(diǎn)1：吊籃（125°）法蘭下，測(cè)點(diǎn)2：流量分配板（135°），測(cè)點(diǎn)3：二次支承柱連接板（180°）。

圖6 脈動(dòng)壓力的監(jiān)測(cè)點(diǎn)示意圖Fig.6 Schematic diagram of monitoring points for pressure fluctuation

圖7 給出了利用LES 模型計(jì)算出的測(cè)點(diǎn)1、測(cè)點(diǎn)2、測(cè)點(diǎn)3 的壓力脈動(dòng)時(shí)程及其功率譜密度曲線（PSD）。從壓力脈動(dòng)的頻譜可以看出，各測(cè)點(diǎn)處的壓力在頻域中沒(méi)有占優(yōu)的頻率，且頻率的分布范圍較寬，PSD 幅值隨頻率增加遞減較快，這種特點(diǎn)符合寬帶隨機(jī)性，作用在堆內(nèi)構(gòu)件上的流體激勵(lì)力是一種寬帶隨機(jī)過(guò)程。

圖7 LES 模型預(yù)測(cè)的各測(cè)點(diǎn)的壓力脈動(dòng)時(shí)程及PSDFig.7 Time history and PSD of pressure fluctuationfor each monitoring point predicted by LES model

圖8 給出了采用SST 模型和k-ε 模型計(jì)算出的壓力脈動(dòng)時(shí)程及其功率譜密度曲線?？梢钥闯觯琒ST 模型和k-ε 模型計(jì)算得到的壓力脈動(dòng)時(shí)程比LES 模型計(jì)算得到的壓力脈動(dòng)時(shí)程平緩，原因在于SST 和k-ε 模型為雷諾時(shí)均模型，在計(jì)算過(guò)程中抹平了湍流的若干細(xì)節(jié)，而LES 模型可以獲得比較精細(xì)的流場(chǎng)壓力。

圖8 SST 模型和k-ε 模型預(yù)測(cè)的測(cè)點(diǎn)1 的壓力脈動(dòng)時(shí)程及PSDFig.8 Time history and PSD of pressure fluctuationfor monitoring point 1 predicted by SST model and k-ε model

2.3 3 種湍流模型結(jié)果對(duì)比

為了更進(jìn)一步對(duì)比這3 種湍流模型對(duì)反應(yīng)堆堆內(nèi)構(gòu)件流體激振力的預(yù)測(cè)情況，將其統(tǒng)計(jì)值（均方根值）列于圖9 中。

圖9 各測(cè)點(diǎn)的壓力脈動(dòng)均方根值Fig.9 RMS(Root mean square)value of pressure fluctuation for each monitoring point

從圖9 可以看出：（1）在預(yù)測(cè)吊籃測(cè)點(diǎn)的壓力脈動(dòng)時(shí)，SST 模型、k-ε 模型的壓力脈動(dòng)均方根值均遠(yuǎn)小于實(shí)驗(yàn)值，嚴(yán)重低估了吊籃上的壓力激勵(lì)水平，而LES 模型預(yù)測(cè)的結(jié)果與實(shí)驗(yàn)值吻合較好。（2）在預(yù)測(cè)流量分配板、二次支承柱連接板的壓力脈動(dòng)時(shí)，SST 模型、k-ε 模型、LES 模型的均方根值均比實(shí)驗(yàn)結(jié)果大，在工程設(shè)計(jì)時(shí)具備保守性。（3）綜合對(duì)比分析各測(cè)點(diǎn)的壓力脈動(dòng)時(shí)程、PSD、均方根值，可以得出，LES 模型可以較好地預(yù)測(cè)反應(yīng)堆堆內(nèi)構(gòu)件的流體激振力，同時(shí)又具有一定的保守性，可以用于堆內(nèi)構(gòu)件流體激振力的計(jì)算。

2.4 流體激振力特性

吊籃是堆內(nèi)構(gòu)件中的關(guān)鍵結(jié)構(gòu)，包括吊籃法蘭、圓柱形筒體、3 個(gè)出口管嘴和堆芯支承板等部件。以吊籃組件為代表結(jié)構(gòu)，分析作用于其上的流體激振力特性。根據(jù)反應(yīng)堆冷卻劑的流動(dòng)特點(diǎn)，流體激振力主要作用于吊籃筒體表面，因此，選擇吊籃組件外壁布置監(jiān)測(cè)點(diǎn)。吊籃組件外壁面共布置460 個(gè)壓力監(jiān)測(cè)點(diǎn)，在z方向分為23 層，每層20 個(gè)點(diǎn)，以圓心角18°均分，如圖10 所示。

圖10 吊籃組件表面上的監(jiān)測(cè)點(diǎn)示意圖Fig.10 Schematic diagram of monitoring points on the core barrel surface

圖11 給出了吊籃組件表面所有測(cè)點(diǎn)脈動(dòng)壓力均方根值的分布情況，從圖中可以看到，在進(jìn)出口位置處（測(cè)點(diǎn)標(biāo)高約-4～-1 m），壓力脈動(dòng)的均方根值變化劇烈且值較大。在遠(yuǎn)離進(jìn)出口以下的位置處（測(cè)點(diǎn)標(biāo)高約-8～-4 m），各測(cè)點(diǎn)的壓力脈動(dòng)均方根值較小且相差不大。

圖11 吊籃組件表面測(cè)點(diǎn)的壓力脈動(dòng)均方根值分布情況Fig.11 The distribution of pressure fluctuation RMS value for monitoring points on the core barrel surface

為了進(jìn)一步詳細(xì)分析各個(gè)標(biāo)高、各個(gè)周向角度處的壓力脈動(dòng)情況，選取位于進(jìn)口處、出口處、進(jìn)口與出口之間不同標(biāo)高的典型測(cè)點(diǎn)以及第10 層、第18 層一周的測(cè)點(diǎn)進(jìn)行分析。

圖12a 給出了第6 個(gè)周向監(jiān)測(cè)點(diǎn)（108°處，位于進(jìn)口處）的壓力脈動(dòng)均方根值隨監(jiān)測(cè)點(diǎn)層數(shù)的分布情況，可以看到，壓力脈動(dòng)均方根值從吊籃組件底部到進(jìn)口處逐漸增大，在靠近進(jìn)口處下方達(dá)到峰值，然后隨著監(jiān)測(cè)點(diǎn)層數(shù)的增加先減小再增大。圖12b 給出了第2 個(gè)周向測(cè)點(diǎn)（36°處，位于出口處）的壓力脈動(dòng)均方根值隨監(jiān)測(cè)點(diǎn)層數(shù)的分布情況，可以看到，壓力脈動(dòng)均方根值從吊籃組件底部到出口標(biāo)高處逐漸增大，在靠近出口標(biāo)高處下方達(dá)到最大值，然后隨著標(biāo)高增加先減小再增大，需要注意的是，由于監(jiān)測(cè)點(diǎn)在層數(shù)18 及21 落入吊籃出口位置，所以該處沒(méi)有監(jiān)測(cè)點(diǎn)。圖12c 給出了第5 個(gè)周向測(cè)點(diǎn)（90°處，位于進(jìn)口與出口之間）的壓力脈動(dòng)均方根值隨監(jiān)測(cè)點(diǎn)層數(shù)的分布情況，流場(chǎng)在該部分區(qū)域發(fā)展的較為充分，受到的擾動(dòng)也小，所以最具有代表性，可以看出，壓力脈動(dòng)從吊籃組件吊籃底部到進(jìn)出口標(biāo)高處逐漸增大，在靠近進(jìn)出口標(biāo)高下方處達(dá)到最大值，然后隨著標(biāo)高增加先減小再增大。

圖12 典型周向測(cè)點(diǎn)的壓力脈動(dòng)均方根值隨監(jiān)測(cè)點(diǎn)層數(shù)的分布情況Fig.12 The distribution of pressure fluctuation RMS value at typical circumferential monitoring points with the layer of monitoring points

圖13a 給出第10 層測(cè)點(diǎn)（吊籃組件下部）的壓力脈動(dòng)均方根值沿周向角度的分布情況，可以看出，第10 層各測(cè)點(diǎn)的壓力脈動(dòng)均方根值是有明顯差別的，是由于各點(diǎn)的周向角位置不同，對(duì)應(yīng)出、入口的位置不同導(dǎo)致的，大致為3 個(gè)波峰和3 個(gè)波谷，3 個(gè)波峰對(duì)應(yīng)出口下方附近，3 個(gè)波谷對(duì)應(yīng)入口下方附近。圖13b給出了第18 層測(cè)點(diǎn)（進(jìn)出口標(biāo)高處）的壓力脈動(dòng)均方根值沿周向角度的分布情況，其分布規(guī)律與第10 層測(cè)點(diǎn)的類似，但由于其標(biāo)高位于進(jìn)出口標(biāo)高處，所以各測(cè)點(diǎn)間的壓力脈動(dòng)值差異較大。

圖13 典型測(cè)點(diǎn)的壓力脈動(dòng)均方根值沿周向角度的分布情況Fig.13 The distribution of pressure fluctuation RMS value with the circumferential angle of monitoring points

3 結(jié)論

（1）根據(jù)堆內(nèi)構(gòu)件的流致振動(dòng)特性建立了堆內(nèi)構(gòu)件流場(chǎng)區(qū)域的幾何模型，在平衡計(jì)算資源和計(jì)算效率的前提下，確定了流場(chǎng)網(wǎng)格的劃分方式和網(wǎng)格尺寸，建立了合理的流場(chǎng)網(wǎng)格模型。

（2）大渦模擬模型可以較好地預(yù)測(cè)反應(yīng)堆堆內(nèi)構(gòu)件的流體激振力，進(jìn)而建立了堆內(nèi)構(gòu)件流致振動(dòng)分析的流場(chǎng)模型，獲得了三維流場(chǎng)特性。

（3）作用在堆內(nèi)構(gòu)件的流體激振力符合寬帶隨機(jī)性，并且計(jì)算得到的壓力脈動(dòng)均方根值與相應(yīng)測(cè)點(diǎn)位置的實(shí)驗(yàn)值在同一量級(jí)，同時(shí)又具有一定的保守性，可以用于反應(yīng)堆堆內(nèi)構(gòu)件的流致振動(dòng)響應(yīng)分析。