Ibrahim Al Asmi

Kai Knobloch

Roman Le Goff Latimier

Thibaut Esence

Kurt Engelbrecht

Hamid Ben Ahmed

Thermocline thermal storage modeling towards its predictive optimal management

Keywords:

The increasing penetration of renewable energies brings into sharp relief the potential of thermal storages, particularly in a multi-energy network. Indeed, they offer low storage cost without environmental concerns. But since physical phenomena during the operation of thermocline thermal storages are complex, associated physical models come with long calculation times. Consequently, their integration within the management of multi-energy networks is challenging at the current time. This work proposes meta modeling based on accurate model's simulation data as an alternative to physical modeling. The developed meta model is compared to three different models with varying complexity. Compared to a typical physical model, the meta model calculation time is reduced by a factor of 1000 while losses as well as internal temperature distribution estimations are provided. In an experimental validation based on four different thermal storage configurations with storage capacities between 450 and 2900 kWh th , a normalized root mean square deviation lower than 3%, between the meta model and the physical model results, is observed. Moreover, physical models often fail to depict real-life behavior of thermal storages, for example due to manufacturing defects or highly dimension-dependent phenomena like flow channeling. In order to overcome this practical challenge, the construction and use of a meta model on the basis of existing experimental data is successfully demonstrated with a normalized root mean square deviation between model and experiments below 4%.

Introduction

In France, 50 to 80 % of electricity consumption has an end use of heat [START_REF] Rte | futurs énergétiques 2050" -consommation et production : les chemins de l'électricité de rte pour la neutralité carbone[END_REF]. Moreover, to achieve a low-carbon economy, heat electrification in residential and industrial sectors is essential and will involve an increase electricity use for Ibrahim AL ASMI et al.: Preprint submitted to Elsevier heating applications [START_REF] Millet | First analysis of the technical potential of electrification of thermal industrial processes by mature technologies -synthesis (6[END_REF]. The coherent production of heat and power, for example in combined heat and power plants, is of increasing interest and paves the way for higher energy generation efficiencies [START_REF] Jimenez-Navarro | Coupling the heating and power sectors: The role of centralised combined heat and power plants and district heat in a European decarbonised power system[END_REF]. Therefore, there is a strong relationship between electrical and thermal carriers. Today, the coupling between these two vectors is rarely taken into account when designing and managing energy networks. While decoupling this problem into several single-vector networks results in a sub-optimal solution, the consideration of a coupled thermal and electrical network could lead to a global optimal solution. This coupling nevertheless comes up against a great technological complexity, from the design of each of the components required for such a so-called multi-energy network, to their modeling and the overall management of the whole [START_REF] Mohammadi | Energy hub: From a model to a concept -a review[END_REF]. The efforts made to overcome the difficulties of a multi-energy network are intensified since the penetration of renewables is increasing and numerous heat resources can be exploited much more widely in the future than today [START_REF] Mancarella | MES (multi-energy systems): An overview of concepts and evaluation models[END_REF]. This is the case for the waste heat generated by industries on the margins of their processes which represents 3500 TWh/year at high temperature (>150 °C) and the equivalent of 2450 Mt CO 2 ∕year in terms of global carbon dioxide emissions [6]. Several renewable sources, such as geothermal energy [START_REF] Moya | Geothermal energy: Power plant technology and direct heat applications[END_REF] and waste incineration [START_REF] Makarichi | The evolution of waste-to-energy incineration: A review[END_REF], also account for a large proportion of available heat resources and make thermal energy storage (TES) an important component in the technological chain of a multienergy network. With economic and environmental costs of 15 $/kWh installed and 1 g CO 2 ,eq ∕kWh installed , respectively, high temperature TES (HTTES) is considered competitive already today [START_REF] Forsberg | Converting excess low-price electricity into hightemperature stored heat for industry and high-value electricity production[END_REF] [START_REF] Marongiu | Modeling of high temperature thermal energy storage in rock beds -experimental comparison and parametric study[END_REF].

State of the art

Thermal storage can be classified according to several criteria [START_REF] Arteconi | State of the art of thermal storage for demand-side management[END_REF], for example i) temperature level (𝑇 <150 °C, 150<𝑇 <300 °C and 𝑇 >300 °C), ii) physical principle (sensible heat, latent heat [START_REF] Xu | Analytical considerations of thermal storage and interface evolution of a pcm with/without porous media[END_REF] and thermochemical [START_REF] Shi | Multi-physics modeling of thermochemical heat storage with enhance heat transfer[END_REF]), iii) storage medium and heat transfer fluid (oil and rock, hot water, molten salt, air and rock, chemical reactants, etc.) and iv) geometry and layout (two tanks storage with physical separation between hot and cold medium or thermocline storage with thermal separation within a single tank). In this work, we focus on sensible HTTES with a packed bed where air is used as transfer fluid and rocks are selected as storage medium, classified as a thermocline thermal storage [START_REF] Mertens | Dynamic simulation of integrated rock-bed thermocline storage for concentrated solar power[END_REF]. The choice of this system among other HTTES systems is for its low cost and advanced research and development achieved during the last years [START_REF] Arteconi | State of the art of thermal storage for demand-side management[END_REF]. In this system, the charge phase consist of blowing hot air through the porous medium. The air circulates in direct contact with the solid and convective heat transfer between the air and the solid is the predominant heat transfer mechanism. As a result, the temperature of the porous medium along the bed increases and a temperature gradient (thermocline) appears as shown in Fig. 1. During discharge, the energy is recovered by cold air circulating in the opposite direction through the porous medium. The air recovers the heat and warms up in contact More details on thermocline systems can be found in [START_REF] Mertens | Dynamic simulation of integrated rock-bed thermocline storage for concentrated solar power[END_REF].

with the hot solid. Typically, in order to charge the system, the air flow rate and entry temperature are regulated to achieve the desired thermal power 𝑃 𝑠𝑡𝑜 . For the discharge phase, the air outlet temperature is used to calculate 𝑃 𝑠𝑡𝑜 . The storage energy is called 𝐸 𝑠𝑡𝑜 and can be estimated by applying a thermal energy balance on the storage materials. The losses 𝑃 𝑙𝑜𝑠𝑠 of the system are composed of heat losses at the walls and during the charging phase. Indeed, during the charge phase, the thermocline moves towards the outlet section. Consequently, the outlet temperature increases when the thermocline approaches the storage end and a part of the charge power is lost as heat in exhaust air as shown in Fig. 1. In a multi-energy network, the HTTES system needs to be operated by a predictive controller that determines the optimal charging and discharging power, 𝑃 𝑠𝑡𝑜 , at each time step by solving an optimization problem along the time horizon. A model is required for the controller to estimate the evolution of the storage energy, 𝐸 𝑠𝑡𝑜 , and heat losses, 𝑃 𝑙𝑜𝑠𝑠 , which are essential to solve the optimization problem. In fact, calculating a feasible optimized command along a time horizon requires anticipation of future events. The controller model is necessary for estimating the system dynamic which links current state (storage energy) to its future value, with respect to the physical constraints and according to command variables. Additionally, the model can also provide cost function terms such as economic storage losses. Finally, the optimization model needs to be accurate but fast in order reach the best performance in terms of cost function score and calculation time. In the literature, a large variety of HTTES models is discussed and increasing complexity can be observed. Among other criteria, most time-dependent models can be classified according to their dimensions which range from 0D (no spatial dependency, only time dependent) to 2D models [START_REF] Marongiu | Modeling of high temperature thermal energy storage in rock beds -experimental comparison and parametric study[END_REF] [START_REF] Xu | Analytical considerations of local thermal non-equilibrium conditions for thermal transport in metal foams[END_REF][16], derived from the numerical solution of partial differential equations (PDE) in two space dimensions, which reflect energy conservation within the system. Among the PDE models, 1D models solving PDE according to the axial direction of the storage are most commonly used in the literature. These models are accurate but require significant calculation times. In general, they are used for numerical simulation purposes in order to accurately estimate the technical performance of the storage for a few operating cycles. For example, in a study of three hot water storage models of different complexity level in [START_REF] Campos Celador | Implications of the modelling of stratified hot water storage tanks in the simulation of CHP plants[END_REF], the authors conclude that the modeling of heat losses and stratification phenomena has a strong impact on the economic feasibility of projects. However, no comparative optimization work was carried out to study the impact of those models on the optimal management controller performance. In [START_REF] Hoffmann | A thermocline thermal energy storage system with filler materials for concentrated solar power plants: Experimental data and numerical model sensitivity to different experimental tank scales[END_REF], the author developed 1D models solving the heat equations on the phases constituting the storage, namely the solid (rock), liquid (oil) and solid (wall) phases. The evaluation of the state of the storage and the heat losses via these models is relatively time-consuming. These models were compared and validated using a prototype with a storage capacity of a few kWh th . The same model has been adapted for a rock bed HTTES where it is assumed that the air is homogeneously distributed via the diffuser without radial effect [START_REF] Hoffmann | Experimental and numerical investigation of a thermocline thermal energy storage tank[END_REF]. [START_REF] Esence | Étude et modélisation des systèmes de stockage thermique de type régénératif solide/fluide[END_REF] suggested an analytical model capable of estimating the state of thermal storage, without going through the resolution of the PDE. This model significantly reduces the computation time compared to other 1D models. However, the author concludes that there are strong limitations to the analytical modeling of such a system since it is limited to partial cycle operation without extraction of the thermocline. Moreover, such a model is not able to estimate heat losses during the load phase and neglects mixing and diffusion processes during rest periods. For the optimization of thermal storage, two categories are identified in the literature, namely management optimization and sizing. Most of these works use simplified 0D models in order to get rid of the computation time limitation. [START_REF] Cox | Real time optimal control of district cooling system with thermal energy storage using neural networks[END_REF], [START_REF] Jørgensen | Economic mpc for a linear stochastic system of energy units[END_REF] and [START_REF] Verrilli | Model predictive control-based optimal operations of district heating system with thermal energy storage and flexible loads[END_REF] suggest a Model Predictive Control (MPC) solving the stochastic optimization problem in order to calculate the management of a thermal storage coupled to a thermal source and heat network. The storage is represented by a 0D model including a linear term for the estimation of heat losses to the ambient. In general, the authors do not differentiate the model used for optimization from the model used for performance simulation. They conclude that this type of model is suitable for the optimization of large-scale systems. Moreover, the design optimization of a multi-energy system at the scale of a region was the subject of the studies [START_REF] Limpens | Energyscope td: A novel open-source model for regional energy systems[END_REF] and [START_REF] Gabrielli | Optimal design of multi-energy systems with seasonal storage[END_REF]. The models involved in this tool are of type 0D in order to guarantee, here again, reasonable calculation times.

Scope of this work

The goal of this work is to propose a thermal storage meta model with a suitable trade-off between accuracy and calculation time for future real-time management and design optimization of thermocline storage systems. In contrast to the literature [START_REF] Esence | Étude et modélisation des systèmes de stockage thermique de type régénératif solide/fluide[END_REF][26] [START_REF] Bayón | Study of thermocline tank performance in dynamic processes and stand-by periods with an analytical function[END_REF], this approach has been applied in a novel way including an extensive validation with real-life installations. Such a model can be used inside a controller to optimize the commands of the storage inside a complex environment such as a multi-energy network with strong renewable energy penetration. This type of model is less accurate than physical models proposed by [START_REF] Marongiu | Modeling of high temperature thermal energy storage in rock beds -experimental comparison and parametric study[END_REF][15] [START_REF] Xu | Thermal transport in microchannels partially filled with micro-porous media involving flow inertia, flow/thermal slips, thermal non-equilibrium and thermal asymmetry[END_REF], which solve partial differential equations and are used for accurate simulation, but gives local information within acceptable calculation time for use inside a storage controller. Indeed, the review of the state of the art demonstrates that the modeling of a thermal storage is essentially on two levels: simplistic models with short calculation times, and complex numerical models with high accuracy and higher dimensions representation. Those two criteria are essential, in our opinion, for co-optimization of an energy system with storage. It is then necessary to study the impact of the optimization model on time and accuracy. Therefore, the purpose of this contribution is to validate an intermediate level of storage modeling adapted to co-optimization studies. In fact, a new approach for thermal storage meta modeling is proposed and compared to existing models in the literature. Moreover, four experimental configurations with different materials, geometries and operations are selected for the meta model validation. For each installation, the meta model simulation results are compared to the physical model as well as the experimental results. Here, material temperature evolution, the state of charge and the outlet temperature are evaluated during charge and discharge phases. Finally, it is shown that a meta model based on experimental results instead of accurate simulation's data is able to describe the dynamics of a storage system with unique characteristics, for example caused by a manufacturing defect.

Sec. 2 presents the main equations governing the models selected in the literature as well as the main contribution of this paper which is the elaboration of a meta model. Sec. 3 describes the studied test facilities and the experimental protocol as well as the instrumentation. In Sec. 4, the comparison between simulation and experimental results are shown and analysed. This includes a proposed alternative for meta modeling based on experimental results instead of physical model simulation results. Finally, the conclusions and perspectives of this work are presented in Sec. 5.

Modeling

The precise operation management of a thermal storage requires the estimation of the temporal evolution of the energy state 𝐸 𝑠𝑡𝑜 and the thermal losses 𝑃 𝑙𝑜𝑠𝑠 . These losses are mainly due to heat losses to the ambient through the walls as well as the loss of heat contained in the warm exhaust air while charging. The estimation of the temperature evolution within the storage is necessary in order to avoid rough assumptions on the storage temperature (hypothesis of 0D models). This Sec. presents two 0D models, one 1D model based on partial differential equations and one 1D meta model thoroughly developed in this work. 2D models are excluded from the scope of this work because of their high time costs.

Ideal Storage: 0D-Ideal

This model assumes an ideal storage without losses at the walls and with a perfect separation between the hot and cold zones. In this case, the temperature gradient is supposed to be equal to infinity, meaning a thermocline of zero thickness. The dynamics of the system according to this model can be expressed as follows:

𝐸 𝑡+Δ𝑡 𝑠𝑡𝑜 = 𝐸 𝑡 𝑠𝑡𝑜 + Δ𝑡 ⋅ 𝑃 𝑡 𝑠𝑡𝑜 (1)
Storage heat losses are zero at all times except when the system is charged even though it is already fully charged. In that case, losses are constant and equal to the charging power. This model underestimates the losses in general but has a fast calculation time, quantitatively lower than 1 ms per iteration.

Uniform temperature storage: 0D-Uniform

In this model, storage is considered, in its entirety, at uniform temperature 𝑇 𝑚𝑒𝑎𝑛 and the temperature gradient is equal to zero, meaning that the thermocline is completely horizontal and its thickness is equal to the storage length. The evolution of the storage energy as well as the power losses during the charge phase can be calculated with Eq. 2a and 2b, where 𝑇 𝑚𝑒𝑎𝑛 is updated at each time step using Eq. 2c. As explained in the introduction Sec. 1, the losses of such thermal storage is composed of two aspects. Here, power loss 𝑃 𝑙𝑜𝑠𝑠 is the power which is not transferred from the air to the storage materials during the charging phase. This power is lost in the air leaving the storage tank to the outside. Therefore, 𝑃 𝑙𝑜𝑠𝑠 is proportional to the air mass flow ṁ and the storage outlet temperature (here equal to 𝑇 𝑡 𝑚𝑒𝑎𝑛 ). This model does not consider the second aspect of loss, the loss at walls, unlike other models such as the model of Sec. 2.3. It is to be noted that the 0D-Uniform model overestimates storage loss at the outlet section. The computation time of 0D models is low, which explains their frequent use in energy system optimization work. However, these models poorly estimate the temperature within the storage (for example the storage outlet temperature for calculating losses, efficiency when coupled to conversion systems in multi-energy networks or aging of equipment). Therefore, 1D models are needed when system optimization depends on local information for finer objective functions.

𝐸 𝑡+Δ𝑡 𝑠𝑡𝑜 = 𝐸 𝑡 𝑠𝑡𝑜 + Δ𝑡 ⋅ 𝑃 𝑡 𝑠𝑡𝑜 -Δ𝑡 ⋅ 𝑃 𝑙𝑜𝑠𝑠 (𝑇 𝑡 𝑚𝑒𝑎𝑛 , ṁ) (2a) with 𝑃 𝑡 𝑙𝑜𝑠𝑠 = ṁ𝑐 𝑝 (𝑇 𝑡 𝑚𝑒𝑎𝑛 -𝑇 𝑒𝑥𝑡 ) (2b)
and

𝑇 𝑡 𝑚𝑒𝑎𝑛 = 𝐸 𝑡 𝑠𝑡𝑜 𝑐 𝑝 𝑚 + 𝑇 𝑒𝑥𝑡 (2c)
where 𝑇 𝑚𝑒𝑎𝑛 in the storage mean temperature, 𝑇 𝑒𝑥𝑡 is the ambient temperature, ṁ is the air mass flow, 𝑚 is the mass of storage materials and 𝑐 𝑝 is the mass specific heat .

Accurate Model: 1D-PDE

This model solves the heat equation on the three phases constituting the storage as shown below by considering only one spatial dimension, the axial variation of the temperature. Unlike 2D models, this model supposes homogeneous radial air distribution without considering radial effects like mixing.

• Fluid :

𝜀 ( 𝜌𝑐 𝑝 ) 𝑓 ( 𝜕𝑇 𝑓 𝜕𝑡 + 𝑢 𝜕𝑇 𝑓 𝜕𝑥 ) = (3)
𝑘 eff 𝑓 𝜕 2 𝑇 𝑓 𝜕𝑥 2 + ℎ 𝑣 (𝑇 𝑠 -𝑇 𝑓 ) + ℎ 𝑤 𝐴 𝑓 ↔𝑤 𝑉 𝑓 + 𝑉 𝑠 (𝑇 𝑤 -𝑇 𝑓 )
• Solid :

(1 -𝜀) ( 𝜌𝑐 𝑝 ) 𝑠 𝜕𝑇 𝑠 𝜕𝑡 = ( 4
)
𝑘 eff 𝑠 𝜕 2 𝑇 𝑠 𝜕𝑥 2 + ℎ 𝑣 (𝑇 𝑓 -𝑇 𝑠 ) + ℎ 𝑤 𝐴 𝑠↔𝑤 𝑉 𝑓 + 𝑉 𝑠 (𝑇 𝑤 -𝑇 𝑠 )
• Wall :

( 𝜌𝑐 𝑝 ) 𝑤 𝜕𝑇 𝑤 𝜕𝑡 = 𝑘 𝑤 𝜕 2 𝑇 𝑤 𝜕𝑥 2 + ℎ ext 𝐴 𝑤↔ext 𝑉 𝑤 (𝑇 ext -𝑇 𝑤 )+ ℎ 𝑤 ( 𝐴 𝑓 ↔𝑤 𝑉 𝑤 (𝑇 𝑓 -𝑇 𝑤 ) + 𝐴 𝑠↔𝑤 𝑉 𝑤 (𝑇 𝑠 -𝑇 𝑤 ) ) (5)
where 𝜀 is the void fraction, 𝜌 is the mass density, 𝑢 is the local fluid velocity, 𝑘 𝑒𝑓 𝑓 is the effective conductivity, 𝑘 𝑤 is the wall effective conductivity, ℎ 𝑣 is the volume heat transfer coefficient, ℎ 𝑤 is the wall heat transfer coefficient, ℎ 𝑒𝑥𝑡 is the heat transfer coefficient with the ambient, 𝐴 is the area of exchange between two media, and 𝑉 is the volume of the medium. The numerical resolution of these equations, by the finite difference method, gives the spatio-temporal evolution of the temperature. All empirical parameters from Eq. 3, 4 and 5 are calculated as expressed in reference [START_REF] Hoffmann | A thermocline thermal energy storage system with filler materials for concentrated solar power plants: Experimental data and numerical model sensitivity to different experimental tank scales[END_REF]. With the help of an energy balance, it is possible to calculate the storage energy as shown in Eq. 6, the thermal losses at the walls and the load losses (the latter are proportional to the storage outlet temperature).

𝐸 sto = (1 -𝜀)𝑆 bed 𝜌 s ∫ 𝐿 0 𝑐 p (𝑥) ( 𝑇 s (𝑥) -𝑇 ext ) d𝑥 (6)
where 𝑆 bed is the section of the rock bed, 𝜌 s the density of rocks and 𝑐 p is the specific heat capacity of the rock. The computation time of the model is relatively long, quantitatively 1-10 s to simulate 1 h of operation. Therefore, the model is adapted for short time horizon simulation purposes and not for optimization. More details on considered models calculation time will be provided in Sec. 2.4.3.

Meta model: 1D-MM

This subsection will present the main contribution of this work which is the meta modeling approach based on a 1D-PDE model. This model consists of an analytical approximation of the temporal-spatial temperature distribution 𝑇 𝑡 (𝑥). It takes two input values which are the command 𝑃 𝑠𝑡𝑜 and the initial storage temperature curve 𝑇 𝑡 (𝑥) in order to give the final storage temperature curve 𝑇 𝑡+Δ𝑡 (𝑥). Afterwards, both the state of the storage and the heat losses can be calculated by an energy balance. The operation of the meta model is illustrated in Fig. 2. Its use, once built, is fast because it avoids solving costly PDE during the simulation phase. Consequently, calculation times are low, quantitatively 0.001 s to simulate 1 h of storage operation. As a result, the model is suitable for use in optimization and simulation algorithms over large time horizons.

1D MM

Operation of the meta model

This model operates in three steps:

• generally, in 1D models, the spatial temperature distribution is a vector of size 𝑛 𝑥 (𝑛 𝑥 = 100 in [START_REF] Hoffmann | A thermocline thermal energy storage system with filler materials for concentrated solar power plants: Experimental data and numerical model sensitivity to different experimental tank scales[END_REF]). In order to reduce the size of the problem, we use the logistic function along with a least square algorithm in order to approximate the initial temperature curve 𝑇 𝑡 (𝑥) with only four parameters (see Eq. 7). It is therefore possible to calculate a quadruplet of 𝑇 𝑚𝑖𝑛 , 𝑇 𝑚𝑎𝑥 , 𝑧 𝑐 and 𝑠 for each temperature curve as shown in Fig. 3

• the initial state 𝑋 𝑡 along with the power command 𝑃 𝑠𝑡𝑜 is then used to estimate the final state 𝑋 𝑡+Δ𝑡 by linear interpolation using the interpolation matrix built in advance. This matrix maps the final state evolution according to different initial states and power command values and its size depends on the desired precision. Therefore, the matrix construction can be time-costly and is performed as explained in 2.4.2. This step output is the final state 𝑋 𝑡+Δ𝑡 ;

• the use of the logistic function again allows to calculate the final temperature distribution 𝑇 𝑡+Δ𝑡 (𝑥) based on the final state 𝑋 𝑡+Δ𝑡 . This temperature curve is then used to estimate the state of the storage and the heat losses by a simple energy balance.

Storage bed length [m] Storage temperature [°C] 1D-PDE LF + LS LF + LS 1D-PDE

Construction of the meta model

As explained in the previous Sec., the meta model second operating step uses an interpolation matrix that should be constructed in advance. The aim of this subsection is to explain the algorithm used to build such a matrix. The role of the interpolation matrix is to map the final state 𝑋 𝑡+Δ𝑡 evolution according to the initial state 𝑋 𝑡 and power command 𝑃 𝑡 𝑠𝑡𝑜 . To do this, the state 𝑋, which depends on 𝑇 𝑚𝑖𝑛 , 𝑇 𝑚𝑎𝑥 , 𝑧 𝑐 , 𝑠 is discretized on a four dimensional grid composed of the vectors 𝑻 𝒎𝒊𝒏 , 𝑻 𝒎𝒂𝒙 , 𝒛 𝒄 , 𝒔. Each vector has size 𝑖 so 𝑋 has 𝑖 4 possible values. The power command 𝑃 𝑠𝑡𝑜 is also discretized to give 𝑷 𝒔𝒕𝒐 of size 𝑗. In this study, 𝑖 = 10 and 𝑗 = 20 are selected, which gives a model based on a total of 210 000 simulations. The size of discretization vectors determines the matrix size and more discretized vectors means higher interpolation quality. It is therefore important to choose the appropriate matrix size for best model precision and calculation time ratio. It is to be mentioned that in this work, the discretization is done on an equidistant grid.

Once the discretization has been achieved, for each initial state and power command, the logistic function is used to calculate an initial temperature distribution 𝑇 𝑡 . A 1D-PDE model simulation, with the investigated time period Δ𝑡, is then performed. The choice of Δ𝑡 for the model time step simulation is determined by the seeked controller time resolution (here Δ𝑡 = 1 h). Later on, the final temperature distribution 𝑇 𝑡+Δ𝑡 is obtained and approximated to 𝑋 𝑡+Δ𝑡 using the logistic function and a least square algorithm. The whole results are then saved into a matrix which will be used later on for linear interpolation. The algorithm 1 resumes the procedure for the matrix building.

Calculation time of the meta model

As discussed, the model needs to be both fast and accurate in order to achieve a good performance of the optimal controller. In this work, two 0D models have been discussed and showed limitations when estimating local variables such as storage axial temperature. This temperature can be essential for certain optimization problems, which is why 0D models are not used in this work. The time costs of all discussed models are summarized in Tab. 1 Intel ® Core i7-6820HQ CPU @ 2.70GHz

System description

The experimental validation in Sec. 4 is based on three HTTES facilities shown in Fig. 4a, 5a and 6a of which one facility was tested with two different storage material sizes. An overview of the resulting four different set ups, referred to as configurations in the following, is provided in Tab. 2. The choice of different HTTES configurations makes it possible to validate models presented in Sec. 2 with a large variety of conceptual as well as operational parameters. Being classified as unpressurized gas/solid packed-bed storages, all systems use atmospheric air as the heat transfer fluid (HTF) and solids as heat storage material (HSM). While configurations 1, 3 and 4 are characterized by a cuboid geometry and a horizontal air flow direction, configuration 2 was designed with a droplet-like geometry for a vertical air flow direction. The maximum storage energy capacity 𝐸 max ranges from 450 kWh th in configuration 4 to 2900 kWh th in configuration 1, calculated for an ambient temperature of 0 °C. It should be noted that the storage capacity increases for higher heater temperatures or lower ambient temperatures and does not include the contribution of air. Further differences between the configurations highly relevant for the validation in Sec. 4 are the different particle diameter 𝐷 P of HSM ranging from 8 to 30 mm as well as the charge and discharge power 𝑃 CH/DC between 27.5 and 320 kW.

Configuration 1: Eco-Stock ®

This horizontal HTTES, destined for industrial applications, has been developed by the company Eco-Tech Ceram. It includes an electrical heater shown on the left side of the process flow diagram in Fig. 4a. The storage in red is composed of a container that encapsulates a part of the pipes where air circulates into and out from the storage medium. In addition, the container includes the storage porous medium composed of bauxite materials and the insulation layers to limit the heat losses through walls. The Eco-Stock ® is equipped with a total of 33 temperature sensors of which 9 are along the central axis of the unit, as shown in Fig 4b . This study will be based on the values recorded by the 9 sensors located on the central axis in order to avoid wall effects which are especially important for small-scale [START_REF] Marongiu | Modeling of high temperature thermal energy storage in rock beds -experimental comparison and parametric study[END_REF], [START_REF] Touzo | Experimental and numerical analysis of a packed-bed thermal energy storage system designed to recover high temperature waste heat: an industrial scale up[END_REF] and [START_REF] Knobloch | Degradation of a Rock Bed Thermal Energy Storage System, under review[END_REF]. Please note that 4 % of rock mass in configuration 3 is characterized by 𝐷 p = 16-22 mm. systems. Further information for configuration 1 can be found in [START_REF] Touzo | Experimental and numerical analysis of a packed-bed thermal energy storage system designed to recover high temperature waste heat: an industrial scale up[END_REF]. on material-and system-level is also analyzed by Knobloch et al. [START_REF] Knobloch | Degradation of a Rock Bed Thermal Energy Storage System, under review[END_REF]. The combination of a cut-off conical top and hemispherical bottom housing gives the HTTES a dropletlike shape, hence the name Droplet. One novel approach is that the heaters as well as inlet and outlet pipes are located on top of the storage to simplify maintenance and allow underground installation of the rocks but with all equipment above the storage material. By employing two fans, air enters and exits the rock bed in the vertical direction during charge and discharge. This flow scheme uses natural stratification to its advantage, as the hottest rocks are located at the top of the storage. A pipe inside the rock bed makes it possible to reverse the flow direction for charge and discharge in order to obtain a flat thermocline. This inner pipe acts as an outlet during charge and inlet during discharge, in both cases for relatively cold air. The Droplet rock bed has a volume of 3.2 m 3 and is filled with a total rock mass of 5394 kg. Whereas the rocks below the inner pipe (219 kg) are characterized by sizes between 16 and 22 mm, rocks between 8 and 11 mm are used for the main part of the rock bed (5175 kg). Particle sizes are determined by sieving through corresponding meshes, meaning that rocks passed through 22 mm and 11 mm openings, respectively, but can be longer in some direction due to their irregular shape. Five different hemispherical and conical housings, consisting of steel, Superwool, Rockwool, bricks and concrete, cover the rock material. These insulations not only reduce the heat losses to the ambient but also support the HTTES structure by balancing mechanical forces. To log experimental data, a data acquisition is installed, presented in Fig. 5b. In total, the main system components including fans and heaters are complemented by 53 temperature sensors, 8 strain gauges, 2 flow meters, 3 pressure sensors and 1 energy meter. For this study, the length of the bed is defined from the storage top to storage bottom and one temperature per vertical length step is estimated by averaging the 5 temperature sensors in central radial position (see Fig. 5b) from 3 directions (E, SW, NW) for each corresponding height. Consequently, 15 temperature sensors in the center of the air flow are considered, as listed in Tab. 2.

Configuration 2: Droplet

Configurations 3 & 4: Shoebox

While Soprani et al. published design and experimental results of configuration 3 in 2019 [START_REF] Soprani | Design and testing of a horizontal rock bed for high temperature thermal energy storage[END_REF], Marongiu et al. investigated, among others, the rock size in his experimental comparison and parametric study for configuration 3 and 4 in 2019 [START_REF] Marongiu | Modeling of high temperature thermal energy storage in rock beds -experimental comparison and parametric study[END_REF]. The only difference between configuration 3 and 4 in this work is the particle diameter of the HSM. While configuration 4 is characterized by particle diameters of 8-11 mm, configuration 3 contains particles with a diameter between 16-22 mm. The used facility, with its process flow diagram in Fig. 6, includes the 1.5 m 3 rock bed and was designed as a simple construction pilot for conceptual testing (rock sizes, layers to suppress natural convection) as well as operational (heater temperature, flow rate, etc.). The cuboid rock bed is filled with 2495 kg Swedish diabase leading to an estimated porosity of 0.45. It is surrounded by an insulation layer and then enclosed in an outer shell of 4 mm thick steel, which guaranteed an airtight inner environment thanks to a high-temperature sealing. A mix of high-temperature resistant bricks and Superwool are used for thermal insulation. In order to structurally support the rock bed, hard insulation is required on the bottom and the sides. By mounting Superwool inside the storage directly in contact with the rocks, porosity at walls and hence flow channeling due to preferential airflow channels are reduced. The top of the rock bed was insulated with 200 mm of soft thermal insulation. An additional 100 mm layer of soft insulation was installed on the outside surface of the steel housing. Two vertical stainless steel grids with 6 mm holes contained the rocks at the hot and cold ends and separated them from entrance and exit regions for the distribution and collection of the airflow. The flow areas before and after the rock bed were also insulated with two 100 mm insulation layers. In order to test different flow settings, two inlet and two outlet ports are employed. Insulated steel pipes with a diameter of 200 mm connect the rock bed to the fan and heaters (two 15 kW heaters arranged in parallel, each connected to a PID temperature controller). Besides the energy meter and temperature sensors at several locations outside the rock bed for energy balances, temperate sensors are also installed inside the rock bed as shown in Fig. 6b. As explained by the legend, temperature sensors aligned in the middle of the rock bed are marked with dots, while triangles represent temperature sensors close to the wall with a distance of 50 mm from it, all on the same side. However, no temperature sensors were installed directly at the walls or in the corners of the flow cross section to ensure that no flow channeling occurrs. Similar to the approach used in configuration 1 and 2, one central temperature per length step is estimated for comparison purposes. Here, row 3 in Fig. 6b is selected as the most suitable selection of temperature sensors since it is the only row which does not include sensors installed close to the storage walls and it is close to the central axis of the horizontal air flow direction. Therefore, all available sensors from row 3 are used for this work, meaning that 3 temperature sensors are considered for configuration 3 (column 1, 3 and 4 due to a defect in column 2) and 4 temperature sensors are considered for configuration 4 (column 1 to 4). This approach includes uncertainties due to the fact that the selected sensor row is not exactly in central position but still allows a comparison to the other configurations since sensors affected by significant wall effects are avoided and inevitable uncertainties stemming from the (dis-) charging with two pipes are considered as more influential. However, Sec. 4.2 investigates a meta model approach which is promising in terms of satisfactorily depicting such challenging HTTES installations without downturns in calculation time.

Meta model validation

The meta model 1D-MM is compared to the physical model 1D-PDE and validated on the basis of four TES configurations (see Tab. 2). The experimental results consist of a charging and discharging cycle where data has been recorded each second by means of the temperature sensors as described in Sec. 3. Please note that just one temperature per bed length step is evaluated in all investigated TES configurations and a homogeneous radial air distribution is assumed.

Meta modeling based on 1D-PDE model

For each configuration, the 1D-PDE has been calibrated taking into account the bed geometry, storage medium and other physical properties of the concerned configuration as mentioned in Tab. 2. Later on, the 1D-MM is constructed based on 210 000 simulations from the 1D-PDE, covering all initial states and power commands as described in Sec. 2.4. In order to calculate the thermal storage energy and other physical quantities essential for its management optimisation, storage material temperature needs to be estimated correctly.

Fig. 7 shows storage temperature evolution during a charging cycle for the 4 configurations described in Tab. 2. Results illustrate the experimental measurements, the 1D-PDE and the 1D-MM results. The y-axis represents the temperature of the solid material (storage medium) 𝑇 𝑠 along the bed length (x-axis). Each curve represents the storage temperature measurements for a given time. The charging phase lasts 15 hours for the Shoebox as well as the Droplet and 6 hours for the Eco-Stock ® facility. The curves are plotted each 3 hours for the Shoebox and Droplet while hourly values are given for the Eco-Stock ® facility. In most optimisation problems addressing thermal storage management and sizing, storage energy and outlet temperature evolution are of primary importance in order to estimate the storage dynamic and losses. Fig. 8 shows both storage energy (black) and outlet temperature (blue) evolution during the charging phase for the four studied configurations. First y-axis represents the storage energy 𝐸 𝑠𝑡𝑜 evolution as a function of charging time (x-axis), while the second yaxis represents the storage outlet temperature progress as a function of charging time (x-axis).

For the storage energy 𝐸 𝑠𝑡𝑜 estimation, the storage temperature 𝑇 𝑠 is being used as expressed in Eq. 6. Tab. 3 presents the mean absolute error (MAE) and the normalized root mean square deviation (NRMSD) of the storage temperature, the outlet temperature and the storage energy, based on experimental results stemming from the four investigated configurations during charging and discharging phases. MAE is expressed in • C for temperatures and in kWh for storage energy errors while the NRMSD is expressed in %. The MAE for all studied criteria is calculated as the absolute difference between the two compared results averaged on the whole simulation period as following:

MAE 𝑇 𝑠 = 𝑛 ∑ 𝑖=0 | | | 𝑻 𝒔,𝒊 𝑚𝑜𝑑𝑒𝑙 -𝑻 𝒔,𝒊 𝑒𝑥𝑝 | | | 𝑛 (8a)
. Table 3: Error evaluation of studied models based on experimental results issued from the four investigated facilities. * First value of error corresponds to the mean absolute error calculated according to Eq. 8. * * Second value of error corresponds to the normalized root mean square deviation calculated according to Eq.9

Models

MAE 𝐸 𝑠𝑡𝑜 = 𝑛 ∑ 𝑖=0 | | | 𝐸 𝑚𝑜𝑑𝑒𝑙 𝑠𝑡𝑜,𝑖 -𝐸 𝑒𝑥𝑝 𝑠𝑡𝑜,𝑖 | | | 𝑛 (8b)
where 𝑻 𝒔 𝑚𝑜𝑑𝑒𝑙 and 𝑻 𝒔 𝑒𝑥𝑝 are the storage temperature vectors issued from the model and from the experiments, respectively, as shown in Fig. 7. 𝐸 𝑚𝑜𝑑𝑒𝑙 𝑠𝑡𝑜 and 𝐸 𝑒𝑥𝑝 𝑠𝑡𝑜 are the storage energy estimations issued from model and experiments, respectively, as shown in black in Fig. 8. 𝑛 is the number of iterations used to calculate the mean error. For the NRMSD calculation, Eq.9 is used as defined in [33][34]. This metric has been chosen in order to avoid scaling problems, that may occur for low energy and temperature values, with mean relative error (MRE) metric (see the Appendix).

√ √ √ √ 𝑛 ∑ 𝑖=0 (𝑥 𝑒𝑥𝑝 𝑖 -𝑥 𝑚𝑜𝑑𝑒𝑙 𝑖 ) 2 𝑛 ( 9b
)
where 𝑥 𝑚𝑎𝑥 and 𝑥 𝑚𝑖𝑛 are the maximum and minimum in the investigated experimental values respectively. In general, results show good agreement between the two investigated models for all TES configurations. Indeed, the NRMSD, in terms of storage temperature, storage energy and outlet temperature, remains lower than 4% for all configurations when comparing the meta model 1D-MM with the physical model 1D-PDE. Consequently, the 1D-MM model faithfully represents the 1D-PDE model that is used for its construction. When comparing the model and experimental results, higher deviations are observed, especially in configuration 4 (Shoebox with 𝐷 p = 11 mm). As summarized in Tab.

3, experiments from configuration 4 differ strongly from 1D-MM and 1D-PDE for both charge and discharge phase, which is mainly because the general model assumption of homogeneous radial air distribution is highly questionable here due to i) the fact that used sensors are not exactly from the central axis (see Sec. 3.3), ii) the (dis-) charge with two instead of one pipe as in configuration 1 and 2, iii) potentially larger insulation degradation due to longer operation time before data acquisition as well as iv) significant wall effects since the rock bed is directly embedded in soft insulation. Furthermore, the configuration 2 (Droplet) is characterized by a vertical air flow direction which uses the buoyancy effects to its advantage, meaning that the buoyancy force does not affect the radial distribution since it works in (discharge) or 180°against (charge) air flow direction. Therefore, higher outlet temperature and storage energy errors are observed when comparing the 1D-PDE model to experiments for configuration 4. For example, storage temperature and outlet temperature MAE for the configuration 4 exceeds 30 °C (>5%) during charging phase. Also outlet temperature MAE exceeds 50 °C (>10%) during the discharging phase in the same facility. In fact, it can be suspected that this is due to the presence of preferential air passage zones with lower density thanks to higher temperature. This results in an increasing non-homogeneity of the air distribution through the cross section during the operational phase. Consequently, the temperature evolution does not respect the assumption of homogeneous radial air distribution used by the 1D-PDE, which explains the difference between the experimental and numerical results. The same observation can be made for the meta model 1D-MM which follows the trend of the physical model 1D-PDE. It can be explained by the fact that 1D-MM is based on 1D-PDE. Therefore, it commits errors where the physical model 1D-PDE itself is not accurate. To address this limitation, a meta modeling approach based on experiments (without solving physical equations) is being explained and discussed in the next subsection.

Meta modeling based on experimental data

The goal of this Sec. is to propose an alternative for PDEbased models. Indeed, storage system dynamics may differ from hypothesis assumed by 1D-PDE models in many cases. Especially the experimental results from configuration 4 (Shoebox 𝐷 p = 11 mm) are characterized by a significant deviation from 1D-PDE and 1D-MM model results, as shown in the previous subsection. In this case, two solutions are possible; 2D-PDE models or meta models based on experiments. Here, meta modeling based on experiments approach is chosen due to the low calculation time and the availability of suitable data sets.

Impact from experimental protocol

This experiment-based approach can be complex and time-consuming when conducted tests are numerous and require further infrastructure. For this reason, this Sec. will focus on the impact of the number of conducted tests and registered measurements frequency to limit the construction phase duration and complexity. First, the TES configuration should be equipped with a minimum of four temperature sensors in order to estimate the unknown parameters of Eq. 7. The relationship between the number of tests to be conducted 𝑁 𝑡𝑒𝑠𝑡𝑠 and the tests power commands can be defined as shown below:

𝑃 𝑖 = 𝑃 𝑚𝑎𝑥 𝑁 𝑡𝑒𝑠𝑡𝑠 𝑖, 𝑖 ∈ 1, ...𝑁 𝑡𝑒𝑠𝑡
with 𝑃 𝑚𝑎𝑥 = 𝐸 𝑚𝑎𝑥 Δ𝑡 here Δ𝑡 = 1h [START_REF] Marongiu | Modeling of high temperature thermal energy storage in rock beds -experimental comparison and parametric study[END_REF] The number of measurements to be collected during a test is governed by the chosen time step between two measurements Δ𝑡 𝑚 . Fig. 9a describes a typical conducted test, showing the evolution of the power command 𝑃 𝑠𝑡𝑜 and the normalized storage energy Ê𝑠𝑡𝑜 according to the test time. Typical tests start with an empty storage which is being charged with the related power command 𝑃 𝑖 until storage saturation ( Ê𝑠𝑡𝑜 = 100 %). Once the storage is saturated, the charging process is considered finished. Afterwards, the discharging phase begins with a power command equal to (b) Meta model storage energy NRMSD evolution depending on measurements frequency, defined by the time period between two measurement Δ𝑡 𝑚 and number of test 𝑁 𝑡𝑒𝑠𝑡𝑠 used for its construction. The NRMSD is calculated by using a Monte Carlo evaluation process, where for each point in the matrix 10000 iterations with randomly chosen initial state and power commands have been operated. -𝑃 𝑖 until the storage is empty again ( Ê𝑠𝑡𝑜 = 0 %). During the whole test, measurements are collected each Δ𝑡 𝑚 , which represents the time step between two measurements.

Based on collected data, the 1D-MM can be constructed with HTTES experimental datasets using the algorithm from 1. In fact, the collected data represents spatial temperature measurements at different times of the conducted test. Each measurement can be considered as a temperature distribution 𝑇 (𝑥) which can be transformed to a state vector 𝑋 by using Eq. 7 and a least square routine. For each measurement, it is possible to associate a final state vector 𝑋 𝑡+Δ𝑡 to an initial state vector 𝑋 𝑡 separated by a time step Δ𝑡 corresponding to the desired model time step. Later on, a multi-dimensional linear interpolation on non regular grid is carried out to build the model [35]. This model is able to estimate final state 𝑋 𝑡+Δ𝑡 depending on a randomly chosen initial state 𝑋 𝑡 and a power command 𝑃 𝑠𝑡𝑜 . As an example, configuration 4 is studied in order to evaluate the performance of the developed model. Due to the lack of experimental datasets for a wide range of power commands, the complexity and the time cost to conduct such tests, the 1D-PDE model is be used to simulate experimental data. These datasets are used to determine the design of the experiments so that an acceptable model error with a minimum of conducted experimental tests is ensured. Once the model is constructed, its accuracy needs to be evaluated over a range of operating points. To do so, a set of operating points is randomly drawn, in a Monte Carlo process. The accuracy of the model is then evaluated by calculating the storage energy NRMSD of the meta model over the whole set of drawn operating points. This error is calculated by applying Eq. 8 and 9 as explained in the last Sec.. In other words, an initial state and a power command, within the tolerable power range of the storage system (here 𝑃 𝑚𝑎𝑥 = 450 kW), will be selected randomly for the performance evaluation. Each simulation of the Monte Carlo process has a fixed time step (here Δ𝑡 = 1h). The number of Monte Carlo evaluation iterations is selected based on the convergence of an error criteria.

Here the error criterion is the NRMSD of 1D-MM estimated storage energy 𝐸 𝑠𝑡𝑜 . The results show that 10000 samples is a suitable number for the error criterion to converge in this case study. Fig. 9b shows the evolution of the storage energy NRMSD for different measurements frequencies Δ𝑡 𝑚 and experimental tests numbers 𝑁 𝑡𝑒𝑠𝑡𝑠 to study the sensibility of 1D-MM regarding the experimental protocol. Axis x describes the number of tests 𝑁 𝑡𝑒𝑠𝑡𝑠 (each test has a unique power command value as explained in Eq. 10) used for the model construction. Axis y describes the time period value separating two collected measurements Δ𝑡 𝑚 (the lower time period between two collected measurements Δ𝑡 𝑚 , the more data used for model construction). Fig. 9b illustrates the non sensitivity of model error regarding the chosen time step between two selected measurements Δ𝑡 𝑚 within the studied range. 8 minutes as a time step for data measurement during conducted tests seems to be a suited value. In other words, data being collected with a frequency higher than 8 minutes does not increase the performance of the studied model. The number of tests necessary for the model construction 𝑁 𝑡𝑒𝑠𝑡𝑠 has a higher impact on the model precision according to Fig. 9b. In fact, for a model constructed on only one experimental test, the model storage energy NRMSD is 15%. On the other hand, starting from three conducted tests, the storage energy NRMSD is lower than 5% and can reach 4% for higher number of experimental tests.

Experiments-based 1D-MM for configuration 4

As said before, physical equation-based modeling approaches are limited in the case of the configuration 4 (Shoebox (𝐷 𝑝 = 11mm)). To address this issue, a meta model was constructed based on a series of experimental measurements using the meta modeling approach based on experience. An example of this case study with real experimental data is shown in Fig. 10. As few experimental data are available, one test 𝑁 𝑡𝑒𝑠𝑡𝑠 = 1 has been selected for model construction and evaluation. Indeed, the same test dataset has been separated into two datasets, one for model construction and another for the validation process. The time step of collected measurements Δ𝑡 𝑚 for model construction in this case is of 8 minutes. The model is validated on the basis of a full charging cycle, where model results will be compared to experiments each three hours. Fig. 10a shows the evolution of the axial storage temperature according to experimental data and results from the proposed meta model. Indeed, the proposed meta model based on experimental data is able to depict the performance of configuration 4. This results in a storage temperature MAE of 1 °C (<1%) when comparing the new 1D-MM results to experimental data. Fig. 10b shows the evolution of the storage energy and outlet temperature according to experimental data and the proposed meta model based on experimental data. A storage energy and outlet temperature NRMSD lower than 2% is observed when comparing the results from the proposed model to the experimental data. Those errors are lower than expected errors according to Fig. 9b. This is explained by the fact that only one experimental test dataset has been used for feeding the model and for testing its performance due to the lack of experimental data issued from this same facility. However, meta modeling based on experimental data shows good adequacy with experimental data in this exemplary study case.

Conclusion

Thermal storage will play an essential role in future coupled energy networks. In order to design and manage such storage systems in an optimal way, models are required to estimate local variables within a reasonable calculation time. Therefore, this work investigates an intermediate level of storage modeling adapted to co-optimization studies: With a NRMSD of 3% compared to the physical model, the meta model shows similar performance as the physical model itself but with a calculation time up to a factor of 1000 smaller. The results obtained with the meta model based on the physical model are successfully validated based on experimental data from four TES configurations, as long as the physical models are suited for the studied configuration. Model limitations that occur for the physical model also apply to the corresponding meta model. This work also demonstrates that if physical models are not suitable for certain TES configurations, for example because of manufacturing defects, limited data acquisition or highly dimension dependent phenomena in specific storage geometries, it is possible to build a meta model on the basis of experimental data. It is indicated that, under a reasoned choice of temperature sensors, number of tests and measurement frequency, a high global performance characterized by a NRMSD between model and experiments below 4% can be achieved despite a low calculation time. Further work is required to investigate the relevance of mentioned models used in model predictive control, where the impact of storage models on optimal energy system management and design should be studied. Appendix A. Error evolution with different metrics as a function of the storage energy (a storage of 450 kWh capacity) estimated by two models. The results stem from a Monte Carlo process with 10 000 samples for a meta model with 𝑁 𝑡𝑒𝑠𝑡𝑠 equal to 5 tests and Δ𝑡 𝑚 equal to 300 s (as shown in Fig. 9b). For each Monte Carlo sample, an initial storage energy and a power command are being selected randomly. Two models are then used to estimate the storage energy evolution and results error are calculated thereafter. The scale for the absolute error represented by gray dots is on the left, while for relative errors (NRMSD and MRE) the scale is on the right.

Figure 1 :

 1 Figure 1: HTTES temperature axial distribution (thermocline) for an initially half charged storage 𝑇 𝑡 0 𝑖𝑛𝑖𝑡𝑖𝑎𝑙 (gray), 1h of 1.5 MW charge 𝑇 𝑡 0 +1ℎ 𝑐ℎ𝑎𝑟𝑔𝑒 (in red), 1h of 1.5 MW discharge 𝑇 𝑡 0 +1ℎ 𝑑𝑖𝑠𝑐ℎ𝑎𝑟𝑔𝑒 (blue).More details on thermocline systems can be found in[START_REF] Mertens | Dynamic simulation of integrated rock-bed thermocline storage for concentrated solar power[END_REF].

Figure 2 :

 2 Figure 2: Meta model Operating Diagram, illustrating the main steps in the calculation of system state evolution and losses. 𝐿𝐹 : Logistic Function, 𝐿𝑆 : Least Square algorithm.

Figure 3 :

 3 Figure 3: Evolution of the storage temperature for 1 MW load during 1 h estimated by the 1D-PDE model and approximated by LF (logistic function) and LS (least square algorithm).

  (a) Process flow diagram during the charging phase (green for cold air and red for hot air) [29]. (b) Location of temperature sensors in side cut view. Sensors considered in this work are circled in red (9 sensors, see Tab. 2). All values are given in mm [29].

Figure 4 :

 4 Figure 4: Configuration 1: Eco-Stock ® -Horizontal HTTES from Eco-Tech Ceram [29].

Figure 5 :

 5 Figure 5: Configuration 2: Droplet -Vertical HTTES from DTU Energy [31].

Figure 6 :

 6 Figure 6: Configuration 3 + 4: Shoebox -Horizontal HTTES from DTU Energy [32].

Fig

  Fig 5a shows the process flow diagram of the vertical HTTES system, here called configuration 2. Knobloch et al. presented the design, a further development of a previous cuboid HTTES (configuration 3 and 4), and testing of configuration 2 [31] in 2022. The degradation of configuration 2

3 ( 4 (

 34 Configuration 1 (Eco-Stock ® ): Charge with 320 kW th . Shoebox with 𝐷 p = 22 mm): Charge with 27.5 kW th . Shoebox with 𝐷 p = 11 mm): Charge with 28.8 kW th .

Figure 7 :

 7 Figure 7: Storage material temperature during charging process for the four investigated TES configurations. The meta model 1D-MM and the physical model 1D-PDE are compared to experimental measurements. Please note that the 1D-MM in this Fig. is based on the physical model 1D-PDE.

3 ( 4 (

 34 Configuration 1 (Eco-Stock ® ): Charge with 320 kW th . Shoebox with 𝐷 p = 22 mm): Charge with 27.5 kW th . Shoebox with 𝐷 p = 11 mm): Charge with 28.8 kW th .

Figure 8 :

 8 Figure 8: Storage energy and outlet temperature during charging processes for the four investigated TES configurations. The meta model 1D-MM and the physical model 1D-PDE are compared to experimental measurements. Please note that the 1D-MM in this Fig. is based on the physical model 1D-PDE.

  Description of the test process needed for meta model based on experimental data. 𝑃 𝑠𝑡𝑜 is the power command used for charging and discharging a storage facility during a test. Ê𝑠𝑡𝑜 is the normalized storage energy, here shown for a full charging and discharging cycle. Black dots represent collected measurements used for meta model construction, the frequency of measurements is defined by Δ𝑡 𝑚 . Time between two measurements[s]

Figure 9 :

 9 Figure 9: Description of the experiment-based meta model test protocol, and the sensibility of the model performance regarding the number of conducted tests and measurements frequency. Here, configuration 4 (Shoebox with 𝐷 p = 11 mm) is studied.

  Storage energy and outlet temperature.

Figure 10 :

 10 Figure 10: Results from the meta model 1D-MM based on experimental data (instead of a physical model) compared to experimental data during a charging process at 28.8 kW th . Here, only configuration 4 (Shoebox with 𝐷 𝑝 = 11 mm) is studied as an example.

Algorithm 1 :

 1 Procedure for the construction of the meta model identification matrix.

	Data: 𝑇 𝑡 𝑚𝑖𝑛 ∈ 𝑻 𝒎𝒊𝒏 , 𝑇 𝑡 𝑚𝑎𝑥 ∈ 𝑻 𝒎𝒂𝒙 , 𝑧 𝑡 𝑐 ∈ 𝒛 𝒄 , 𝑠 𝑡 ∈ 𝒔, 𝑃 𝑡 𝑠𝑡𝑜 ∈ 𝑷 𝒔𝒕𝒐
	Result: 𝑓 ∶ 𝑋 𝑡 ; 𝑃 𝑡 ⟼ 𝑋 𝑡+Δ𝑡
	with 𝑋 𝑡+Δ𝑡 = [𝑇 𝑡+Δ𝑡 𝑚𝑖𝑛 , 𝑇 𝑡+Δ𝑡 𝑚𝑎𝑥 , 𝑧 𝑡+Δ𝑡 𝑐 and 𝑋 𝑡 = [𝑇 𝑡 𝑚𝑖𝑛 , 𝑇 𝑡 𝑚𝑎𝑥 , 𝑧 𝑡 𝑐 , 𝑠 𝑡 ]	, 𝑠 𝑡+Δ𝑡 ]
	initialization;
	foreach 𝑇 𝑡 𝑚𝑖𝑛 , 𝑇 𝑡 𝑚𝑎𝑥 , 𝑧 𝑡 𝑐 , 𝑠 𝑡 do foreach 𝑃 𝑡 𝑠𝑡𝑜 do calculate 𝑇 𝑡 (𝑥) with Eq. 7 ;
	calculate 𝑇 𝑡+Δ𝑡 (𝑥) by solving Eq. 3 4 5 ;
	calculate 𝑇 𝑡+Δ𝑡 𝑚𝑖𝑛 , 𝑇 𝑡+Δ𝑡 𝑚𝑎𝑥 , 𝑧 𝑡+Δ𝑡 𝑐 Eq. 7 and least square algorithm , 𝑠 𝑡+Δ𝑡 with
	end
	end

Table 1 :

 1 [START_REF] Rte | futurs énergétiques 2050" -consommation et production : les chemins de l'électricité de rte pour la neutralité carbone[END_REF], where 𝑁 𝑐𝑜𝑛𝑠 is the number of required simulations and 𝑡 𝑐𝑜𝑛𝑠 is the simulations time cost in order to construct the meta model. 𝑡 𝑠𝑖𝑚 is the model time cost required to simulate 1h of storage operation. 𝑡 𝑜𝑝𝑡 is the time cost to solve a typical thermal storage optimization problem over 24 hours, hour by hour, using a Model Predictive Control algorithm (MPC)[START_REF] Jørgensen | Economic mpc for a linear stochastic system of energy units[END_REF][START_REF] Al Asmi | Impact of coupling thermal and electrical carriers on the optimal management of a multi-energy network[END_REF]. The 1D-MM allows to reduce the time costs by a factor of 1000 or 150 for simulation or control purpose, respectively, compared to the 1D-PDE model. Summary of time costs of studied models.

	Model	𝑁 𝑐𝑜𝑛𝑠	𝑡 𝑐𝑜𝑛𝑠 [s]	𝑡 𝑠𝑖𝑚 [s] 1	𝑡 𝑜𝑝𝑡 [s]
	0D-Ideal	-	-	< 0.001	7
	0D-Uniform	-	-	< 0.001	7
	1D-PDE	-	-	1 -10	≈ 1 000 000
	1D-MM	210 000	207 000	≈ 0.001	6500

Table 2 :

 2 Summary of HTTES configurations used in this work. Physical properties of the HSMs are given in

	Configuration	1	2	3	4
	Facility	Eco-Stock ®	Droplet	Shoebox	Shoebox
	Flow	Horizontal	Vertical	Horizontal Horizontal
	HTF	Air	Air	Air	Air
	HSM	Bauxite	Diabase	Diabase	Diabase
	Diameter (m)	1.7	1.65 (ø)	1	1
	Length (m)	3.08	1.5	1.5	1.5
	𝑉 bed (𝑚 3 ) 𝐷 p (mm)	8.9 30	3.2 8-11	1.5 16-22	1.5 8-11
	𝐸 max (kWh th ) 𝑇 storage,max (°C)	2900 525	1007 600	450 560	450 560
	𝑃 CH (kW) 𝑃 DC (kW) 𝜀 (-)	320 320 0.4	42 58 0.45	27.5 27.5 0.45	28.8 28.8 0.45
	Sensors used	9	15	3	4

Acknowledgements

The authors would like to express their gratitude to the Association Nationale de Recherche et de Technologie (ANRT France 2018/1534) and the Danish Energy Technology Development and Demonstration Program (EUDP 64016-0027) for funding this work.

Declaration of interests

The authors have no competing interests to declare.