Abstract
Ice slurry is an advanced secondary refrigerant that has been attracting considerable attention for the past decade due to the growing concerns regarding energy shortage and environmental protection. To stimulate the potential applications of ice slurry, the corresponding pressure drop of this refrigerant must be comprehensively investigated. The flow of ice slurry is a complex phenomenon that is affected by various parameters, including flow velocity, ice particle size, and ice mass fraction. To predict the pressure drop of ice slurry flow in pipes, a mixture computational fluid dynamic model was adopted to simulate a two-phase flow without considering ice melting. The numerical calculations were performed on a wide range of six ice particle sizes (0.1, 0.3, 0.5, 0.75, 1, and 1.2 mm) and ice mass fraction ranging within 5%–20% in the laminar range of ice slurry flow. The numerical model was validated using experimental data. Results showed that the ice volumetric loading and flow velocity have a direct effect on pressure drop; it increases with the increase in volumetric concentration and flow velocity. The findings also confirmed that for constant ice mass fraction and flow velocity, the pressure drop is directly and inversely related to the particle and pipe diameters, respectively. Moreover, the rise in pressure drop is more significant for large ice particle diameter in comparison to smaller size ice particles at high values of ice concentration and flow velocity.
1. Introduction
The acceleration of global warming exerts a significant effect on environmental sustainability and is directly linked with increased greenhouse gas emissions and power consumption [1]. The building sector, which is one of the main energy consumption sectors, accounts for approximately 40% of the total energy consumption and carbon dioxide emission. A large portion of this percentage can be attributed to heating, ventilation, air conditioning, and refrigeration (HVAC and R) systems, which are designed for human comfort [2, 3]. Therefore, the development of sustainable and environmental friendly materials is crucial to improve the energy efficiency of HVAC and R systems [4]. Phase change materials (PCMs) are potential substitutes for such systems due to their enhanced energy transport density and cold energy storage medium [5–7].
Ice slurry, a type of PCM, is a mixture of ice particles dispersed in aqueous solution. This material is considered as an environmental friendly secondary refrigerant [8, 9]. The size of ice particles ranges from 1 mm to 1.2 mm in diameter, so it can be pumped through pipes or stored in tanks [10]. Owing to the combined effect of latent and sensible heats, ice slurry possesses high energy density, excellent heat transfer capabilities, and fast cooling rates [11, 12]. Consequently, this refrigerant has become a promising alternative to the traditional single-phase fluid used in cold energy storages and cold chains.
However, the implementation of ice slurry in engineering applications is not as straightforward as anticipated because of the numerous complexities involved [13, 14]. One of the main problems is transporting the ice slurry in the cooling network to fulfil the cold requirements. Ice concentration and particle size are the critical parameters that significantly affect the ability to transport ice slurry through piping without particle aggregation, pipe plugging, and large pressure drop. If the values of above mentioned parameters are optimised, ice slurry can be transported as a single-phase flow (i.e., without any complication) in a wide range of conditions.
The pressure drop in the two-phase flow of ice slurry is an important consideration in engineering applications that has been extensively studied by many researchers [15–18]. For ice slurry flow in pipes, the pressure drop initially increases up to 10% due to flow velocity and ice mass fraction and eventually increases because of the ice mass fraction [19]. Bel et al. [20] obtained the same inference and stated that the frictional pressure drop for narrow tubes increases with the increase in ice concentration and the rate of augmentation is high in the low Reynolds number (Re) range [20]. Pressure drop increases due to the increase in flow velocity and ice concentration [18, 21–23], and the sensitivity of such drop to the variation in both factors is strongly influenced by the flow regime [24]. On the contrary, Knodal et al. [17] and Lieu et al. [25] reported that the pressure drop in low velocity regions decreased with the increase in ice concentration.
Related experiments rarely provided comprehensive information on the velocity distribution, ice concentration distribution, and particle size of ice slurry in pipes due to technical difficulties. Given the advancement in the numerical techniques and the complexity involved in experimental setups, computational fluid dynamics (CFD) emerge as an excellent tool for obtaining extensive knowledge on two-phase flows, particularly ice slurry flow and accurate pressure drop estimation in pipes. Some researchers developed numerical modelling approaches to simulate ice slurry flow in pipes [26–29].
Wang et al. [30] adopted the Eulerian–Eulerian model to evaluate the pressure drop of an isothermal ice slurry flow under various pipe configurations (e.g., vertical, horizontal, and elbow pipes). The results confirmed that the pressure drop in all configuration increases with ice concentration and flow velocity. Similarly, Tian et al. [31] investigated the effect of ice concentration, pipe, and ice particle size on safe transportation of ice slurry by using the Eulerian–Eulerian multiphase model. The findings revealed that the safe transportation of ice slurry in cooling systems can be achieved at an ice concentration of up to 20% and ice particle size of 100 μm [31].
Niezgoda-Zelasko and Zalewski [32] combined the mixture multiphase and Eulerian–Eulerian models to investigate the pressure drop in the laminar flow range of a horizontal pipe with a pipe diameter of 16 mm. The Eulerian–Eulerian model underestimates the pressure drop, whereas the former provides a reasonably correct description of the pressure drop with respect to experimental results. The maximum percentage error of the numerical predictions obtained through the mixture multiphase model of pressure drop is approximately 15% of the experimental results [33].
Within the available literature, however, there is not much work that addresses the application of a numerical model to describe the ice slurry flow [18]. Therefore, if ice slurry is categorized as solid-liquid flow, the mixture multiphase CFD model is appropriate [34]. It is reported that the mixture method has been extensively applied to ice slurry horizontal slurry pipe flows and to water and kerosene mixtures, which have density ratio nearly identical to ice slurry for a wide range of particle sizes (ranging from 0.1 mm to 2 mm) [32, 35–38]. The mixture model is recommended for multiphase flow with small differences in density between the phases [34]. Given that the particle relaxation time is short, the mixture model formulation can also be successfully applied to the two-phase flow with a significant density difference [39].
Zhang and Shi [7] used the Eulerian–Eulerian model to investigate the isothermal ice slurry flow under various ice particle diameters (0.1 mm and 0.27 mm). However, they only considered the ice particle velocity and solid particle distribution and disregarded the particle-dependent pressure drop characteristics. Moreover, they did not investigate the ice particle size-dependent flow characteristics of ice slurry in different pipe diameters. Clarifying the influence of ice particle diameter on pressure drop and flow structure is important, especially for choosing a suitable ice generation technique.
On the basis of the discussion above, the studies on the particle size-dependent pressure drop characteristics are lacking. Moreover, previous studies only considered limited particle sizes and disregarded the particle-based pressure drop characteristics of ice slurry flow. To demonstrate the effect of ice particle diameter on pressure drop, the present study aims to analyse the pressure drop estimation of isothermal ice slurry flow based on pipe diameter and particle size using an aqueous solution of 10.3% ethanol in horizontal pipes with diameters of 9 mm and 23 mm. To the best of the authors’ knowledge, no previous study has investigated the particle size-dependent flow characteristics of isothermal ice slurry within such dimeter range. The mixture CFD model is applied to describe the pressure drop characteristics of isothermal ice slurry flow without considering ice melting. Simulations are performed considering the laminar range of the ice slurry flow for six different particle diameters (0.1-1.2 mm). The ice mass fraction for each particle diameter ranges within 5%–20%.
2. Problem Description
This study aims to investigate the pressure drop and flow characteristics of isothermal ice slurry flowing through horizontal pipes with diameters of 9 mm and 23 mm by using a mixture multiphase model. Ice slurry is a mixture of ice and aqueous solution of the carrier liquid (10.3% ethanol–water). Simulations were performed for six different particle diameters (0.1-1.2 mm) in the laminar regime. An ice mass fraction of 5%–20% was set for each particle diameter. The effects of various parameters, including flow velocity, ice mass fraction, and ice particle diameter, on the pressure drop characteristics of the isothermal ice slurry were analysed using numerical calculations. The mathematical modelling and property calculation correlations for ice slurry are discussed in the following subsections.
2.1. Mathematical Modelling
The mixture multiphase model is adopted to simulate the isothermal ice slurry flow. This model is a simplified model that utilises mixture properties to capture the two-phase flows, in which the liquid and solid phases move at different velocities but assume local equilibrium over short spatial length scales. An algebraic slip formulation is used to describe the relative velocity between the liquid and solid phases, and the numerical solution is reduced to a set of equations for the mixture phase, which comprises continuity, momentum, and particle volume fraction equations for the secondary phase. Given the small number of equations involved in the mixture model, the stability and convergence of the numerical solution are encouraging which greatly reduces the computational efforts. The fundamental governing equations for the mixture model are given below. The continuity equation for the mixture multiphase model is expressed aswhere denotes the density of mixture, and represents the local velocity. The momentum equation for the mixture model is written aswhere subscript denotes the phase ( represents the liquid phase and solid phase, respectively), represents mixture viscosity, is the ice volume fraction, and is the drift velocity defined as the ith phase velocity relative to the mixture phase. The drift velocity is defined by the relationship proposed by Taivassalo and Kallio [34].where is the slip velocity, which is defined using an algebraic slip formulation introduced by Schiller [40].
The drag function proposed by Schiller [40] is expressed aswhere denotes the particle diameter, and represents the Re of the particle.
The viscosity and mixture density in the mixture multiphase model are respectively defined as
From the continuity equation, the volume fraction equation for dispersed phase is expressed as
The boundary conditions for the inlet, outlet, and pipe walls in mixture and dispersed phases are as follows: (1) uniform velocity is applied to the pipe inlet in the mixture phase, and ice volume fraction is adopted for the dispersed phase, (2) atmospheric boundary condition is adopted for the outlet in mixture and dispersed phases, and (3) a no-slip boundary condition is applied to the pipe wall in mixture and dispersed phases. The influence of temperature change on the flow is excluded due to the isothermal conditions (Figure 1).

2.2. Physical Properties of the Ice Slurry
The corelations for the determination of the physical properties of the ice slurry are obtained from the handbook of ice slurry [41].where and , respectively, denote the ice volume fraction and density of the ice particles, and are the densities of the carrier liquid, and and represent the viscosity of the ice slurry and carrier liquid (10.3% ethanol–water solution), respectively. Table 1 presents the physical properties of the carrier liquid and ice particles calculated at a phase equilibrium temperature of 268.65 K.
3. Numerical Details
The commercial CFD software ANSYS Fluent (version 13) was used to investigate the pressure drop of the isothermal ice slurry in a horizontal pipe with diameters of 9 mm and 23 mm. The pipe length was set according to the relation L ≥ 100 D to attain the fully developed flow. Variable density-based meshes were generated using the ICEM software; meshes near the wall are finer than the far ones. A mesh independence study was performed to validate the high precision of the results using minimum computational resources by considering three different levels of mesh refinement and varying these levels by a factor of three (total number of elements for D = 23 mm: = 67,155, = 188,889, and = 284,138; for D = 9 mm: = 109,020, = 301,500, and = 603,980). The marginal difference in the numerical calculations of the pressure drop was evident for all grid sizes. Therefore, grid size was adopted in the subsequent simulations to save computation time (Figures 2 and 3).


All governing equations were discretised using the finite volume approach, where the second-order upwind and SIMPLEC schemes were used for convective diffusive terms and pressure velocity coupling, respectively. The convergence requirement for all variables was fixed to 10−4.
4. Results and Discussion
To confirm the accuracy of the mixture CFD model, the numerical results are compared with the experimental results obtained by Grozdek et al. [21]. The validation is conducted under a wide range of flow velocities and ice mass fractions. The validation is conducted over laminar range of ice slurry flow . For a two-phase flow of ice slurry, the flow pattern transition from laminar to turbulent can be observed within the range of Reynolds number from 2000 to 2300 depending upon the pipe diameter and ice mass fraction [21, 32]. It should be noted that for constant flow velocity and pipe diameter, the Reynolds number decreases as the ice mass fraction increases. Since viscosity increases and mixture density decreases with the increase in ice mass fraction. The flow behaviour was evaluated on the basis of Reynolds number given bywhere is the Reynolds number of mixture, and the mixture density is , and the mixture viscosity was calculated from equations (9) and (10).
A pipe diameter of 9 mm is considered to match the experimental conditions in the study of Grozdek et al. [21]. As shown in Figures 4(a)–4(c), the mixture CFD model provides a reasonable estimation of the pressure drop of isothermal ice slurry flow. All numerical predictions by the CFD model are in close agreement with the experimental data, and the relative errors for all considered cases are restricted to 15%. However, at high ice mass fraction of 10% and specifically 15%, the general concavity of the graph is different than those of at ice mass fraction of 5%. This could be well explained by the phenomenon of laminarization. The high concentration of ice delays the transition from laminar to turbulent flow. As the concentration of ice rises, the critical number of Reynolds where transition happens also increases and can reach even higher values than = 2300. These results are consistent with the other findings [32, 42].

(a)

(b)

(c)
4.1. Pressure Drop
The pressure drop for ice slurry flow is an important characteristic because it reflects the pumping power requirement. In this section, the influence of ice particle diameter on the pressure drop characteristics of the isothermal ice slurry flowing through a horizontal pipe with diameters of 9 mm and 23 mm is discussed. The results have been discussed in the laminar range of flow for both pipes evaluated on the basis of bulk Reynold number. The maximum Reynolds number is 1798 and 1998 for pipe diameter 23 mm and 9 mm, respectively.
The results are presented in terms of pressure gradient as a function of variable ice particle size for a wide range of flow velocity and ice mass fraction. For a given ice concentration, the pressure drop per unit length increases with the increase in ice particle diameter for all flow velocities (Figure 4). However, at a fixed ice concentration, the increase in pressure drop is high for large ice particle diameters at high flow velocity. For low velocities (velocities approximately below 0.3592 m/s), the effect of ice particle diameter on the pressure drop is small.
For small particle diameters (0.1–0.3 mm), the rise in the pressure drop at low velocities is marginal. As the ice particle size increases from 0.3 mm to 1.2 mm, the effect of particle diameter on pressure drop becomes increasingly pronounced (Figures 5(a)–5(d)). For instance, at fixed ice mass fraction (5%) and flow velocity (0.1796 m/s), the percentage increase in pressure drop as the particle diameter increases from 0.1 mm to 0.3 mm is 10% (Figure 5(a)). When the particle diameter is 1.2 mm, the percentage increase in pressure drop is 37%, which is approximately 3.5% of the original value. For the same ice concentration and flow velocity of 0.5389 m/s, the percentage increase in pressure drop when the particle diameter increases from 0.1 mm to 0.3 mm is 13%; the increment when the particle diameter is 1.2 mm is approximately 76%, which is almost 5.5% of the pressure drop at 0.1 mm. Large ice particle size leads to a larger drag force which causes larger relative velocity of ice particles. Because of this, larger lift force will act on the ice particle which will increase its flotation tendency. Consequently, the heterogeneity of ice particles distribution will increase with the increase in the size of the particle, thereby resulting in additional frictional losses. Furthermore, the effect of ice particle reinforcement becomes increasingly evident when the particle diameter is large [30].

(a)

(b)

(c)

(d)
At a fixed velocity, the pressure drop of slurry with high volume fraction is large. For instance, at a constant velocity of 0.1796 m/s, the pressure drop increases by an average of 2%, 3%, and 4% for each additional 5% of ice concentration as the particle size increases from 0.1 mm to 1.2 mm; the increment at an ice concentration of 5% is 1% (Figures 5(a)–5(d)). The reason for this phenomenon is that the addition of more ice particles, i.e., increasing ice mass fraction, results in an increased viscosity of ice slurry, which alters the density of the mixture as well as the number of Reynolds (results in decrease of the number of Reynolds, which leads to the accumulation of ice particles in the walls, and raises the potential of ice blockage, which increases resistance to flow) [15, 43]. This affects the safe transport of ice slurry in pipes at high ice concentrations.
The pressure drop per unit length is illustrated as a function of ice particle diameter (pipe diameter = 23 mm) in Figures 6(a)–6(d). Figure 6 demonstrates that the pressure gradient increases with the increase in ice particle size. However, the pressure drop at a pipe diameter of 23 mm is lower than that at smaller ones. For a given velocity, the pressure drop increases with the decreases in pipe diameter and vice versa. Liu et al. [44] concluded that the ice particles in pipes with large diameters provide a coherent structure to the fluid and reduce the amount of energy dissipated by the viscous frictions. On the basis of the results, the pressure drop associated with ice slurry pipe flow is a function of particle size. Ice slurry with large particle size that flows through pipes with small diameters yields unwanted increments in pressure drop and maximises the pumping power penalties.

(a)

(b)

(c)

(d)
4.2. Velocity Contours
Figures 7–9 show the contours of the velocity distribution of the solid phase at the perpendicular cross-section in the outlet of the 9 mm pipe. The contours correspond to varying inlet velocities (0.1796, 0.3592, and 0.5389 m/s) at different ice particle diameters (0.1, 0.75, and 1.2 mm) and constant ice mass fraction of 10%.

(a)

(b)

(c)

(a)

(b)

(c)

(a)

(b)

(c)
The velocity contours indicate that the local velocity along the pipe wall is minimal and almost stagnant because viscous forces are dominating the vicinity of the wall. Moreover, when the ice slurry flows in the low velocity region with high ice concentration, the ice particles cluster together at the walls, and the value of the velocity increases towards the centre of the pipe. The Figures 7–9 depict that given a constant inlet velocity, the local velocity value increases with the increase in ice particle size towards the centre of the pipe. In addition, the velocity distribution becomes increasingly homogeneous as the inlet velocity increases when the volume fraction of ice is constant [43]. Furthermore, no difference is observed in the velocity distribution between the solid and liquid phases of all cases.
4.3. Velocity and Concentration Profile
One way to assess the flow patterns of ice slurry flow is to examine the velocity distribution in the cross-section of the flow. The influence of ice particle diameter on the mixture velocity distribution profiles in the vertical diameter direction along the horizontal cross-section of the pipe outlet is illustrated in Figures 10 and 11. The mixture velocity distributions of three different particle diameters (0.1, 0.75, and 1.2 mm) at ice mass fractions of 5% and 10% are plotted under different inlet velocities (pipe diameter = 9 mm).

(a)

(b)

(a)

(b)
Figures 10 and 11 show mixture velocity profile distribution in the cross-section of the flow at distance L = 0.5 (mid-length) and at the outlet cross-section of the pipe. It has been observed that at around 0.5 m distance from the entrance of the pipe, the velocity profile is still developing. However, after this length, the ice slurry flow is fully developed, and the velocity profile is parabolic in nature. Hence, the entrance length is around 0.5 m for the present case, and the ice slurry flow is fully developed after this distance (data after L = 0.5 m showed the fully developed flow and are not included to avoid complexity). The velocity fields for all particle size subjected to the present study are symmetrical with respect to the vertical diameter, and the extremes values of velocity occur on the centre of diameter.
Similarly, the ice concentration distribution along the vertical diameter for pipe diameter = 9 mm at constant ice mass fraction of 10%, for different inlet velocities, and ice particle diameter of 0.1 mm, 0.75 mm, and 1.2 mm are presented in Figure 12. It is clear from the figure that for constant ice volume fraction, heterogeneity of ice slurry flow increases with the increase in ice particle diameter and flow velocity. The larger particles exert larger drag force, which in turn results in the larger relative velocity and increases the heterogeneity of ice particles.

(a)

(b)
5. Conclusion
The study presented the CFD modelling of the two-phase flow of isothermal ice slurry through horizontal pipes with diameters of 9 mm and 23 mm within the laminar range. The pressure drop characteristics are evaluated under six ice particle sizes (0.1, 0.3, 0.5, 0.75, 1, and 1.2 mm) at an ice mass fraction range of 5%–20%. The results of the present study are validated by the experimental findings. The pressure drop predictions are in excellent agreement with the experimental results obtained by Grozdek et al. [21] for a wide range of inlet velocity and ice mass fraction under the laminar regime. The following conclusions are derived from the numerical calculations.(i)The pressure drop per unit length significantly increases with the increase in ice mass fraction and flow velocity, and the pressure loss rapidly rises at high ice concentration and flow velocity.(ii)Pressure drop per unit length increases with the increase in ice particle size and decreases with the increase in pipe diameter.(iii)The increase in the pressure drop per unit length of ice particles with large diameters at high flow velocity and ice concentration is significant. The rate of increase for the particle with a diameter of 1.2 mm is approximately 3–5 times higher than that of the particle with a diameter of 0.1 mm as the flow velocity increases from 0.1796 m/s to 0.5381 m/s. Moreover, the pressure drop increases by an average of 2%, 3%, and 4% for every 5% increment in ice concentration; the corresponding drop is only 1% when the ice concentration is 5%.(iv)The local velocity near the pipe wall is minimal and almost stagnant due to the dominance of the viscous forces near the wall region. This value, however, increases towards the centre of the pipe.
Results suggest that pipes with small diameters are not recommended to transport ice slurry because of their higher pressure drop than those with large diameters. In addition, given that a large particle size results in undesirable increments in pressure drop, the size of the ice particle must be kept small to reduce the pumping power penalties. In the future, the effect of ice particle size on the performance of ice slurry should be explored using various performance evaluation criteria.
Nomenclature
: | Pipe length, m |
: | Reynolds number |
: | Particle diameter, m |
: | Pipe diameter, m |
: | Inlet velocity, m/s |
: | Gravitational acceleration, m/s2 |
: | Temperature, K |
: | Drag force |
: | Drift velocity, m/s |
: | Slip velocity, m/s |
: | Phase index |
: | Pressure, Pa |
: | Viscosity, Pa·s |
: | Density, kg/m3 |
: | Volume fraction |
: | Liquid phase, aqueous solution |
: | Solid phase, ice particles |
: | Mixture, ice slurry |
: | Interaction between the solid and liquid phases |
: | Particle. |
Data Availability
The data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This study was supported by the National Research Foundation of Korea and funded by the Korean government (MSIP, 2020R1A2B5B02002512 and 2020R1A4A1018652).