A thorough investigation of thermochemical heat storage system from particle to bed scale

Salt hydrates are promising candidates for the long-term thermochemical heat storage (TCES) in the building environment. In such storage systems, the surplus of energy will be exploited in an endothermic reaction to dehydrate the salt hydrates. Once it is demanded, the stored energy will be released through an exothermic reaction by hydrating the salt, which results in an increase in the mass and temperature of salt particles as well as changes in the species of material. In order to construct an improved storage system, it is very important to deeply investigate the details of the (de)hydration processes in salt hydrates. Poor heat and mass transfer is the bottle neck in this technology. Therefore, the main objective of this work is to investigate how heat and mass transfer influence the (de)hydration in a closed TCESsystem. The novelty of this work is to provide a high degree of detailed information about (de)hydration of TCM (thermochemical material) in a bed by calculating transport phenomena for each single particle while considering their interactions with each other. This is achieved by applying and developing the Extended Discrete Element Method (XDEM) as a numerical modeling tool and Thermogravimetric Analysis (TGA) measurements. Comparisons are carried out for the results of the hydration and dehydration process in a single particle with the measurements which shows a very good agreement. Moreover, impact of particle size on the hydration process is also studied. Further, simulations for the hydration process in a chain of six potassium carbonate (K2CO3) particles are performed in order to understand the mechanism of heat and mass transfer inside the packed beds. 2021 The Author(s). Published by Elsevier Ltd. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).


Introduction
The loss-free shifting of renewable energy over longer periods of time in general and to be able to supply the peak energy demand can make a completely sustainable energy supply possible. To answer the question of how much renewable energy can contribute to the decarbonization of the energy sector, one needs to answer how much energy storage can be improved. In order to boost this energy transition, there is an urgent need for the development of heat storage systems. Thermal energy storage is the solution for a key bottleneck against having 100% solar heat supply for space heating and domestic hot water (DHW) in the building environment. Among the available technologies, thermochemical energy storage (Sunku Prasad et al., 2019;Müller et al., 2019;Gravogl et al., 2019) (TCES) creates opportunities for loss-free seasonal storage of heat. The working principle is based on a reversible (de)sorption reaction. Among different categories of thermochemical materials (TCM), salt hydrates are promising candidates with relatively high storage density (greater than 2 GJ/m 3 ) to be used in the building environment (Donkers et al., , 2016. In such a system, the surplus of energy will be exploited in an endothermic reaction to dehydrate the TCM (Han et al., 2018;Schmidt et al., 2017). In winter, this stored energy will be released through an exothermic reaction by hydrating the TCM (Solé et al., 2013;Lahmidi et al., 2006;Bales et al., 2005;Cot-Gores et al., 2012;N'Tsoukpoe et al., 2014).
The reactor is an important part of the TCES system, where (de) hydration of salt hydrates happen. The (de)hydration involves a coupled heat and mass transfer which determines the performance of the system like specific power, energy density, thermal efficiency, and stability. The heat transfer in the reactor bed is related to the thermal conductivity of particles and of the gas phase, heat transfer coefficient and the thermal contact between the particles. One issue related to salt hydrates is their low instinct thermal conductivity which causes low heat transport. In the case of a closed system , there is an additional thermal resistance caused by the lack of area contact between the wall of heat exchanger and particles. Material macrostructure can also limit mass transport (Scapino et al., 2017). The heat and mass transfer in such packed beds has been a subject of great research (Michel et al., 2016;Mamani et al., 2018), however, the influence of applying different sizes and shapes of the active material on heat and mass transfer is at present unclear. In order to optimize the configurations of TCMs inside heat exchangers of thermochemical storage reactors, the fundamental heat and mass transfer properties of active materials should be known. Often it is very difficult to investigate the processes inside the reactor experimentally because of the limited access inside the bed and also due to the very costly experimental setups. In contrast, computational modeling is a more reliable approach in order to deeply investigate and understand the complex physical and chemical processes within the packed beds.
Detailed reviews on different approaches that have been used for modeling and solving heat transfer in granular materials are presented in (Solé et al., 2015;Balasubramanian et al., 2010). For modeling complex thermochemical energy storage systems con-siderable achievements have been gained (Solé et al., 2015;Balasubramanian et al., 2010). Balasubramanian et al. (Balasubramanian et al., 2010) developed a mathematical model to study the energy storage ability of salt hydrates. Michel et al. (Michel et al., 2012) developed another mathematical model for the thermochemical reaction of solid and gaseous phases by considering moist airflow through the porous bed of Strontium bromide. Chen et al. (Chen et al., 2019) presented a two-dimensional pseudo-homogeneous model to simulate reaction kinetics and thermodynamics in conical as well as cylindrical ammonia dissociation reactors. Takasu et al. (Takasu et al., 2019) evaluates the performance of a TCES system based on Li 4 SiO 4 /zeolite/CO 2 for thermal energy storage. Li et al. (Li et al., 2016) described and analyzed the operating principle and working performance of the thermochemical multilevel sorption thermal battery for energy storage. Michel et al. (Michel et al., 2016) reported the local operation and reactive bed behavior of such systems. Pal et al. (Pal et al., 2014) reported study on thermal energy storage in porous materials with adsorption and desorption of moisture. Mette et al. (Mette et al., 2014) investigates the water vapor adsorption and kinetics of zeolite both experimentally and numerically. Wu et al. (Wu et al., 2009) performed numerical analysis of an open-type thermal storage system using composite sorbents.  studied the dehydration/hydration of granular beds for thermal storage applications. Narayanan et al. (Narayanan et al., 2014) developed a detailed computational model based on a governing adsorption dynamics in a single adsorption layer and pellet. Xia et al. (Xia et al., 2019) explored mass transfer in porous thermochemical heat storage materials. Chen et al. (Chen et al., 2020) performed topology optimization for heat transfer enhancement in thermochemical heat storage. Robino and Boer (Boer, R.d., Rubino, A., 2012) presented a model for open sorption reactors. Ernewein and Lorente (Malley-Ernewein and Lorente, 2019) studied the geometrical features of an open reactor considering pressure drop and temperature distribution. They found increased number of parallel layers of salt will increase the ratio between the stored heat and power. In the recent past Gaeini et al. (Gaeini et al., 2017) present a 2D model for studying moisture and heat transport in a packed bed reactor. Most recently, Risthaus et al.  (Risthaus et al., 2020) numerically analyzed the hydration of calcium oxide in a fixed bed reactor and found that the reaction kinetics at low pressures (8.7-50 kPa) is very sensitive towards pressure and temperature. Also, they concluded that the thermal losses have a significant influence and thus have to be accounted for in the models. Stengler et al. (Stengler et al., 2021) investigate the performance of thermochemical energy storage system based on strontium bromide using finite element method and discovered that heat transfer coefficient between the reactive bulk phase and the heat exchanger is the main contributor for storing maximum thermal power. Kant et al. (Kant et al., 2021) developed a 3dimensional numerical model and analyzed the performance of a K 2 CO 3 -based thermochemical energy storage system by considering a honeycomb structured heat exchanger. In a parametric study of geometrical parameters of the honeycomb heat exchanger they found that the reduction of cell size of the honeycomb provides better heat transport as well as an improved reaction rate. However, many critical phenomena such as intra particle temperature gradient, diffusion of gas into the particle and particle-particle and particle-wall interaction models have been neglected in the available literature, while they have significant effects on the heat and mass transfer and therefore the system performance of TCESsystems. This will be even more crucial when using salt hydrates as they might experience an aging effect during multi-cycling performance (Knoll et al., 2017;Barreneche et al., 2015). For that reason, it is important to resolve transport phenomena at different scales.
In contrast to some of the previous works, we pursue a 3dimensional computational model which is capable of resolving the aforementioned critical phenomena. The purpose of this study is to address, how heat and mass transfer influence the (de)hydration in a closed TCES-system. For this purpose, both experimental and numerical investigations by considering potassium carbonate (K 2 CO 3 ) (Gaeini et al., 2019) as thermochemical heat storage material are performed. The focus will be hence on modeling heat and mass transfer in a single or chain of particles by resolving the governing equations at particle scale. To have an acceptable prediction of (de)hydration of a TCM reactor, available models in the literature are very dependent on certain inputs such as porosity distribution in the bed and effective thermal conductivity of the bed. These inputs can be provided with correlations extracted from welldefined measurements. Considering the effect of particle size, shape and packing method on these inputs, large number of measurements are required to provide suitable correlations for a comprehensive study. The novel numerical model proposed in this work is a reliable approach that avoid the abovementioned issue by using the basic properties of the material (e.g., particle size, particle shape, thermal properties) for the calculations. This is achieved by using a CFD-DEM (Tsuji et al., 1992(Tsuji et al., , 1993 approach in which each salt grain (solid particle) is modeled with a Lagrangian approach while the surrounding gas phase is modeled with a Eulerian approach. The CFD-DEM method is now fully developed and widely applied in granular flows and fluidized beds. It combines computational fluid dynamics for the continuous phase and the discrete element method for the particle phase and makes it possible to investigate the complex processes in a Eulerian-Lagrangian framework. This gives the possibility to study complex physics with high degree of details and resolution in both particle and bed scales.
The obtained results have been compared and validated with the thermal gravimetric analysis (TGA) measurements. Also, the impact of particle size on the hydration process is investigated and a detailed mechanism study of heat transfer between particle-particle and particle-wall is performed to interrogate the inherent physics.

Material
The studied material is K 2 CO 3 (CAS number 584-0-7) produced by Evonik. The material is produced in 83-85% hydrated form, which means the initial material is almost completely a sesquihydrate. The grains seem like single grains and are sieved in a fraction between 1 and 1.4 mm diameter. No further preparation is performed before the hydration/dehydration experiments.

TGA-DSC setup
The hydration/dehydration kinetics are studied in two setups: a home-build vacuum TGA and TGA/DSC of Seteram. The home-build vacuum TGA is a vacuum oven (Binder: VD023-230 V) combined with a balance. The temperature in this setup is controlled within 1 K. The balance is a FUTEK 200-200 g sensor. On top of the balance a water heater/cooled plate is connected, which keeps the temperature of the sample plate constant over time. As water flows through the sample plate a high noise level is present in the measurement. The error in the mass is there for 50 mg. Sample sizes of minimal 10 g should be used to have enough accuracy (2.5% error at 10 g hydrated sample). This setup is used to perform hydration experiments under controlled conditions at low pressure (pure water). For accurate kinetic hydration experiments the sample is pre-cycled for 10 times. At that moment, the hydration kinetics in this setup was constant. The hydration conditions are: 40°C and 12 mbar water vapor pressure. The dehydration conditions are: 40°C with constant evacuation of the setup. After the 10 pre-cycles different hydration experiments are performed with varying temperature (35-50°C) and different water vapor pressures (12-19.9 mbar). During these experiments the oven and the sample plate temperature are kept constant. In between these experiments the samples are dehydrated at 40°C under a constant evacuation of the oven. Before switching hydration/dehydration conditions the full conversion was reached. The water vapor pressure in the setup is controlled by an external evaporator, whereby the temperature of the evaporator is controlled by a thermostatic bath. The saturation vapor pressure of pure water is used to control the absolute humidity in the system.
As the vacuum TGA cannot be heated above 70°C, the dehydration is performed in a Seteram TGA/DSC Sensys EVO. The temperature of the sample is controlled within 0.1 K. Due to the geometry of the sample holder, the diffusion of the water vapor into the pile of particles in the sample holder is limited. This causes a small pressure-drop and consequently makes it impossible to perform kinetics hydration experiments in this setup. On the other hand, the dehydration experiment has sufficient accuracy as the pressure differences in this kind of experiments are much larger. In this setup only a limited amount of TCM (~10-40 mg) can be tested but with high accuracy (<2 lg error, <0.02% uncertainty on a sample of 10 mg). Also, in this setup the material is pre-cycled for 10 times, while the material is dehydrated at 120°C and 12 mbar vapor pressure and hydrated at 40°C and 12 mbar vapor pressure. After the pre-cycles the temperature during dehydration is varied from 90 to 130°C. Before switching between hydration/dehydration conditions the full conversion was reached.

Numerical modeling
The packed bed reactors are widely used in industrial applications and involve intra-particle transport processes in the presence of non-isothermal reactive particulate flows. The complexity involved in the modeling of intra-particle transport processes and reactions has attracted the attention of many researchers and engineers. In the present work, the thermal interaction between particles (conduction, radiation) and the particle's interaction with its surroundings (conduction, convection) is resolved with the extended discrete element method (XDEM). The XDEM is an advanced numerical simulation tool for modeling and solving multiphysics problems, which was proposed and developed by Peters (Peters and Peters, 2003) and will be further extended in this work for the application in thermochemical heat storage; XDEM presents the coupling between computational fluids dynamics (CFD) and discrete element method (DEM) for modeling the continues and discrete phase simultaneously. Such a coupled numerical method allows you to explore numerous mechanisms and processes in fixed or fluidized beds and is very effective for addressing the challenges in reactor engineering. This method treats thermochemical conversions and motion of particles by considering each particle as an individual entity. In the case of porous particles, diffusion of gas in the pores volume is accounted for. A particle can exchange heat and mass with its surroundings subject to the applied boundary conditions at the interface. The temperature distribution and the species are accounted for by a system of onedimensional transient conservation equations (Peters and Peters, 2003;Ha and Choi, 1994). In this framework, investigations and predictions in a de-coupled mode i.e. considering just granular phase without CFD part are also possible (Peters et al., 2010).
The equation of conservation of mass within the pore volume of a porous particle is written as follows: @ @t The term on the right-hand side of the Eq. (1) deals with the mass exchange between solid and gas phases due to chemical reaction. Where q f is the fluid phase density and m ! f is the advective velocity and e f represents the porosity of the particle. The fluid density q f is the sum of the densities of all the species present in the gaseous phase i.e.
The species equation can be written as @ @t It is considered that the transport of gaseous species obeys Darcy's law (Darcy, 1856) inside the pore space of particles. For momentum conservation Darcy's law states as Since Darcy's law is derived for the flows with low Reynolds numbers hence it is applicable for certain flow regimes.
In order to obtain the energy balance equation, the temperature of the different interacting phases inside the particle i.e. gas, liquid or solid is assumed to be the same at interphase or at the point of interaction of the two phases. The thermal mass of the gas phase qc p is negligible compared to the mass of fluid and solid phase, so in the heat transport through bulk motion, the diffused gaseous species within the pore spaces are neglected. Therefore, the energy equation for a porous medium based on the homogeneous model as reported by Faghri (Faghri and Zhang, 2006) is written as The source term on the right-hand side of Eq. (5) represents the exemption or consumption of heat due to chemical reaction k and H k is the reaction enthalpy, where k eff is the effective thermal conductivity. Further, this formulation allows choosing the geometry of the domain to be either an infinite plate, infinite cylinder or sphere by selecting the value of n equal to 0; 1 or 2 respectively. Thermodynamic equilibrium for intra-particle fluid is assumed and considered it as an ideal gas.
In order to complete the model for a single particle, the symmetric boundary condition is applied at the center of the particle because the particle shape is considered as either sphere, infinite cylinder of the infinite plate. Àk eff @hTi @r For the heat and mass transfer, the following boundary conditions are applied at the particle surface Àk eff @hTi @r where T 1 denotes the ambient gas temperature and q i;1 is the density of ambient gas species i. The heat flux _ q 00 rad represents the potential radiative heat exchange with the surroundings and the heat flux _ q 00 cond represents conductive heat exchange through physical contact with other objects. D, b and h are the diffusion, mass transfer and heat transfer coefficients respectively. (See Fig. 1)

Reaction kinetics
The focus of this numerical work is on the validation of the experimental results for the (de)hydration of the potassium carbonate (K 2 CO 3 ). Fig. 2 shows the interaction of particles with each other and with their surroundings. The basic reversible reaction behind the hydration and dehydration process of K 2 CO 3 particle can be expressed as: Generally, it is assumed that the rate of the reaction can be parameterized in terms of temperature T, partial pressure P and conversion X (Vyazovkin et al., 2011). So, for a single-step mechanism, the following kinetic model must hold here K defines the rate coefficient and hðP; P eq Þaccounts for the pressure dependence. The functionf ðXÞ describes the reaction mechanism. The conversionX should be defined in such a way that the following expression must be satisfied Based on the TGA measurements (explained in the section 2), the reaction kinetics was calculated as bellow In case of hydration reaction K ¼ 2:7 Â 10 À9 ½1=s; E ¼ À34828 ½J=mol and q ¼ 0:7, where in the dehydration case K ¼ 225 ½1=s; E ¼ 43382 ½J=mol and q ¼ 0:8.

Contact model
To deal with solid-solid interactions a contact model is required that can handle thermal heat exchange between solids. In the thermal heat transfer process through particulate or granular media, the solid particles exchange heat fluxes with other neighboring particles or with solid walls through conduction and radiation. For this additional surface source or sink terms for solid phase are introduced. The neighboring particles of each particle are detected by standard DEM algorithms and the exchange of thermal heat fluxes are calculated accordingly. By considering the boundary condition for a single particle pas given by Eq. (8), i.e. heat fluxes the particle p is receiving form its N number of neighboring particles by conduction and radiation are calculated by using the _ q 00 p; cond ¼ where F p!j is the view factor or area of contact between the particle p and its neighbor particle j and k denotes the thermal conductivity of the particles.

Validation
Firstly, test cases are performed for the hydration and dehydration of K 2 CO 3 particles in order to validate our mathematical model. Particles are modeled as perfect sphere, as it is the closest shape to the real irregular shape of particles. Moreover, It is assumed that particle size does not change during hydration and dehydration. The material properties are listed in Table 1. Due to the uniform conditions for all particles in the measurement (i.e. each particle is in contact with its neighbor particles on the left and right and the heat exchanger wall underneath), for the sake of validation, we only calculate the hydration and dehydration of one single particle. The obtained numerically calculated and predicted results are compared with the experimental findings. The goal of this comparison study is to provide a confirmation that the developed model is implemented accurately in the present method.

Hydration of a single particle layer
The hydrated K 2 CO 3 particle is considered as spherical in shape with a diameter D p ¼ 1 mm surrounded with a constant pressure of P atm ¼ 12mbar. The simulations of the hydration process are performed at four different values of temperature as applied in the corresponding experiments i.e. T ¼ 35 C; 40 C; 45 C and 50 C.
The conversion X is given by the expression in Eq. (15) where m hyd ; m intÀhyd and m f Àhyd are the masses of the K 2 CO 3 particle at hydrate, initial hydrate, and final hydrate state, respectively. The comparison of numerical predictions with experimental results for the conversion of the particle during the hydration process at different values of applied temperature with respect to the time is shown in Fig. 3 (a). Overall, the numerical results show good agreement with the experimental results for the conversion of a particle with time at all the different temperature conditions applied around the particle. It is noticed that the conversion rate of the particle at low temperature i.e. T ¼ 35 C is high and this rate  of conversion decreases as the temperature increases. At lower temperatures, the pressure difference between the equilibrium pressure (associated with local particle temperature) and local pressure of water vapor on the particle is the driving force for the reaction. According to the P-T diagram reported by Sögütoglu et al. (Sögütoglu et al., 2018), this driving force is higher at a lower temperature. That is why hydration is faster at T ¼ 35 C. Also, It can be seen that the reaction kinetics is not constant with time and it decreases as the reaction time progress. This kind of behavior is well known for such kinds of reactions between solid and gas . Further, the change in temperature of the particle surface during the reaction is calculated. The plots in Fig. 3 (b) illustrate the comparisons of the numerically calculated results and the obtained experimental results for the change in the particle's temperature with respect to time during the reaction. It can be seen in Fig. 3 (b) that although the temperature of the particle is slightly different from the experimental values during the reaction, the results of our mathematical model agree very well with the experimental data in the beginning and at the ending time of the reaction. The temperature difference between measurements and model can be explained by following two reasons uncertainty on the location of the thermocouple: initially it is tried to make a good contact between thermocouple and surface of salt particles. However, during the hydrationdehydration cycle, the particle size varies slightly. This can affect the contact between thermocouple and particle surface and consequently affect the temperature measurements. irregular shape of particles: the real particle shape is irregular (but close to sphere). This means each particle might have a slightly different contact surface with the heat exchanger wall. This certainly affects the particle temperature. However, in the model spherical particle has been used which has a known contact surface with the wall.
Next, we considered that the temperature applied to the particle is fixed at T ¼ 45 C, where the value of the applied pressure around the particle varies. Simulations are performed by considering four different values of the pressure according to the relative experiments and we select P ¼ 12mbar; 14:4mbar; 17mbar and 19:9mbar. Fig. 4 shows that the conversion of a particle with respect to time at constant temperature and varying applied pressure values. The experimental data shows a slight decrease in the first seconds due to the Buoyance effect, wherefor no solid correction was possible for the given experimental setup. The conversion data is reliable after 200 s after switching from dehydration to hydration. The numerically predicted results are showing a good agreement overall with the experimental results. From the obtained results, it is found that the higher the applied pressure higher is the conversion rate of the particle, as it can be seen that at P ¼ 19:9 mbar the conversion of the particle completed earlier compared to the rest of the cases with lower applied pressure. This is due to the higher driving force for the reaction at a higher pressure.

Dehydration of a single layer of particles
The setup of this case is similar to the case of the hydration of a single layer of particles of diameter D p ¼ 1 mm, where the applied vapor pressure is 12mbar. The simulations of the hydration process are performed at four different values of temperature as applied in the corresponding experiments i.e., T ¼ 90 C; 100 C; 110 C; 120 C and 130 C. The comparison of numerical results with experimental results for the conversion of the particle during the dehydration process at different values of applied temperature with respect to the time is shown in Fig. 5. The numerical results are in perfect agreement with the experimental results for the conversion of a particle with time at all the different temperature conditions applied around the particle. It is observed that the dehydration reaction within the particle at high temperatures completed earlier as compared to the cases with low temperatures. This is due to the higher driving force at a higher temperature compared to a lower temperature which is closer to the equilibrium temperature at 12mbar. Further, it is noticed that there is a big difference between the reaction completion times when T ¼ 90 C and T ¼ 100 C, but this difference significantly decreases at higher values of applied temperature.
The comparison of the obtained numerical results with the experimental results is in good agreement for both hydration and dehydration cases. This provides us a confirmation that our mathematical model is implemented accurately which gives us the confidence to further explore the inherent physics of this reversible reaction and perform an in-depth mechanism study of the heat and mass transport in a packed bed of TCM. . Comparison with experimental data (a) Conversion of a bed of single K 2 CO 3 particle with respect to time at different temperatures, the error in conversion is max. 2.5% (b) Temperature at the surface of a single K 2 CO 3 particle with respect to time at different applied temperatures.

Particle size effects on hydration
Here we study the impact of particle size on the hydration process. The setup of this case is similar to the case for the validation of hydration of a single particle i.e., a single particle is placed at the top of a plate with a constant temperature T ¼ 35 C and vapor pressure P ¼ 12 mbar is applied on the particle. The particles of different diameters are used ranges from 1 mmto20 mm:The time taken by the particles for the complete conversion with respect to different particle sizes is shown in Fig. 6. As it can be seen that when the diameter of the particle increases the time required for the hydration to be completed also increases. This is due to the fact  that mass transfer is in a larger particle is not as good as a small particle. Moreover, the smaller particle has a larger surface to volume fraction therefore heat transfer from the particle to the surrounding is better, hence hydration will be faster in smaller particles. It is observed that the time required for full hydration will be about 12.6 times more for the particle with a diameter of 20 mm compare to the particle with 1 mm diameter.
Here we look closely into the particle to study the detailed physics of the hydration process. For this purpose, we look at the particle with a diameter of 12 mm placed at the top of a plate with a constant temperature and at constant vapor pressure. The gradient of the mass fraction (conversion) of hydrated salt over time within the particle of 12 mm diameter is shown in Fig. 7 (a). Here the radius is normalized, it means 0 is the center of the particle and 1 is the surface of the particle. It can be seen that the hydration starts at the surface of the particle while inside the particle the reaction starts with a delay. We can see there is a considerable gradient of conversion within the particle, as the hydration is completed at the surface after 10,000 s while at the center it takes 16,500 s to finish the hydration process. Fig. 7 (b) shows the tem-perature gradient inside the particle with respect to time along the radius of the particle. It is found that there is a slight temperature gradient along the radius of the particle and the temperature profile is almost straight all the time along the radius. This small temperature gradient is one of the main causes of a long time needed for the hydration in a large particle. This will be discussed later in this section. Fig. 8 (b) shows the water vapor pressure along the radius of the particle with respect to time. During the hydration the water vapor is consumed, therefore the pressure decreased with time. Due to the pressure difference between inside and outside of the particle, the vapor constantly enters into the particle. The vapor is consumed on the way inside the particle due to the hydration. This leads to a pressure drop inside the particle from the surface to the center. It can be seen in Fig. 8 (b) that near the end of the reaction the pressure gradient is higher because the reaction rate is high at that time (see Fig. 8 (a)). The reaction rate inside the particle is shown in Fig. 8 (a). In the beginning, the reaction rate is high close to the surface of the particle and almost zero at the center, this means that there is no considerable hydration taking place at the center. The reaction front (ration layer) moves towards the center over time and this move can be seen with the increase in the reaction rate. Regardless of the location, reaction rate increases over time, this is mainly due to the lower particle temperature when the hydration approaches the end, see Fig. 7 (b). The reaction rate becomes maximum when the hydration is taking place at the center of the particle.
In Fig. 9 you can see the conversion and the vapor pressure at three different locations in the particle i.e. at the center when r ¼ 0:07, in the middle when r ¼ 0:52 and close to the surface of particle where r ¼ 0:92. It is noticed that pressure near the center reduces over time due to vapor consumption at upstream and it intensifies once the hydration rate increases close to the center. Looking at r ¼ 0:52, the pressure profile can be divided into three sections with an almost constant slope: (1st)~500 s to~6000 s: the pressure drop is mainly due to vapor consumption upstream as the hydration rate at this location, i.e., at (r ¼ 0:52), is relatively low. (2nd) From~6000 s to~11,500 s: during this period the slope is steeper, it is mainly due to hydration around r ¼ 0:52. An increase in the slope of conversion also proves that (3rd) from~11,500 to 14,000 s: pressure increases since the conversion Fig. 6. Particle conversion time with respect to particle diameters. Fig. 7. The gradient of (a) conversion and (b) temperature of the particle over time along its radius.
A. Mahmoudi, Pim A.J. Donkers, K. Walayat et al. Chemical Engineering Science 246 (2021) 116877 at r ¼ 0:52 approaches to the end and less vapor is consumed, therefore vapor pressure approaches to 12 mbar (outside pressure). A similar kind of behavior can be seen at r ¼ 0:07 but with a steeper reduction and rise in the pressure and the conversion profile.
The pressure distribution inside the particle at different time instants is shown in Fig. 10 (a). As you can see there is always a pressure gradient inside the particle at different time instants. The gradient increases over time and once the hydration is finished it reduces again (compare pressure profile at time 14,000 s and 15,000 s). This has been also seen in the 3D graph as well in Fig. 8 (b). Fig. 10 (b) is showing the temperature distribution inside the particle at different time instants and you can see at any time there is a slight temperature gradient inside the particle. It can be seen that, although there is a pressure drop in the particle due to the consumption of water vapor by hydration, there is always vapor available with considerable pressure. However, near the center of particle hydration, this is almost negligible for a period of time (see Fig. 7 (a)).
At the time t ¼ 1000s, the vapor pressure at the center is 9:4 mbar and the associated equilibrium temperature to this pressure is 56:53 C, where the temperature at the center of the particle is 56:43 as shown in Table 2. Although, the temperature is slightly below equilibrium temperature, but the difference is very small, Fig. 8. The (a) reaction rate inside the particle and (b) vapor pressure along the radius of particle with respect to time. Fig. 9. Conversion and vapor pressure at three different locations in the particle.
A. Mahmoudi, Pim A.J. Donkers, K. Walayat et al. Chemical Engineering Science 246 (2021) 116877 therefore the driving force for hydration is almost negligible. This is the reason why the hydration does not proceed close to the center of the particle. Table 2 shows that the vapor pressure and temperature at 8000s and 14000s as well. It can be seen that even at 8000s the temperature at the center is very close to the equilibrium temperature but at 14000s there is a considerable difference, that is why the hydration is relatively fast at this time at the center of the particle (see Fig. 7 (a)).

Hydration and dehydration of a chain of particles
Here we want to study the behavior of particles during the hydration and dehydration once they are placed in a packed bed. In this case, the particles exchange energy among other particles as well as with the heat exchanger. The schematic diagram of a packed bed of the TCM in a closed system is shown in Fig. 11. The temperature of the heat exchangers along the three side walls is T ¼ 35 C and constant water vapor pressure P ¼ 12 mbar is sup-plied to the bed from the left side. This pressure continuously reduces along with the x À direction due to the reaction. However, for the first column of the particle, it can be assumed that vapor pressure for all particles is the same as supplied from the inlet. For simplicity, we just model one column of particles to see their behavior during hydration and dehydration processes. As the vapor pressure and the temperature of the incoming vapor are known for the first column of particles, therefore, solving gas flow around the particles are not necessary. This model can be a representative of the entire bed and its study can describe the behavior of the whole bed. Due to the symmetry, we only modeled half of this column and the particle at the bottom is in contact with heat exchanger and the rest of the particles are arranged on its top such that it makes a chain of particles.
In this case, we have six particles of dehydrated K 2 CO 3 arranged in such a way that the particles make a chain-like structure as shown in Fig. 12. The particles are considered as spherical in shape and the diameter of each particle D p ¼ 1 mm. Constant pressure Fig. 10. Conversion and vapor pressure at three different locations in the particle.

Table 2
The vapor pressure and temperature distribution inside the particle. P ¼ 12 mbar is applied around the particles. The 1st particle is in contact with the heat exchanger and the temperature of this wall is set as T wall ¼ 35 C for the first 10 h and then it is set as T wall ¼ 90 C while the other particles are connected with each other. The plots in Fig. 12 shows the vapor pressure in the reactor and the heat exchanger temperature over time. According to the temperature profile of the heat exchanger, the hydration reaction is expected for the first 10 h and dehydration reaction for the last 10 h, which represents a complete reaction cycle for the charging and discharging of the TCM. Fig. 13 shows the conversion of the particles with respect to time during the hydration and the dehydration process. The conversion rate of the first particle is very fast and is completed first then in the second particle, then in the third and so on. The furtherer the particles are from the heat exchanger, the slower is their conversion. It can be seen that the conversion in the particle 5 and 6 reached to only 20 percent and 5 percent after 10 h respectively because the amount of heat transferred from these particles through conduction is less as compared to the particles which are closer to the heat exchanger. Further, when the reverse reaction starts for dehydration, it is again seen that the conversion of the first particles is fast compared to the others and with increasing distance from the heat exchanger the conversion rate of the particles drops due to the low heat conduction.
Moreover, the temperature of the particles with respect to the time during one complete hydration-dehydration cycle is shown in Fig. 14. It shows that the first particle has better heat transfer with the heat exchange compared to other particles, therefore, the heat is generated due to the hydration that can be taken away from the first particle better and its temperature remains low. It is found that the larger the difference between particle temperature and equilibrium temperature, it provides a larger driving force for the reaction which means the reaction will be faster. The particles 5 and 6, in the beginning, undergo hydration and their temperature increases to the equilibrium temperature associated with 12 mbar which is T ¼ 60:03 C. Because of the low heat transfer through conduction from these particles, they cannot release the heat to the environment with the same rate at which it is produced by the reaction. Therefore, their temperature remains at the equilibrium temperature, and hence, hydration does not proceed anymore. The reaction does not start until the temperature of their neighboring particle drops down from the equilibrium temperature value. For particle 5, this occurs after about 3 h when the temperature of the particle 4 is considerably below the equilibrium temperature (look at the zoomed figure). The same behavior can be seen for particle 6 when time is 8 h but with a lower rate, because the temperature difference between particles 6 and 5 is less than that of particles 5 and 4.
At 10 h, when the temperature of the heat exchanger shifts and jumps from below to above the equilibrium temperature, dehydration starts in particles 1, 2 and 3 with a small delay. However, particle 4 takes more time to sense the temperature rise as compared to particles 1, 2, and 3. Therefore, the hydration is continued in this particle for about 0.2 h as its temperature is below equilibrium temperature. After this delay, when its temperature goes above the equilibrium, the dehydration process starts in it.
If you look at the temperature profile of particle 1, as the conversion proceeds and approaches to the end, since the reaction rate decreases and there is less salt to be hydrated, therefore less heat is   produced. This leads to a temperature drop with a higher gradient towards the heat exchange temperature. When the hydration in the first particle is completed (at about 1.5 h) no heat is generated in this particle, however, it receives heat from particle 2 by conduction (as a result of the heat generated by the hydration in the particle 2). Therefore, the temperature decreases with a lower rate and the slope of the temperature profile changes. If there is only one layer of particles on the heat exchanger, the temperature reaches to heat exchanger temperature see Fig. 15. However, when there is more than one layer of particles, the temperature profile is slightly different. As can be seen in Figs. 14 and 15, the slope of temperature profile decreases after the hydration finished in particle 1 or it is continued with a lower rate. This is due to the heat generated by the particle 2 which is transferred by conduction to the first particle. Similar behavior can be observed at 4 h (when the hydration is completed in particle 2) in the temperature profile of particles 1 and 2 due to hydration of particle 3. To explain this behavior in a better way, it is useful if we compare it with the case of a single particle on a heat exchanger. Fig. 15 shows the comparison of the conversion and temperature change with respect to time of a single particle when it is in contact with the heat exchanger alone and when it is in contact with both heat   15. Comparison of the conversion and temperature of a single particle in contact with heat exchanger alone and with heat exchanger as well as with neighbor particle.
A. Mahmoudi, Pim A.J. Donkers, K. Walayat et al. Chemical Engineering Science 246 (2021) 116877 exchanger as well as with the neighboring particle in a chain of particles during hydration reaction. The conversion rate in the case of the one-layer particle on the heat exchanger (single particle) is higher because it has only heat transfer with heat exchanger. But in the case of the multi-layer bed, particle in contact with the heat exchanger as well as with other particles, it also receives heat from its neighbor particle. Therefore, its temperature is higher compared to the single particle and this leads to a lower driving force for the hydration reaction. Further, the temperature and the net heat conduction between particle-particle and between particle-wall of the heat exchanger with respect to time for the first three particles in the chain of particles are shown in Fig. 16. The negative values for the net heat conduction mean that the heat is dissipated out of the particle. Since in our case the packed bed is in the vacuum so there is no thermal convection involved and the heat transfer mechanism is only based on conduction and radiation. The heat generated in the first particle is transferred to the heat exchanger while at the same time it (particle 1) receives heat from its neighbor particle (particle 2). Comparing the temperature difference between particles 1 and 2 and the wall, it is obvious that the net heat conduction must be outward i.e. from particle 1 to heat exchanger.
It is noticed that, as the time progress and conversion proceeds in particle 1, the temperature of particle 1 decreases, therefore the temperature difference between particle 1 and heat exchanger decreases which reduces the net heat conduction, but it is still outward. A similar kind of behavior is also seen in particles 2 and 3; the net heat conduction is reduced after their conversion is completed. However, the net heat conduction is still negative, that is why their temperatures keep decreasing over time. It also shows that the flow of heat from particles during hydration is towards the heat exchanger.
The net conduction heat transfer from particle 1 is larger than particles 2 and 3, even after complete conversion. The main reason behind this is the high thermal conductivity of the heat exchanger wall. The temperature difference between particle 1 and the wall is less than the temperature difference between particles 1 and 2. Even with a lesser temperature difference between particle 1 and wall, it causes negative heat conduction from particle 1 and larger in amount than that of the net heat conduction of particle 2.

Conclusion
This paper presents a detailed mathematical model of the coupled heat and mass transfer phenomena in K 2 CO 3 salt hydrate in the Extended Discrete Element Method (XDEM) framework for the thermochemical heat storage process. K 2 CO 3 is considered as one of the most promising materials in terms of cyclability and energy density for the thermochemical heat storage and based on the TGA results we modeled the reaction kinetics of the reversible (de)hydration reaction. The reaction kinetics of the hydration and dehydration processes of K 2 CO 3 particles are studied with different applied conditions. First, we simulate numerical test cases for the validation of hydration and dehydration processes in a K 2 CO 3 particle. For the sake of simplicity of the problem and because of the uniform conditions around the particle, we consider a single spherical K 2 CO 3 particle for the validation of hydration and dehydration processes. The obtained numerical results and the TGA measurements shows a perfect agreement with each other which gives us a confidence to enhance our investigations and to further explore the inherent physics of the reversable thermochemical reaction. It is found that during hydration process the reaction rate is high at low temperature when constant pressure is applied on particle and due to high reaction rate the particle's temperature decreases much faster. In the case of constant temperature and varying pressure, the reaction rate of hydration reaction is high at higher pressure because of the higher driving force for the reaction. Where in dehydration process with constant pressure the reaction rate is high at higher temperature. So, we conclude that the hydration and dehydration reaction rates are directly and indirectly proportional to the vapor pressure and temperature respectively. It is also found that the particle's size has an impact Fig. 16. Temperature and the net heat conduction between particles and particle-wall.
A. Mahmoudi, Pim A.J. Donkers, K. Walayat et al. Chemical Engineering Science 246 (2021) 116877 on the hydration reaction and the conversion time of the particle increases with increase in the size of particle. Further, it is seen that the reaction in the core of particle starts and ends with a time delay compared to its surface, also there is always a small temperature gradient along the radius. The consumption of water vapor inside the particle during hydration leads to pressure drop inside the particle from surface and the reaction rate becomes maximum when reaction front reaches the center of particle. Finally, we simulate a hydration dehydration cycle in a chain of six particles which represents an entire packed bed. It is found that smaller the distance between particles and heat exchanger higher is the reaction rate. The net conduction heat transfer in the particles which are closer to the heat exchanger is high that leads to a higher reaction rate in particle 1 (which in contact with heat exchanger) than in particle 2 and higher reaction rate in particle 2 than in particle 3 and so on. Furthermore, it is observed that in multi-layer bed of particles the rate of hydration rection is low compared to a single particle because in case of multi-layer bed a particle also receives heat from its neighboring particles.
CRediT authorship contribution statement

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.