- Email: [email protected]
Studies of the deteriorated turbulent heat transfer regime for the gas-cooled fast reactor decay heat removal system
Nuclear Engineering and Design 237 (2007) 1033–1045
Studies of the deteriorated turbulent heat transfer regime for the gas-cooled fast reactor decay heat removal system Jeong Ik Lee ∗ , Pavel Hejzlar, Pradip Saha, Mujid S. Kazimi Center for Advanced Nuclear Energy Systems, Department of Nuclear Science and Engineering, Massachusetts Institute of Technology, Cambridge, MA 02139, USA Received 21 November 2006; received in revised form 11 January 2007; accepted 11 January 2007
Abstract Increased reliance on passive emergency cooling using natural circulation of gas at elevated pressure is one of the major goals for the gas-cooled fast reactor (GFR). Since GFR cores have high power density and low thermal inertia, the decay heat removal (DHR) in depressurization accidents is a key challenge. Furthermore, due to its high surface heat flux and low velocities under natural circulation in any post-LOCA scenario, three effects impair the capability of turbulent gas flow to remove heat from the GFR core, namely: (1) acceleration effect, (2) buoyancy effect and (3) property variation. This paper reviews previous work on heat transfer mechanisms and flow characteristics of the deteriorated turbulent heat transfer (DTHT) regime and decomposes governing non-dimensional groups to provide an insight for the GFR DHR system design. It is shown that by applying the developed methodology, the GFR’s DHR system has a potential for operating in the DTHT regime and the gas DTHT regime is different from the liquid or super-critical fluid’s DTHT regime. A description of an experimental facility designed and built to investigate the DTHT regime is provided together with the initial test results. The initial runs were performed in the forced convection regime to verify facility operation against well-established forced convection correlations. The results of the three runs at Reynolds numbers of 6700, 8000 and 12,800 show good agreement with the Gnielinski correlation for heat transfer, which is considered the best available correlation for forced convection over a wide range of Reynolds and Prandtl numbers. However, even in the forced convection regime, the effect of variation of the fluid properties was found to be significant. © 2007 Elsevier B.V. All rights reserved.
1. Introduction The gas-cooled fast reactor (GFR) concept for Generation IV has generated considerable interest and is under development in the U.S., France, and Japan. One of the key candidates is a blockcore configuration first proposed by MIT (Hejzlar, 2001), which was shown in (Williams et al., 2004), to have the potential to operate in deteriorated turbulent heat transfer (DTHT) regime or in the transition between the DTHT and normal forced turbulent or laminar convection regime during post-loss-of-coolant accident (LOCA) conditions. The definition and the physical phenomena of the DTHT regime will be explained in the next
∗
Corresponding author. Tel.: +1 617 253 6188. E-mail addresses: [email protected] (J.I. Lee), [email protected] (P. Hejzlar), [email protected] (P. Saha), [email protected] (M.S. Kazimi). 0029-5493/$ – see front matter © 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.nucengdes.2007.01.003
section. This is contrary to most industrial applications where operation is in a well-defined and well-known turbulent forced convection regime. As a result, an important need emerged to develop heat transfer correlations that make possible rigorous and accurate predictions of decay heat removal (DHR) during post-LOCA in these regimes. Therefore, a collaborative INL/MIT project was initiated, to acquire experimental data (to the extent possible) for heat transfer correlations for all regimes including the transition regions in such a manner that there would be no discontinuities at the boundaries. The paper is organized as follows: first, a short literature review on the DTHT regimes will be provided followed by a simple analysis that shows the GFR DHR system has a strong potential to operate in this regime. Next section will give a simple comparison between gas, liquid and super-critical fluid to show that existing experimental data of liquid and super-critical fluid in the DTHT regime may not be appropriate for understanding gas DTHT. An overview of the loop construction with a detailed
1034
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Nomenclature A Bo* cp cts D f g G Gr* h k kts Kν Nu Pr Pw Re qb+ qw Q t T U x
area buoyancy parameter = Gr* /Re3.425 Pr0.8 specific heat at constant pressure (fluid) specific heat of test section material pipe diameter friction factor gravitational acceleration mass flux D4 /kν2 Grashof number with heat flux = gβqw heat transfer coefficient thermal conductivity of gas thermal conductivity of test section material acceleration parameter = (ν/Ub2 )(dUb /dx) ≈ + (4qb /Re) Nusselt number = hD/k Prandtl number = ν/α wetted perimeter Reynolds number =Ub D/ν /Gc T non-dimensional heat flux = qw p b wall heat flux volumetric flow rate time temperature flow velocity flow direction
Greek letters Λ volume α thermal diffusivity = k/ρcp β thermal expansion coefficient = −(1/ρ)(∂ρ/∂T)p μ viscosity ν kinematic viscosity ρ density Subscripts b bulk i inner ts test section w wall
description of each major loop component will be presented. Finally, the forced circulation data obtained from the loop will be presented and compared to the established correlations to validate the loop performance. 2. Review and decomposition of the deteriorated heat transfer regime controlling parameters In McEligot and Jackson (2004), the conditions of onset of DTHT regime based on three different physical effects: (1) flow acceleration, (2) buoyancy effect and (3) thermo-physical properties variation are well summarized. The reason it is called deteriorated turbulent heat transfer regime is because normal turbulent convective heat transfer is hindered by these phenomena
and the fluid’s heat transfer capability drops significantly in this region. This regime was first found in the super-critical fluid heat transfer and, later on, it was also found in general gas and liquid heat transfer systems, where the ratio of heating to flow rate is high. Even though the onset of this regime due to acceleration and buoyancy effects is well defined, the heat transfer coefficient and the friction factor for this regime have not gained universal agreement among the investigators. In contrast, thermo-physical property variation has a wide range of experimental data and is well established as reviewed in McEligot and Jackson (2004). 2.1. Acceleration effect The acceleration effect is also known as a “laminarization due to a favorable pressure gradient”. The decrease in turbulence transport occurs whenever a flow directional acceleration exceeds a certain value. The stream wise acceleration can be quantified by using the definition of acceleration parameter, Kν = (ν/Ub2 )(dUb /dx). In a heated system the acceleration occurs from gas expansion. This effect is different from the buoyancy effect, since the effect is not related to the gravitational force. When a flow condition satisfies the onset of laminarization due to acceleration, the heat transfer coefficient of the fluid quickly drops to the laminar heat transfer coefficient value, even though the Reynolds number is well above the critical Reynolds number. Applying energy balance and continuity equation with perfect gas and constant cross section assumptions to the original definition of Kν , an alternative definition of an acceleration parameter can be obtained for the heated case in the form of a non-dimensional number Kν ≈ 4qb+ /Re (McEligot et al., 1969). Based on McEligot and Jackson (2004), we can normalize the non-dimensional number that represents acceleration effect with the onset criterion (NKν = Kν /3 × 10−6 ). Thus, whenever the normalized value exceeds one, it indicates that the heat transfer regime is changing from normal forced convection regime to the DTHT regime due to the acceleration effect. It is also possible to separate the normalized acceleration parameter into two groups. One is the controlled group1 (CKν ), which includes geometry, heat flux and the volumetric flow rate. The other is the fluid properties group (PKν ), which is dependent on the system conditions such as pressure and operation temperature, and is not directly controlled. It should be noted that all the non-dimensional groups presented in this paper are evaluated at the fluid bulk temperature. Later in the paper, analysis will be focused on the properties group to determine which system pressure and operating temperature will likely force the fluid flow into the DTHT regime: 4qb+ Kν ≈ 3 × 10−6 3 × 10−6 Re qw APw μβ = = CKν PKν Q2 3 × 10−6 cp ρ2
NKν =
(1)
1 This group was designated “controlled”, because it consists of parameter that can be controlled by designer.
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
1035
2.2. Buoyancy effect The buoyancy effect originates from the density gradient due to heating, and reduces the turbulent heat transport when the flow direction is the same as the buoyancy force. Hall and Jackson (1969) explained the decrease in turbulence by the shear stress re-distribution in the flow. On the other hand Petukhov and Polyakov (1988) reported that the phenomenon is governed by two terms that are competing in the turbulent energy equation. One is the velocity gradient and the other is the turbulence work that needs to be provided to overcome a stable density gradient. Nevertheless, the criterion for onset of buoyancy induced DTHT is similar in both cases. In this paper we will follow definition of McEligot and Jackson (2004), which is an extended version of Hall and Jackson (1969) discussion. Similar to the acceleration effect, the buoyancy group will be normalized with the onset criterion (NBo* = Bo* /6 × 10−7 ). The buoyancy effect can also be viewed as a product of a control group (CBo* ) and a properties group (PBo* ) (Eq. (2)). The buoyancy properties group will be utilized as a measure of how the system pressure and operating temperature combination shifts the system to operate in the DTHT regime: Bo∗ Gr ∗ = 6 × 10−7 6 × 10−7 Re3.425 Pr 0.8 A4 gμ0.625 β qw = Q3.425 Pw 0.575 2.704 × 10−7 cp 0.8 k0.2 ρ1.425
NBo∗ =
= CBo∗ PBo∗
(2)
2.3. Thermo-physical properties variation When strong heating is provided to a fluid, thermo-physical properties vary across the flow area due to the radial temperature gradient of the fluid. Therefore, the heat transfer coefficient based on the assumption of constant fluid properties does not hold anymore and needs to be modified. Since liquid and gas properties changes with the temperature are different, the modification function is also different. However, thermo-physical property variation does not have a strong impact on the heat transfer coefficient like the acceleration effect or buoyancy effect. This effect has to be considered seriously if the ratio between wall temperature and bulk temperature exceeds 1.1, which causes 5% decrease in the Nusselt number McEligot and Jackson (2004).
Fig. 1. Acceleration properties group.
in the analysis, since it is relatively well known compared to the other two effects. It is clear from the figures that as temperature increases and pressure decreases, helium and carbon dioxide property groups increase. This indicates that the normal gas turbulent heat transfer has stronger tendency to fall into DTHT regime when the system depressurizes and the operating temperature increases. This situation occurs in the GFR system during the loss-ofcoolant accident. Therefore, both acceleration and buoyancy DTHT criteria need to be checked when designing the GFR DHR system. The figures also show that for the same controlled variable group, helium is more susceptible to the DTHT regime than carbon dioxide. Thus, more attention is needed when designing a helium-cooled system. Fig. 3 shows one of the possible operating ranges calculated with LOCA-COLA, an in-house code (Williams et al., 2004), of a GFR DHR system (Cochran et al., 2004). The Reynolds number is normalized to 4000 where fully developed adiabatic turbulent flow regime starts, and normalized DTHT parameters are NKν and NBo* . Curves in Fig. 3 show the axial variation of NKν and NBo* for given system pressure and heated diameter of the hot channel with helium and carbon dioxide as an operating fluid.
3. Simple analysis of GFR DHR system design Based on Eqs. (1) and (2) we can perform a simple parametric study with the gas properties. Since, the GFR design is considering helium and carbon dioxide as candidates for the coolant, both gases are evaluated. Figs. 1 and 2 show both helium and carbon dioxide acceleration properties groups and buoyancy properties groups for different pressure and temperature, respectively. It is noted that thermo-physical property variation effect is neglected
Fig. 2. Buoyancy properties group.
1036
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 3. GFR DHR operation range.
It can be observed that there is a possibility for GFR DHR system to operate not only in the DTHT regime but also in the laminar to turbulent transition regime. The laminar to turbulent transition criteria are not fully understood even in the adiabatic flow situation and a complete set of study for heated flow is rarely found (Lee, 2005). 4. DTHT regime for different fluids In this section, behavior of the buoyancy and acceleration parameters during the experimental runs will be discussed. Since the MIT experimental facility is a single flow channel, which will be further described in details in the next section, Eqs. (1) and (2) need to be rearranged to fit the single channel analysis. This is because along the channel the mass flow rate is constant not the volumetric flow rate. Eqs. (3) and (4) are rearranged forms of Eqs. (1) and (2): μβ qw APw NKν = = Cm˙ Kν Pm˙ Kν (3) m ˙2 3 × 10−6 cp ∗
NBo =
A4 qw m ˙ 3.425 Pw 0.575
= Cm˙ Bo∗ Pm˙ Bo∗
ρ2 gμ0.625 β 2.704 × 10−7 cp 0.8 k0.2
Fig. 4. Nitrogen acceleration properties group.
(4)
Figs. 4 and 5 plot the acceleration properties group and buoyancy properties group behavior of nitrogen at three different pressures with varying temperature. Since most of the gas thermal properties behave similarly with the temperature, nitrogen is chosen as an example. First observation from the comparison of Fig. 1 against Fig. 4 and Fig. 2 against Fig. 5 is the trend reversal. This is because the density term in the properties group moved to the controlled group side (since m ˙ = ρQ). Since the bulk temperature of the fluid will increase as it flows toward the downstream due to the heating, it can be predicted from the figures that both the
Fig. 5. Nitrogen buoyancy properties group.
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 6. Liquid water acceleration properties group.
buoyancy parameter and the acceleration parameter of the gas flow will decrease along the downstream. Therefore, it can be concluded that the maximum values of the buoyancy and acceleration numbers for the gas flow will be at the channel inlet when the control group is fixed as a constant. The control group will be approximately a constant for the MIT facility along the channel, since the test section has circular tube shape with a constant diameter and is designed to operate with uniform heat flux boundary condition. Another interesting observation is that the buoyancy parameter properties group of nitrogen increases with the pressure while the acceleration parameter properties group is pressureindependent. This indicates that by increasing the operating pressure of the facility, the buoyancy effect will increase relative to the acceleration effect. Various correlations were developed from the experimental data collected with working fluid of water and super-critical phase fluid (see Celata et al., 1998; Jackson et al., 1989; Kakac¸ et al., 1987). Thus, it is of interest to check how these fluids behave compared to gases. Figs. 6–9 show the properties group trends for liquid water and super-critical CO2 . If we compare the acceleration number properties group of water and nitrogen (Figs. 4 and 6), their behavior is also different. The nitrogen case shows steady decrease with increasing temperature, but water shows parabolic behavior with temperature. However, in contrast to the buoyancy parameter (Figs. 5 and 7), the order of magnitude of the acceleration property group is the same for both water and nitrogen cases. From the liquid water buoyancy properties group behavior (Fig. 7), one can see that the trend is totally different from that of the gas. The liquid water buoyancy parameter will increase or stay the same downstream the channel, since the properties group increases with the fluid bulk temperature. This is because water is incompressible fluid and the thermal properties variation due to the temperature change in the low-pressure range is insignificant compared to gas. One can also conclude from the figures that the order of magnitude for the water buoyancy properties group is four times higher than that of the low-pressure gas.
Fig. 7. Liquid water buoyancy properties group.
Fig. 8. Super-critical CO2 acceleration properties group.
Fig. 9. Super-critical CO2 buoyancy properties group.
1037
1038
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Therefore, it can be predicted that water experiments can reach buoyancy induced DTHT regime relatively easier than gas experiments, if the control group is the same. In addition, the channelwise behavior due to the buoyancy force will also be different. When we compare the super-critical carbon dioxide to nitrogen, the situation is different from the comparison of water to nitrogen. The trend with temperature is the same for both fluids but the order of magnitude change with temperature is different. This is because of a dramatic jump in both property groups for the super-critical CO2 near the critical point (Figs. 8 and 9), where the properties undergo a steep change. Therefore, it can be predicted that an experiment using a super-critical fluid will enter the DHTH regime, induced by either buoyancy effect or acceleration effect, where the temperature is near the critical point. From the comparison of all fluids, the following observations can be made: 1. The channel behavior of the buoyancy number will be very different among liquid, gas and super-critical fluids. 2. For gases, both maximum buoyancy number and maximum acceleration number will occur at the inlet of the heated channel.
3. For water, the maximum buoyancy number will be at the outlet of the heated channel and the maximum acceleration number will depend on the bulk fluid temperature. 4. For the super-critical carbon dioxide, both the maximum buoyancy number and maximum acceleration number will be at the point where the bulk temperature is near the critical point. 5. By pressurizing the gas, the buoyancy effect can be made more pronounced while the acceleration effect is immune from the pressure change for the same control group. The simple analysis presented in this section shows that the gas DTHT and the liquid or super-critical fluid DTHT regime are going to have different characteristics along the single channel and this prohibits us from utilizing the liquid or super-critical fluid experimental data for the gas system analysis. Therefore, the general heat transfer coefficient and friction factor in the gas DTHT regime and all the transitions from normal turbulent flow to the gas DTHT regime are not well defined yet, experimental studies in these heat transfer regimes are necessary to develop reliable correlations for the design of GFR DHR system.
Fig. 10. Schematic diagram of experimental facility.
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
5. Description of the experimental facility Fig. 10 is a simplified diagram, with dimensions, of the experimental facility (Cochran et al., 2004). The facility can be viewed as three different systems, as indicated by boxes in Fig. 10. The first system is the main loop, where most of the measurements are made. The main loop can be further divided into three sections—the test section, the chimney section and the downcomer section. The second system involves the compressor compartment located at the bottom part of Fig. 10. This system is necessary to operate the facility in a high flow rate range and to initiate a natural circulation flow. The last system is the gas charging and the vacuum system. This system is on the top right side in Fig. 10. The vacuum system is needed to pump out the air from the main loop before the test gases such as pure nitrogen, carbon dioxide and helium are charged into the main loop through the gas charging system. Each system and sub-system will be described in more detail in the following subsections. 5.1. Main loop With the exception of the test section, the main loop is constructed of 1 in. Type 316 stainless steel tube. The test section and flow development region, ahead of the test section, are made of a 0.75 in. Type 316 tubing. The test section wall thickness is 0.065 in. and the inner diameter of the test section is 0.62 in. (15.7 mm). The 1 in. stainless steel tube has a wall thickness of 0.049 in. and an inner diameter of 0.902 in. (22.9 mm). Connecting the test section and flow development section to the larger tubing of the loop is a simple reducing union at the entrance and exit of the test section. Both reducing sections change the tube sizing from 1 to 0.75 in. All the tube connections are done with compression type stainless steel fittings purchased from SWAGELOK. The test section is 1.976 m long and the chimney section is 3.797 m long. Power taps at the top and bottom of the test section connect a dc power supply providing up to 1500 A, at 10 V to the test section. This power supply provides for either voltage or current control of the power applied to the loop. Remote control of the power supply by the test control software is provided by a DA converter installed in a HP 3852A Data Acquisition unit. Pressure taps for measuring differential pressure across the test section have been welded to the tubing adjacent to the power taps. As can be seen from the diagram (Fig. 10) the tap at the entrance to the test section is also used for measuring system pressure. Differential pressure taps have also been installed on the downcomer section. All pressure taps have been connected to solenoid valves for selecting measurements from either the test section or the downcomer section. Teflon tubing has been used to connect pressure tap output to a differential pressure transducer. This transducer is a MKS Type 120A high accuracy pressure transducer. The range of this pressure transducer can be selected from 1 to 25,000 mmHg (133.3 Pa–3.33 MPa) with resolution of 10−6 of the full-scale range. Accuracy is ±0.12% of reading. The unit has its own readout equipment. Teflon tubing has been used for providing electrical isolation of the differential pressure transducer. Solenoid valves are con-
1039
trolled with relays in the same HP3852. Output of the system pressure transducer is monitored by a National Instrument Data Acquisition (NI-DAQ) system. Two flow thermocouples have been installed adjacent to the flow development and test sections. The first is 5.188 in. below the heated (test) section with the other installed 4.750 in. above the outlet of the heated section. This layout was necessary since the test section is heated using a direct Joule heating method, and installing flow thermocouples in the test section would create non-uniform heating in the vicinity of the thermocouple. Moreover, it was impractical to install them near the power taps. The flow thermocouples are connected to the NI-DAQ. 5.2. Instrumentation Temperatures of the loop were to be measured by Type K thermocouples spot-welded directly to the tubing. Initial testing showed significant errors in temperature readings during operation of the dc power supply. Analysis of the readings determined the cause of this problem to be electrical noise from the power supply. The problem was resolved by installing ungrounded sheathed thermocouples on the test section. Ungrounded thermocouples provide isolation from the electrical noise. Sheathed thermocouples required a different method of attachment to the test section. Each thermocouple was brazed to a small piece of shim stock, which was then spot-welded to the test section tubing. These sheathed thermocouples are used to monitor loop temperatures while dc heating current is applied to the loop. Due to concerns about slight errors with the sheathed thermocouples, temperatures used for determining heat transfer parameters are measured with thermocouples spot-welded directly to the tubing. Readings are taken during brief (2 s) interruptions of dc heating current. Sheathed thermocouples are connected to the NI-DAQ and the welded on thermocouples are connected to the HP3852. Twenty welded thermocouples and four sheathed thermocouples are installed in the test section and nine welded thermocouples and three sheathed thermocouples are installed in the chimney section (see Fig. 11). All the thermocouples are K-type and manufactured from OMEGA. To prove that the test section wall temperature does not drop appreciably during heating current interruptions, a simple idealized problem was solved. Eq. (5) is the solution to a lumped parameter model of the tube wall, selected to illustrate temperature response of the wall with time: hAi T (t) − Tb = exp − t (5) T (0) − Tb ρts cts ts Using Eq. (5), the time that is required for a test section temperature to drop 1 K from the steady state temperature, given the heat transfer coefficient of a gas to be 100 W/m2 K, and the temperature difference between the bulk gas and test section wall, at steady state, to be 200 K (conservatively high values for this facility), can be calculated to be 227 s. Therefore, even though the problem setup is simplified, it confirms that wall temperature changes during the short 2-s power-off period are negligible. This is mainly because stainless steel test section total
1040
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 11. Instrumentations and components in test section and chimney section.
heat capacity is much higher than the gas heat transfer capability. This was further confirmed during the experimental runs, where no visible temperature drops could be observed during the offpower period when the data acquisition system was collecting the data. The test section and the chimney section are under a layer of insulation. The insulation material is mineral wool. The test section insulation is equipped with thermocouples installed at six axially different positions and three azimuthally different positions (total 18) to measure the thermal resistance of the insulation and the heat losses during the experiment. On the chimney section, thermocouples are installed on the insulation at nine axially different positions and three azimuthally different positions, but at some places only one azimuthal thermocouple is installed (total 21). These insulation thermocouples will be used for not only measuring the heat losses but also to control the guard heater operation.
Fig. 11 shows the instrumentation and the components that form the test section and the chimney section. All the dimensions of the distances between thermocouples and their positions are indicated. The downcomer section includes a heat exchanger, pressure taps and electrical break. The gas charging and vacuum systems are connected to the downcomer section. One inch stainless steel tube is used for the downcomer section and the same compression type stainless steel fittings are used for the riser section. The length of the heat exchanger is 2 m and it is helically twisted with 30 cm coil diameter (6 and 2/3 turns). The building cooling water provides relatively constant temperature water (approximately 283 K). Two flow thermocouples (shown in Fig. 12) are installed near the heat exchanger; one at the inlet and the other at the outlet of heat exchanger. The two flow thermocouples are connected to the NI-DAQ. Building cooling water
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
1041
Fig. 12. Instrumentations and components in downcomer section.
temperature is also monitored by two thermocouples welded on to the heat exchanger. These thermocouples are also connected to the NI-DAQ. Pressure taps are installed through the connections on the downcomer, as indicated in Fig. 12. The pressure taps are controlled with the solenoid valves described in the previous section. The differential pressure transducer is very sensitive and will be used together with the hot-wire in the downcomer for the determination of flow rate in very low Reynolds number flow regimes. An electrical break is essential to eliminate a parallel electrical path around the test section. This insulating break was manufactured and installed in the downcomer section. Since pressure drop in the downcomer section is to be measured, great care was taken to achieve maximum smoothness of the inner surface of the break and smooth transition between the break and the tube. Four welded thermocouples and two sheathed thermocouples were installed on the downcomer section with the same
reasoning as explained in the previous section. The sheathed thermocouples are wired to the NI-DAQ and welded on thermocouples are connected to the HP3852. The downcomer section is also equipped with the same insulation material used on the test section. The insulation is needed since it is desirable to maintain a constant gas temperature along the downcomer for reliable pressure drop measurements and to prevent personnel contact with the downcomer tube wall in case that downcomer temperature exceeds 40 ◦ C in some runs. 5.3. Blower section and gas charging system A blower system is required for forced circulation at high Reynolds number testing and to initiate natural circulation. A three cylinders compressor unit supplies the input necessary to circulate the gas around the loop. Gas from the downcomer section of the loop is routed to an accumulator tank at the end of the downcomer section. The accumulator provides smooth-
1042
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 13. Instrumentations and components in compressor.
ing of the suction flow to the compressor. The compressor has been re-plumbed to connect all cylinders in parallel. After passing through the cylinders the gas is fed into a second larger accumulator tank to smooth the outgoing flow to the loop. No compression of the gas is done as it passes through the unit during operation. When the loop is vacuumed before the new gas is charged, the compressor unit is isolated from the remainder of the loop, by six valves installed in the inlet/outlet lines immediately adjacent to the compressor cylinders, to prevent the lubrication oil in the compressor from seeping out and contaminating the test section (Fig. 13). Three additional ball valves installed in the loop, near the compressor unit provide the isolation of the compressor from the loop for natural circulation testing. A variable frequency drive unit (VFD) provides speed control of the compressor unit. A relief valve is provided at the outlet to the compressor to prevent possible over-pressure event. After the tank outlet, a rotameter is installed to measure the flow rate in the experiments. The vacuum and gas charging systems control the gas purity in the main loop and the gas pressure. The vacuum system will first remove the air from the main loop except for a small volume in the compressor by closing six valves. After the loop is vacuumed, gas charging system is connected to the main loop and provides the test gases such as nitrogen, carbon dioxide and helium. The gas charging system is controlled with a pressure regulator. During the experiment, the gas charging system is going to be connected to the system to maintain a constant pressure and make up for any minor leakage in the main loop. 6. Forced convection results Three different experimental conditions were tested as a first check of the facility performance. These conditions were selected to perform experiments at normal forced convection turbulent heat transfer regime, with a slight effect from the thermo-physical properties variation of the fluid. Table 1 shows
the conditions for different experimental sets. All the tests were performed with nitrogen gas. For the run designation, the first three digits indicate the Reynolds number (for example 067 stands for Re at the inlet of about 6700) and last two digits indicate the qb+ (e.g., 13 stands for qb+ = 0.0013). qb+ is the ratio of wall heat flux to inlet flow enthalpy flux, which is a key indicator for heating to flow rate ratio. The three experimental sets cover the range of local Reynolds numbers from 4700 to 12,800, which are usually considered within the turbulent flow regime in the pipe flow. The correlations selected for comparison are Dittus–Boelter, and Gnielinski. Dittus–Boelter correlation was selected because it is the most widely used correlation for forced convection turbulent heat transfer. However, the Dittus–Boelter equation used in the comparison differs from the widely known form, which uses the value of 0.023 as a leading constant instead of 0.021 presented in Eq. (6). This is because the widely known form is more appropriate for liquids. Gnielinski’s correlation (Eq. (7)) is a newer correlation, which exhibits better accuracy than Dittus–Boelter and is valid for a larger range of Reynolds and Prandtl numbers.
Table 1 Experimental conditions
Mean power (W) Mass flow rate (kg/s) Mean pressure (MPa) Inlet Reynolds number Outlet Reynolds number qb+ inlet Inlet bulk temperature (K) Outlet bulk temperature (K) Minimum wall temperature (K) Maximum wall temperature (K) Wall to bulk temperature ratio
06713 run
07915 run
12823 run
397.5 0.0015 0.1380 6738 4624 0.0013 298.6 494.4 330.5 584.9 1.09–1.30
515.0 0.0018 0.1533 7992 5304 0.0015 299.0 521.0 337.2 624.1 1.11–1.35
1187.1 0.0028 0.2279 12871 7438 0.0023 299.4 622.1 363.9 783.3 1.17–1.52
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Modified Dittus–Boelter correlation (McEligot et al., 1966) Nu = 0.021Re0.8 Pr 0.4
(6)
Gnielinski correlation (Gnielinski, 1976) Nu =
(f/8)(Re − 1000)Pr √ 1 + 12.7 f/8(Pr 2/3 − 1)
(7)
The friction factor that is used in the Gnielinski correlation is given in Eq. (8). Colebrook correlation (valid for 4000 < Re < 107 ) (Kaka¸c et al., 1987) 2 16 f = (8) 1.5635 ln(Re/7) Before presenting the results, it is important to note that the Gnielinski correlation given in Kakac¸ et al. (1987), omitted a correction factor that accounts for the effect of property variation in the film for cases where the temperature difference between the wall and bulk is large. This correction factor is defined in the original Gnielinski’s publication (1976) and it is given in Eq. (9). Similar correction factor exists for the Dittus–Boelter correlation—see Eq. (10) (McEligot et al., 1966): Tb 0.45 Correction factor (CF) = , Twi NuGnielinski modified = NuGnielinski Eq. (7) × CF CF =
Tb Twi
Fig. 14. Nusselt number ratio vs. axial position (06713 run).
(9)
0.5 ,
NuDittus–Boelter modified = NuDittus–Boelter Eq. (6) × CF
(10)
The correction factor is necessary for the gas since the wall to bulk temperature ratio range is 1.1–1.5 in all experimental runs affecting significantly properties in the film, and thus heat transfer coefficient. Table 1 shows the range of wall to bulk temperature ratio for three cases. Figs. 14–16 plot the ratio between the experimental Nusselt numbers and the Nusselt numbers calculated from correlations. Since the correlations shown above are for fully developed region, Figs. 14–16 are magnified with the thermal entrance region data points left out for the reader’s convenience. In the fully developed region, experimental runs 06713 and 07915 match the correlations better than the 12823 run. Experiments 06713 and 07915 have nearly the same ratio between 0.95 and 1.0 with respect to the Gnielinski correlation with the correction factor in the fully developed region. Since the 06173 and 07915 runs have smaller errors than the 12823 run in terms of energy, the larger discrepancy of 12823 run with the correlation can be explained. Both correlations, the Dittus–Boelter and the Gnielinski, show good agreement with the experimental data after the correction factor is incorporated. This indicates that the property variation due to the large temperature gradient in the fluid has
Fig. 15. Nusselt number ratio vs. axial position (07915 run).
Fig. 16. Nusselt number ratio vs. axial position (12823 run).
1043
1044
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 17. Nusselt number comparison with error bars (06713 run).
Fig. 19. Nusselt number comparison with error bars (12823 run).
a significant effect on the heat transfer even for the forced convection regime. It can be also concluded that both the Dittus–Boelter and Gnielinski correlations predict the experimental data at low Reynolds number regime (4000–13,000) satisfactorily. The last point of experimental Nusselt number in each plot rises above 1.0, potentially due to the under-estimation of the axial heat losses at the top of the heated section. This may require installation of a thermocouple right after the heated section to reduce the axial heat loss at the outlet. The error analysis will be discussed further in the following section. Figs. 17–19 show the experimental Nusselt number and the Nusselt number from the Gnielinski correlation that includes the correction factor. The Gnielinski correlation has an error of ±20% (Gnielinski, 1976) and it is represented in the figures representing error bars with the uncertainty of the experimental Nusselt number.
Figs. 17–19 show the 06713 and 07915 runs to have very good agreement with the Gnielinski correlation with the property temperature dependence factor. The 12823 run gives slightly higher Nusselt number than the correlation. This was due to the larger error in estimating the thermal resistance of the insulation data for reduction of the 12823 run, since the wall temperature of the 12823 run is significantly above the temperature at which the thermal resistance was measured. Fig. 20 compares the experimental Nusselt number versus the Reynolds number against the Gnielinski correlation with the correction factor for all three runs. The entrance and the end points, where the uncertainty is the highest were excluded. It is expected that further reduction in the uncertainty can be achieved after: (1) a rotameter is replaced with a turbine meter, (2) guard heaters are installed on top of the primary insulation
Fig. 18. Nusselt number comparison with error bars (07915 run).
Fig. 20. Nusselt number comparison (Re vs. Nu).
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
and (3) more thermocouples on the insulation to determine the heat loss with higher accuracy. 7. Summary and conclusions The work presented in this paper can be summarized as follows: 1. It was shown that gas-cooled fast reactor (GFR) decay heat removal (DHR) system is likely to operate after LOCA in the deteriorated turbulent heat transfer (DTHT) regime. Since the heat transfer coefficient and friction factor correlations are not well established in this regime, new experiments were necessary to develop reliable correlations. 2. The non-dimensional groups reflecting the buoyancy and acceleration effects have been decomposed into two parts. One representing the design controlled parameters, such as geometry and heating rate. The other representing the fluid properties, which are not directly controlled in the design of a facility. 3. An experimental facility was constructed to investigate the DTHT regime. However, due to the power supply noise, a different temperature measurement procedure was developed, which combines welded thermocouples and sheathed thermocouples. It was shown by a simple analysis that the new measurement procedure gives accurate data, which was also confirmed by experimental runs. 4. Three forced circulation experiments were performed and showed very good agreement with the Gnielinski correlation that incorporates the effect of temperature dependent properties in the film, confirming that there is no significant defect on measurement system and the experimental loop. However, a few minor adjustments were identified to be done, in order to reduce measurement uncertainty. 5. The initial results already confirm that the Nusselt number correlation in forced convection has to incorporate temperature dependence of properties in the film to correctly predict gas heat transfer in the GFR channels since the wall-to-bulk temperature ratios exceed the value of 1.1. This is important as it leads to reduced heat transfer coefficient and thus to higher wall temperatures.
1045
Acknowledgments The authors would like to acknowledge the financial support of the Idaho National Laboratory. They thank P. Cochran of MIT for providing preliminary analysis of the experimental apparatus, P. Stahle of MIT for his valuable effort in the construction of the facility and Dr. Glenn E. McCreery and Dr. D.M. McEligot of INL for valuable comments on the facility design. References Celata, G.P., D’Annibale, F., Chiaradia, A., Cumo, M., 1998. Upflow turbulent mixed convection heat transfer in vertical pipes. Int. J. Heat Mass Transfer 41, 4037–4054. Cochran, P., Saha, P., Hejzlar, P., McEligot, D.M., McCreery, G.E., Schultz, R.R., 2004. Scaling Analysis and Selection of Test Facility for Fundamental Thermal-Hydraulic Studies related to Advanced Gas-cooled Reactor, INL/EXT-05-00158. Idaho National Laboratory and Massachusetts Institute of Technology, Department of Nuclear Science and Engineering. Gnielinski, V., 1976. New equations for heat and mass transfer in turbulent pipe and channel flow. Int. Chem. Eng. 16 (2), 359–387. Hall, W.B., Jackson, J.D., 1969. Laminarization of a turbulent pipe flow by buoyancy forces, ASME 69-HT-55. Hejzlar, P., 2001. A modular, gas turbine fast reactor concept (MFGR-GT). In: Proceedings of Transactions of the American Nuclear Society, vol. 84, Milwaukee, Wisconsin, June 17–21. Jackson, J.D., Cotton, M.A., Axcell, B.P., 1989. Studies of mixed convection in vertical tubes. Int. J. Heat Fluid Flow 10, 2–15. Kakac¸, S., Shah, R.K., Aung, W., 1987. Handbook of Single-phase Convective Heat Transfer. John Wiley & Sons (Chapter 15). Lee, J.I., 2005. The Flow Structure under Mixed Convection in a Uniformly Heated Vertical Pipe. Master Thesis, Nuclear Science and Engineering Department, Massachusetts Institute of Technology. McEligot, D.M., Coon, C.W., Perkins, H.C., 1969. Relaminarization in tubes. Int. J. Heat Mass Transfer 13, 431–433 (Shorter Communications). McEligot, D.M., Ormand, L.W., Perkins, H.C., 1966. Internal low Reynolds number turbulent and transitional gas flow with heat transfer. J. Heat Transfer 88, 239–245. McEligot, D.M., Jackson, J.D., 2004. Deterioration criteria for convective heat transfer in gas flow through non-circular ducts. Nucl. Eng. Design 232, 327–333 (Short Communication). Petukhov, B.S., Polyakov, A.F., 1988. Heat Transfer in Turbulent Mixed Convection. Hemisphere Publishing Corporation. Williams, W.C., Hejzlar, P., Driscoll, M.J., 2004. Decay heat removal from GFR core by natural convection. In: International Congress on Advances in Nuclear Power Plants ICAPP ’04, Paper 4166, Pittsburgh, USA, June 13–17.
Studies of the deteriorated turbulent heat transfer regime for the gas-cooled fast reactor decay heat removal system Jeong Ik Lee ∗ , Pavel Hejzlar, Pradip Saha, Mujid S. Kazimi Center for Advanced Nuclear Energy Systems, Department of Nuclear Science and Engineering, Massachusetts Institute of Technology, Cambridge, MA 02139, USA Received 21 November 2006; received in revised form 11 January 2007; accepted 11 January 2007
Abstract Increased reliance on passive emergency cooling using natural circulation of gas at elevated pressure is one of the major goals for the gas-cooled fast reactor (GFR). Since GFR cores have high power density and low thermal inertia, the decay heat removal (DHR) in depressurization accidents is a key challenge. Furthermore, due to its high surface heat flux and low velocities under natural circulation in any post-LOCA scenario, three effects impair the capability of turbulent gas flow to remove heat from the GFR core, namely: (1) acceleration effect, (2) buoyancy effect and (3) property variation. This paper reviews previous work on heat transfer mechanisms and flow characteristics of the deteriorated turbulent heat transfer (DTHT) regime and decomposes governing non-dimensional groups to provide an insight for the GFR DHR system design. It is shown that by applying the developed methodology, the GFR’s DHR system has a potential for operating in the DTHT regime and the gas DTHT regime is different from the liquid or super-critical fluid’s DTHT regime. A description of an experimental facility designed and built to investigate the DTHT regime is provided together with the initial test results. The initial runs were performed in the forced convection regime to verify facility operation against well-established forced convection correlations. The results of the three runs at Reynolds numbers of 6700, 8000 and 12,800 show good agreement with the Gnielinski correlation for heat transfer, which is considered the best available correlation for forced convection over a wide range of Reynolds and Prandtl numbers. However, even in the forced convection regime, the effect of variation of the fluid properties was found to be significant. © 2007 Elsevier B.V. All rights reserved.
1. Introduction The gas-cooled fast reactor (GFR) concept for Generation IV has generated considerable interest and is under development in the U.S., France, and Japan. One of the key candidates is a blockcore configuration first proposed by MIT (Hejzlar, 2001), which was shown in (Williams et al., 2004), to have the potential to operate in deteriorated turbulent heat transfer (DTHT) regime or in the transition between the DTHT and normal forced turbulent or laminar convection regime during post-loss-of-coolant accident (LOCA) conditions. The definition and the physical phenomena of the DTHT regime will be explained in the next
∗
Corresponding author. Tel.: +1 617 253 6188. E-mail addresses: [email protected] (J.I. Lee), [email protected] (P. Hejzlar), [email protected] (P. Saha), [email protected] (M.S. Kazimi). 0029-5493/$ – see front matter © 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.nucengdes.2007.01.003
section. This is contrary to most industrial applications where operation is in a well-defined and well-known turbulent forced convection regime. As a result, an important need emerged to develop heat transfer correlations that make possible rigorous and accurate predictions of decay heat removal (DHR) during post-LOCA in these regimes. Therefore, a collaborative INL/MIT project was initiated, to acquire experimental data (to the extent possible) for heat transfer correlations for all regimes including the transition regions in such a manner that there would be no discontinuities at the boundaries. The paper is organized as follows: first, a short literature review on the DTHT regimes will be provided followed by a simple analysis that shows the GFR DHR system has a strong potential to operate in this regime. Next section will give a simple comparison between gas, liquid and super-critical fluid to show that existing experimental data of liquid and super-critical fluid in the DTHT regime may not be appropriate for understanding gas DTHT. An overview of the loop construction with a detailed
1034
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Nomenclature A Bo* cp cts D f g G Gr* h k kts Kν Nu Pr Pw Re qb+ qw Q t T U x
area buoyancy parameter = Gr* /Re3.425 Pr0.8 specific heat at constant pressure (fluid) specific heat of test section material pipe diameter friction factor gravitational acceleration mass flux D4 /kν2 Grashof number with heat flux = gβqw heat transfer coefficient thermal conductivity of gas thermal conductivity of test section material acceleration parameter = (ν/Ub2 )(dUb /dx) ≈ + (4qb /Re) Nusselt number = hD/k Prandtl number = ν/α wetted perimeter Reynolds number =Ub D/ν /Gc T non-dimensional heat flux = qw p b wall heat flux volumetric flow rate time temperature flow velocity flow direction
Greek letters Λ volume α thermal diffusivity = k/ρcp β thermal expansion coefficient = −(1/ρ)(∂ρ/∂T)p μ viscosity ν kinematic viscosity ρ density Subscripts b bulk i inner ts test section w wall
description of each major loop component will be presented. Finally, the forced circulation data obtained from the loop will be presented and compared to the established correlations to validate the loop performance. 2. Review and decomposition of the deteriorated heat transfer regime controlling parameters In McEligot and Jackson (2004), the conditions of onset of DTHT regime based on three different physical effects: (1) flow acceleration, (2) buoyancy effect and (3) thermo-physical properties variation are well summarized. The reason it is called deteriorated turbulent heat transfer regime is because normal turbulent convective heat transfer is hindered by these phenomena
and the fluid’s heat transfer capability drops significantly in this region. This regime was first found in the super-critical fluid heat transfer and, later on, it was also found in general gas and liquid heat transfer systems, where the ratio of heating to flow rate is high. Even though the onset of this regime due to acceleration and buoyancy effects is well defined, the heat transfer coefficient and the friction factor for this regime have not gained universal agreement among the investigators. In contrast, thermo-physical property variation has a wide range of experimental data and is well established as reviewed in McEligot and Jackson (2004). 2.1. Acceleration effect The acceleration effect is also known as a “laminarization due to a favorable pressure gradient”. The decrease in turbulence transport occurs whenever a flow directional acceleration exceeds a certain value. The stream wise acceleration can be quantified by using the definition of acceleration parameter, Kν = (ν/Ub2 )(dUb /dx). In a heated system the acceleration occurs from gas expansion. This effect is different from the buoyancy effect, since the effect is not related to the gravitational force. When a flow condition satisfies the onset of laminarization due to acceleration, the heat transfer coefficient of the fluid quickly drops to the laminar heat transfer coefficient value, even though the Reynolds number is well above the critical Reynolds number. Applying energy balance and continuity equation with perfect gas and constant cross section assumptions to the original definition of Kν , an alternative definition of an acceleration parameter can be obtained for the heated case in the form of a non-dimensional number Kν ≈ 4qb+ /Re (McEligot et al., 1969). Based on McEligot and Jackson (2004), we can normalize the non-dimensional number that represents acceleration effect with the onset criterion (NKν = Kν /3 × 10−6 ). Thus, whenever the normalized value exceeds one, it indicates that the heat transfer regime is changing from normal forced convection regime to the DTHT regime due to the acceleration effect. It is also possible to separate the normalized acceleration parameter into two groups. One is the controlled group1 (CKν ), which includes geometry, heat flux and the volumetric flow rate. The other is the fluid properties group (PKν ), which is dependent on the system conditions such as pressure and operation temperature, and is not directly controlled. It should be noted that all the non-dimensional groups presented in this paper are evaluated at the fluid bulk temperature. Later in the paper, analysis will be focused on the properties group to determine which system pressure and operating temperature will likely force the fluid flow into the DTHT regime: 4qb+ Kν ≈ 3 × 10−6 3 × 10−6 Re qw APw μβ = = CKν PKν Q2 3 × 10−6 cp ρ2
NKν =
(1)
1 This group was designated “controlled”, because it consists of parameter that can be controlled by designer.
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
1035
2.2. Buoyancy effect The buoyancy effect originates from the density gradient due to heating, and reduces the turbulent heat transport when the flow direction is the same as the buoyancy force. Hall and Jackson (1969) explained the decrease in turbulence by the shear stress re-distribution in the flow. On the other hand Petukhov and Polyakov (1988) reported that the phenomenon is governed by two terms that are competing in the turbulent energy equation. One is the velocity gradient and the other is the turbulence work that needs to be provided to overcome a stable density gradient. Nevertheless, the criterion for onset of buoyancy induced DTHT is similar in both cases. In this paper we will follow definition of McEligot and Jackson (2004), which is an extended version of Hall and Jackson (1969) discussion. Similar to the acceleration effect, the buoyancy group will be normalized with the onset criterion (NBo* = Bo* /6 × 10−7 ). The buoyancy effect can also be viewed as a product of a control group (CBo* ) and a properties group (PBo* ) (Eq. (2)). The buoyancy properties group will be utilized as a measure of how the system pressure and operating temperature combination shifts the system to operate in the DTHT regime: Bo∗ Gr ∗ = 6 × 10−7 6 × 10−7 Re3.425 Pr 0.8 A4 gμ0.625 β qw = Q3.425 Pw 0.575 2.704 × 10−7 cp 0.8 k0.2 ρ1.425
NBo∗ =
= CBo∗ PBo∗
(2)
2.3. Thermo-physical properties variation When strong heating is provided to a fluid, thermo-physical properties vary across the flow area due to the radial temperature gradient of the fluid. Therefore, the heat transfer coefficient based on the assumption of constant fluid properties does not hold anymore and needs to be modified. Since liquid and gas properties changes with the temperature are different, the modification function is also different. However, thermo-physical property variation does not have a strong impact on the heat transfer coefficient like the acceleration effect or buoyancy effect. This effect has to be considered seriously if the ratio between wall temperature and bulk temperature exceeds 1.1, which causes 5% decrease in the Nusselt number McEligot and Jackson (2004).
Fig. 1. Acceleration properties group.
in the analysis, since it is relatively well known compared to the other two effects. It is clear from the figures that as temperature increases and pressure decreases, helium and carbon dioxide property groups increase. This indicates that the normal gas turbulent heat transfer has stronger tendency to fall into DTHT regime when the system depressurizes and the operating temperature increases. This situation occurs in the GFR system during the loss-ofcoolant accident. Therefore, both acceleration and buoyancy DTHT criteria need to be checked when designing the GFR DHR system. The figures also show that for the same controlled variable group, helium is more susceptible to the DTHT regime than carbon dioxide. Thus, more attention is needed when designing a helium-cooled system. Fig. 3 shows one of the possible operating ranges calculated with LOCA-COLA, an in-house code (Williams et al., 2004), of a GFR DHR system (Cochran et al., 2004). The Reynolds number is normalized to 4000 where fully developed adiabatic turbulent flow regime starts, and normalized DTHT parameters are NKν and NBo* . Curves in Fig. 3 show the axial variation of NKν and NBo* for given system pressure and heated diameter of the hot channel with helium and carbon dioxide as an operating fluid.
3. Simple analysis of GFR DHR system design Based on Eqs. (1) and (2) we can perform a simple parametric study with the gas properties. Since, the GFR design is considering helium and carbon dioxide as candidates for the coolant, both gases are evaluated. Figs. 1 and 2 show both helium and carbon dioxide acceleration properties groups and buoyancy properties groups for different pressure and temperature, respectively. It is noted that thermo-physical property variation effect is neglected
Fig. 2. Buoyancy properties group.
1036
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 3. GFR DHR operation range.
It can be observed that there is a possibility for GFR DHR system to operate not only in the DTHT regime but also in the laminar to turbulent transition regime. The laminar to turbulent transition criteria are not fully understood even in the adiabatic flow situation and a complete set of study for heated flow is rarely found (Lee, 2005). 4. DTHT regime for different fluids In this section, behavior of the buoyancy and acceleration parameters during the experimental runs will be discussed. Since the MIT experimental facility is a single flow channel, which will be further described in details in the next section, Eqs. (1) and (2) need to be rearranged to fit the single channel analysis. This is because along the channel the mass flow rate is constant not the volumetric flow rate. Eqs. (3) and (4) are rearranged forms of Eqs. (1) and (2): μβ qw APw NKν = = Cm˙ Kν Pm˙ Kν (3) m ˙2 3 × 10−6 cp ∗
NBo =
A4 qw m ˙ 3.425 Pw 0.575
= Cm˙ Bo∗ Pm˙ Bo∗
ρ2 gμ0.625 β 2.704 × 10−7 cp 0.8 k0.2
Fig. 4. Nitrogen acceleration properties group.
(4)
Figs. 4 and 5 plot the acceleration properties group and buoyancy properties group behavior of nitrogen at three different pressures with varying temperature. Since most of the gas thermal properties behave similarly with the temperature, nitrogen is chosen as an example. First observation from the comparison of Fig. 1 against Fig. 4 and Fig. 2 against Fig. 5 is the trend reversal. This is because the density term in the properties group moved to the controlled group side (since m ˙ = ρQ). Since the bulk temperature of the fluid will increase as it flows toward the downstream due to the heating, it can be predicted from the figures that both the
Fig. 5. Nitrogen buoyancy properties group.
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 6. Liquid water acceleration properties group.
buoyancy parameter and the acceleration parameter of the gas flow will decrease along the downstream. Therefore, it can be concluded that the maximum values of the buoyancy and acceleration numbers for the gas flow will be at the channel inlet when the control group is fixed as a constant. The control group will be approximately a constant for the MIT facility along the channel, since the test section has circular tube shape with a constant diameter and is designed to operate with uniform heat flux boundary condition. Another interesting observation is that the buoyancy parameter properties group of nitrogen increases with the pressure while the acceleration parameter properties group is pressureindependent. This indicates that by increasing the operating pressure of the facility, the buoyancy effect will increase relative to the acceleration effect. Various correlations were developed from the experimental data collected with working fluid of water and super-critical phase fluid (see Celata et al., 1998; Jackson et al., 1989; Kakac¸ et al., 1987). Thus, it is of interest to check how these fluids behave compared to gases. Figs. 6–9 show the properties group trends for liquid water and super-critical CO2 . If we compare the acceleration number properties group of water and nitrogen (Figs. 4 and 6), their behavior is also different. The nitrogen case shows steady decrease with increasing temperature, but water shows parabolic behavior with temperature. However, in contrast to the buoyancy parameter (Figs. 5 and 7), the order of magnitude of the acceleration property group is the same for both water and nitrogen cases. From the liquid water buoyancy properties group behavior (Fig. 7), one can see that the trend is totally different from that of the gas. The liquid water buoyancy parameter will increase or stay the same downstream the channel, since the properties group increases with the fluid bulk temperature. This is because water is incompressible fluid and the thermal properties variation due to the temperature change in the low-pressure range is insignificant compared to gas. One can also conclude from the figures that the order of magnitude for the water buoyancy properties group is four times higher than that of the low-pressure gas.
Fig. 7. Liquid water buoyancy properties group.
Fig. 8. Super-critical CO2 acceleration properties group.
Fig. 9. Super-critical CO2 buoyancy properties group.
1037
1038
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Therefore, it can be predicted that water experiments can reach buoyancy induced DTHT regime relatively easier than gas experiments, if the control group is the same. In addition, the channelwise behavior due to the buoyancy force will also be different. When we compare the super-critical carbon dioxide to nitrogen, the situation is different from the comparison of water to nitrogen. The trend with temperature is the same for both fluids but the order of magnitude change with temperature is different. This is because of a dramatic jump in both property groups for the super-critical CO2 near the critical point (Figs. 8 and 9), where the properties undergo a steep change. Therefore, it can be predicted that an experiment using a super-critical fluid will enter the DHTH regime, induced by either buoyancy effect or acceleration effect, where the temperature is near the critical point. From the comparison of all fluids, the following observations can be made: 1. The channel behavior of the buoyancy number will be very different among liquid, gas and super-critical fluids. 2. For gases, both maximum buoyancy number and maximum acceleration number will occur at the inlet of the heated channel.
3. For water, the maximum buoyancy number will be at the outlet of the heated channel and the maximum acceleration number will depend on the bulk fluid temperature. 4. For the super-critical carbon dioxide, both the maximum buoyancy number and maximum acceleration number will be at the point where the bulk temperature is near the critical point. 5. By pressurizing the gas, the buoyancy effect can be made more pronounced while the acceleration effect is immune from the pressure change for the same control group. The simple analysis presented in this section shows that the gas DTHT and the liquid or super-critical fluid DTHT regime are going to have different characteristics along the single channel and this prohibits us from utilizing the liquid or super-critical fluid experimental data for the gas system analysis. Therefore, the general heat transfer coefficient and friction factor in the gas DTHT regime and all the transitions from normal turbulent flow to the gas DTHT regime are not well defined yet, experimental studies in these heat transfer regimes are necessary to develop reliable correlations for the design of GFR DHR system.
Fig. 10. Schematic diagram of experimental facility.
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
5. Description of the experimental facility Fig. 10 is a simplified diagram, with dimensions, of the experimental facility (Cochran et al., 2004). The facility can be viewed as three different systems, as indicated by boxes in Fig. 10. The first system is the main loop, where most of the measurements are made. The main loop can be further divided into three sections—the test section, the chimney section and the downcomer section. The second system involves the compressor compartment located at the bottom part of Fig. 10. This system is necessary to operate the facility in a high flow rate range and to initiate a natural circulation flow. The last system is the gas charging and the vacuum system. This system is on the top right side in Fig. 10. The vacuum system is needed to pump out the air from the main loop before the test gases such as pure nitrogen, carbon dioxide and helium are charged into the main loop through the gas charging system. Each system and sub-system will be described in more detail in the following subsections. 5.1. Main loop With the exception of the test section, the main loop is constructed of 1 in. Type 316 stainless steel tube. The test section and flow development region, ahead of the test section, are made of a 0.75 in. Type 316 tubing. The test section wall thickness is 0.065 in. and the inner diameter of the test section is 0.62 in. (15.7 mm). The 1 in. stainless steel tube has a wall thickness of 0.049 in. and an inner diameter of 0.902 in. (22.9 mm). Connecting the test section and flow development section to the larger tubing of the loop is a simple reducing union at the entrance and exit of the test section. Both reducing sections change the tube sizing from 1 to 0.75 in. All the tube connections are done with compression type stainless steel fittings purchased from SWAGELOK. The test section is 1.976 m long and the chimney section is 3.797 m long. Power taps at the top and bottom of the test section connect a dc power supply providing up to 1500 A, at 10 V to the test section. This power supply provides for either voltage or current control of the power applied to the loop. Remote control of the power supply by the test control software is provided by a DA converter installed in a HP 3852A Data Acquisition unit. Pressure taps for measuring differential pressure across the test section have been welded to the tubing adjacent to the power taps. As can be seen from the diagram (Fig. 10) the tap at the entrance to the test section is also used for measuring system pressure. Differential pressure taps have also been installed on the downcomer section. All pressure taps have been connected to solenoid valves for selecting measurements from either the test section or the downcomer section. Teflon tubing has been used to connect pressure tap output to a differential pressure transducer. This transducer is a MKS Type 120A high accuracy pressure transducer. The range of this pressure transducer can be selected from 1 to 25,000 mmHg (133.3 Pa–3.33 MPa) with resolution of 10−6 of the full-scale range. Accuracy is ±0.12% of reading. The unit has its own readout equipment. Teflon tubing has been used for providing electrical isolation of the differential pressure transducer. Solenoid valves are con-
1039
trolled with relays in the same HP3852. Output of the system pressure transducer is monitored by a National Instrument Data Acquisition (NI-DAQ) system. Two flow thermocouples have been installed adjacent to the flow development and test sections. The first is 5.188 in. below the heated (test) section with the other installed 4.750 in. above the outlet of the heated section. This layout was necessary since the test section is heated using a direct Joule heating method, and installing flow thermocouples in the test section would create non-uniform heating in the vicinity of the thermocouple. Moreover, it was impractical to install them near the power taps. The flow thermocouples are connected to the NI-DAQ. 5.2. Instrumentation Temperatures of the loop were to be measured by Type K thermocouples spot-welded directly to the tubing. Initial testing showed significant errors in temperature readings during operation of the dc power supply. Analysis of the readings determined the cause of this problem to be electrical noise from the power supply. The problem was resolved by installing ungrounded sheathed thermocouples on the test section. Ungrounded thermocouples provide isolation from the electrical noise. Sheathed thermocouples required a different method of attachment to the test section. Each thermocouple was brazed to a small piece of shim stock, which was then spot-welded to the test section tubing. These sheathed thermocouples are used to monitor loop temperatures while dc heating current is applied to the loop. Due to concerns about slight errors with the sheathed thermocouples, temperatures used for determining heat transfer parameters are measured with thermocouples spot-welded directly to the tubing. Readings are taken during brief (2 s) interruptions of dc heating current. Sheathed thermocouples are connected to the NI-DAQ and the welded on thermocouples are connected to the HP3852. Twenty welded thermocouples and four sheathed thermocouples are installed in the test section and nine welded thermocouples and three sheathed thermocouples are installed in the chimney section (see Fig. 11). All the thermocouples are K-type and manufactured from OMEGA. To prove that the test section wall temperature does not drop appreciably during heating current interruptions, a simple idealized problem was solved. Eq. (5) is the solution to a lumped parameter model of the tube wall, selected to illustrate temperature response of the wall with time: hAi T (t) − Tb = exp − t (5) T (0) − Tb ρts cts ts Using Eq. (5), the time that is required for a test section temperature to drop 1 K from the steady state temperature, given the heat transfer coefficient of a gas to be 100 W/m2 K, and the temperature difference between the bulk gas and test section wall, at steady state, to be 200 K (conservatively high values for this facility), can be calculated to be 227 s. Therefore, even though the problem setup is simplified, it confirms that wall temperature changes during the short 2-s power-off period are negligible. This is mainly because stainless steel test section total
1040
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 11. Instrumentations and components in test section and chimney section.
heat capacity is much higher than the gas heat transfer capability. This was further confirmed during the experimental runs, where no visible temperature drops could be observed during the offpower period when the data acquisition system was collecting the data. The test section and the chimney section are under a layer of insulation. The insulation material is mineral wool. The test section insulation is equipped with thermocouples installed at six axially different positions and three azimuthally different positions (total 18) to measure the thermal resistance of the insulation and the heat losses during the experiment. On the chimney section, thermocouples are installed on the insulation at nine axially different positions and three azimuthally different positions, but at some places only one azimuthal thermocouple is installed (total 21). These insulation thermocouples will be used for not only measuring the heat losses but also to control the guard heater operation.
Fig. 11 shows the instrumentation and the components that form the test section and the chimney section. All the dimensions of the distances between thermocouples and their positions are indicated. The downcomer section includes a heat exchanger, pressure taps and electrical break. The gas charging and vacuum systems are connected to the downcomer section. One inch stainless steel tube is used for the downcomer section and the same compression type stainless steel fittings are used for the riser section. The length of the heat exchanger is 2 m and it is helically twisted with 30 cm coil diameter (6 and 2/3 turns). The building cooling water provides relatively constant temperature water (approximately 283 K). Two flow thermocouples (shown in Fig. 12) are installed near the heat exchanger; one at the inlet and the other at the outlet of heat exchanger. The two flow thermocouples are connected to the NI-DAQ. Building cooling water
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
1041
Fig. 12. Instrumentations and components in downcomer section.
temperature is also monitored by two thermocouples welded on to the heat exchanger. These thermocouples are also connected to the NI-DAQ. Pressure taps are installed through the connections on the downcomer, as indicated in Fig. 12. The pressure taps are controlled with the solenoid valves described in the previous section. The differential pressure transducer is very sensitive and will be used together with the hot-wire in the downcomer for the determination of flow rate in very low Reynolds number flow regimes. An electrical break is essential to eliminate a parallel electrical path around the test section. This insulating break was manufactured and installed in the downcomer section. Since pressure drop in the downcomer section is to be measured, great care was taken to achieve maximum smoothness of the inner surface of the break and smooth transition between the break and the tube. Four welded thermocouples and two sheathed thermocouples were installed on the downcomer section with the same
reasoning as explained in the previous section. The sheathed thermocouples are wired to the NI-DAQ and welded on thermocouples are connected to the HP3852. The downcomer section is also equipped with the same insulation material used on the test section. The insulation is needed since it is desirable to maintain a constant gas temperature along the downcomer for reliable pressure drop measurements and to prevent personnel contact with the downcomer tube wall in case that downcomer temperature exceeds 40 ◦ C in some runs. 5.3. Blower section and gas charging system A blower system is required for forced circulation at high Reynolds number testing and to initiate natural circulation. A three cylinders compressor unit supplies the input necessary to circulate the gas around the loop. Gas from the downcomer section of the loop is routed to an accumulator tank at the end of the downcomer section. The accumulator provides smooth-
1042
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 13. Instrumentations and components in compressor.
ing of the suction flow to the compressor. The compressor has been re-plumbed to connect all cylinders in parallel. After passing through the cylinders the gas is fed into a second larger accumulator tank to smooth the outgoing flow to the loop. No compression of the gas is done as it passes through the unit during operation. When the loop is vacuumed before the new gas is charged, the compressor unit is isolated from the remainder of the loop, by six valves installed in the inlet/outlet lines immediately adjacent to the compressor cylinders, to prevent the lubrication oil in the compressor from seeping out and contaminating the test section (Fig. 13). Three additional ball valves installed in the loop, near the compressor unit provide the isolation of the compressor from the loop for natural circulation testing. A variable frequency drive unit (VFD) provides speed control of the compressor unit. A relief valve is provided at the outlet to the compressor to prevent possible over-pressure event. After the tank outlet, a rotameter is installed to measure the flow rate in the experiments. The vacuum and gas charging systems control the gas purity in the main loop and the gas pressure. The vacuum system will first remove the air from the main loop except for a small volume in the compressor by closing six valves. After the loop is vacuumed, gas charging system is connected to the main loop and provides the test gases such as nitrogen, carbon dioxide and helium. The gas charging system is controlled with a pressure regulator. During the experiment, the gas charging system is going to be connected to the system to maintain a constant pressure and make up for any minor leakage in the main loop. 6. Forced convection results Three different experimental conditions were tested as a first check of the facility performance. These conditions were selected to perform experiments at normal forced convection turbulent heat transfer regime, with a slight effect from the thermo-physical properties variation of the fluid. Table 1 shows
the conditions for different experimental sets. All the tests were performed with nitrogen gas. For the run designation, the first three digits indicate the Reynolds number (for example 067 stands for Re at the inlet of about 6700) and last two digits indicate the qb+ (e.g., 13 stands for qb+ = 0.0013). qb+ is the ratio of wall heat flux to inlet flow enthalpy flux, which is a key indicator for heating to flow rate ratio. The three experimental sets cover the range of local Reynolds numbers from 4700 to 12,800, which are usually considered within the turbulent flow regime in the pipe flow. The correlations selected for comparison are Dittus–Boelter, and Gnielinski. Dittus–Boelter correlation was selected because it is the most widely used correlation for forced convection turbulent heat transfer. However, the Dittus–Boelter equation used in the comparison differs from the widely known form, which uses the value of 0.023 as a leading constant instead of 0.021 presented in Eq. (6). This is because the widely known form is more appropriate for liquids. Gnielinski’s correlation (Eq. (7)) is a newer correlation, which exhibits better accuracy than Dittus–Boelter and is valid for a larger range of Reynolds and Prandtl numbers.
Table 1 Experimental conditions
Mean power (W) Mass flow rate (kg/s) Mean pressure (MPa) Inlet Reynolds number Outlet Reynolds number qb+ inlet Inlet bulk temperature (K) Outlet bulk temperature (K) Minimum wall temperature (K) Maximum wall temperature (K) Wall to bulk temperature ratio
06713 run
07915 run
12823 run
397.5 0.0015 0.1380 6738 4624 0.0013 298.6 494.4 330.5 584.9 1.09–1.30
515.0 0.0018 0.1533 7992 5304 0.0015 299.0 521.0 337.2 624.1 1.11–1.35
1187.1 0.0028 0.2279 12871 7438 0.0023 299.4 622.1 363.9 783.3 1.17–1.52
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Modified Dittus–Boelter correlation (McEligot et al., 1966) Nu = 0.021Re0.8 Pr 0.4
(6)
Gnielinski correlation (Gnielinski, 1976) Nu =
(f/8)(Re − 1000)Pr √ 1 + 12.7 f/8(Pr 2/3 − 1)
(7)
The friction factor that is used in the Gnielinski correlation is given in Eq. (8). Colebrook correlation (valid for 4000 < Re < 107 ) (Kaka¸c et al., 1987) 2 16 f = (8) 1.5635 ln(Re/7) Before presenting the results, it is important to note that the Gnielinski correlation given in Kakac¸ et al. (1987), omitted a correction factor that accounts for the effect of property variation in the film for cases where the temperature difference between the wall and bulk is large. This correction factor is defined in the original Gnielinski’s publication (1976) and it is given in Eq. (9). Similar correction factor exists for the Dittus–Boelter correlation—see Eq. (10) (McEligot et al., 1966): Tb 0.45 Correction factor (CF) = , Twi NuGnielinski modified = NuGnielinski Eq. (7) × CF CF =
Tb Twi
Fig. 14. Nusselt number ratio vs. axial position (06713 run).
(9)
0.5 ,
NuDittus–Boelter modified = NuDittus–Boelter Eq. (6) × CF
(10)
The correction factor is necessary for the gas since the wall to bulk temperature ratio range is 1.1–1.5 in all experimental runs affecting significantly properties in the film, and thus heat transfer coefficient. Table 1 shows the range of wall to bulk temperature ratio for three cases. Figs. 14–16 plot the ratio between the experimental Nusselt numbers and the Nusselt numbers calculated from correlations. Since the correlations shown above are for fully developed region, Figs. 14–16 are magnified with the thermal entrance region data points left out for the reader’s convenience. In the fully developed region, experimental runs 06713 and 07915 match the correlations better than the 12823 run. Experiments 06713 and 07915 have nearly the same ratio between 0.95 and 1.0 with respect to the Gnielinski correlation with the correction factor in the fully developed region. Since the 06173 and 07915 runs have smaller errors than the 12823 run in terms of energy, the larger discrepancy of 12823 run with the correlation can be explained. Both correlations, the Dittus–Boelter and the Gnielinski, show good agreement with the experimental data after the correction factor is incorporated. This indicates that the property variation due to the large temperature gradient in the fluid has
Fig. 15. Nusselt number ratio vs. axial position (07915 run).
Fig. 16. Nusselt number ratio vs. axial position (12823 run).
1043
1044
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
Fig. 17. Nusselt number comparison with error bars (06713 run).
Fig. 19. Nusselt number comparison with error bars (12823 run).
a significant effect on the heat transfer even for the forced convection regime. It can be also concluded that both the Dittus–Boelter and Gnielinski correlations predict the experimental data at low Reynolds number regime (4000–13,000) satisfactorily. The last point of experimental Nusselt number in each plot rises above 1.0, potentially due to the under-estimation of the axial heat losses at the top of the heated section. This may require installation of a thermocouple right after the heated section to reduce the axial heat loss at the outlet. The error analysis will be discussed further in the following section. Figs. 17–19 show the experimental Nusselt number and the Nusselt number from the Gnielinski correlation that includes the correction factor. The Gnielinski correlation has an error of ±20% (Gnielinski, 1976) and it is represented in the figures representing error bars with the uncertainty of the experimental Nusselt number.
Figs. 17–19 show the 06713 and 07915 runs to have very good agreement with the Gnielinski correlation with the property temperature dependence factor. The 12823 run gives slightly higher Nusselt number than the correlation. This was due to the larger error in estimating the thermal resistance of the insulation data for reduction of the 12823 run, since the wall temperature of the 12823 run is significantly above the temperature at which the thermal resistance was measured. Fig. 20 compares the experimental Nusselt number versus the Reynolds number against the Gnielinski correlation with the correction factor for all three runs. The entrance and the end points, where the uncertainty is the highest were excluded. It is expected that further reduction in the uncertainty can be achieved after: (1) a rotameter is replaced with a turbine meter, (2) guard heaters are installed on top of the primary insulation
Fig. 18. Nusselt number comparison with error bars (07915 run).
Fig. 20. Nusselt number comparison (Re vs. Nu).
J.I. Lee et al. / Nuclear Engineering and Design 237 (2007) 1033–1045
and (3) more thermocouples on the insulation to determine the heat loss with higher accuracy. 7. Summary and conclusions The work presented in this paper can be summarized as follows: 1. It was shown that gas-cooled fast reactor (GFR) decay heat removal (DHR) system is likely to operate after LOCA in the deteriorated turbulent heat transfer (DTHT) regime. Since the heat transfer coefficient and friction factor correlations are not well established in this regime, new experiments were necessary to develop reliable correlations. 2. The non-dimensional groups reflecting the buoyancy and acceleration effects have been decomposed into two parts. One representing the design controlled parameters, such as geometry and heating rate. The other representing the fluid properties, which are not directly controlled in the design of a facility. 3. An experimental facility was constructed to investigate the DTHT regime. However, due to the power supply noise, a different temperature measurement procedure was developed, which combines welded thermocouples and sheathed thermocouples. It was shown by a simple analysis that the new measurement procedure gives accurate data, which was also confirmed by experimental runs. 4. Three forced circulation experiments were performed and showed very good agreement with the Gnielinski correlation that incorporates the effect of temperature dependent properties in the film, confirming that there is no significant defect on measurement system and the experimental loop. However, a few minor adjustments were identified to be done, in order to reduce measurement uncertainty. 5. The initial results already confirm that the Nusselt number correlation in forced convection has to incorporate temperature dependence of properties in the film to correctly predict gas heat transfer in the GFR channels since the wall-to-bulk temperature ratios exceed the value of 1.1. This is important as it leads to reduced heat transfer coefficient and thus to higher wall temperatures.
1045
Acknowledgments The authors would like to acknowledge the financial support of the Idaho National Laboratory. They thank P. Cochran of MIT for providing preliminary analysis of the experimental apparatus, P. Stahle of MIT for his valuable effort in the construction of the facility and Dr. Glenn E. McCreery and Dr. D.M. McEligot of INL for valuable comments on the facility design. References Celata, G.P., D’Annibale, F., Chiaradia, A., Cumo, M., 1998. Upflow turbulent mixed convection heat transfer in vertical pipes. Int. J. Heat Mass Transfer 41, 4037–4054. Cochran, P., Saha, P., Hejzlar, P., McEligot, D.M., McCreery, G.E., Schultz, R.R., 2004. Scaling Analysis and Selection of Test Facility for Fundamental Thermal-Hydraulic Studies related to Advanced Gas-cooled Reactor, INL/EXT-05-00158. Idaho National Laboratory and Massachusetts Institute of Technology, Department of Nuclear Science and Engineering. Gnielinski, V., 1976. New equations for heat and mass transfer in turbulent pipe and channel flow. Int. Chem. Eng. 16 (2), 359–387. Hall, W.B., Jackson, J.D., 1969. Laminarization of a turbulent pipe flow by buoyancy forces, ASME 69-HT-55. Hejzlar, P., 2001. A modular, gas turbine fast reactor concept (MFGR-GT). In: Proceedings of Transactions of the American Nuclear Society, vol. 84, Milwaukee, Wisconsin, June 17–21. Jackson, J.D., Cotton, M.A., Axcell, B.P., 1989. Studies of mixed convection in vertical tubes. Int. J. Heat Fluid Flow 10, 2–15. Kakac¸, S., Shah, R.K., Aung, W., 1987. Handbook of Single-phase Convective Heat Transfer. John Wiley & Sons (Chapter 15). Lee, J.I., 2005. The Flow Structure under Mixed Convection in a Uniformly Heated Vertical Pipe. Master Thesis, Nuclear Science and Engineering Department, Massachusetts Institute of Technology. McEligot, D.M., Coon, C.W., Perkins, H.C., 1969. Relaminarization in tubes. Int. J. Heat Mass Transfer 13, 431–433 (Shorter Communications). McEligot, D.M., Ormand, L.W., Perkins, H.C., 1966. Internal low Reynolds number turbulent and transitional gas flow with heat transfer. J. Heat Transfer 88, 239–245. McEligot, D.M., Jackson, J.D., 2004. Deterioration criteria for convective heat transfer in gas flow through non-circular ducts. Nucl. Eng. Design 232, 327–333 (Short Communication). Petukhov, B.S., Polyakov, A.F., 1988. Heat Transfer in Turbulent Mixed Convection. Hemisphere Publishing Corporation. Williams, W.C., Hejzlar, P., Driscoll, M.J., 2004. Decay heat removal from GFR core by natural convection. In: International Congress on Advances in Nuclear Power Plants ICAPP ’04, Paper 4166, Pittsburgh, USA, June 13–17.