Thermo-Mechanical Investigations of Packed Beds for High Temperature Heat Storage: Continuum Modeling

Abstract

Thermal energy storage (TES) systems are central elements for various types of new power plant concepts, whereat packed beds represent a promising storage inventory option. Due to thermal expansion and shrinking of the packed bed’s particles during cyclic thermal charging and discharging operation, high technical risks arise, and possibly lead to material failure. In order to accurately design the heat storage system, suitable tools for calculating induced forces and stresses are mandatory. Continuum models offer time efficient simulation results, but are in need of effective packed bed parameters. This paper introduces a methodology for applying a simplified continuum model and presents first results for an exemplarily large-scale application.

1. Introduction

The application of thermal energy storage (TES) is a key technology for power plants by offering improved operating flexibility for conventional power plants and extended operating time for power plants driven by renewable energy sources [1,2]. Other power generating applications for TES are adiabatic compressed air energy storage systems [3,4]. Further TES concepts exist in industrial heat management, e.g., in the steel industries (in the form of “Cowper” storage), in glass producing industries and air purification systems [5,6]. Here, the aim is a decrease in operational costs and an increase in system efficiency and flexibility.

When air is used as a heat transfer medium, regenerator-type heat storage with a packed bed inventory is a straightforward design option, offering cost-effective solutions and flexibility in terms of arrangement and material choices for the storage medium. Possible choices for packed bed materials are ceramic pebbles or simply natural stone [7]. One of the major challenges for such a packed bed TES system is its vulnerability to mechanical failures caused by the punctiform contacts between particles and the encircling container wall [8]. During the thermal charging and discharging processes of the heat storage system, the pebbles or particles of the packed bed start to expand and shrink, respectively. As a consequence, thermo-mechanical loads are induced which lead to high technical risks for the storage inventory itself and the encircling containment insulation of the heat storage system. Due to the operating mode consisting of repeated cyclic thermal loading and unloading, the induced thermo-mechanical loads potentially lower the life expectancy of the storage system. Therefore, knowledge of the occurring stresses and loads are required in order to adequately adjust the packed bed installation and container design if necessary. Hence, appropriate numerical models and simulation tools for calculating these loads are needed. Conventional approaches for mechanical packed bed modeling are based on either particle discrete models [9,10] or continuum models [11]. As discussed in [12], both of these approaches offer advantages and disadvantages: While the particle discrete models are able to model packed beds with high accuracy at the expense of high computing times, continuum models generally offer faster simulation results, but are in need of complex model parametrizations.

The work at hand aims to provide a methodology for a simplified approach to approximate the thermo-mechanical stresses of a large-scale regenerator type TES with a packed bed as the storage medium. The present work is a continuation of a previously published paper [12] and builds upon its findings. Therefore, it is recommended to first study the findings in [12] in order to obtain additional background information and deeper insight of the discussed subjects in the present work. Subsequently, the results of an exemplarily large-scale regenerator type TES, as described in [4], is used as a basis for a validation of the developed simulation approach.

2. Mechanical Packed Bed Modelling

This paragraph covers the theoretical background of the present study. It will give a brief overview of the mechanical modelling approaches for packed beds. These modeling techniques are the basis for the implementation of a simplified approach for calculating thermo-mechanical stresses of the large-scale TES application [4].

An early approach to calculate packed bed stresses was derived by Janssen [13] in the context of silo technology. By analyzing a theoretical disc element of a silo packed bed, he derived the following analytical correlation of vertical stresses (${\sigma }_v$) acting on the container wall:

$${\sigma }_v=\frac{g·{\rho }_b·A}{\lambda ·\tan \left({\phi }_x\right)·U}·\left(1-e^{\frac{-\lambda ·\tan \left({\phi }_x\right)·U}{A}}\right).$$

(1)

This correlation is dependent on the bulk density (${\rho }_b$), cross-section ($A$), horizontal load ratio ($\lambda$), wall friction angle (${\phi }_x$), the containment circumference ($U$), and gravity ($g$). The required bulk properties can be determined with various soil mechanical tests. More information on this approach can be gathered from multiple sources, e.g., [13]. While this model enables the calculation of stresses in stationary and thermally uncharged packed beds, it is not employable for applications with additional thermal loads.

For a more detailed modelling of the packed bed, two fundamentally different approaches exist: the first is a particle discrete model, where each particle is modelled as independent physical object; the second is a continuum model, where the packed bed is perceived as a macroscopic, cohesive material.

Particle discrete models can be implemented using the Discrete Element Method (DEM). The DEM enables the computation of forces between interacting particles. Here, each particle is represented as a friction–spring–damper system. Thus, particle–particle and particle–wall contacts will result in spring, friction, and damping forces. The algorithm used in DEM enables the determination of time dependent movement and rotation of particles. More details on the governing equations and algorithms are given in [12]. Previous work has shown that the DEM produces highly accurate results and allows valuable insights into the basic phenomena [14,15]. The particle discrete approach is well suited for packed beds with a relatively small number of particles or for a short simulation time. However, due to the high computing time required to solve the multitude of underlying numerical operations for a large number of particles, the DEM is unfit to simulate large-scale applications.

As discussed in the previous work [12], a continuum based finite element model approach for packed beds offers a significant speed up in calculation time in comparison to discrete models. The major challenge of the application of a continuum model is the adequate parametrization of the packed bed properties. In the case of continuous packed bed modeling, “effective” material properties are needed, instead of bulk material properties. Especially in the area of soil mechanics and nuclear fusion technology, different constitutive equations have been developed to describe the continuous characteristics of packed beds. In the context of the investigation of pebble beds for the use in ceramic breeder blankets, Reimann [16] developed a power law function for the effective Young’s Modulus, dependent on stress and temperature, by conducting uniaxial compression tests. The function can be converted to the following form:

$$E_{eff}=C_1·\left(1+C_2·T^3\right)·{\sigma }^{C_3}$$

(2)

A similar approach was developed by An et al [17]. The effective Young’s Modulus is a key parameter when describing a packed bed as a continuum and provides insight into the mechanical behavior of the packed bed. Further information about continuum models and relevant parameters can be gathered from multiple sources [11,18,19].

Here, the parameters $C_i$ can be determined by evaluating the stress–strain relations of a selected packed bed during uniaxial compression. Other effective parameters like effective creep, effective thermal expansion and effective thermal conduction were investigated in [20,21,22]. However, the results in this paper are limited to the identification of the effective Young’s Modulus and the corresponding Reimann parameters. Therefore, the presented continuum model is a simplified model. Further effective parameters still need to be determined in future works in order to enhance the accuracy of the continuum model.

3. Methodology

In order to simulate the thermomechanical stresses during thermal cycling of a packed bed, a methodology was developed aiming for time efficient calculations. This methodology is based on the utilization of a particle discrete model for identifying mechanical properties (the effective Young’s modulus) of a specific packed bed (detailed in [12]), as well as thermal profiles of an exemplarily large-scale application (detailed in [4]). This information is implemented into a continuum model, which enables time efficient calculations of the thermo-mechanical induced stresses of the packed bed. Figure 1 displays a flow chart of the involved steps of this methodology and the sequence in which those steps are carried out.

Based on a selected packed bed of ceramic pebbles, a particle discrete model for simulating uniaxial compression tests was implemented into a Discrete Element Method tool. The mechanical properties of the pebbles, e.g., the Young’s modulus and the particle diameter, are used as input for the particle discrete model. During uniaxial compression, the selected packed bed is confined in all directions, while being compressed and sequentially relaxed with a top cover plate. This procedure of compression/decompression is repeated for multiple cycles. The resulting stress/strain relations of the packed bed are recorded, which enable the determination of mechanical packed bed properties, e.g., the effective Young’s Modulus. The model was validated with the experimental results of a custom-built test rig for uniaxial compression tests. Further information on these steps can be gathered from the preceding work [12], in which the modelling approach and experimental setup are described in detail.

Based on the stress/strain results of the uniaxial compression simulations, it is possible to calculate specific continuum parameters [16]. Those parameters are used as central elements for the implementation of a continuum model in a Finite Element Model (FEM) tool. The work at hand is based on a simplified approach, in which only the effective Young’s modulus is implemented. The implementation of the continuum model is presented in the next section. In addition to the mechanical aspects, thermally induced forces need to be taken into consideration as well. For this purpose, the thermal profile shown in Figure 2 is applied. The selected profile shows the temperature spread between the end of thermal charging and discharging process of the large-scale TES described in [4]. This specific profile was chosen since it represents the maximum temperature spread during TES operation in [4] and thus can be considered as the “worst-case” profile, leading to maximum thermo-mechanical loads.

With the combined mechanical input from the particle discrete simulations and the thermal profiles derived from the external source [4], it is possible to utilize the continuum model for calculating the contact stresses between the selected packed bed and a container wall of the large scale application described in [4]. Due to the lack of large scale packed bed TES in commercial applications, there currently are no literature data available regarding the thermo-mechanical behavior of large scale TES. Therefore, the data available in [4] is the only source for validation of the model at the moment.

4. FEM Model

The effective Young’s Modulus was implemented in the FEM software ANSYS. Figure 3 illustrates the applied algorithm of the model and its iterative sequence. It relies on two external inputs: a thermal profile for setting the thermal boundary conditions and the Reimann coefficients derived from uniaxial compression tests for the specific packed bed. This quasi stationary model is able to calculate the acting stresses on the container wall of a packed bed TES for one specific time. For a time-resolved calculation, the model needs to be initialized each time from new.

At first, the packed bed is confined by geometrical boundary conditions and material properties. The ceramic packed bed (see [12] for material properties) was modeled with a packed bed diameter of 10 m and packed bed height of 40 m, according to [4]. The continuum was subdivided into a grid of 68 × 12 elements. This low number of elements was chosen for time-efficient calculations while still obtaining reliable results. Initially, a uniform effective Young’s Modulus is assigned for all elements of the grid. In the next step, temperatures corresponding to the temperature profile in Figure 2 are assigned to the elements. The implemented temperature profile is shown in Figure 2. Here, the temperature distribution is considered uniform in radial direction but changes in axial direction over bed height.

Based on the thermally induced forces, due to thermal expansion, and mechanically induced forces, due to self-weight of the bulk, the axial stresses are calculated for each element. The resulting stresses, together with the previously determined Reimann parameters, serve as input for the calculation of an updated effective Young’s Modulus. The updated Young’s Moduli are assigned to the elements of the grid, and in an iterative sequence, the solving process is repeated until the stresses reach an equilibrium state. The assessment of the contact pressures between the packed bed and the confining wall is based on the converged results.

5. Results

The stress–strain relation of a ceramic packed bed has been determined via experimental and simulation results as presented in [12] and are used as a basis for the present work. The parameters are fitted to approximate the experimental and simulation data following Reimann’s approach.

The effective Young’s Modulus, which is calculated from the stress/strain ratio of uniaxial compression tests, is plotted as a function of stress in Figure 4. Here, the experimental results are shown as well as the results from the DEM simulations. Again, details on the procedure of those tests are described in [12].

For all four different cases shown in Figure 4a–d, the effective Young’s Modulus increases at a higher rate for low stresses compared to higher stress levels. This means that the particles still have the possibility to move and to rearrange themselves at low stress states due to the low compaction of the packed bed. The increase rate slows down since particle movement is increasingly restrained with continuing compression.

Comparing the compression (Figure 4a,b) and decompression (Figure 4c,d) separately, the effective Young’s Modulus at 600 °C is slightly higher than at 20 °C at any given stress state, which is the result of the increasing Young’s Modulus of the ceramic bulk material at higher temperatures. The differences between the gradients of compression and decompression can be explained by energy dissipation effects. The maximum effective Young’s Modulus is significantly higher for the 600 °C packed bed is about 15% to 20% higher compared to the 20 °C packed bed. Again, this is a result of the temperature dependency of the ceramic bulk material’s Young’s Modulus.

Compared to the bulk material, the effective Young’s Modulus for the packed bed at a given stress range is over two orders of magnitude smaller, which matches with the results of [23].

As displayed, the effective Young’s Modulus of the packed bed depends on temperature and stress state and thus can be expressed as a continuous function. In this paper, Reimann’s approach (see Equation (1)) is applied. The parameters ${\mathrm{C}}_{\mathrm{i}}$ were fitted to match the simulation data. Since the compression cycle yield higher stress and stiffness compared to the decompression cycle, the following results for the Reimann parameters are only given for the compression cycle since this represents the more critical case. The resulting parameters are listed for a range of different packing factors (investigated in [12]) of the packed bed in

Table 1. Reimann parameters fitted for the simulation data.

Packing Factor [%]	C1	C2	C3
56.35	104.29	8.03 × 10−10	0.206
56.64	154.11	10.4 × 10−10	0.238
56.96	231.28	5.68 × 10−10	0.307
57.82	325.01	5.66 × 10−10	0.373
58.70	359.97	5.27 × 10−10	0.371

.

The parameters above are implemented in the simplified continuum model described in the previous paragraph. The results for the contact pressure between packed bed and container wall are displayed in Figure 5. Additionally, results of the Janssen approach [13] for the “cold” state of the packed bed is plotted for comparison. The relevant material properties for calculating the Janssen approach are listed in

Table 2. Material properties for the Janssen approach.

$\mathbf{Bulk}\mathbf{Density}{\mathsf{\rho }}_{\mathbf{b}}[\mathbf{kgm}$−3]	Horizontal Load Ratio [-]	Wall Friction Angle [°]	Height [m]	Diameter [m]
1440	0.4	48.7	40	10

.

For both cases (cold and hot state), DEM results from [4] for the large-scale application are also plotted. These results were obtained using a 2D-DEM model to reduce the number of simulated particles and thus computation time. While this decreases the overall accuracy of the simulation results compared to a 3D analysis, the results are still expected to be of satisfying quality.

Figure 5a shows the results for the packed bed in a cold state, whereas Figure 5b displays the results for the thermally charged packed bed.

For the cold packed bed, the results of the continuum model with packing factor 56.35% and 58.7% are in good agreement with the DEM data and Janssen’s approach. Near the bottom of the packed bed, there is a significant drop in pressure. The reasons for this are numerical inaccuracies due to the low number of elements and the chosen boundary conditions of the model. Furthermore, it can be observed that the two results for different packing factors deviate within up to 20% of each other. The more rigid packed bed with the higher packing factor exhibits slightly lower pressures. This result is unexpected at first, but can be explained by the limitations of the model: In the cold state of the packed bed, the occurring pressure is purely based on mechanical, gravity driven effects. Due to the higher effective Young’s Modulus of the packed bed with the denser packing factor, the strain is lower, both vertically and horizontally. Thus, the exerted stresses on the container wall are lower for the less densely packed bed configuration.

In comparison, the dependency of the packing factor on the results for the heated packed bed is more distinctive. In general, the results yield higher stresses for increasing packing factors. In the thermally charged state, the occurring stresses are primarily a result of thermal expansion. Due to the rigid walls of the model, the thermal expansion of the packed bed is converted directly into pressure on the wall. The pressures for a thermally charged packed bed simulated with the continuum model reach their peak values at around 15 m of bed height. The reason behind this is the applied thermal profile (see Figure 2), which has its highest temperature differences at that height. Compared to the DEM result of the thermally cycled packed bed, the continuum model tends to significantly overestimate the stresses of the packed bed for the upper half of the packed bed. In contrast, the stresses on the lower part of the packed bed fall of notably for the continuum model compared to the DEM results. Considering stochastical variance of the DEM results and the limitations of the continuum model due to simplifications, the agreement between those two approaches still is sufficient. In order to improve the accuracy of the continuum model, the changing condition of the packed bed, e.g., the densification and packing factor during thermal cycling, must be considered.

6. Conclusions and Outlook

In order to assess thermo-mechanical stresses in high temperature thermal energy storage, a methodology was developed for time efficient calculations based on a continuum model. In this context, the required effective Young’s Modulus for a selected packed bed was determined by evaluating stress–strain relations during uniaxial compression. For an exemplary application, the resulting pressures of a thermally charged packed bed have been simulated with the continuum model and compared to results from an external source [4]. The overall magnitude of the occurring stresses was found to be of satisfying agreement to the results of [4]. While providing an important first estimate of the occurring stresses, the results also indicate that a more refined model is needed for a more detailed analysis. This can be achieved by utilizing elasto-plastic models, which also captures the plastic characteristic of the packed bed, e.g., its deformation by densification. However, the experimental and simulation effort for this approach is considered significantly higher than the methodology presented in the present study.