Annals of Nuclear Energy 92 (2016) 175–185
Contents lists available at ScienceDirect
Annals of Nuclear Energy journal homepage: www.elsevier.com/locate/anucene
CFD aided approach to design printed circuit heat exchangers for supercritical CO2 Brayton cycle application Seong Gu Kim 1, Youho Lee 1, Yoonhan Ahn 1, Jeong Ik Lee ⇑ Department of Nuclear and Quantum Engineering, Korea Advanced Institute of Science and Technology, 373-1 Guseong-dong, Yuseong-gu, Daejeon 305-701, Republic of Korea
a r t i c l e
i n f o
Article history: Received 17 June 2015 Received in revised form 15 October 2015 Accepted 10 January 2016
Keywords: Supercritical CO2 S-CO2 Brayton cycle Printed circuit heat exchanger S-CO2 PCHE CFD Supercritical fluids
a b s t r a c t While most conventional PCHE designs for working fluid of supercritical CO2 require an extension of valid Reynolds number limits of experimentally obtained correlations, Computational Fluid Dynamics (CFD) code ANSYS CFX was used to explore validity of existing correlations beyond their tested Reynolds number ranges. For heat transfer coefficient correlations, an appropriate piece-wising with Ishizuka’s and Hesselgreaves’s correlation is found to enable an extension of Reynolds numbers. For friction factors, no single existing correlation is found to capture different temperature and angular dependencies for a wide Reynolds number range. Based on the comparison of CFD results with the experimentally obtained correlations, a new CFD-aided correlation covering an extended range of Reynolds number 2000–58,000 for Nusselt number and friction factor is proposed to facilitate PCHE designs for the supercritical CO2 Brayton cycle application. Ó 2016 Elsevier Ltd. All rights reserved.
1. Introduction The supercritical carbon dioxide (S-CO2) Brayton cycle is a promising power cycle for it has high thermal efficiency in relatively low temperature ranges. The S-CO2 cycle takes advantages of both the Rankine cycle and the Brayton cycle. The Rankine cycle achieves high efficiency thanks to low pumping power and the Brayton cycle achieves high efficiency primarily because of high turbine inlet temperature. The S-CO2 cycle employs both advantages because (i) it requires low CO2 compressing work near the critical point, and (ii) it has a higher turbine inlet temperature than that of the steam Rankine cycle (Ahn and Lee, 2014). Many researchers have proposed various applications of the supercritical CO2 Brayton cycle. To design an efficient power conversion system for a sodium-cooled fast reactor (SFR), several gas mixtures (CO2 mixed with H2, N2, O2, He, Ar, Kr, Xe) were evaluated as working fluids for the supercritical Brayton cycle (Jeong et al., 2011). Preliminary design of S-CO2 Brayton cycle and radial type turbomachinery coupled to System-integrated Modular Advanced Reactor (SMART) were investigated in the previous studies (Yoon et al., 2012; Lee et al., 2014). Various S-CO2 cycles for a power conversion system coupled to Molten Carbonate Fuel Cell (MCFC) hybrid
⇑ Corresponding author. Tel.: +82 42 350 3829; fax: +82 42 350 3810. E-mail addresses:
[email protected] (S.G. Kim),
[email protected] (Y. Lee),
[email protected] (Y. Ahn),
[email protected] (J.I. Lee). 1 Tel.: +82 42 350 3829; fax: +82 42 350 3810. http://dx.doi.org/10.1016/j.anucene.2016.01.019 0306-4549/Ó 2016 Elsevier Ltd. All rights reserved.
system were studied. Several S-CO2 cycles were examined in terms of thermal efficiency, net electricity, and the volume for MCFC hybrid system. Previous studies demonstrated that the S-CO2 Brayton cycle with a simple recuperated layout can be an efficient bottoming cycle for the MCFC hybrid system than the gas turbine power generation system (Bae et al., 2014). The small modular high temperature gas-cooled reactor (SM-HTGR) system coupled to helium Brayton cycle and S-CO2 Brayton cycle was considered as a distributed power generation system (Bae et al., 2015). The S-CO2 cycle components such as turbomachineries and heat exchangers can be designed in more compact forms than those of the Rankine cycle (Yoon et al., 2012). The fluids with moderate heat transfer capacity, usually in either gaseous or supercritical state, are the typical working fluids for a power conversion cycle that employs PCHEs as intermediate heat exchangers. The large heat transferring area in a relatively small volume of PCHE is rooted in its manufacturing procedure: the PCHE is made by stacking multiple etched plates, followed by a diffusion bonding process (Lee and Lee, 2014). As a consequence of the manufacturing process, narrowly spaced semi-circular micro channels (on the order of 103 m), and serving as conduits for heat transferring fluids, are formed. The zigzag fluid channel shape provides large advantages by producing vortices at the channel corners which improves Nusselt number at the expense of an increase in friction factor (Pra et al., 2008). Heat transfer and pressure drop performances are critical factors for designing a heat exchanger. In experiments, it is difficult
176
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185
Nomenclature CP Nu Re Q _ m hout hin A Twall Tbulk h
specific heat capacity at constant pressure [kJ/kg-K] Nusselt number Reynolds number heat load [W] mass flow rate [kg/s] enthalpy at the outlet [kJ/kg] enthalpy at the inlet [kJ/kg] heat transfer area [m2] wall temperature [°C] bulk temperature [°C] local heat transfer coefficient [W/m2-K]
to measure and obtain local flow parameters in micro channels of PCHEs. Therefore, Computational Fluid Dynamics (CFD) analysis can be a useful approach to obtain local flow information in PCHE flow channels. In this paper, CFD analysis of PCHE with S-CO2 as a working fluid is discussed. The Nusselt number and pressure drop correlations proposed in previous experimental studies have a limited validity due to restricted experimental conditions. In the course of designing real PCHEs or evaluating system performance subject to unexpected flow transients, the developed correlations may be subject to extrapolating beyond the explored flow conditions and geometries. Yet, the validity of extrapolating developed correlations to experimentally unexplored flow conditions and geometries should be questioned. In this study, we computationally explore heat transfer and pressure drop behavior beyond the tested flow conditions and geometries, and assess the validity of existing correlations for extrapolation. Such a computational support for extrapolation of existing correlations can advance the design and analysis of PCHE. In this study, heat transfer performance and pressure drop of PCHE with wide Reynolds number ranges were investigated with CFD analysis and the results were compared to the existing correlations to examine if the existing correlations can be used in beyond the experimented Reynolds number ranges. Furthermore, with the support of CFD analysis, we propose a new set of correlation for the Nusselt number and friction factor in the extended Reynolds number range. 2. Previous studies on PCHE with CO2 as a working fluid In previous studies, empirical correlations were derived from combining experimental and numerical results. Ishizuka et al. (2005) conducted an experimental study to examine heat transfer and pressure drop performance of PCHE in Tokyo Institute of Technology (TIT) supercritical CO2 loop. In the study, empirical correlations for pressure loss coefficient, local heat transfer coefficient, and overall heat transfer coefficient were proposed. Since measurement systems were installed at the inlet and the outlet of the PCHE, local values of pressure and temperature inside the flow channels were found by assuming their linear distributions from the inlet to the outlet. With the assumed local values, the local heat transfer and pressure loss coefficients were obtained. The ranges of Reynolds number are 2400–6000 and 5000–13,000 in the hot side and the cold side, respectively. Nikitin et al. (2006) proposed an empirical heat transfer correlation for the PCHE in the TIT supercritical CO2 loop. In order to obtain accurate local heat transfer coefficients, the amount of heat loss in the PCHE was found by a 2-D CFD analysis. Ngo et al. (2007) conducted an experimental study for a micro-channel heat exchanger (MCHE) with S-shaped and zigzag fins. A new MCHE with S-shaped fin was developed
h Nu Nu f Pin Pk
qk
Vk Dh Lk f
channel averaged local heat transfer coefficient [W/m2-K] Nusselt number channel averaged Nusselt number friction factor pressure drop at the kth plane [Pa] density at the kth plane [kg/m3] velocity at the kth plane [m/s] hydraulic diameter [m] active length from the inlet [m] channel averaged friction factor
for a recuperator of S-CO2 cycle for nuclear applications by using the commercial CFD code, FLUENT (Ngo et al., 2007). The pressure drop correlation and the Nusselt number correlation for local convective heat transfer were proposed. In the study, the range of Reynolds number was 3500–22,000. Hesselgreaves (2001) proposed a correlation of wavy (corrugated or herringbone) fins taking into account the geometric factors – the wave pitch and width. This correlation is based on the experimental studies with wavy fins conducted by Kays and London (1984) and Oyakawa et al. (1989). Kim et al. (2008) performed a 3-D numerical study of the S-CO2 PCHE with zigzag fins, and they proposed an airfoil shaped fin which has much lower pressure drop without sacrificing heat transfer performance. The numerical result of the zigzag channel was compared to the experimental results from Ishizuka et al. and the CFD-results showed a good agreement with the experimental results for the outlet temperature and pressure. Lee and Kim (2014) performed a 3-D numerical analysis of PCHE to examine effectiveness, Nusselt number, and the Colburn j factor. In this study, the authors proposed a new design of PCHE with a different thickness and the reversed plate configuration to enhance heat transfer capability. The correlations developed in the discussed studies are summarized in Table 1. The empirical correlations from the experimental studies have limited Reynolds number ranges since the facilities were made in a small lab scale. The Colburn j factor correlations proposed by Hesselgreaves (2001) have a sufficiently high Reynolds number range up to 100,000. Yet, since the correlation was originally proposed for folded wavy fins between flat plates, a validity-check for the PCHE geometry is required.
3. In-house heat exchanger design code The heat exchanger design code developed by KAIST research team, KAIST-HXD is utilizing the correlation developed by Ngo et al. (2007). The code estimates pressure drop, heat transfer performance, and, size of a PCHE. The input parameters are pressure, temperature and mass flow rates at the inlet of both hot side and cold side. This code calculates effectiveness of the heat exchanger and outlet temperature of both sides, and returns whether the design achieves target performance. The volume and mass of the heat exchanger can be estimated with dimensions of channel variables: diameter, pitch, plate thickness, length, and number of channels. This code uses properties of CO2 from the NIST REFPROP 8.0 database (Lemmon et al., 2007). However, due to the restricted Reynolds number ranges of the correlation, velocity in the flow channel is set to have an upper limit. This may result in unrealistic PCHE dimensions with excessive number of channels to meet the flow rate limit set by the correlation. If we can use an extended
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185
177
Table 1 Empirical correlations developed in the previous studies of PCHE.
Ishizuka et al. (2005)
Nikitin et al. (2006)
Ngo et al. (2007)
Valid range
Heat transfer coefficient (HTC) or Nusselt number
Friction factor
Hot side 2400 6 Re 6 6000 Cold side 5000 6 Re 6 13000
Local HTC h ¼ 0:2104 Re þ 44:16 Overall HTC U ¼ 0:1106 Re þ 15:943
Hot side
Hot side 2800 6 Re 6 5800 Cold side 6200 6 Re 6 12100
Local HTC
Hot side
hhot ¼ 2:52 Re0:681
f h ¼ ð1:402 106 0:087 106 Þ Re þ ð0:04495 0:00038Þ Cold side
3500 6 Re 6 22000 0:75 6 Pr 6 2:2
f h ¼ 2:0 106 Re þ 0:0467 Cold side f c ¼ 2:0 106 Re þ 0:1023
hcold ¼ 5:49 Re0:625 Nu ¼ ð0:1696 0:0144Þ Re
f c ¼ ð1:545 106 0:099 106 Þ Re þ ð0:09318 0:00090Þ 0:6290:009
Pr
0:3170:014
Hot side f hot ¼ ð0:3390 0:0285Þ Re0:1580:009 Cold side f cold ¼ ð0:3749 0:1293Þ Re0:1540:036
Hesselgreaves (2001)
600 < Re < 3000 104 < Re < 105
600 < Re < 3000: 0:75 j ¼ 0:4 Re0:4 2b K 104 < Re < 105 : j ¼ 0:4 Re
0:36 2b 0:75
ðKÞ
600 < Re < 3000: f ¼ 5j 104 < Re < 105 : 1:5 f ¼ 4:8 Re0:36 2b K
where j ¼ Re Nu Pr1=3 b: channel width K: channel pitch
Table 2 Operating condition of a recuperator. Hot side inlet condition Cold side inlet condition Mass flow rate Channel diameter Minimum effectiveness Pressure drop
432.7 °C, 7.58 MPa 149.9 °C, 20.0 MPa 175.34 kg/s 1.9 mm 93% 1.5% of system pressure
correlation in higher Reynolds number range for such a case, the heat exchanger can be designed with more realistic dimensions. We compared differences in heat exchanger dimensions designed with the different correlations listed in Table 1; recuperators – each subject to the different correlations – were designed to meet the identical performance metrics under the same environment. The operating condition of the recuperator was determined from the cycle design of a reference supercritical CO2-cooled micro modular reactor (MMR) (Kim et al., 2013) (see Table 2). The operating conditions and cycle layout for MMR are described in Fig. 1. It is worth noting that since the MMR concept assumes containing the reactor and power conversion system in a module, the volume of heat exchanger is a key factor to the feasibility of the concept. The average Reynolds number of channel was set to be in the valid range of the correlation. The target effectiveness was set to be higher than 93%. The maximum pressure drop inside the heat exchanger was set to be 1.5% of the maximum system pressure (300 kPa). Design results of PCHE with the different correlations are shown in Table 3. In design result I, the designed heat exchanger has a pan-cake external shape with a large frontal area and a short length. Because the flow speed in the channel has a limit due to the low Reynolds number limit (Re < 6000), the number of channel was increased to match the total mass flow rate. Such a pan-cake PCHE shape is practically hard to be realized because of the difficulty in designing a large size PCHE header that properly distributes flow into individual channels. In addition, the volume of design I is 4.7 times bigger than that of PCHE design III. In design result II, the length of heat exchanger is 0.639 m and the frontal area is still large. Assuming the frontal area is square, the side of frontal section is longer than the length of heat exchanger. Design
Fig. 1. Schematic figure of simple recuperated S-CO2 Brayton cycle for MMR.
result III shows small frontal area, compact size and practical external dimensions. As shown in the design results, when the Reynolds number limit is low, the PCHE design requires an excessive number of channels with an unrealistically short length and flat external shape. From this result, it is inferable that realistic PCHE designs would run with Reynolds number higher than the highest valid Reynolds numbers of most of the available correlations – a key reason why the extension of Reynolds number are needed. We performed extensions of the existing correlations and compared them with results of CFD analysis. 4. PCHE analysis model 4.1. Analysis model for PCHE channel The geometry and configuration of CFD-simulated PCHE were made to be the same as those of the experimental study performed
178
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185
Table 3 Design results of PCHE recuperator by using KAIST-HXD code.
Reynolds number range Average Reynolds number Heat load Relative volume (volume divided by volume III) Pure mass Length Relative frontal area (frontal area divided by frontal area III) Pressure drop
Design result I with T. Ishizuka correlation
Design result II with T.L. Ngo correlation
Design result III with J.E. Hesselgreaves correlation
2400 < Re < 6000 Hot side: 5921 Cold side: 5171 67.1 MW 4.738 11,558 kg 0.414 m 8.280
3500 < Re < 22,000 Hot side: 21,542 Cold side: 18,850 67.1 MW 2.009 4901 kg 0.639 m 2.275
10,000 < Re < 100,000 Hot side: 56,336 Cold side: 39,380 67.1 MW 1.000 1597 kg 0.724 m 1.000
10.12 kPa (hot) 3.36 kPa (cold)
187.2 kPa (hot) 61.1 kPa (cold)
295.3 kPa (hot) 98.4 kPa (cold)
External shape
Fig. 2. Geometry and boundary condition of PCHE analysis model.
by Ishizuka et al. (2005). A channel configuration consisting of two hot channels and one middle cold channel with the vertical and horizontal periodic conditions was used as the reference PCHE channel unit as shown in Fig. 2. The original cross sectional shape of the PCHE specimen is illustrated in Figs. 3 and 5. Since one cold side plate is placed between two hot side plates, the reference PCHE channel unit is composed with two hot channels and one cold channel to represent the original configuration. The hot and cold channels are semi-circular, whose sizes are 1.9 mm and 1.8 mm in diameter, respectively. The structural material of the PCHE plates is AISI316, which has a thermal conductivity of 16.2 W/m-K. The PCHE facility in the experiment of Ishizuka et al. (2005) has different fin angles for the hot (32.5°) and cold side (40.0°), making actual length of the cold side 10% longer than that of the hot side. To facilitate periodic boundary conditions on the sides of the tested PCHE unit, the fin angle of analysis model was assumed to be the same for both the cold and hot channels at either 32.5° or 40.0° (see Fig. 4). Thus, two analysis models were generated and compared to the selected correlations. Detail dimensions of PCHE specimen is listed in Table 4.
Fig. 3. Cross-sectional shape of PCHE.
4.2. Grid generation Volume meshing elements were generated by using ANSYS ICEM CFD software. 18 sheets of prism layers are produced near the wall, and y+ of the first node adjacent the wall is close to 1.0. To confirm the grid independence of the PCHE analysis model, grid sensitivity study was carried out and results are shown in Table 5. As the number of nodes exceeds 1,204,052, channel averaged heat transfer coefficients at the both sides converged. Hence, the grid system with 1,204,052 nodes and 3,772,579 elements was used as the reference grid scheme for the CFD analysis. Tetrahedral, and hexahedral meshes were used. The typical meshing scheme used in CFD analysis is illustrated in Fig. 6. 4.3. Fluid property implementation Real gas property (RGP) format table was used for CO2 property implementation. The properties of CO2 are taken from the NIST REFPROP 8.0 fluid property database (Lemmon et al., 2007). The tables can be generated in prescribed pressure and temperature
179
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185
Fig. 4. Fin angle of PCHE specimen (Ishizuka et al., 2005) and analysis model 1 (32.5°) and analysis model 2 (40.0°).
Table 5 Result of grid sensitivity study. Number of nodes
Number of elements
Average heat transfer coefficient (W/m2-K)
Error (%)
303,998
850,164
678,175
2,036,925
1,204,052 (Reference grid) 2,283,447
3,772,579
1850 4016 1994 4042 2122 4253 2166 4237
14.6 5.2 7.9 4.6 2.0 +0.4 +0.0 +0.0
7,643,631
(hot) (cold) (hot) (cold) (hot) (cold) (hot) (cold)
Fig. 5. PCHE flow channel configuration and dimension (Ishizuka et al., 2005).
Table 4 Dimensions of PCHE (Ishizuka et al., 2005).
Plate material Plate thickness Horizontal pitch Wall width Channel width Channel depth Hydraulic diameter Heat transfer area Channel active length Channel length
Hot channel
Cold channel
SS316L 1.63 mm 4.50 mm 0.60 mm 1.90 mm 0.90 mm 1.15 mm 0.697 m2 1000 mm 896 mm
SS316L 1.63 mm 3.62 mm 0.70 mm 1.80 mm 0.90 mm 1.15 mm 0.356 m2 1100 mm 896 mm Fig. 6. Mesh of the PCHE analysis model.
ranges with an in-house MATLAB program. Table resolution can be prescribed by the user. The CFX solver utilizes the lookup table method to refer the CO2 properties from the RGP table and calculates output values by using bilinear interpolation. The RGP table, pressure ranging from 1 MPa to 20 MPa, and temperature ranging
from 290 K to 650 K, with resolution of 500 by 500 was used. Although the table resolution can be made up to the 2000 by 2000, operating condition of the PCHE recuperator is far from the critical point of CO2 under little property variations, so the tables are not required to have such a high resolution (Kim et al., 2014).
180
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185
Table 6 Property errors of RGP table in the hot side and the cold side. Enthalpy [kJ/kg] Hot side (3.196 MPa, 243.0 °C) CFX 703.26 NIST 702.59 Error (%) 0.09 Cold side (10.49 MPa, CFX NIST Error (%)
137.4 °C) 553.84 553.19 0.12
CP [J/kg-K]
Viscosity [Pa-s]
Thermal conductivity [W/m-K]
Density [kg/m3]
1065.2 1064.31 0.08
2.4929e05 2.4898e05 0.13
0.03577 0.03571 0.19
33.787 33.4152 1.11
1315.1 1313.9 0.09
2.2774e05 2.2730e05 0.19
0.03252 0.03244 0.25
163.05 162.67 0.23
To estimate the error between real properties obtained from the NIST REFPROP and calculated properties within the CFX code, a comparison was made at the mid-point of the hot and cold sides as summarized in Table 6. With the 500 by 500 table, density showed 1.1% of property error while the other properties showed error less than 0.3%. This degree of property error is considered acceptable for the purpose of CFD analysis. Hence, the table with 500 by 500 resolution was mainly used for the CO2 property implementation. 4.4. Problem setup A commercial code, ANSYS CFX 14.5 was used for the CFD analysis with k–x Shear Stress Transport (SST) turbulence model (ANSYS, 2011). The k–x SST is chosen because standard k–e model tends to underestimate pressure drop according to the previous numerical study performed by Val Abel et al. (2011). In the previous study (Val Abel et al., 2011), the results with k–e models estimated 30% less pressure drop than both an experimentally measured pressure drop and a prediction of k–x SST model. A similar underestimation for pressure drop was observed in the result of analysis model 1 with the standard k–e and RNG k–e turbulence models. The total energy equation was used for energy conservation, and conjugate heat transfer mode in the fluid and solid domain was employed in this analysis. 5% turbulent intensity and automatic calculated length scale were applied. The convergence criteria was met when the RMS residual of mass, momentum and energy is less than 104 and the imbalance of mass, momentum and energy is less than 0.5%. Boundary conditions were determined to be the same as the experimental conditions of Ishizuka et al. (2005). In this experiment, inlet conditions of the hot side and the cold side were 280 °C, 3.2 MPa and 108 °C, 10.5 MPa, respectively. However, the channel length of the analysis model used in this study is much shorter than that of the experiment. Therefore, the inlet temperature of cold side which is 54 mm far from the hot side inlet was referred from the previous CFD result with 896 mm channel length performed by Kim et al. (2008). In this numerical study, the calculated results of pressure and temperature at both outlets showed a good agreement with the experimental results. Temperature distributions at the hot and cold sides are shown in Fig. 7. Based on this result, temperature at the cold side inlet was set to be 221.8 °C. The outlet conditions were controlled to mass flow rate from 30 kg/h to 400 kg/h with an increment of 10–20 kg/h. Note that the simulated mass flow rate range includes that of the experiment in order to cover flow conditions beyond the valid range of correlation. 4.5. Post-processing For post-processing, numbers of planes were created in the fluid domain of each analysis model. All the fluid properties and flow variables were read from those planes, and the bulk and wall temperatures were calculated by using a CFX CEL expression. Heat load
Fig. 7. Result of temperature profile at the hot side and the cold side of PCHE conducted by Kim et al. (2008).
Q of the each channel was found from the enthalpy differences from the inlet to the outlet. Local heat transfer coefficient was obtained at each plane, and channel averaged local heat transfer coefficients were obtained by using the local values. After that, channel averaged Reynolds number and Prandtl number were used to calculate Nusselt number and friction factor values of existing correlations for a comparison with the CFD results. The correlation for the local heat transfer coefficient developed by Ishizuka et al. was converted to a correlation for Nusselt number by multiplying l , where l is the characteristic length of channel and k is the chank nel averaged value of thermal conductivity. The Colburn j factor correlation proposed by Hesselgreaves was converted to a Nusselt number correlation with geometric parameters of analysis model 1 (32.5° model) and 2 (40° model). The heat load was obtained using following expression
_ out hin Þ Q ¼ mðh
ð1Þ
_ mass flow rate, hout: enthalpy at the outlet, where Q: heat load, m: hin: enthalpy at the inlet. The local heat transfer coefficient was obtained by using the following expression
h¼
Q =A T wall T bulk
¼1 h n
Xk
h i¼1 i
ð2Þ
where A: heat transfer area, Twall: wall temperature, Tbulk: bulk tem channel averaged local heat transfer coefficient. perature; h: The bulk temperature was obtained by using following expression
R q h V dA T bulk ¼ R q C P V dA
ð3Þ
where q: density, h: static enthalpy, V: velocity, CP: specific heat capacity at constant pressure.
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185
181
The Nusselt number was obtained by using following expression
Nu ¼
hlocal l k
Nu ¼
1 Xk Nui i¼1 n
ð4Þ
where hlocal: local heat transfer coefficient, Nu: channel averaged Nusselt number. The Fanning friction factor f was obtained by using following expression
f ¼
Xk ðP in Pk Þ Dh f ¼ 1 f 2 i¼1 i n Lk 2qk V k
ð5Þ
where f: friction factor, Pin Pk: pressure drop at the k-th plane, qk: density at the k-th plane, Vk: velocity at the k-th plane, Dh: hydraulic diameter, Lk: active length from the inlet, f : channel averaged friction factor. 5. Results 5.1. Results of Nusselt number
Fig. 9. Result of analysis model 1 (32.5°) – Nusselt number in the cold side compared to the existing correlations.
Nusselt numbers obtained from CFD analysis in the hot and cold channels were compared to those calculated from the existing correlations. The ranges of the existing correlations were extended up to Reynolds number of 58,000. In the extended Reynolds number ranges – outside the originally proposed valid range – the correlations are shown as dotted or dashed line in Figs. 8–11. The channel averaged results from CFD analysis are compared to the calculated values from correlations, and relative differences are summarized in Tables 7 and 8. As can be seen in Figs. 8–11, the CFD-obtained Nusselt numbers tend to be enveloped by the Ishizuka’s and Ngo’s correlation. It is worth noting that the Hesselgreaves’s correlation – found for a general plate-fin type heat exchanger – gives, in general, a poor agreement with CFD-results. Analysis model 1 (32.5°) – hot side (Fig. 8): the fin angle of the hot side (32.5°) is the same as the experiment of Ishizuka et al. (2005). As anticipated, the CFD-results give a good agreement with Ishizuka’s prediction for the Reynolds number range (2400 < Re < 6000) provided by the correlation. Ishikuka’s prediction is also shown to give an accurate prediction for CFDobtained results below the explored minimum Reynolds number (Re < 2400). However, above the tested Reynolds number
Fig. 8. Results of analysis model 1 (32.5°) – Nusselt number in the hot side compared to the correlations.
Fig. 10. Result of analysis model 2 (40.0°) – Nusselt number in the hot side compared to the existing correlations.
Fig. 11. Result of analysis model 2 (40.0°) – Nusselt number in the cold side compared to the existing correlations.
182
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185
Table 7 Relative differences (channel-averaged differences between CFD-result and correlation) of Nusselt number results of analysis model 1 in comparison to the correlations. Analysis model 1 (32.5°)
Ishizuka’s correlation (Hot side) (%)
Ishizuka’s correlation (Cold side) (%)
Ngo’s correlation (Hot side) (%)
Ngo’s correlation (Cold side) (%)
Hesselgreaves’s correlation (Hot side) (%)
Hesselgreaves’s correlation (Cold side) (%)
Relative difference in the valid range Relative difference in the extended range
+5.18
+8.04
11.18
14.15
+26.85
+20.77
+23.66
+28.77
15.25
33.12
+35.67
+31.71
Table 8 Relative differences (channel-averaged differences between CFD-result and correlation) of Nusselt number results of analysis model 2 in comparison to the correlations. Analysis model 2 (40.0°)
Ishizuka’s correlation (Hot side) (%)
Ishizuka’s correlation (Cold side) (%)
Ngo’s correlation (Hot side) (%)
Ngo’s correlation (Cold side) (%)
Hesselgreaves’s correlation (Hot side) (%)
Hesselgreaves’s correlation (Cold side) (%)
Relative difference in the valid range Relative difference in the extended range
+5.58
+8.35
8.84
29.01
+43.52
+21.33
+22.06
12.44
+10.16
61.84
+43.77
+35.11
range by Ishizuka (Re > 6000), the Ngo’s correlation gives a better prediction for Nusselt number whereas Ishizuka’s correlation gives an increasing difference with CFD-results. Analysis model 1 (32.5°) – cold side (Fig. 9): while Ishizuka’s prediction gives an acceptable agreement with CFD-obtained results within the valid correlation range and (5000 < Re < 13,000), we recommend to use Hesselgreaves’s or Ngo’s correlation for the Nusselt number with Reynolds number greater than 10,000. It is worth noting that Ishizuka’s correlation was developed with the cold fin angle of 40°, which is different from the CFD cold fin angle 32.5°. Such a difference in the fin angle is anticipated to be a cause for the observed departure from the CFD results in the valid correlation range. Analysis model 2 (40.0°) – hot side (Fig. 10): one can note that predictability of the correlations for hot channel CFD results does not noticeably change from changing the fin angle 32.5° (see Fig. 8) to 40° (see Fig. 10). Hence, Ishizuka’s correlation is recommended for the Reynolds number below 6000 – the maximum of the experimentally tested range. Above Reynolds number 6000, the Ngo’s correlations give a good agreement with CFD results.
Analysis model 2 (40.0°) – cold side (Fig. 11): fin angle of the cold side (40.0°) is same as the experiment of Ishizuka et al. (2005). The CFD-results give the best agreement with Ishizuka’s prediction throughout the simulated Reynolds number range (4000–58,000). That is, Ishizuka’s correlation can be safely extrapolated beyond the tested Reynolds number range (5000–13,000). It is noteworthy that in comparison with the CFD results, the Ngo’s correlation underestimates Nusselt numbers throughout the simulated Reynolds number range. Hesselgreaves’s correlations starts to exhibit an agreement from the CFD-obtained results above Reynolds number around 45,000.
5.2. Result of friction factor CFD-obtained channel-averaged friction factors are compared with those obtained from the correlations and their relative differences are summarized in Tables 9 and 10. The correlations of Ishizuka were obtained at fin angles of 32.5°, and 40° for the hot, and cold channels, respectively. The Ngo’s correlation is based on the experimental results with the fin angle of 52°. Hesselgreaves’s
Table 9 Relative differences of friction factors results of analysis model 1 in comparison to the correlations. Analysis model 1 (32.5°)
Ishizuka’s correlation (Hot side) (%)
Ishizuka’s correlation (Cold side) (%)
Ngo’s correlation (Hot side) (%)
Ngo’s correlation (Cold side) (%)
Hesselgreaves’s correlation (Hot side) (%)
Hesselgreaves’s correlation (Cold side) (%)
Relative difference in the valid range Relative difference in the extended range
22.04
+55.43
+51.01
+60.09
+30.92
+28.10
More than 100
More than +100
+50.79
+60.17
+43.63
+46.33
Table 10 Relative differences of friction factor results of analysis model 2 in comparison to the correlations. Analysis model 2 (40.0°)
Ishizuka’s correlation (Hot side) (%)
Ishizuka’s correlation (Cold side) (%)
Ngo’s correlation (Hot side) (%)
Ngo’s correlation (Cold side) (%)
Hesselgreaves’s correlation (Hot side) (%)
Hesselgreaves’s correlation (Cold side) (%)
Relative difference in the valid range Relative difference in the extended range
More than 100
1.86
16.04
+8.91
+19.35
+18.9
More than 100
More than 100
21.29
+3.83
17.21
+15.5
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185
183
Analysis model 1 (32.5°) – hot side (Fig. 12): the Ishizuka’s correlation (32.5°) gives the closest agreement with the CFD results in the valid Reynolds number range. However, the correlation assumes linear behavior of friction factor with respect to Reynolds number, and this inherently limits an extension of the correlation beyond the tested range. That is, one can clearly see in Fig. 12 with the CFD-results and the other correlations that the marginal decrease in friction factor with increasing Reynolds number decreases, which exhibits a significant departure from the linear behavior. An exponent type correlation would be a proper form when considering the wide range of Reynolds number for the friction factor correlations. Analysis model 1 (32.5°) – cold side (Fig. 13): CFD results of cold side shows a good accordance with the Hesselgreaves’s correlation in the valid range. Relative difference is 28.1% in this valid Reynolds number range. The CFD results demonstrate a large difference with the Ishizuka’s and Ngo’s correlation. The difference comes from the disagreement of the angles between the correlations (40° for Ishizuka, and 52° for Ngo) and the simulation (32.5°). Ishizuka’s correlation gives almost a factor of two
difference in friction factor for the angle 32.5° and 40°. Such sensitivity of fin angles on the friction factor is later confirmed by CFD analysis. Hence, the loss of predictability with Ishizuka’s and Ngo’s correlation for the simulation at fin angle 32.5° predominantly pertains to the disagreement of the fin angle. Analysis model 2 (40.0°) – hot side (Fig. 14): although the fin angle for the Ngo’s correlation (52°) is far greater than that of the model (40°), the CFD results give only little higher friction factors than those obtained by the Ngo’s correlation. This analysis model 2 employs the fin angle of 40.0° and it significantly increases the friction factor from the simulation with the fin angle of 32.5°. The cold channel gives slightly lower friction factors than those of the hot channel because of larger density. The Hesselgreaves’s correlation becomes similar value to that of the Ngo’s after 15,000 Reynolds number range. This correlation shows a good agreement in the valid Reynolds number range. Analysis model 2 (40.0°) – cold side (Fig. 15): CFD results showed little difference with the Ishizuka’s correlation in the valid Reynolds number range. This is, in part, because CFDmodel has the same the cold channel fin angle (40°) and boundary conditions with the experiment Ishizuka conducted to obtain the correlation. Again, the assumed linear behavior of Ishizuka’s friction factor with Reynolds number limits an
Fig. 12. Result of analysis model 1 (32.5°) – friction factor in the hot side compared to the existing correlations.
Fig. 14. Result of analysis model 2 (40.0°) – friction factor in the hot side compared to the existing correlations.
Fig. 13. Result of analysis model 1 (32.5°) – friction factor in the cold side compared to the existing correlations.
Fig. 15. Result of analysis model 2 (40.0°) – friction factor in the cold side compared to the existing correlations.
correlation can deal with various fin angles as it requires the angle as an input.
184
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185 Table 11 Recommended combination of existing correlations for PCHE as S-CO2 power cycle recuperator design.
Hot side – fin angle 32.5° Hot side – fin angle 40.0°
Cold side – fin angle 32.5°
Cold side – fin angle 40.0°
Nusselt number
Friction factor
2400 < Re < 6000: Ishizuka’s correlation 6000 < Re < 30,000: Ngo’s correlation 5000 < Re < 10,000: Ishizuka’s correlation 10,000 < Re < 58,000: Ngo’s correlation 3500 < Re < 58,000: Ishizuka’s correlation
Fin angle 32.5: Hesselgreaves’s correlation Fin angle 40.0: Ngo’s correlation
Table 12 New correlation derived from the CFD results.
Fin angle 32.5°
Fin angle 40.0°
Nusselt number
Friction factor
2000 < Re < 58,000, (0.7 < Pr < 1.0)
2000 < Re < 58,000 (0.7 < Pr < 1.0)
Nu ¼ ð0:0292 0:0015Þ Re0:81380:0050
f ¼ ð0:2515 0:0097Þ Re0:20310:0041
R2 ¼ 0:999
R2 ¼ 0:983
2000 < Re < 55,000 (0.7 < Pr < 1.0)
2000 < Re < 55,000 (0.7 < Pr < 1.0)
Nu ¼ ð0:0188 0:0032Þ Re0:87420:0162
f ¼ ð0:2881 0:0212Þ Re0:13220:0079
R2 ¼ 0:995
R2 ¼ 0:908
Table 13 Design result of PCHE recuperator by using KAIST-HXD code with new correlation.
Reynolds number range Average Reynolds number Heat load Relative volume (volume divided by volume III) Pure mass Length Relative frontal area (frontal area divided by frontal area III) Pressure drop
Design result II with T.L. Ngo correlation
Design result IV with new correlation
3500 < Re < 22,000 Hot side: 21,542 Cold side: 18,850 67.1 MW 2.009 4901 kg 0.639 m 2.275 Hot side: 187.2 kPa Cold side: 61.1 kPa
2000 < Re < 58,000 Hot side: 40,107 Cold side: 35,687 67.1 MW 1.011 2398 kg 0.650 m 1.126 Hot side: 299.5 kPa Cold side: 98.2 kPa
extrapolation beyond the experimented Reynolds numbers. Throughout the Reynolds number ranges in Fig. 15, the Ngo’s correlation, in general, gives the best agreement with the CFD results among the compared correlations. 5.3. PCHE design with new correlation Consequently, we are able to identify which correlation best fits CFD results for a range of Reynolds numbers. That is, we can propose a new correlation for a wide Reynolds number range with a piece-wise function of appropriately extended existing correlations as shown in Table 11. The proposed correlation mainly combines Ishizuka’s correlation in the low Reynolds number range and Ngo’s correlation in the high Reynolds number range. Since friction factor mainly depends on the fin angle of a channel, it was categorized according to the fin angles. Based on the result of CFD analysis, new sets of correlation were obtained using non-linear curve fitting. The results of both the hot and cold sides were used to derive a correlation that covers a wide range of Reynolds numbers. However, the results were restricted in small range of Prandtl number. This is because the operating range of a recuperator in the supercritical CO2 power cycle behaves very close to normal gas even though it is in the supercritical state. Thus, the new correlation only depends on Reynolds numbers in restricted range of Prandtl number (Pr near 0.7 only). New sets of correlation are shown in Table 12. Also, standard error and
correlation coefficient for each correlation were found. The previously discussed recuperator in Section 3 was redesigned with this new correlation in Table 11. The same pressure drop restriction (1.5% of the maximum system pressure) and minimum target effectiveness (93%) applied for design I–III were used. The newly designed PCHE was compared to design result II as shown in Table 13. The external dimensions underwent a noticeable change in more practical, yet manufacturable, way thanks to the new correlation that deals with a wide range of Reynolds number. However, pressure drop in both sides are increased as much as 1.6 times compared to the design result II due to the long channel length and the high flow velocity. The volume and frontal area of design result IV are slightly larger than those of the design result III, because the new correlation has a lower Reynolds number limit (Re < 58,000) than the limit of the Hesselgreaves’s correlation. 6. Conclusions CFD is a potent tool to explore flow characteristics that are sometime too costly to be examined by experiments; an extension of existing heat transfer and friction factor correlations for an increasing Reynolds number is a kind. CFD analyses were performed to examine heat transfer and pressure drop performance of a PCHE for the S-CO2 power cycle application. The CFD results for Nusselt number and friction factors with a wide range of Reynolds number enables to explore if existing correlations can
S.G. Kim et al. / Annals of Nuclear Energy 92 (2016) 175–185
be extended beyond the originally proposed valid Reynolds number range. For heat transfer coefficients, an appropriate piecewising with Ishizuka’s, Ngo’s, and Hesselgreaves’s correlation is found to cover a wide range of Reynolds numbers. For friction factors, no single existing correlation is found to represent different temperature and angular dependencies for a wide Reynolds number range. Yet, a proper piece-wising with Hesselgreaves’s, Ngo’s correlations can well cover a wide range of Reynolds numbers. PCHE design is subject to a significant change with the new correlation that can cover a wider range of Reynolds number. This illuminates that an appropriate combination of experimental work and computational simulation should facilitate understanding physical phenomena, and advance component and system designs. Acknowledgements Authors gratefully acknowledge that this research is supported by the National Research Foundation of Korea (NRF2013M2A8A1041508) and funded by the Korean Ministry of Science, ICT and Future Planning. References Ahn, Yoonhan, Lee, Jeong Ik, 2014. Study of various Brayton cycle designs for small modular sodium-cooled fast reactor. Nucl. Eng. Des. 276, 128–141. ANSYS, 2011. ANSYS CFX-Solver Theory Guide, Release 14.0. Tech., rep., November. Bae, Seong Jun, Ahn, Yoonhan, Lee, Jekyoung, Lee, Jeong Ik, 2014. Various supercritical carbon dioxide cycle layouts study for molten carbonate fuel cell application. J. Power Sources 270, 608–618. Bae, Seong Jun, Lee, Jekyoung, Ahn, Yoonhan, Lee, Jeong Ik, 2015. Preliminary studies of compact Brayton cycle performance for small modular high temperature gas-cooled reactor system. Ann. Nucl. Energy 75, 11–19. John E. Hesselgreaves, 2001. Compact Heat Exchangers – Selection, Design and Operation, Pergamon. Takao, Ishizuka, Yasuyoshi, Kato, Yasushi, Muto, Konstantin, Nikitin, Ngo Lam, Tri, Hiroyuki, Hashimoto, 2005. Thermal–hydraulic characteristic of a printed circuit heat exchanger in a supercritical CO2 loop. In: The 11th International Topical Meeting on Nuclear Reactor Thermal–Hydraulics (NURETH-11), Avignon, France, October 2–6.
185
Jeong, Woo Seok, Lee, Jeong Ik, Jeong, Yong Hoon, 2011. Potential improvements of supercritical recompression CO2 Brayton cycle by mixing other gases for power conversion system of a SFR. Nucl. Eng. Des. 241, 2128–2137. Kays, London, 1984. Compact Heat Exchangers, third ed. Krieger. Kim, Dong Eok, Kim, Moo Hwan, Cha, Jae Eun, Kim, Seong O., 2008. Numerical investigation on thermal–hydraulic performance of new printed circuit heat exchanger model. Nucl. Eng. Des. 32 (238), 3269–3276. Seong Gu, Kim, Min Gil, Kim, Seong Jun, Bae, Jeong Ik, Lee, 2013. Preliminary design of S-CO2 Brayton cycle for KAIST micro modular reactor. In: Transactions of the Korean Nuclear Society Autumn Meeting, Gyeongju, Korea, October 24–25. Kim, Seong Gu, Lee, Jekyoung, Ahn, Yoonhan, Lee, Jeong Ik, Addad, Yacine, Ko, Bockseong, 2014. CFD investigation of a centrifugal compressor derived from pump technology for supercritical carbon dioxide as a working fluid. J. Supercrit. Fluids 86, 160–161. Lee, Sang-Moon, Kim, Kwang-Yong, 2014. Thermal performance of a double-faced printed circuit heat exchanger with thin plates. J. Thermophys. Heat Transfer 28, 251–257. Lee, Youho, Lee, Jeong Ik, 2014. Structural assessment of intermediate printed circuit heat exchanger for sodium-cooled fast reactor with supercritical CO2 cycle. Ann. Nucl. Energy 73, 84–95. Lee, Jekyoung, Lee, Jeong Ik, Yoon, Ho Joon, Cha, Jae Eun, 2014. Supercritical carbon dioxide turbomachinery design for water-cooled small modular reactor application. Nucl. Eng. Des. 270, 76–89. Lemmon, E.W., Huber, M.L., McLinden, M.O., 2007. NIST reference fluid thermodynamic and transport properties – REFPROP. Nat. Inst. Stand. Technol. Ngo, Tri Lam, Kato, Yasuyoshi, Nikitin, Konstantin, Ishizuka, Takao, 2007. Heat transfer and pressure drop correlations of microchannel heat exchangers with S-shaped and zigzag fins for carbon dioxide cycles. Exp. Therm. Fluid Sci. 32, 560–570. Nikitin, Konstantin, Kato, Yasuyoshi, Ngo, Lam, 2006. Printed circuit heat exchanger thermal–hydraulic performance in supercritical CO2 experimental loop. Int. J. Refrig. 29, 807–814. Oyakawa, Kenyu, Shinzato, Takao, Mabuchi, Ikuo, 1989. The effect of the channel width on heat-transfer augmentation in a sinusoidal wave channel. Jpn. Soc. Mech. Eng. 32 (3), 403–410, Series II. Pra, Franck, Tochon, Patrice, Mauget, Christian, Fokkens, Jan, Willemsen, Sander, 2008. Promising designs of compact heat exchangers for modular HTRs using the Brayton cycle. Nucl. Eng. Des. 238, 3160–3173. Eric N., Val Abel, Mark H., Anderson, Michael L., Corradini, 2011. Numerical investigation of pressure drop and local heat transfer of supercritical CO2 in printed circuit heat exchangers. In: Supercritical CO2 Power Cycle Symposium, Boulder, Colorado, May 24–25. Yoon, Ho Joon, Ahn, Yoonhan, Lee, Jeong Ik, Addad, Yacine, 2012. Potential advantages of coupling supercritical CO2 Brayton cycle to water cooled small and medium size reactor. Nucl. Eng. Des. 245, 223–232.