Thermal-hydraulic performance of printed circuit heat exchangers with zigzag flow channels

International Journal of Heat and Mass Transfer 130 (2019) 356–367

Contents lists available at ScienceDirect

International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt

Thermal-hydraulic performance of printed circuit heat exchangers with zigzag flow channels Minghui Chen a,⇑, Xiaodong Sun a,⇑, Richard N. Christensen b a b

Department of Nuclear Engineering and Radiological Sciences, University of Michigan, Ann Arbor, MI 48109, United States Nuclear Engineering Program, University of Idaho, Idaho Falls, ID 83401, United States

a r t i c l e

i n f o

Article history: Received 23 August 2018 Received in revised form 21 September 2018 Accepted 7 October 2018

Keywords: HTGR PCHE Thermal hydraulics Thermal boundary conditions CFD

a b s t r a c t Printed circuit heat exchanger (PCHE) is one of the leading candidates to be employed in advanced nuclear reactors and next generation concentrated solar power applications due to its compactness and capability for high-temperature, high-pressure applications with high effectiveness. In the current study, thermal-hydraulic performance of a zigzag-channel PCHE with high-pressure, high-temperature helium on both the hot and cold sides was simulated using a computational fluid dynamics (CFD) software package STAR-CCM + . Comparisons between the experimental data and CFD simulation results showed good agreement with some discrepancies in the pressure drop and heat transfer results. Local thermal-hydraulic performance analyses indicated that a fully-developed flow condition was not observed in the PCHE, mainly due to the nature of the zigzag channels, leading to periodic flow disturbance at each of the zigzag bends. It was also found that the local and global heat transfer coefficients were considerably different in the PCHE. Furthermore, thermal boundary conditions showed that the fluid temperatures and heat fluxes were not uniform along the azimuthal direction of a cross section of the flow channel and that the helium temperature distribution for each segment along the flow direction presented a wavy profile. However, the distribution of the helium bulk temperature along the flow direction was approximately linear. For the heat flux distributions, although they were significantly different at different segments, the trend of the heat flux for each segment along the fluid flow direction was similar. Finally, effects of several parameters on the thermal-hydraulic performance of the PCHE were investigated, including the fluid and solid thermophysical properties, radius of curvature at zigzag bends, channel configuration, channel length pitch in the flow direction, and zigzag pitch angle. No considerable enhancement in the Nusselt numbers was observed when the zigzag pitch angles were greater than 30°. Ó 2018 Elsevier Ltd. All rights reserved.

1. Introduction Printed circuit heat exchanger (PCHE) is one of the leading candidates to be employed in advanced nuclear reactors and next generation concentrated solar power applications due to its compactness and capability for high-temperature, high-pressure applications with high effectiveness [1]. The Very-HighTemperature Reactor (VHTR) is one promising candidate for advanced reactors due to its capability of generating electricity with high efficiency and providing high-temperature process heat for industrial applications. The helium temperature at the reactor core outlet is reduced to 800 °C in phase I of the reactor development stage and is planned to be increased to 1000 °C or above in the future development [2]. The VHTR can support a number of ⇑ Corresponding authors. E-mail addresses: [email protected] (M. Chen), [email protected] (X. Sun). https://doi.org/10.1016/j.ijheatmasstransfer.2018.10.031 0017-9310/Ó 2018 Elsevier Ltd. All rights reserved.

applications that require such high coolant temperatures. The electric power generation may use a Rankine cycle with high-pressure steam generators, or a direct Brayton cycle with the primary helium coolant as the working fluid, or an indirect Brayton cycle with a secondary helium as the working fluid. Potential process heat applications may include hydrogen production, coal gasification, petroleum refining, bio-fuels production, and production of chemical feedstocks for use in the fertilizer and chemical industries [3]. The efficiencies of electricity generation and process heat applications of VHTRs are significantly dependent upon the performance of intermediate heat exchangers (IHXs). The IHX in a VHTR is a key component for transferring thermal energy from the primary coolant to a secondary coolant, which could be helium, supercritical carbon dioxide (sCO2), or molten salt. In addition, the IHX serves as an interface that isolates the VHTR primary system from the plants for electricity generation and/or process heat

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

applications. Therefore, it must be sufficiently robust to maintain the entire system integrity during normal and off-normal conditions. Since helium typically has a low heat transfer capability due to its low volumetric thermal capacity and low thermal conductivity, a compact heat exchanger with a high surface area to volume ratio, generally at 700 m2/m3 or higher [4], is preferable to be employed as an IHX in VHTRs. PCHEs stand out from several heat exchanger candidates due to their high effectiveness, compactness, robustness, and their ability to withstand high pressures and temperatures [5]. The first PCHE was developed in 1980s in Australia and a company called Heatric was founded based on the developed PCHE fabrication technology [6]. PCHEs are platetype compact heat exchangers in which flow channels (typically, channels with a small hydraulic diameter) are etched into flat metal plates using a photochemical machining process. There are several PCHE designs with respect to flow channel geometries, such as straight, zigzag, S-shape finned, and airfoil finned channels [7]. The etched metal plates are grinded and lapped to remove scratches on the plate surfaces followed by making the etched plates flat and parallel. At the end, the etched plates are stacked together with a prescribed arrangement configuration and diffusion bonded to create a high-integrity solid block before flow distribution headers are attached to the heat exchanger block [8]. Several studies have been performed on the steady-state thermal-hydraulic performance of PCHEs over the past several years. Kim et al. [9] proposed a mathematical expression for predicting the thermal performance of cross flow, parallel flow, and counter flow PCHEs with semicircular straight channels based on an extensive numerical study for sCO2. A similar study was later performed by the same research group for wavy-channel PCHEs and the effects of waviness were analyzed [10]. Chu et al. [11] carried out an experimental study on a PCHE under water-to-sCO2 heat transfer condition and the obtained results showed that sCO2 has better heat transfer capability than water when the mass flow rates on both sides were the same. This confirms that sCO2 has superior heat transfer performance than its gas phase. A numerical study on PCHEs with new airfoil fins under sCO2-to-sCO2 condition was performed by Cui et al. [12] in a recent study and the authors concluded that PCHEs with airfoil fins have better thermalhydraulic performance than PCHEs with zigzag channels. However, they did not conduct a study on zigzag-channel PCHEs. These studies mainly focused on the thermal-hydraulic performance of PCHEs with a working fluid of sCO2. There are also many studies on PCHEs under helium fluid conditions. Mylavarapu [3] designed and constructed a high-temperature helium test facility (HTHF) that can facilitate compact heat exchanger thermalhydraulic testing using helium as the working fluid at temperatures and pressures up to 800 °C and 3.0 MPa with helium mass flow rates ranging from 15 to 49 kg/h. Chen et al. [13,14] presented the fabrication procedure of a counter-flow PCHE with zigzag channels that can be used in high-temperature, high-pressure conditions and performed a thermal-hydraulic investigation of the fabricated PCHE for steady-state and transient conditions at high temperatures. Bartel et al. [15] performed a comparative study of PCHE designs for a typical VHTR design and concluded that the zigzag-channel PCHE with an inclination angle of 15° provided the best performance in heat transfer while incurring a moderate pressure drop penalty as compared to the straight-channel PCHEs. Yoon et al. [16] conducted a comprehensive numerical study on the flow and heat transfer in zigzag-channel PCHEs with various geometric parameters. General explicit correlations for both the friction factor and Nusselt number as a function of geometric parameters and Reynolds number were developed through numerical simulations. The correlations are applicable to a wide range of geometric parameters and operation conditions. It was indicated

357

that the friction factor for zigzag-channel PCHEs was mainly influenced by the zigzag channel geometry while the Nusselt number was influenced by the overall heat exchanger design, including the plenum sections. Aneesh et al. [17] also carried out a numerical study on straight-channel PCHEs with helium as the working fluid and Alloy 617 as the heat exchanger structural material. They investigated the effects of wavy channel configurations and concluded that the trapezoidal wavy channel PCHE has the best performance. Although a number of studies have been performed on PCHEs, a limited number of numerical and experimental research has been performed on zigzag-channel PCHEs with helium as the working fluid at high temperatures. The local thermal-hydraulic performance of zigzag-channel PCHEs has not been thoroughly investigated with several unknowns for different effects on their performance. The aim of the current research is therefore to investigate the thermal-hydraulic performance of zigzag-channel PCHEs under high-temperature helium conditions. In this paper, a computational fluid dynamics (CFD) code, i.e., STAR-CCM+, was used to simulate detailed thermal-hydraulic performance of a fabricated PCHE with a simplified two-channel geometry model. The thermal boundary conditions, and effects of thermophysical properties and geometrical parameters on the thermal-hydraulic performance of PCHEs with zigzag flow channels are numerically examined. 2. Printed circuit heat exchanger and test facility 2.1. Printed circuit heat exchanger A laboratory-scale zigzag-channel PCHE was fabricated. A total of eight hot and eight cold plates, made of Alloy 617, with 11 zigzag channels in each of the 1.6-mm thick plates were first stacked in an alternating manner with one 12.7-mm thick Alloy 617 plate added each at the top and bottom of the stack and then diffusion bonded together to form a monolithic metal block, i.e., the heat exchanger core. Fig. 1 shows a section of a hot-side or cold-side plate. The cross section of the fluid passages is approximately semicircular with a nominal diameter of 2.0 mm and a nominal pitch of 2.5 mm in the span-wise direction. The longitudinal configuration of the flow passage is zigzag with roundness at the tip of each bend. The radius of curvature at the bends is 4.0 mm. The angle between the flow direction in the zigzag channels and the edge of the block, i.e., the direction along the heat exchanger length, is 15°. The PCHE was designed in such a way that each side can withstand a maximum pressure of 3.0 MPa at 800 °C temperature in the HTHF. The overall dimensions of the PCHE block are 339.1 mm in length, 126.0 mm in width, and 50.8 mm in height. The detailed fabrication procedures and heat exchanger geometries can be found in Ref. [8]. 2.2. High-temperature helium test facility Experiments were performed to examine the thermal-hydraulic characteristics of the fabricated zigzag-channel PCHE under various steady-state conditions. Fig. 2 is a schematic of the HTHF showing the helium flow path for heat exchanger testing, where one of the two previously fabricated straight-channel PCHEs was replaced by the current zigzag-channel PCHE (i.e., the PCHE being tested and studied) [13]. The HTHF was designed and constructed to facilitate thermal-hydraulic performance testing of compact heat exchangers, including diffusion-bonded heat exchangers, at temperatures and pressures up to 800 °C and 3.0 MPa, respectively. The readers can refer to Ref. [3] for more details about the HTHF design.

358

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

Fig. 1. Images of a section of a hot-side or cold-side plate [13].

Fig. 2. Layout of the high-temperature helium test facility [13].

3. Numerical modeling 3.1. Numerical simulation model A three-dimensional CFD study was carried out using STARCCM+. Modeling a full-size laboratory-scale PCHE is highly time

consuming and can lead to prohibitively high computational costs. Therefore, a simplified two-channel geometry model, consisting of one hot-side flow channel, one cold-side flow channel, and the associated solid structure metal, as shown in Fig. 3, was modeled in the current simulations. The cross-sectional shape of both the cold- and hot-side flow channels is semicircular with a nominal

359

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

  /  /2  : er ¼  1 /1 

Fig. 3. Simplified geometry model with its geometrical parameters and boundary conditions [20].

diameter of 2.0 mm. The geometrical parameters of the computational domain were created to match the nominal values of the core region of the fabricated PCHE. The plate thickness is 1.6 mm and the modeled block width is 2.5 mm. The shape of the channels along the fluid flow direction is zigzag with a 4.0-mm-radius roundness at the tip of each bend. The geometry model was comprised of a total of 16 zigzag bends per flow channel. The physical properties of the plate material Alloy 617 were found from Special Metals [18]. The specific heat was formulated as a temperaturedependent function and the thermal conductivity values were tabulated as a table file. Both the specific heat function and thermal conductivity table were imported into STAR-CCM+. The density of Alloy 617 was considered as a constant of 8360 kg/m3. The temperature-dependent thermophysical properties of helium were exported from NIST chemistry webbook [19] in a table file and the file was implemented into STAR-CCM+. Note that the dependence of the helium density on its pressure was neglected since the pressure drop in this study along the PCHE channels was relatively small (<1%) as compared to the operating pressure of 3.0 MPa. Three boundary conditions were used in the model: fluid inlet, fluid outlet, and wall. The fluid inlets were specified as the mass flow rate inlets on both the hot and cold sides. The fluid outlets were specified as the pressure outlets on both sides. The left and right walls, as shown in Fig. 3, were coupled and set as a periodic boundary type. Similarly, a thermally-coupled periodic boundary type for the top and bottom walls was used. The front and back surfaces were specified as adiabatic walls. In this numerical study, the helium inlet temperatures were 800 and 350 °C on the hot side and cold side, respectively. The helium outlet absolute pressures on both sides were set to 3.0 MPa. To perform 3-D conjugate heat transfer simulations, steady state, laminar flow, and constant helium inlet temperatures were assumed while the thermal radiation effects and viscous dissipation were neglected in the simulations.

N 1X ðDV k Þ N k¼1

#1=3 ;

ð1Þ

where DV k is the volume of the kth cell and N is the total number of cells used for the entire computational domain. Three sets of grid sizes were selected with a refinement factor greater than 1.3 as suggested in [21]. The refinement factor is defined as,

r¼

hcoarse : hfine

The relative error can be obtained as,

3.3. Numerical data reduction and methodology The global pressure loss (equivalent Fanning friction) factor for each fluid, f g , is defined as,

fg ¼


A2c Dpf dh q ; _2 2lm

ð4Þ

where dh , Ac , and l denote the channel hydraulic diameter, the area of a channel cross section in the span-wise direction, and the actual _ is the mass flow rate. The mean fluid flow length, respectively. m
, is calculated using the arithmetic mean temperafluid density, q ture based on the helium inlet and outlet temperatures on each side. Dpf is the pressure loss, which can be expressed as,

Dpf ¼ Dpt  Dpa ;

ð5Þ

where Dpt and Dpa are the total pressure drop and the pressure drop due to flow acceleration or deceleration, respectively. There was no hydrostatic pressure drop since the heat exchanger was horizontally orientated. For a small segment along the span-wise direction, the local pressure loss factor is defined as,

fi ¼


i A2c Dpf ;i dh q ; _2 2Dli m

ð6Þ

where the subscript, i, denotes the i segment. The mean pressure loss factor over the entire flow channel can be computed as,

To examine the spatial convergence of the CFD simulations, a grid independence study was first conducted. The cell or grid size h is defined as,

"

The STAR-CCM+ software package provides several meshing strategies that can be used for different applications. In this study, parts-based meshing was selected with surface and polyhedral meshers. Since the largest helium inlet flow rate gives the most complex flow scenario that needs the finest mesh size to resolve the small flow structures, a grid independence study was carried out for the simplified simulation model with the maximum helium mass flow rate of 1.1238  104 kg/s per channel on both the hot and cold sides. Three grids with different numbers of cells were generated, as listed in Table 1. The helium outlet temperatures and pressure drop across each side were monitored. The helium outlet temperatures on both sides obtained from the coarse and medium mesh cases only had small differences compared to the fine mesh case. The relative errors for the hot-side and cold-side helium outlet temperatures between the medium and fine mesh results were 0.30 and 0.26%, respectively. For the pressure drops, both relative errors between the medium and fine mesh results were 3.8% for the hot and cold sides. In the current study, the fine mesh with 7.21 million cells was considered sufficient and was therefore used in the simulations. An image of the fine mesh is shown in Fig. 4.

th

3.2. Grid independence study

h¼

ð3Þ

ð2Þ


f ¼ 1 n

 Pn

i¼1 Dpf ;i

2

 Pn  q
i A2c dh  i¼1 ; _2 i¼1 Dli m

Pn

ð7Þ

where n is the total number of segments along the flow channel. Three heat transfer coefficient definitions were used in this numerical study, which are the local heat transfer coefficient, hj ;
and global heat transfer channel mean heat transfer coefficient, h; coefficient, hg . The local heat transfer coefficient is defined as:

hj ¼

q_ 00j jT f ;j  T w;j j

:

ð8Þ th

where the subscript, j, denotes the j segment in the span-wise direction; q_ 00 is the surface-averaged heat flux; T w is the surfaceaveraged temperature on the wall that is in contact with the

360

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

Table 1 Results of grid independency study. Mesh type

No. of cells, Million

Grid size, 106 m

Dph , Pa

Dpc , Pa

T h;o , °C

T c;o , °C

Fine Medium Coarse Relative error

7.21 2.39 0.73

61.5 89.0 132.0

27939.8 26890.3 24535.0 3.8%

23139.0 22256.8 20171.6 3.8%

511.5 513.0 508.6 0.30%

642.7 641.0 645.2 0.26%

4. Results and discussions 4.1. Numerical and experimental comparisons

Fig. 4. An image of the adopted fine mesh.

surrounding fluid; and T f is the reference fluid temperature. For the local heat transfer coefficient, the fluid bulk (mean) temperature is used for T f . The fluid mean temperature is given as [22],

R

qucp TdA T f ;c ¼ R : qucp dA Ac Ac

ð9Þ

The mean heat transfer coefficient over the entire channel can be computed as,


¼P h m

High-temperature experiments were carried to obtained the thermal-hydraulic performance of the fabricated hightemperature PCHE with zigzag channels. The pressure loss factors for the zigzag flow channels with the uncertainty bars representing the calculated uncertainties at the 95% confidence level for both the hot and cold sides were obtained. The relative uncertainties for the range of the Reynolds numbers tested in the experiments varied from 3.7 to 7.2% in those experimental data. The uncertainties of the pressure loss factors from the current experimental data were within ±6% in the laminar flow regime [7]. The pressure loss factors obtained from the experimental data were approximately 17% less than the global pressure loss factors obtained from the numerical simulations, as shown in Fig. 5. One of the reasons for the discrepancies could be due to potential distortions in the cross-sectional shape and dimensions of the flow channels after diffusion bonding. The local heat transfer coefficients could not be computed from the experimental data since no internal wall and fluid temperature measurements were conducted in the experiments. Instead, the global heat transfer coefficients were obtained from the steady-state experimental testing. Fig. 6 shows the comparison of the Nusselt numbers obtained from STAR-CCM+ and the experimental values computed from Eq. (11). Although the numerically obtained Nusselt numbers were slightly higher than the experimental results, they had very similar trends and most of the simulation results were within the experimental uncertainty bars. 4.2. Local thermal-hydraulic performance analysis It is challenging to obtain the local helium temperatures and pressures inside PCHEs from experiments due to the small size of channels. One advantage of the numerical simulations is that the

Pm

_ 00 j¼1 qj As;j

j¼1 ðAs;j jT f ;j

 T w;j jÞ

;

ð10Þ

where m is the total number of segments along the span-wise direction. The global heat transfer coefficient is often used in experimental data reduction and it is given by,

hg ¼

jT f ;g

q_ 00g ;  T w;g j

ð11Þ

where T f ;g is defined as the arithmetic mean of the fluid inlet and outlet temperatures on either the hot or cold side and T w;g is the arithmetic mean of the four terminal temperatures of a heat exchanger. They can be written as,

T f ;g ¼

ðT in þ T out Þ ; 2

ð12Þ

T w;g ¼

ðT h;in þ T h;out þ T c;in þ T c;out Þ : 4

ð13Þ

Fig. 5. Comparison of the experimental and numerical global pressure loss factors.

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

Fig. 6. Comparison of the experimental and numerical global Nusselt numbers.

local information, such as the temperatures, pressures, and velocities inside the heat exchanger, can be readily obtained. The base case studied here is for both the hot-side and cold-side inlets having the same mass flow rate of 5.108  105 kg/s, which gives a balanced heat capacity on both sides. Fig. 7 shows the temperature distributions in the fluid and solid metal inside the heat exchanger and Fig. 8 displays the gauge pressure (i.e., relative to the system pressure of 3.0 MPa) distributions on both the hot and cold sides along the fluid flow directions of the heat exchanger. It was observed that the fluid bulk mean temperature distribution on each side was approximately linear, which was expected due to the equal heat capacities on both sides. In addition, the pressure drop on the hot side was slightly larger than that on the cold side, as shown in Fig. 8. The helium flow on the cold side accelerated from the inlet to the outlet, while the helium velocity decreased from the inlet to the outlet on the hot side, as shown in Fig. 9. The changes in the helium temperature along the helium flow directions resulted in a decrease in the helium density on the cold side while an increase on the hot side of the heat exchanger. It can also be seen that the velocity on the hot side was higher than that on the cold side inside the heat exchanger, which is due to the smaller helium density on the hot side. Fig. 10 shows the contours of the local helium pressure, velocity magnitude, and temperature on several cross sections inside the flow channel. The velocity profile inside the heat exchanger changed along the fluid flow direction, which indicated that a fully-developed flow condition was not reached in the current simulations. This is attributed to the wavy nature of the zigzag flow channels as well as the large temperature differences along the flow channels, resulting in large fluid property variations. Fig. 11 shows the local heat transfer coefficients on both the hot and cold sides. Compared to the cold side, the higher velocities and larger thermal conductivities on the hot side led to higher average heat transfer coefficients. The helium flow was interrupted from becoming fully-developed at each zigzag bend and thereby restarting the boundary layers, which increased heat transfer coefficient. As shown in Fig. 11, each zigzag bend led to a higher heat transfer coefficient. This may result in a slightly wavy temperature distribution in the solid plates along the heat exchanger. It should be noted that the peak of the heat transfer coefficient, due to the complex flow phenomena, is not necessary right after a bend. It was observed that the highest heat transfer coefficients were at the inlets on both the hot and cold sides of the heat exchanger. This effect was also shown in Fig. 7 that smaller temperature

361

Fig. 7. Local temperature distributions along the fluid flow directions from the STAR-CCM+ simulations.

Fig. 8. Gauge pressure inside the heat exchanger from the STAR-CCM+ simulations.

Fig. 9. Velocity distributions inside the heat exchanger from the STAR-CCM+ simulations.

362

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

Fig. 11. Local heat transfer coefficient distributions inside the heat exchanger from the STAR-CCM+ simulations.

Fig. 12. Local heat flux distributions inside the heat exchanger from STAR-CCM+ simulations.

Fig. 10. Fluid contour plots of (a) pressure, (b) velocity magnitude, and (c) temperature.

differences between the fluid and metal walls were presented near the inlet section of each side of the heat exchanger. As the helium temperature decreased along the hot-side fluid flow direction, the overall trend of the heat flux slightly decreased, except for the hot-side outlet due to the entrance effect of the coldside inlet. It can be seen from Fig. 12 that the heat fluxes on both heat exchanger terminals had the largest values since the thickness of the thermal boundary layers in the entrance regions was small.

The effect of the bends on the heat flux was also displayed. It was shown that the heat flux had a substantially larger value in each zigzag bend. The change in the helium flow direction, as well as the local flow separation and enhanced flow mixing, led to the higher heat transfer rate. Similarly, the peak was not necessary right after a bend. In this flow condition and channel geometry, the peak of the heat flux was found to be near the middle of the straight section. Figs. 13 and 14 present the local pressure loss factors and local Nusselt numbers for both the cold and hot sides. The global pressure loss factors calculated by Eq. (4) were 0.03696 and 0.03531 for the hot and cold sides, respectively. The mean pressure loss factors calculated by Eq. (7) were 0.03581 and 0.03665 for the hot and cold sides, respectively. The mean pressure loss factor on the hot side was slightly smaller than that on the cold side by approximately 3%. It can be seen that the global values could predict the mean values over the entire channel since the discrepancies between the mean and global values were about 3.5%. As can be seen from Fig. 13, the pressure loss factor peaks on the hot and cold sides were different, which can be attributed to the location difference of the entrance region for a heat exchanger with a countercurrent flow configuration.

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

Fig. 13. Local pressure loss factors from the STAR-CCM+ simulations.

363

[4]. Therefore, it is essential to examine thermal boundary conditions when the fluid flow in a heat exchanger is in the laminar flow regime. Three important thermal boundary conditions for heat exchangers are: (1) constant wall temperature, both axially and peripherally throughout the passage length; (2) constant wall heat transfer rate in the axial direction and constant wall temperature in the azimuthal direction at any axial location; and (3) constant wall heat transfer rate in the axial direction as well as in the azimuthal direction. The simplified two-channel heat exchanger simulation model was divided into eight segments from each crosssectional area of the flow channels to identify the thermal boundary conditions for the zigzag-channel PCHE. A schematic of the CFD model with the cross-sectional area divisions and the associated notations is shown in Fig. 15. The temperature and heat flux distributions along the helium flow direction for each segment can be obtained from the CFD simulations. Since the geometry is axial symmetric with respect to the vertical plane passing through the two circle centers, the results for only the first four segments (i.e., h1, h2, h3, h4, c1, c2, c3, and c4) are presented. It can be seen from Figs. 16 and 17, the wall temperature distribution was not uniform in the azimuthal/peripheral direction of a channel crosssection. The temperatures on the top flat wall were smaller than those on the bottom semicircular wall on the hot side. On the cold side, the temperatures on the top flat walls are however larger than those on the bottom semicircular wall on the cold side, opposite to the hot side. The heat flux boundary conditions are shown in Figs. 18 and 19. The heat flux distributions on different segments were significantly different. However, the trend of heat flux for each segment along the helium flow direction was similar. The largest values were observed at every other bend for the first four segments while the smallest values on the straight sections of the zigzag bends. It is expected that the largest values occur at the other bends for the rest four heat exchanger segments, resulting in higher heat flux values after each of the bends in the zigzag channels. The actual thermal boundary conditions did not correspond to any of the

Fig. 14. Local Nusselt numbers from the STAR-CCM+ simulations.

Similarly, the mean Nusselt numbers calculated from the heat transfer coefficients by Eq. (10) were 11.66 and 12.34 for the hot and cold sides, respectively. The global Nusselt numbers obtained from the heat transfer coefficients by Eq. (11) were 9.18 and 9.94 on the hot and cold sides, respectively. The mean Nusselt number on the cold side was greater than that on the hot side by approximately 6%. However, the Nusselt numbers obtained from the global heat transfer coefficients underestimated the mean values by 21.3 and 19.4% for the hot side and cold side, respectively. For laminar flows in a heat exchanger, the heat transfer correlation is a strong function of not only the flow passage geometry but also the thermal boundary condition. Also, different reference temperatures would give the different Nusselt number values. Therefore, it should be cautious about the definition of the heat transfer coefficient as well as the operation conditions when it comes to using the heat transfer correlations for the heat exchanger thermal designs. 4.3. Various effects on thermal-hydraulic performance of PCHEs 4.3.1. Thermal boundary condition [23] The Nusselt number is strongly dependent on the thermal boundary condition and flow passage geometry for laminar flows

Fig. 15. Schematic of the divided geometry model and the associated notations.

364

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

Fig. 16. Hot-side boundary temperature distributions on the first four segment walls.

Fig. 19. Cold-side surface-averaged heat flux distributions on the first four segment walls.

simulation results may be caused by, in addition to the potential variations in the geometrical parameters due to chemically etching and diffusion bonding processes, the different thermal boundary conditions since the thermal boundary conditions of the zigzag channels in the experiments are much more complex than those in the numerical simulations, considering the thermal insulation and flow fluctuations. Figs. 20 and 21 present the helium temperature distributions for each segment on the hot and cold sides, respectively. As can be seen, the helium temperature distribution for each segment presented a wavy shape. However, the fluid bulk temperature distributions along the entire channel length were approximately linear. Fig. 17. Cold-side boundary temperature distributions on the first four segment walls.

Fig. 18. Hot-side surface-averaged heat flux distributions on the first four segment walls.

three idealized boundary conditions described above. Hence, the existing correlations in the literature may not truly reflect the actual heat transfer conditions in our experiments. Therefore, the discrepancies between the experimental data and the CFD

4.3.2. Effect of the variations of the fluid and solid thermophysical properties The effect of the variations of the fluid and solid thermophysical properties was studied through CFD simulations. For the first case, a simulation with constant fluid and solid properties was carried out using the same boundary conditions from the previous simulations. Table 2 shows the comparisons between these two simulation cases. As illustrated, the variation in the fluid properties could affect the pressure drop and heat transfer in a heat exchanger. Heat is transferred from the hot helium to the metal wall on the hot side, which sets up a temperature gradient in the helium in the radial direction at each cross section. Helium on the hot side has higher temperatures near the channel centerline and lower temperatures near the wall. Since the gas viscosity decreases with decreasing temperature, helium near the wall has a smaller viscosity than that near the channel centerline. Compared to the simulation with constant properties, helium near the wall shall have a larger velocity to satisfy the mass conservation. The decrease in the gas viscosity near the wall yields a lower friction factor and therefore a lower pressure loss. The increased velocity near the wall may lead to a more efficient convective heat transfer process, resulting in a higher heat transfer coefficient. Vice versa, a higher friction factor, a higher pressure loss, and a lower heat transfer coefficient are expected on the cold side. However, it should be emphasized that the variation in the helium velocity may also depend on the channel cross-sectional shape, surface interruptions, developing flows, and thermal boundary conditions. As can be seen from Table 2, compared with the simulation with temperature-dependent fluid properties, the pressure drop factors in the simulation with constant fluid and solid properties were

365

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367 Table 3 Comparisons of the simulation results using helium and nitrogen gases.

Num h Num c m fh m fc

Fig. 20. Hot-side helium temperature distributions for each heat exchanger segment along the PCHE length.

Fig. 21. Cold-side helium temperature distributions for each heat exchanger segment along the PCHE length. Table 2 Comparisons of the simulation results between the temperature-dependent and constant fluid and solid properties.

Num h Num c m

fh

m

fc

Temperature-dependent properties

Constant properties

Relative difference (%)

11.66 12.34 0.0358

8.908 9.142 0.0373

23.6 25.9 4.2

0.0366

0.0351

4.0

slightly lower on the cold side but higher on the hot side. However, the mean Nusselt numbers for the simulation with constant properties were 23.6 and 25.9% lower on the hot side and cold side, respectively. This was mainly attributed to the complex thermal boundary conditions as described in Section 4.3.1 and the channel shape that led to multiple developing flows inside the heat exchanger. The second case is to study the thermal-hydraulic performance with a different gas. The simulation was carried out using nitrogen gas with the same Reynolds numbers as the helium gas at both the hot-side and cold-side inlets. The comparison is listed in Table 3 and the discrepancies between the two simulations were less than 1.5%. It was expected since the Prandtl numbers of the two gases were similar.

Helium gas

Nitrogen gas

Relative difference (%)

11.66 12.34 0.0358

11.70 12.51 0.0360

0.36 1.33 0.50

0.0366

0.0364

0.68

4.3.3. Effect of the radius of curvature at zigzag bends A roundness was introduced at the tip of each bend in the zigzag flow channels of the fabricated PCHE. CFD simulations were conducted to study the performance of a simplified PCHE model with sharp bends (i.e. with a zero radius of curvature for the sharp bends) as well as a radius of curvature of 8.0 mm in the flow channels. The base case is the simplified two-channel model with a radius of curvature of 4.0 mm. Table 4 shows comparisons of the mean pressure loss factor and the mean Nusselt number among the three zigzag channels with the rounded and sharp bends. e1 and e2 indicated the percentages of deviation from the base case. Both the mean pressure loss factors and mean Nusselt numbers in the zigzag channel with the sharp bends were larger than those in the zigzag channels with rounded bends. This was expected since the flow in the channels with the sharp bends was more chaotic than that in the channels with rounded bends, resulting in enhancing the heat transfer but increasing the pressure loss. Increasing the radius of curvature from 4.0 to 8.0 mm reduced the pressure loss factor by about 5.6 and 5.3% on the hot side and cold side, respectively. However, the Nusselt numbers on the hot side and cold side were also reduced by about 5.0 and 4.7%, respectively. Whether it is worth introducing a radius of curvature to the bends is based on multiple factors, such as the capital cost, heat exchanger effectiveness, and the entire system efficiency, and this discussion is out of the current research scope. For the case studied, it is noted that the heat exchanger performance was not varied considerably by increasing the radius of curvature from 4.0 to 8.0 mm. 4.3.4. Effect of different channel configuration The effect of different channel configurations was examined over five different channel arrangements that are shown in Fig. 22. Type a is the base case that was used for the simulations in the previous sections. The flow length of the heat exchanger was reduced to half of the previous simulation model while the inlet helium temperatures were kept the same. Simulations were conducted for three mass flow rates in the laminar flow regime. Table 5 only compares the heat exchanger effectiveness for the five configurations. It was observed that Types a and e had the smallest effectiveness values while Types c and d, in general, had the largest values. However, fabrication of Types c and d heat exchangers is challenging. 4.3.5. Effect of different channel pitch length and zigzag angle A comprehensive numerical study on the friction factor and heat transfer for zigzag-channel PCHEs with various geometric

Table 4 Comparisons of simulation results with the different radii of curvature for the zigzagchannel bends.

Num h Num c m

fh

m

fc

Base case R = 4.0 (mm)

R = 8.0 (mm)

e1 (%)

R = 0 (mm)

e2 (%)

11.66 12.34 0.0358

11.08 11.76 0.0338

5.0 4.7 5.6

12.09 12.77 0.038

3.7 3.5 6.3

0.0366

0.0347

5.3

0.039

5.9

366

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

Fig. 22. Schematic of the five different channel configurations studied.

Table 5 Comparison of effectiveness’s for the five different channel configurations. Mass flow rate (kg/s)

Effectiveness (%)

5

4.00  10 6.00  105 8.00  105

Type a

Type b

Type c

Type d

Type e

70.73 65.08 60.76

71.11 65.84 61.80

71.04 66.09 62.13

71.02 66.04 62.11

70.84 65.80 61.76

parameters was carried out at Idaho National Laboratory [16]. General explicit correlations, Eqs. (14)–(16), for both the friction factor and Nusselt number as a function of the geometric parameters were developed from extensive CFD analyses:

f ¼

 2:3833aþ0:26648 15:78 6:7268 lR þ expð6:6705aÞ Re 1000 Dh 4:3551a  1:0814 þ 100

Nuh ¼ ð0:71a



lR þ 0:289Þ Dh Nuc ¼ ð0:18a



lR þ 0:457Þ Dh

Fig. 23. Effect of the channel pitch lengths on the friction factors and Nusselt numbers.

 

0:087

0:038

ð14Þ

0:11ða0:55Þ2 0:004

Re

lR Dh

aþ0:54

Pr0:56

ð15Þ

Pr0:58

ð16Þ

  0:23ða0:74Þ2 0:004

Re

lR Dh

aþ0:56

with:

200 6 Re 6 2000 5 6 a 6 45

:

4:09 6 DlR 6 12:27 h

where lR is the length of single straight section. These correlations developed based on mean pressure loss factors and mean Nusselt numbers allow for accurate determination of the thermal-hydraulic performance and heat exchanger effectiveness. The effects of the channel pitch lengths and zigzag angles were analyzed based on the obtained numerical correlations. As shown in Fig. 23, both the friction factor and Nusselt number

Fig. 24. Effect of the channel zigzag angles on the friction factors and Nusselt numbers.

M. Chen et al. / International Journal of Heat and Mass Transfer 130 (2019) 356–367

decreased as the channel pitch lengths increased when the zigzag angle was fixed at 15°. This was expected since the extent of entrance effects per unit length was reduced as the channel pitch lengths increased, resulting in smaller friction factors and Nusselt numbers. Also, as zigzag angles increased, both the friction factors and Nusselt numbers increased. However, it can be seen from Fig. 24 that there was not significant gain in the Nusselt numbers when the zigzag pitch angle was greater than 30°. Therefore, further increasing the zigzag angle beyond 30° is not recommended. 5. Conclusions In this paper, a three-dimensional CFD study of a simplified two-channel geometry model for a fabricated zigzag-channel PCHE was carried out. Comparisons of the numerical results with the experimental data showed good agreement in the channel Nusselt numbers and pressure loss factors. The relatively minor discrepancies could be attributed to the lack of the detailed geometrical measurements of the flow channels in the PCHE after the diffusion bonding process and different thermal boundary conditions used in the simulations. It was observed that the wall temperatures were not uniform along the azimuthal direction of a channel cross section. In addition, the helium temperature distribution for each segment of a cross-sectional area of the channel presented a wavy shape along the flow direction. However, the global helium temperature distributions along the entire channels were approximately linear. For the heat flux distributions, although they were significantly different on the channel walls of different segments, the heat flux for each segment along the fluid flow direction was similar. Local thermal-hydraulic performance analyses also indicated that a fully-developed flow condition was never observed in the zigzag channels, due to the nature of zigzag channels, leading to periodic flow disturbance at each bend. Several effects on the thermal-hydraulic performance of the PCHE were studied, including the fluid and solid thermophysical properties, radius of curvature at zigzag bends, channel configuration, channel pitch length, and zigzag angle. It was confirmed that the differences between the results obtained from two different gases of similar Prandtl numbers were relatively small. It was also confirmed that the mean pressure loss and mean Nusselt number in the zigzag channel with the sharp bends were larger than those in the zigzag channels with rounded bends. The effects of channel pitch and zigzag pitch angle were analyzed based on the available correlations in the literature. Both the friction factor and Nusselt number decreased as the channel pitch lengths increased when the zigzag angle was fixed. It was observed that little gain was obtained in the Nusselt numbers when the zigzag pitch angles were greater than 30°. Therefore, further increasing the zigzag angle beyond 30° is not recommended for high-temperature helium-to-helium applications. Conflict of interest The authors declare that there is no conflict of interest. Acknowledgements This material is based upon work supported by the U.S. Department of Energy Office of Nuclear Energy’s Nuclear Energy University Program under Award Number DE-NE0008714. Disclaimer This paper was prepared as an account of work sponsored by an agency of the United States Government. Neither the United States

367

Government nor any agency thereof, nor any of their employees, makes any warranty, express or implied, or assumes any legal liability or responsibility for the accuracy, completeness, or usefulness of any information, apparatus, product, or process disclosed, or represents that its use would not infringe privately owned rights. Reference herein to any specific commercial product, process, or service by trade name, trademark, manufacturer, or otherwise does not necessarily constitute or imply its endorsement, recommendation, or favoring by the United States Government or any agency thereof. The views and opinions of authors expressed herein do not necessarily state or reflect those of the United States Government or any agency thereof. References [1] H. Li, Y. Zhang, L. Zhang, M. Yao, A. Kruizenga, M. Anderson, PDF-based modeling on the turbulent convection heat transfer of supercritical CO2 in the printed circuit heat exchangers for the supercritical CO2 Brayton cycle, Int. J. Heat Mass Transf. 98 (2016) 204–218. [2] J. Nestell, T.-L. Sham, ASME Code Considerations for the Compact Heat Exchanger, Oak Ridge National Laboratory, ORNL/TM-2015/401, 2015. [3] S.K. Mylavarapu, Design, Fabrication, Performance Testing, and Modeling of Diffusion Bonded Compact Heat Exchangers in a High-Temperature Helium Test Facility Ph.D. Dissertation, The Ohio State University, 2011. [4] R.K. Shah, D.P. Sekulic, Fundamentals of Heat Exchanger Design, John Wiley & Sons, 2003. [5] S.K. Mylavarapu, X. Sun, R.N. Christensen, R.R. Unocic, R.E. Glosup, M.W. Patterson, Fabrication and design aspects of high-temperature diffusion bonded compact heat exchangers, Nucl. Eng. Des. 249 (2012) 49–56. [6] Heatric, https://www.heatric.com/about.html (accessed in September 2017). [7] M. Chen, Performance Testing and Modeling of Printed Circuit Heat Exchangers for Advanced Nuclear Reactor Applications Ph.D. Dissertation, University of Michigan, 2018. [8] M. Chen, Design, Fabrication, Testing, and Modeling of a High-Temperature Printed Circuit Heat Exchanger Master Thesis, The Ohio State University, 2015. [9] W. Kim, Y.J. Baik, S. Jeon, D. Jeon, C. Byon, A mathematical correlation for predicting the thermal performance of cross, parallel, and counterflow PCHEs, Int. J. Heat Mass Transf. 106 (2017) 1294–1302. [10] Y.J. Baik, S. Jeon, B. Kim, D. Jeon, C. Byon, Heat transfer performance of wavychanneled PCHEs and the effects of waviness factors, Int. J. Heat Mass Transf. 114 (2017) 809–815. [11] W. Chu, X. Li, T. Ma, Y. Chen, Q. Wang, Experimental investigation on sCO2water heat transfer characteristics in a printed circuit heat exchanger with straight channels, Int. J. Heat Mass Transf. 113 (2017) 184–194. [12] X. Cui, J. Guo, X. Huai, K. Cheng, H. Zhang, M. Xiang, Numerical study on novel airfoil fins for printed circuit heat exchanger using supercritical CO2, Int. J. Heat Mass Transf. 121 (2018) 354–366. [13] M. Chen, X. Sun, R.N. Christensen, I. Skavdahl, V. Utgikar, P. Sabharwall, Pressure drop and heat transfer characteristics of a high-temperature printed circuit heat exchanger, Appl. Therm. Eng. 108 (2016) 1409–1417. [14] M. Chen, X. Sun, R.N. Christensen, I. Skavdahl, V. Utgikar, P. Sabharwall, Dynamic behavior of a high-temperature printed circuit heat exchanger: numerical modeling and experimental investigation, Appl. Therm. Eng. 135 (2018) 246–256. [15] N. Bartel, M. Chen, V.P. Utgikar, X. Sun, I.H. Kim, R.N. Christensen, P. Sabharwall, Comparative analysis of compact heat exchangers for application as the intermediate heat exchanger for advanced nuclear reactors, Ann. Nucl. Eng. 81 (2015) 143–149. [16] S.J. Yoon, J. O’Brien, M. Chen, P. Sabharwall, X. Sun, Development and validation of Nusselt number and friction factor correlations for laminar flow in semi-circular zigzag channel of printed circuit heat exchanger, Appl. Therm. Eng. 123 (2017) 1327–1344. [17] A.M. Aneesh, A. Sharma, A. Srivastava, P. Chaudhury, Effects of wavy channel configurations on thermal-hydraulic characteristics of printed circuit heat exchanger (PCHE), Int. J. Heat Mass Transf. 118 (2018) 304–315. [18] Special Metals, www.specialmetals.com/documents/Inconel%20alloy%20617. pdf (accessed on July 15, 2017). [19] NIST Chemistry WebBook, http://webbook.nist.gov/chemistry/ (accessed on June 15, 2017). [20] M. Chen, J.E. O’Brien, S.J. Yoon, R.N. Christensen, P. Sabharwall, X. Sun, Thermal-hydraulic performance of a high-temperature zigzag-channel printed circuit heat exchanger, in: Proceedings of the 8th International Topical Meeting on High Temperature Reactor Technology, November 6–10, 2016. [21] I. Celik, Procedure for Estimation and Reporting of Discretization Error in CFD Applications, Report, West Virginia University, 2000. [22] F.P. Incropera, D.P. Dewitt, T.L. Bergma, A.S. Lavine, Introduction to Heat Transfer, John Wiley & Sons, 2007. [23] M. Chen, X. Sun, R.N. Christense, Numerical investigation of thermal boundary conditions of a high-temperature PCHE with zigzag flow channels, in: Proceedings of the 2018 International Congress on Advances in Nuclear Power Plants, April 8–11, 2018.