Analytical model of solar energy storage using non—Newtonian Fluid in a saturated porous media in fully developed region: carboxymethyl cellulose (CMC) and graphite model

Abbreviations: A: Constant = $ {\left[\frac{-K}{\mu }\frac{\partial p}{\partial z}\right]}^{\frac{1}{n}}\frac{d{T}_{w}}{dz} $ ; CMC: Carboxymethyl cellulose; FDM: Finite-difference methods; $ {{C}_{p}}_{f} $ : Specific heat of the fluid at constant pressure, $ J⁄(Kg.K) $ ; $ {C}_{t} $ : The tortuosity factor given by Christopher and Middleman which equal 25/12; $ h $ : Convective heat transfer coefficient, $ W⁄({m}^{2}.K) $ ; $ \mathrm{H} $ : Consistency index, $ Pa.{S}^{n} $ ; $ K $ : Permeability of porous media, $ {m}^{2} $ ; $ {K}^{*} $ : Modified permeability, $ {m}^{1+n} $ ; $ k $ : Thermal conductivity, $ W/(m.K) $ ; $ {k}_{f} $ : Thermal conductivity of fluid, $ W/(m.K) $ ; $ {k}_{s} $ : Thermal conductivity of solid, $ W/(m.K) $ ; $ {k}_{eff} $ : Effective thermal conductivity, $ W/(m.K) $ ; $ \dot{m} $ : Mass flow rate, $ kg/s $ ; $ n $ : Power law index; $ \mathrm{N}\mathrm{u} $ : Nusselt number; $ {\mathrm{N}\mathrm{u}}_{D} $ : Local Nusselt number; $ {\mathrm{N}\mathrm{u}}_{m} $ : Mean Nusselt number; $ Pr $ : Prandtl number; $ ∆P/L $ : Pressure gradient, $ Pa/m $ ; $ Q $ : Amount of heat stored, Joule; $ {q}^{\text{'}\text{'}} $ : Heat flux per unit area, $ W/{m}^{2} $ ; $ Q $ : The amount of heat stored, J; $ r $ : Radius of the conduit, $ m $ ; $ {r}^{*} $ : Dimensionless radius of the conduit, $ m $ ; $ R $ : Outer radius of the conduit, $ m $ ; $ {Re}_{D} $ : Reynolds number based on Darcy velocity; $ {Re}_{K} $ : Reynolds number based on Permeability; $ T $ : Temperature, $ K $ ; $ {T}_{e} $ : Temperature at the entrance, $ K $ ; $ {T}_{w} $ : Temperature at the wall of the conduit, $ K $ ; $ {T}_{m} $ : The mean temperature, $ K $ ; $ {T}^{*} $ : Dimensionless temperature; $ u $ : Velocity in z-direction, $ m/s $ ; $ {u}_{D} $ : Darcy velocity, $ m/s $ ; $ {u}_{m} $ : The mean velocity, $ m/s $ ; $ z $ : Length of conduit, $ m $ ; $ {z}^{*} $ : Dimensionless length; $ \mu $ : Dynamic viscosity, $ Pa.s $ ; $ {\mu }_{eff} $ : Effective viscosity, $ Pa.{S}^{n} $ ; $ \varepsilon $ : Porosity of porous media; $ {\alpha }_{eff} $ : Effective thermal diffusivity, $ {m}^{2}/s $ ; $ {\rho }_{f} $ : Fluid density, $ kg/{m}^{3} $ ; $ {\rho }_{s} $ : Solid density, $ kg/{m}^{3} $ ; $ \nu $ : Kinematic viscosity, $ {m}^{2}/s $ ; $ \upsilon $ : The velocity at radial (r) direction, $ m/s $

1.   Introduction

The development of efficient energy storage systems is crucial because the ongoing progressive shift away from traditional baseload energy sources and renewables gain a greater foothold in the energy system. Solar energy is the cheapest form of renewable energy, as well as the cleanest, delivering exceptional benefits for both the planet and for businesses. Yet, in order to maximize the return on solar requires an efficient way to store that energy and flexibility for its future use.

The thermal energy storage method which uses the absorption and discharge effects of heat caused by the temperature change of materials is called sensible heat thermal energy storage (SHTES).

Thermal energy conservation and management is a key issue for energy savings in any field with intensive thermal demands.

The storage of sensible heat is one of the best-known and most widespread technologies, with the domestic hot water tank as an example. The storage medium may be a liquid such as water or thermo-oil, or a solid such as concrete or the ground. Thermal energy is stored solely through a change of temperature of the storage medium. The capacity of a storage system is defined by the specific heat capacity and the mass of the medium used ^[1]. Porous materials are of considerable scientific and technological interest for energy storage and conversion ^[2]. Porous material can improve the efficiency of solar thermal energy storage systems. Because it has a high dispersal so that enhances the heat convection. Another benefit for using porous media is the irregular movement of the flow around the beads or through pores which mixes the fluid more effectively.

Dhifaoui et al. (2007) ^[3] studied the storage of thermal energy in a porous bed in a vertical channel, filled with glass beads, heated on one of the vertical walls by a constant heat flux. As a result, the decreases in discharge time led to increase in the system efficiency (the ratio of the stored energy over the energy given) for a fixed channel width and a fixed charge duration and increases in the volume of storage. A parametric study to investigate the effect of the bed dimensions, fluid flow rate, particle material and diameter on the thermal behavior of the packed bed up to 550 ℃ was performed by Hänchen et al. (2007) ^[4].

The behavior of non-Newtonian fluid as flow and heat transfer inside a saturated porous media have been investigated for centuries in order to improve the heat transfer rate, because of their wide applications in medicine, mechanical, industrial, chemical, and civil engineering, such as combustion technology, fuel cells, vehicle's radiator, heat exchanger in building, oil and gas flow in reservoirs, Solar energy collections etc. ^[5].

Wu and Zhao (2011) ^[6] established an experimental study of using porous materials such as metal foams and expanded graphite to enhance the heat transfer capability of PCMs in high- temperature thermal energy storage system. NaNO₃ was selected with 99% purity as the PCM for the experimental investigation.

The result has been reached that the heat transfer rate can be enhanced through the addition of the porous materials by 2.5 times compared to that of pure NaNO₃ in the heating process from 250–300 ℃.

Chamkha et al. (2014) ^[7] analyzed the effect of uniform lateral mass flux on non-Darcy natural convection of non-Newtonian fluid along a vertical cone embedded in a porous medium filled with a nanofluid. Dimensional governing equations are transformed into a non-similar form and then solved numerically by Keller box finite-difference method.

The effect of power-law index parameter n on the dimensionless velocity, temperature and, volume fraction profiles are displayed. The result that has been reached as the power-law index n is increased, both the volume fraction and temperature profiles are increased, but the slip velocity at the cone surface is decreased. Also, the non-Newtonian fluid with a higher power-law index has a greater thermal boundary layer thickness and a greater volume fraction boundary layer thickness.

In this study a theoretical model is built and solved analytically in order to describe the non-Newtonian fluid, (CMC), and Graphite as porous media.

Graphite with highly advanced properties can be applied in various fields, such as thermal transfer devices, electrical transfer devices and electrochemical super capacitors.

The higher thermal conductivity leads to a shorter time to charge or discharge. Graphite has high thermal conductivity so the time used for heat transfer inside the material will be very short, which means better system performance and that way improving the efficiency of thermal energy storage systems ^[8].

Mohammed & Mohammed, (2009) ^[9] studied Rheological Properties for two types of CMC at 30 ℃. This study used this information. Parameters for (CMC) like power index and consistency coefficient were used from this study are shown in Table 1.

The properties of Graphite were used in calculations of this study are shown in Table 2.

The properties of CMC were used in calculations of this study are shown in Table 3.

Survey points

The forced convection problem inside a saturated porous media is studied under the effects of constant surface heat flux, Darcy model as momentum equation is depicted.

In this study, a theoretical model is built and solved analytically in order to describe the non-Newtonian fluid, (CMC), and Graphite as porous media. Also, the effects of some (power index, porosity…. etc.) parameters on temperature variation were studied to show the differences between water and non-Newtonian fluid in energy storage. The study is shows that the temperature and energy storage in winter in Amman city

2.   Mathematical modeling and analysis

2.1.   Nusselt number calculations

Consider the steady two-dimensional forced convection flow of incompressible non-Newtonian fluid moving in circular conduit saturated with porous media is shown in Figure 1:

The flow is said to be fully developed when the normalized temperature profile remains unchanged and the velocity profile is fully developed and remains unchanged.

Considering steady state, incompressible, laminar, forced fluid flow, and all the physical properties of the fluid are isotropic and homogeneous. The governing equations were written as follow:

The continuity equation:

The momentum equation in z direction:

where $ {u}_{D} $ is Darcy velocity

The momentum equation in r direction:

The averaged velocity profile remains unchanged when the flow is fully developed, and thus: $ \partial u/\partial z = 0 $ , $ u = {u}_{D} $ .

The energy equation:

where K is permeability, p is the fluid pressure, μ is the dynamic viscosity, T is the fluid temperature and $ {\alpha }_{eff} $ is the effective thermal diffusivity (m²/s) which equals $ {\alpha }_{eff} = \frac{{k}_{eff}}{{\left(\rho {C}_{p}\right)}_{f}} $ , where $ {k}_{eff} = \varepsilon {k}_{f}+(1-\varepsilon ){k}_{s} $ , $ {k}_{f} $ is thermal conductivity of the fluid and $ {k}_{s} $ is thermal conductivity of the solid and $ \varepsilon $ is porosity of the porous media.

It can be integrated with respect to these boundary conditions:

Temperature profile in the fully developed region are presented by Oosthuizen & Naylor, 1999a ^[12]:

Substitute temperature profile and Darcy velocity into energy equation:

Rewrite Eq 7 by assuming the left term in constant (A):

Integration of Eq 8 twice with respect of the boundary conditions and substitute $ \frac{\partial T}{\partial r} = \frac{{q}^{\text{'}\text{'}}}{k} $ at $ r = R $ because the model is subjected by constant heat flux gives:

Now to find Nusselt number we should find mean temperature $ {T}_{m} $ to substitute it in the heat transfer coefficient:

The Nusselt number is defined in terms of the difference between the wall and the mean Temperatures $ {(T}_{w}-{T}_{m}) $ .The mean temperature, can be found as: ^[12]

Substitution Eq 9 into 11 and the result is integrated over the interval $ [0, R]. $

Finally, mean temperature represented as follows:

Substitute Eq 12 into 10 to get:

Now Nusselt number is founded by substitute Eq 13 in Nusselt number equation $ \left(\frac{hD}{k}\right) $ , where $ D = 2R $ . Noted that the heat flux is concentrated over the upper half side, so because of that, the result of Nusselt number from the upper wall to the center can be evaluated as: ^[13]

2.1.1.   Solutions validation

In the entrance region or developing region where the boundary layer develops and reaches the tube center, the developing in velocity profile is then essentially fully developed before the heat transfer take place, and because of the presence of porous media, the velocity profile is not changing with z along the conduit which means there is only a thermal entrance region.

The energy equation is written in dimensionless form and solved numerically in case of constant surface heat flux.

The velocity in $ r-direction$ was neglected; i.e., $ \upsilon $ = 0, and there is no heat transfer in 𝑟 direction at the center; i.e., $ \partial T / \partial r = 0 ~~at ~~ r = 0 $ . The diffusion of heat in axial direction will be neglected compared to that in the radial direction.

We can rewrite energy Eq 4 after substitute effective thermal diffusivity as:

The initial point of the conduit has no heat transfer so the temperature is uniform and constant and can be presented as; $ T = {T}_{e} $ at z = 0.

The following dimensionless variables are used in case of constant heat flux at the wall:

where $ {T}^{*} $ is dimensionless temperature, $ {T}_{e} $ is the temperature at the entrance point, $ {u}^{*} $ $ $ is dimensionless velocity, $ {u}_{D} $ is the Darcy velocity in the conduit, $ D = 2R $ is the diameter of the conduit, $ {z}^{*} $ is the dimensionless length and $ {r}^{*} $ is the dimensionless radius of the conduit.

By combining dimensionless variables (15) and energy Eq 14, the energy equation can be written as:

We can eliminate dimensionless velocity $ {u}^{*} $ because the velocity through porous media is constant and here is Darcy velocity which means $ u = {u}_{D} $ (i.e., $ {u}^{*} = 1 $ ).

The equation is subjected to the following dimensionless boundary conditions:

The linear second-order differential Eq 16 concerning boundary conditions (17) was solved by finite difference method using MATLAB program.

We can calculate local Nusselt number by substitute equation $ h = \frac{{q}^{\text{'}\text{'}}}{{T}_{w}-{T}_{c}} $ in $ Nu = \frac{hD}{k} $ ^[14] and get:

From Eq 15, the right term represents of $ \frac{1}{{T}^{*}} $ at wall of the conduit so Eq 18 can be written as follows:

where $ {Nu}_{D} $ the local Nusselt number and the dimensionless temperature at constant heat flux at the wall is $ {\left.{T}^{*}\right|}_{{r}^{*} = 0.5} = \frac{{T}_{w}-{T}_{e}}{{q}^{\text{'}\text{'}}D/k} $ .

The Nusselt number is defined in terms of the difference between the wall and the mean temperatures (𝑇 $ {T_w}- {T_m} $ ).

The mean temperature, 𝑇 $ {T_w} $ can be found as: ^[12]

Temperature T can be obtained from dimensionless form Eq 15:

Substitution Eq 21 into 22 and $ u = {u}_{m} $ we get:

The right-side term represents a non-dimensional mean temperature $ {{T}^{*}}_{m} $ . After substituting the dimensionless form, the result will be:

Equation 23 can be solved by MATLAB program and after that, we use the data of to find Nusselt number by this equation:

The results of these equations have been shown in the Figures 2 and 3.

The Nusselt number starts from z* = 0 exponentially and decreases through developing region along entrance length and become constant along fully developed region; both the friction factor and the heat transfer coefficient remain constant in the fully developed region and there is no dependence on the Reynolds or the Prandtl numbers.

Figure 3 represent the dimensionless temperature T* with axial direction z* at different Reynolds number (Re) and Prandtl number (Pr) respectively at selective values of the dimensionless radius, the temperature increases in the axial direction of the flow; this is due to continuous heating of the fluid in the entrance region. The figures show that increasing in Reynolds and Prandtl numbers lead to decreasing in temperature profile and this because of increasing in Re means increasing in inertia forces and increasing in heat transfer coefficient. Also, the increasing in Prandtl number are found to decrease non-Newtonian fluid temperatures and consequently enhance heat transfer rates.

2.2.   Temperature variation and energy with axial direction calculations in fully developed region

Temperature is not constant in the axial direction through the pipe flow. We can calculate the mean temperature variation in the axial direction by applying the energy balance equation lead to:

where $ {u}_{m} $ the mean velocity of the fluid, $ {c}_{p} $ is specific heat. Because the temperature gradient is independent of z in fully developed flow and thus the shape of the temperature profile does not change along the tube.

Christopher & Middleman, 1965 ^[15] developed a modified Blake-Kozeny equation for a power-law, non-Newtonian fluid with laminar flow through packed porous media and this equation is:

where u is the Darcy's velocity; $ {K}^{*} $ is modified permeability and known as: ^[16]

$ {C}_{t} $ is the tortuosity factor given by Christopher and Middleman ^[15] and it is equal $ \left(25/12\right) $ , $ △P/L $ is the pressure gradient; and $ {\mu }_{eff} $ is effective viscosity and given as:

where $ \varepsilon $ is the porosity, and $ K = \frac{{\varepsilon }^{3}{d}^{2}}{150{(1-\varepsilon )}^{2}} $ . Note that $ {\mu }_{eff} $ does not have the units of viscosity but (Pa.sⁿ).

Heat flux at the wall $ {q}_{w}^{\text{'}\text{'}} $ (W/m²) in equation was found as a function of time (t) (hour) for Amman city and the data from PVGIS website. It is taken as 4^th order to give more accurate results. General equation of heat flux in Amman is given by the following:

where $ A, B, C, D, $ and $ E $ are constants and will be found for each month in the Table 4.

The amount of stored sensible heat in a material depends on its heat capacity (energy density) and the thermal diffusivity (rate at which the heat can be released and extracted) ^[18].

where $ Q $ is the amount of heat stored, $ {T}_{i} $ is the initial temperature, $ {T}_{f} $ is the final temperature, m is the mass of heat storage medium, and $ {C}_{p } $ is the specific heat.

3.   Results and dissections

In this region, many figures represent the effect of porosity, pore diameter, and power-law index, for temperature variation through pipe.

Figures 4 and 5 illustrate a comparison between CMC-1 and CMC-2 solutions and water in porous media. The porosity, pressure and pore diameter for both types and water conduit are constant.

The mean temperature variation with axial coordinate for CMC-1 and water at different heat flux are shown in Figure 4. It is clear that the performance of CMC is better than water at the same heat flux. When the concentration of CMC-1 solution is increased from 40 g/l to 50 g/l to 60 g/l the power index (n) is decreased from 0.724 to 0.7182 to 0.7182 respectively and the temperature variation about the end of conduit increases from 15 ℃ to 30 ℃ to 70 ℃ respectively at $ 200w/{m}^{2} $ . But the maximum temperature variation of water is 40 ℃ at $ 800/{m}^{2} $ .

CMC-2 solution has a higher consistency index at the same concentration which means higher viscosity and less power index than CMC-1. So, Figure 5 shows that at the end of conduit the temperature variation at $ 200w/{m}^{2} $ for 40 g/l (CMC-2) is 20 ℃ and for 50 g/l is 80 ℃.

It can be observed that for the same type of solution when the power index is decreased, temperature variation is increased at the same conditions. The less the power index of fluid, the more the shear stress it has and the less shear strain, which means less velocity and the fluid gain more solar energy and heat.

In Figures 4b, c and 5b the temperature is more than 100 ℃ so the system can be operated at pressures that allow this temperature to be achieved without boiling, the solution or the conditions (porosity, particle diameter…...etc.) can be changed to decrease temperature variation.

Figure 6a, b show the effect of pressure on the temperature at different (n) for CMC-1 and CMC-2 respectively at constant porosity and pore diameter at $ 200w/{m}^{2} $ . It is clear as pressure is decreased the temperature is increased and that due to the decreasing velocity of the fluid. The highest temperature variation reading at 40 kPa for all power index of CMC-1 and CMC-2.

The effect of porosity on the temperature is shown in Figure 7, as shown the temperature increases with decreasing porosity and a small change in porosity make a large change in temperature.

Three values of porosity were taken 0.3, 0.25, and 0.23. The temperature is doubled when the porosity increases 0.05. The temperature of CMC-1 reached to 120 ℃ at $ n = 0.7122 $ at $ \varepsilon = 0.23 $ and 150 ℃ for CMC-2 at $ n = 0.0.57 $ at the porosity.

Figure 8 represents the effect of increasing particles or pore diameter on the temperature, where the pressure and the porosity are constant 40 kPa and 0.25 respectively. As it is clear the size of particles affects inversely of temperature which means at $ d = 0.0003 $ m the temperature is the highest value of temperature, because of decreasing fluid velocity and the amounts of solar energy stored are increased.

The radius of the conduit also affects the temperature of storage, Figure 9 shows this effect. It can be observed that the temperature of CMC-1 and $ n = 0.7122 $ is increased from 45 ℃ to 90 ℃ with decreasing the radius of the conduit to half. The other power index also behaves the same way.

The temperature variation at different heat flux in Amman was found and calculated based on the data from PVGIS web site for all months. Figures 10 and 11 represent temperature during the day t (hr) at different heat flux for different months in Amman.

The temperature variation just reaches to approximately 50 ℃ for $ n = 0.724 $ in summer Figure 10a but for $ n = 0.7182 $ and $ n = 0.7122 $ the temperature reach to more than 100 ℃ Figure 10b, c.

The temperature variation of CMC-2, $ n = 0.599 $ in summer reach to 70 ℃ Figure 11a but reach to more than 200 ℃ for $ n = 0.57 $ Figure 11b.

In January and December, the irradiance is at law value so the heating of the fluid as little as possible. Also, the time, which the temperature is higher as possible, is different between months because the solar noon is different.

Figure 12 illustrate the capacity of the storage in January and December at specific conditions and parameters shown on figures. We can change the parameter to design the system at the wanted temperature. It's clear from the picture that decreasing in power index led to increasing of the temperature.

Figure 13 illustrate a comparison of stored energy between CMC-1 and CMC-2 solutions at different power index and the stored energy per unit mass with axial coordinate for CMC-1 at low heat flux $ 200w/{m}^{2} $ . When concentration of CMC solution is increased the power index (n) is also decreased and stored energy increases at the same conditions in 2 m conduit with $ r = 0.03m $ , $ \varepsilon $ = 0.25, d = 0.00035, and $ p = 50 $ kPa. the stored energy of CMC-1 for n = 0.724, n = 0.7182, and n = 0.7122 is approximately 120 kJ, 300 kJ, and more than $ 600 $ kJ respectively.

We can observe that for the same type of solution when the power index is decreased, temperature variation is increased at the same conditions. The less the power index of fluid, the more the shear stress it has and the less shear strain, which means less velocity and the fluid gain more solar energy and heat.

The effect of porosity on the stored energy along the conduit is shown in Figure 14 as shown the stored energy increases with decreasing porosity and a small change in porosity make a large change in temperature.

Figure 15 represents the effect of increasing particles or pore diameter on the energy stored in the conduit, the size of particles effect inversely of temperature that mean at d = 0.0003 m the temperature is the highest value of temperature, because of decreasing fluid velocity which leads to decrease the amounts of solar energy stored are increased and reached to more than 150 kJ/kg along the conduit.

Figure 16 shows the effect of pressure on the stored energy at different (n) for CMC-2 and. It is clear as pressure is decreased the stored energy is increased and that due to the decreasing velocity of the fluid. The stored energy is reached to 300 kJ/kg at 40 kPa for CMC with n = 0.57.

5.   Conclusions

In this work, a model of non- Newtonian fluid—which is CMC- through graphite porous media was built and solved numerically and analytically, one can conclude that:

1.Porous media enhance the stored energy inside energy storage, especially when using porous media with high thermal conductivity.

2.The non-Newtonian fluid (CMC) shows the best temperature variation than water through porous media.

3.The Nusselt number was calculated analytically and numerically and the result was 4.

4.Effect of porosity: It is found that as porosity decreases the velocity of the fluid decreases, as a result, the temperature and energy along with the duct increases.

5.Effect of power index: For the same fluid, the power index is affected by the viscosity because that increase in power index led to a decrease in temperature and energy along the duct at the same conditions.

6.Effect of particle diameter: the size of particles affects inversely temperature and energy stored, because of decreasing fluid velocity, and the amounts of solar energy stored are increased.

7.Effect of pressure: The temperature variation and energy stored are increased when pressure decreases.

8.In winter the heat flux decreases but non-Newtonian fluid gives a good temperature variation for solar energy storage in Amman.

The following are suggestions for further investigations:

1.Study different types of porous media and non-Newtonian fluid for thermal energy storage.

2.Applying this model experimentally and compare it with theoretical study.

3.Study the non-Darcy flow by using Forchheimer equations for non-Newtonian fluid.

Conflict of interest

All authors declare no conflicts of interest in this paper.