中間熱交換器建模仿真研究

張厚明，段天英，王 剛，劉 勇，熊文彬，劉兆陽

（1．國家核安全局華北核與輻射安全監督站，北京100191；2．中國原子能科學研究院，北京102413；3．國家核安全局核與輻射安全中心，北京100082）

中間熱交換器是液態金屬冷卻快中子增殖堆（以下簡稱“快堆”）熱傳輸系統中隔離放射性的一回路載熱劑與非放射性二回路載熱劑的屏障，同時也是將堆芯釋熱傳遞給蒸汽發生器的重要設備。

對于快堆核電廠熱工水力、儀控系統的設計來講，中間熱交換器的建模及仿真研究可以為之提供設計依據。同時，可以將建模及仿真研究的成果用于設計工程仿真機、全范圍仿真機來為操作員提供培訓。另外，仿真模型在高速計算機上運行，通過實時采集快堆核電廠數據提前計算出各個參數的運行趨勢，從而為操縱員提供技術支持。因此，對于上述三個方面，中間熱交換器建模及仿真研究都具有重要意義。本文給出一種計算迅速、精度較高、可以用于控制系統分析、實時計算的中間熱交換器仿真模型。

1 中間熱交換器建模簡介

早期實驗、原型快堆中間熱交換器類型、結構形式較多［1，2］，例如：BN－350采用豎直U形傳熱管束置于水平矩形筒體內；BR－5及Grand Quevilly French Test Station采用U形傳熱管束置于U形筒體內；Fermi及FFTF采用直傳熱管置于豎直圓柱形筒體內，采用可移動的封頭來補償熱膨脹效應；Hallam Nuclear Power Facility采用直傳熱管，殼側有膨脹節來補償熱膨脹效應。

示范、商用大型快堆，例如Super Phénix、CRBRP（未建造）、BN－800，均采用直換熱管、逆流管殼式中間熱交換器。二回路鈉均從頂部進入，通過中心下降管流到底部下腔室，然后流向反轉，向上流過換熱管束。一回路鈉均從上部進入，流經傳熱管間的殼側，從底部流出，如圖1所示［3］。

因此，本文選取某個直換熱管、逆流殼式中間熱交換器作為研究對象。

描述中間熱交換器三維熱工水力的方程組非常復雜。由于動量、能量方程的相互耦合，并且，這些方程的非線性使得方程組很難求解。然而，即使在自然循環工況下，三維計算結果與一維計算結果也相差不大［4］，因此以往的研究中，大多假設：

圖1 CRBRP中間熱交換器概圖Fig．1 Sketch of CRBRP IHXs'

（a）用一根傳熱管代替傳熱管束；

（b）流體物性參數沿徑向保持不變；

（c）忽略載熱劑、管壁的軸向導熱；

（d）忽略中間熱交換器殼體的散熱。

在上述這些假設下，Gunby［5］采用“節距固定”方法將中間熱交換器的傳熱管劃分為10個等長節段，在每個節段上建立了能量守恒微分方程組，采用不同的差分來求解，研究了流量、入口溫度變化等各種工況下中間熱交換器的溫度的分布及變化。研究結果表明：對一、二回路不平衡流動（10∶1），流量大的一側其計算結果差別較大；并且，即使出口溫度計算結果較為可信，各種差分方法計算出的內部溫度分布也均不正確；在要求出口溫度和內部溫度分布計算結果同時可靠時，大擾動情況下，采用這種“節距固定”的方法，劃分節段數目較少時，各種差分方法均難以獲得正確的計算結果。

Brukx［6］針對Gunby方法中出現的內部溫度分布不正確問題，給出了一種將能量守恒微分方程組在空間上離散的差分計算方法，計算結果較好，但這種方法中無法實時計算溫度變化，在控制系統設計，例如使用Matlab／Simulink軟件對控制系統進行設計中，諸多不便。

段天英［7，8］根據能量守恒方程，在小擾動情況下推導出中間熱交換器出口溫度隨流量、入口溫度等參數變化的傳遞函數，并將傳遞函數進行修正，計算結果較為精確，這種方法可以用于擾動不大的控制系統設計中；但是，這種方法無法計算各種瞬態工況下中間熱交換器內部的溫度變化，且不可以用于擾動較大的控制系統設計。

Tzanos［9］給出兩種中間熱交換器模型。第一種為“節距固定”方法：將傳熱管均分為145個節段，在每個節段上采用能量守恒方程進行描述，這樣獲得描述中間熱交換器模型的常微分方程組。由于早期計算機運算速度較慢，Tzanos給出第二種“節距可變”方法：將中間熱交換器按照等能量劃分為10～15個節段，采用Leibnitz方程將每個節段上的能量守恒方程偏微分方程組轉化為溫度導數與節段長度導數有關的常微分方程組。由于出口溫度變化與長度變化有關，“節距可變”方法使得節段平均溫度為節段出入口溫度的算術平均值這一假設在節段數目較少的情況下比“節距固定”的方法更為合理，于是，采用這種階段長度可變的方法可以大大減少常微分方程組維數和計算機機時。但是，這種“節段長度可變”的方法與“節距固定”方法的計算結果仍有一定差別。這種方法在EBR－Ⅱ的TEST－8A試驗［9，10］中得到驗證，最大計算誤差在6K左右。

綜上所述，在以往的對中間熱交換器建模的各種方法中，采用“節距固定”的方法來建立中間熱交換器模型，若取得較好的計算結果，瞬態過程中中間熱交換器內溫度可信，不會出現數值不穩定，那么節段數目的選取應至少為100個左右［9］。本文給出150節段、適于Matlab／Simulink軟件中計算的“節距固定”模型。本模型建立過程中首先提出兩種穩態計算方法，比較兩種方法計算的穩態結果，然后給出一種瞬態模型，并將穩態計算的結果作為瞬態計算的參數，對二回路鈉流量指數衰減的工況進行了計算，以下分別介紹穩態、瞬態計算模型。

2 中間熱交換器穩態計算方法

文獻［11］在假定穩態情況下整個中間熱交換器內鈉、換熱管材料的比熱容及換熱系數均不變，推導出沿高度方向溫度呈指數分布，這種計算方法較為粗糙，計算結果誤差較大。但本文在程序設計中將這種方法計算出的結果作為粗略估計值，從而可以減少穩態計算的循環次數。

本文給出的兩種穩態工況下的計算方法，將中間熱交換器傳熱管分為n個節段，如圖2所示。

圖2 中間熱交換器第i個節段Fig．2 Segment i in IHXs

穩態工況下，中間熱交換器在節段i上的能量守恒方程為：

式（1）中，

其中，W為質量流量，C為比熱容，U為傳熱系數，l為節段長度，T為溫度，λ為熱導率，a為一回路流道截面積，p為節距；角標p為一回路，m為管壁，s為二回路，in為管內，out為管外。

根據式（1）采用包括：等節段長度、等換熱量兩種方法來劃分節段數目并編制程序，具體算法、循環語句如圖3所示。根據這兩種方法，分別調整α、β、δ、ε的值，可以獲得精度較高的計算結果，經過計算兩種算法出的結果相差很小，如圖4所示。與實際值相比，所計算出的值絕對誤差誤差在3K以內。

圖3 穩態計算方法Fig．3 Steady state algorithms

圖4 穩態工況下的計算結果Fig．4 Calculated result with steady state algorithms

由計算結果可知：管內熱阻不超過總熱阻的15%～20%，殼側的傳熱系數明顯小于管側，其熱阻約為總熱阻的30%～40%，管壁熱阻約為總熱阻的40%～55%。因此在設計中，進一步增加管內流速不會大幅降低總熱阻；增加殼側的流速可以有效降低總熱阻，減小管壁厚度可以大大降低總熱阻。這與文獻［1，2］得出的結論是一致的。

鈉、管壁材料的物性參數來自文獻［1，11］；一、二回路的傳熱系數計算公式來自文獻［9］。

3 中間熱交換器瞬態模型

用來描述中間熱交換器瞬態變化較好的模型應包括一、二回路入口流量、溫度變化帶來中間熱交換器內部、出口溫度變化的數學模型，同時，應選取一種計算迅速、易于使用、誤差小的計算方法。根據上述討論，描述中間熱交換器的能量守恒偏微分方程組為方程（3）的形式。

其中，密度ρ、比熱C、傳熱系數U，均為溫度的函數。

在每個節段上，將方程（3）離散為差分方程（4）：

對方程（4）這種差分方法進行了穩定性分析，將方程（4）展開，并轉化為方程（5）的形式。

設3n－3維向量T＝（Tp，Tm，Ts），也即：

那么差分方程（5）可以寫成方程（6）的形式，簡寫為MT′＝NT，其中A－G均為系數矩陣。

根據差分方程組（4）依次建立每個節段上的模型，如圖5所示。并根據上述文獻及穩態計算結果計算出每個節段上的物性參數與各微分方程組的系數。

圖5 第i個節段上的瞬態計算模型Fig．5 Transient state model of segment i in IHXs

通過本文中的瞬態模型中變量取值，帶入上述矩陣方程分析，在各種工況下增長矩陣M－1N的譜半徑始終小于1，符合Von－Neumann條件，故差分格式是穩定的。

根據計算出來的每個微分方程組的系數，可得方程組的最大剛性比達106左右，對于這種剛性常微分方程組，采用ode15s方法求解，選取相對精度為10－4。

采用這種模型和算法，在其他參數均不變，對于為了便于比較，使輸出量趨于同一方向，輸入量分別選取100%功率工況下：一次側流量階躍降低5%、二次側流量階躍上升5%、一次側入口溫度階躍降低5%及二次側入口溫度階躍降低5%，對一、二次側出口溫度影響分別如圖6、圖7所示。

圖6 一回路出口溫度Fig．6 Temperature at primary outlet

圖7 二回路出口溫度Fig．7 Temperature at secondary outlet

計算結果表明：

（1）流通延遲

一次側入口溫度對其出口溫度影響有約6s的延遲、二次側（中間回路）入口溫度對其出口溫度有約2．9s的延遲，這與文獻［12］中相應的分別有5．8s、2．9s延遲的結論一致的；

（2）一次側出口溫度

二次側入口溫度對一次側出口溫度影響較之其他3個參數的影響大。二次側入口溫度階躍降低后，一次側出口溫度迅速降低，經過約80s穩定在最終值附近。二次側入口溫度降低約29．2K，一次側出口溫度降低約23．5K，這與穩態計算結果相差較小，這是因為，如圖5，在瞬態模型建立中將鈉溫度、流量對傳熱系數的影響加入到相應的程序中。

圖8 一回路流量衰減工況計算結果Fig．8 Simulation result of loss of primary flow transient

一、二次側鈉流量對一次側鈉出口溫度的影響均較小，約5K，且動態響應過程也類似，經過約90s一次側出口鈉溫穩定在最終值附近；一次側入口溫度對其出口溫度影響較之鈉流量的影響稍大，動態響應過程也稍長一些，經過約100s穩定在最終值附近；

（3）二次側（中間回路）出口溫度

一次側入口溫度對二次側（中間回路）出口溫度影響較之其他3個參數的影響大。一次側入口溫度階躍降低后，二次側出口溫度迅速降低，經過約90s穩定在最終值附近。一次側入口溫度降低約39．4K，二次側出口溫度降低約36．5K，這與穩態計算結果相差很小。一、二次側鈉流量和對二次側出口溫度影響均較小，并且過程較為緩慢、平滑，經過約90s穩定在最終值附近。

另外，在一回路鈉流量以時間常數為15s指數衰減的工況下進行計算，給出了一、二次側溫度的25s、50s、75s時的計算結果，如圖8所示，其中各組曲線中溫度較高的為一次側、較低的為二次側。

4 結論及展望

本文使用Matlab／Simulink軟件建立的中間熱交換器模型，兩種穩態計算方法計算出的符合性很好，與實際差別較小。瞬態模型計算迅速、數值穩定。與以往的中間熱交換器建模方法相比，本文給出的兩種穩態計算模型考慮了各物性參數變化對計算結果造成的影響。瞬態計算模型考慮了流量及鈉溫變化對換熱系數造成的影響，適用于大擾動的計算，例如主泵墮轉流量衰減。同時，瞬態模型采用Simulink模塊庫中的模塊建模，從而避開了具體算法程序的編寫，可以方便地對模型進行修改。因此本模型可用于中間熱交換器的熱工水力設計以及快堆核電廠的控制系統設計中。

未來將在中國實驗快堆（CEFR）的試驗中對本模型進行驗證，或者在其他公開文獻中獲得數據對其進行驗證。

同時，將進一步對中間熱交換器建模仿真進行研究，從一維向三維熱工、水力耦合模型擴展，從正常工況向自然循環工況擴展，從而獲取一種更為精確的模型。

［1］ 蘇著亭，葉長源，閻鳳文，等．鈉冷快增殖堆［M］．北京：原子能出版社，1991：332，361－362．

［2］ Tang Y S，Coffield Jr R D，Merkley R A．Thermal Analysis of Liquid Metal Fast Breeder Reactors［M］．ISBN 0894480111．La Grange Park，Illinois，USA：American Nuclear Society，1978．319－332，369－373．

［3］ Devlin R W，Bresnahan J D．FFTF and CRBRP intermediate heat exchanger design，testing，and fabrication［C］．Seminar on LMFBR components，Canoga Park，CA，USA，1977

［4］ CheungA C，Kao T T，Cho S M，et al．A Comparison of Between DEMO 1－D and the COMMIX 3－D Analysis of the CRBRP－IHX Under a Natural Circulation Event［C］．21st ASME／AIChE National Heat Transfer Conference．Seattle，Washington，USA，1983．

［5］ GunbyA L．Intermediate Heat Exchanger Modeling for FFTF Simulation［R］．AEC Research ＆Development Report of Battle Northwest Library．BNWL－1367，1970．

［6］ Brukx J F L M．CSDT Simulation Model for Prediction of IHX Transient Behavior Using S／360CSMP［D］．Delft，Netherlands：Delft University of Technology，1974．

［7］ 段天英．中國實驗快堆控制系統的仿真研究［D］．中國原子能科學研究院博士學位論文，1999．

［8］ 段天英．熱交換器仿真模型的修正［J］．中國試驗快堆年報，2000：5．

［9］ Constantine P，Tzanos．Liquid－Metal Fast Reactor Breeder Beactor Intermediate Heat Exchanger Transient Modeling for Faster than Real－Time Analysis［J］．Nuclear Technology，1987，76（3）：337－351．

［10］ Constantine P，Tzanos．Validation of an Intermediate Heat Exchanger Model for Real Time Analysis［C］．American Nuclear Society and Atomic Industrial Forum Joint Meeting．Washington DC，USA，1986．

［11］ 居懷明，徐元輝，李懷萱．載熱質熱物性計算程序及數據手冊［M］．北京：原子能出版社，1990：92－99．

［12］ 中國實驗快堆最終安全分析報告［R］．中國原子能科學研究院，2008．