**About the Editor**

**Roc´ıo Bay ´on** received an M.Sc. degree in Chemistry at the University of Valladolid and a Ph.D. degree at the University Autonoma of Madrid. She was a postdoctoral fellow at the ´ Hahn-Meitner-Institut of Berlin (Germany) and worked as research engineer at IMRA Europe in Sophia Antipolis (France). At present, she works as senior scientist in the Thermal Energy Storage Unit of CIEMAT-PSA in Madrid. She has been working in different scientific areas related to materials and systems for renewable energy applications such as thin film solar cells, selective absorbers and thermal energy storage. She has participated in many national and international research projects, international conferences and is the co-author of about 40 scientific papers. She collaborates as a teacher in various university masters programs and participates in technological collaboration platforms of the International Energy Agency related to thermal storage.

### *Editorial* **Special Issue "Advanced Phase Change Materials for Thermal Storage"**

**Rocío Bayón**

Thermal Energy Storage Unit, CIEMAT-PSA, Av. Complutense 40, 28040 Madrid, Spain; rocio.bayon@ciemat.es; Tel.: +34-913466048

#### **1. Introduction**

Thermal energy storage using phase change materials (PCMs) is a research topic that has attracted much attention in recent decades. This is mainly because the potential use of PCMs as latent storage media not only covers renewable energy and building efficiency applications, but also the temperature control of electronic devices, batteries and even clothes. Although a number of companies worldwide are producing a variety of PCMs, advanced materials with improved properties and new latent storage concepts are required to better meet the specific requirements of different applications. Moreover, the development of common validation procedures for PCMs is an important issue that should be addressed in order to achieve commercial deployment and implementation of these kinds of materials in latent storage systems.

#### **2. Advanced Phase Change Materials for Thermal Storage**

The key subjects included in this special issue were related not only to materials in terms of new PCM formulations and concepts, validation and assessment procedures, characterization and simulation, but also to PCM applications in terms of implementation and testing in storage prototypes, innovative approaches and the simulation of novel storage modules for latent heat. Despite COVID-19 crises and lockdowns in most countries, there were still six papers submitted to this special issue, and five of them were accepted, which proves the quality of the research and the strong interest in the field of latent heat storage. In the following paragraphs, a summary of these papers with their most relevant contributions is presented.

The first paper included in this issue dealt with a procedure for selecting the appropriate PCM for two kinds of innovative compact energy storage systems implemented in residential buildings: the Mediterranean (MED) concept, intended for space cooling, and the continental (CON) concept, used for heating and domestic water [1]. The selection methodology consists of a qualitative decision matrix, which uses some common features of PCMs to assign an overall score to each material so that different options can be compared. The most important PCM features to be considered in the decision matrix for material selection are the melting enthalpy and temperature range, availability, cost and, in the case of the CON concept, the maximum working temperature range. Apart from the qualitative results, the authors made an experimental characterization of the best candidates, consisting of differential scanning calorimetry (DSC) and thermogravimetric analysis (TGA), before making a final decision. This selection process led to various possible candidates, but two commercial PCMs were selected as the most promising ones: RT4 (Tmelt = 4 ◦C) for the MED cooling system and RT64HC (Tmelt = 64 ◦C) for the CON system providing house heating and domestic hot water.

The main contribution of this paper is that it provides a simple, quick tool for prescreening a PCM before being implemented in any application and selecting at least the most promising candidates to be included in further validation tests.

**Citation:** Bayón, R. Special Issue "Advanced Phase Change Materials for Thermal Storage". *Appl. Sci.* **2021**, *11*, 1390. https://doi.org/10.3390/ app11041390

Received: 25 January 2021 Accepted: 30 January 2021 Published: 4 February 2021


**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the author. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

The second paper presented a methodology that allowed comparing latent heat energy storage (LHES) modules with different designs with respect to their compactness and heat transfer performance [2]. Nowadays, many novel and promising heat exchanger designs and concepts have emerged, aiming to enhance the heat transfer inside LHES devices. However, the wide range of experimental conditions that can be found in the literature for their characterization makes it difficult to compare their performance. In light of this, the methodology described in the paper established just two key performance indicators (KPIs) which were minimally influenced by the experimental conditions: the compactness degree (ΦPCM) and a normalized heat transfer performance coefficient (NHTPC). In the paper, these KPIs were calculated for several LHES units already reported in the literature, allowing a leveled performance comparison with regard to operating conditions at different scales, while characteristics like geometry, structural materials and PCMs remained intrinsic. The robustness of the proposed KPIs was confirmed for units at different size scales by varying the heat transfer fluid mass flow rates and temperature levels. The evaluation procedure was applied to various LHES systems, and the most promising designs for different applications were identified and discussed. The authors clarified that the storage units analyzed were application-oriented and, in most cases, a high heat transfer rate was not a requirement, which led to low values for the KPIs. However optimized versions of the evaluated LHES systems are expected to deliver considerably higher performance indicators. Hence, the conclusions of this analysis should be considered as preliminary pictures of the general potential associated to each technological approach. The main contribution of this paper is that the methodology proposed is expected to open new paths in LHES research by allowing the leveled-ground comparison of technologies among different studies, facilitating the evaluation and selection of the most suitable design or designs for a specific storage application.

The third paper accepted in this special issue introduced the use of PCMs for the thermal management of lithium-ion batteries (LIBs), since temperature is an important factor affecting the working efficiency and service life of these devices [3]. In this work, the authors studied the thermal performance of two commercial batteries (Sony and Sanyo) under different working conditions: extreme conditions (inside a closed aluminum tube), natural convection cooling and PCM cooling with and without heat dissipation fins. The PCM used was a composite of wax with expanded graphite (Tmelt = 52 ◦C), and the experimental results showed that the PCM was able to absorb some of the heat produced during both the charge and discharge processes and, hence, effectively reduce the temperature and keep the battery capacity stable. In fact, the tests performed at different discharge rates showed that the temperature decrease attained under PCM cooling was much higher for the Sanyo LIB (between 4.7 ◦C and 12.8 ◦C) than for the Sony LIB (between 1.1 ◦C and 2 ◦C), in both cases being compared with the natural convection experiments. The temperature reduction impact on the Sanyo LIB was stronger because this battery generated more heat due to its larger storage capacity. As for future developments, the authors suggested that further optimization of LIB thermal management in terms of surface temperature reduction could be achieved if a PCM with a higher latent heat was combined with heat dissipation fins. In my opinion, the most interesting contribution of this paper is that electrical storage and thermal storage working together can improve the performance efficiency of energy storage systems, which is one of the main challenges to be addressed and solved in future energetic scenarios.

The fourth paper resulted from the collaboration of several institutions and was developed within the framework of Annex 33/SHC Task 58 Material and Component Development for Compact Thermal Energy Storage, a joint working group of the Energy Storage (ES) and Solar Heating and Cooling (SHC) Technology Collaboration Programmes of the International Energy Agency (IEA) [4]. It consisted of a survey with a detailed description of the experimental devices present in those institutions and used for investigating the long-term stability and performance of PCMs under application conditions [5]. In fact, an important prerequisite to select a reliable material for thermal energy storage

applications is to investigate its performance under real working conditions. In the case of solid–liquid PCMs, the long-term performance in terms of the melting and solidification processes should be ensured along the lifetime of the storage system, taking into account the conditions of the intended application. In this work, the different institutions presented up to 18 experimental set-ups and devices that allowed for performing thermal tests (cycling and constant temperature) not only for conventional PCMs, but also for the ones with stable supercooling, as well as phase change slurries (PCSs). Moreover, the paper introduced appropriate methods to investigate possible degradation mechanisms of PCMs. Considering the diversity of the devices and the wide range of experimental parameters, further work toward a standardization of PCM stability testing is strongly recommended. The main contribution of this paper is that it puts together many experimental facilities currently in use and the know-how of the corresponding institutions for assessing the long-term performance of PCMs, which certainly is a key issue for the commercial implementation of LHTS systems.

The last paper of this special issue presented a compact model of latent heat thermal storage (LHTS) for its integration in multi-energy systems [6]. In this way, the study developed a new modeling approach to quickly characterize the dynamic behavior of an LHTS unit. The thermal power released or absorbed by an LHTS module was expressed only as a function of the current and the initial state of charge. Moreover, the model also allowed for simulating even the partial charge and discharge processes. In general, the results were quite accurate when compared with a 2D finite volume model, with the advantage of the computational effort being much lower. Due to its simplicity, this model can be used to investigate optimal LHTS control strategies at the system level. In light of this, the authors implemented it in two relevant case studies: the reduction of the morning thermal power peak in district heating systems and the optimal energy supply schedule in multi-energy systems. However, this study describes the functioning of the LHTS unit at the system level only on the basis of numerical results. Hence, future work should also test the LHTS unit in a real case application to better quantify the model uncertainties.

The main contribution of this paper is the development of a simple model for LHTS modules that can be implemented in the simulation of multi-energy systems, although the model uncertainties still remain unquantified since it should be validated with data obtained from real applications.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Conflicts of Interest:** The author declares no conflict of interest.

#### **References**


## *Article* **Compact Model of Latent Heat Thermal Storage for Its Integration in Multi-Energy Systems**

#### **Alessandro Colangelo \*, Elisa Guelpa, Andrea Lanzini, Giulia Mancò and Vittorio Verda**

Department of Energy, Politecnico di Torino, 10129 Torino, Italy; elisa.guelpa@polito.it (E.G.); andrea.lanzini@polito.it (A.L.); giulia.manco@polito.it (G.M.); vittorio.verda@polito.it (V.V.)

**\*** Correspondence: alessandro.colangelo@polito.it

Received: 30 October 2020; Accepted: 14 December 2020; Published: 16 December 2020

#### **Featured Application: Dynamic simulation of Latent Heat Thermal Storage at system level; optimization of control strategies in multi-energy systems; investigation of Demand Side Management strategies.**

**Abstract:** Nowadays, flexibility through energy storage constitutes a key feature for the optimal management of energy systems. Concerning thermal energy, Latent Heat Thermal Storage (LHTS) units are characterized by a significantly higher energy density with respect to sensible storage systems. For this reason, they represent an interesting solution where limited space is available. Nevertheless, their market development is limited by engineering issues and, most importantly, by scarce knowledge about LHTS integration in existing energy systems. This study presents a new modeling approach to quickly characterize the dynamic behavior of an LHTS unit. The thermal power released or absorbed by a LHTS module is expressed only as a function of the current and the initial state of charge. The proposed model allows simulating even partial charge and discharge processes. Results are fairly accurate when compared to a 2D finite volume model, although the computational effort is considerably lower. Summarizing, the proposed model could be used to investigate optimal LHTS control strategies at the system level. In this paper, two relevant case studies are presented: (a) the reduction of the morning thermal power peak in District Heating systems; and (b) the optimal energy supply schedule in multi-energy systems.

**Keywords:** latent heat thermal storage; pcm; 0D dynamic model; multi-energy system; district heating; thermal network

#### **1. Introduction**

The thermal energy storage is a key element for increasing the operational flexibility of energy systems, especially when these are supplied by renewable energy sources. When the space available for storage technologies is limited, such as in urban contexts, Latent Heat Thermal Storage (LHTS) through the melting/solidification of Phase Change Materials (PCMs) can reduce the occupied volume up to five times with respect to sensible thermal storages to reach the same energy content. Therefore, an increasing amount of attention is focused on LHTS systems in order to advance their Technology Readiness Level (TRL) and market development. Currently, LHTS technologies have scarcely been introduced to market, and their TRL is lower than seven [1], which indicates that this technology has generally been demonstrated only in operational environments.

In this framework, numerous LHTS modeling approaches have been proposed in the existing literature. LHTSs are commonly analyzed through detailed numerical models, which, in few cases, are validated by experimental results. In fact, the intrinsic physical complexity of the phase change phenomenon can be easily simulated by numerical models [2]. Most of the models are focused on LHTSs design optimization with the aim of enhancing the heat transfer between the PCM and

the heat transfer fluid (HTF). As reported in [3], several enhancement techniques can be pursued. Among these, numerous works investigated the optimal sizing of HTF tube fins. For instance, Niyas et al. [4] performed a numerical study of a shell-and-tube LHTS prototype filled with PCM salts with a full 3D model. Using commercial software based on a finite element scheme, they determined the optimal number of HTF tubes and longitudinal fins per tube that minimizes the overall discharging time. They also monitored a few performance parameters of the identified configuration, such as the melting fraction, charging/discharging times, outlet HTF temperature, and energy stored. Then, results were validated with experiments [5] showing good agreement. Sciacovelli et al. [6] adopted a shape optimization strategy for Y-shaped fins with one or two bifurcations in order to enhance the performance of a shell-and-tube LHTS. A two-dimensional cross-section model, solving the underlying equations with a commercial finite volume code, was implemented. The fins optimization was performed requiring the complete discharge of the LHTS unit both in a short and long-time interval. A further application of CFD modeling to LHTS fins design optimization was presented by Pizzolato et al. [7]. They extended the investigation of the fins design through the combination of topology optimization techniques and ad hoc finite element two-dimensional and three-dimensional models. Here, only heat conduction is considered. Then, the same authors increased their model complexity by including also natural convection [8].

Esapour et al. [9] developed a detailed 3D in-house code to study the performance of an annular LHTS with multiple tubes. More specifically, they determined the effect of the arrangement and number of tubes on the liquid fraction and melting time. Then, in another work [10], the same model was also adopted to investigate the influence of variable operational parameters such as HTF inlet temperature and flow rate on the LHTS performance. Similarly, Saddegh et al. [11] compared the thermal behavior of a shell-and-tube LHTS using a 2D conduction–convection heat transfer model for two different LHTS orientations (vertical and horizontal).

In an attempt to simplify the simulation of LHTS systems, Neumann et al. [12] coupled a 1D model for the heat exchanger tubes and a reduced 3D model for the PCM and fins. This model can be used to improve the geometric and operational parameters of fin-and-tubes LHTSs with a limited computational effort, although the PCM parameters should be carefully calibrated. A different simplification approach was followed by Parry et al. [13]: a 3D computational model for a shell-and-tube LHTS was calibrated with experimental results. This was successively reduced to a one-dimensional radial model by employing an effective diffusivity technique. According to the authors, this model is suited to predict performance over long-time spans. An original simple model was proposed by Tay et al. [14]. This is based on the effectiveness-number of transfer units, and it was calibrated on a tube-in-tank salt-based PCM system with radial circular fins. The correlation proposed by the authors is semi-empirical, and it can be used for sizing and optimizing LHTS systems. A further interesting simplification methodology was suggested by Johnson et al. [15]. First, they analyzed the heat transfer properties of a tube assembly cross-section through a 2D finite volume model. Then, they produced a 1D radial model in Dymola for the ensemble of PCM and fins able to yield comparable results. Finally, they coupled this latter model with a 1D axial model for the HTF again in Dymola and studied the time evolution of the thermal heat flux for different configurations.

Nevertheless, although a few studies investigated the integration of LHTS in existing heating systems, little attention has been devoted to the analysis of LHTSs operating conditions and partial charge/discharge. For instance, Colella et al. [16] analyzed the behavior of an LHTS unit in a District Heating (DH) substation through a 2D finite volume model. Xu et al. [17] studied the performance of an encapsulated LHTS unit in a residential heating system in Sweden by means of a finite element model. Johnson et al. [15] designed an LHTS for the production of high-temperature steam in a cogeneration plant using the simplified model previously mentioned. However, in order to facilitate LHTS market development, their optimal integration in energy systems is essential. Dynamic models that are fast and accurate could positively support the study of LHTSs at a system level.

In this work, a compact 0D model for a modular shell-and-tube LHTS is proposed. The model is able to predict the thermal power released or absorbed by the assembly as a function of its state of charge. Furthermore, it can represent partial charge and discharge operations. Therefore, this model can be easily deployed for simulations at a system level in order to investigate optimal control strategies or demand side management (DSM). Indeed, the proposed model is here applied to two relevant case studies regarding the integration of LHTS in DH and multi-energy systems.

#### **2. Materials and Methods**

In order to study the optimal integration of different energy technologies at the building level, a reliable representation of each energy asset is essential. In this framework, an LHTS model for the simulation of the operating conditions provides quick results while retaining an acceptable physical description. Considering the physical complexity of the phase change in PCMs and the heat transfer in shell-and-tube configurations, a two-step model simplification was adopted. First, a 2D numerical model was developed in Ansys Fluent (2020 R1). Then, a 0D mathematical description was adopted, fitting the results obtained from the former numerical model with an exponential function. The LHTS system used as a reference is a vertical shell-and-tube type, where finned tubes are crossed by water as heat transfer fluid (HTF) and its shell encloses the PCM.

#### *2.1. D Detailed Model*

Figure 1 shows the computational domain of the initial detailed model. Due to the simultaneous presence of a solid and liquid phase in the PCM, both heat conduction and natural convection are relevant heat transfer phenomena in LHTS systems. Thus, a full 3D numerical domain is generally recommended. However, the characteristics of the configuration analyzed in this study allow a substantial domain reduction. In fact, as claimed by Vogel et al. [18], compact LHTS designs hinder the development of natural convection phenomena, which are the main aspects responsible for problem asymmetries. According to [18], the Rayleigh number (*RaW*) and LHTS aspect ratio (A) can be used to determine the relevance of the wall heat flux due to natural convection in the liquid PCM. This correlation asserts that the flow is dominated by natural convection only if *RaW <sup>A</sup>* ≥ 500. Regarding the LHTS configuration proposed in this analysis, the above-mentioned correlation yields *RaW <sup>A</sup>* = 374. Therefore, the heat flux due to natural convection is likely to be negligible both during charging and discharging phases if compared to conduction. Hence, the computational domain is reduced to the horizontal cross-section plane. Here, the temperature gradient is expected to be larger compared to the axial direction [19], and heat conduction is dominant. The numerical domain is also further reduced, taking advantage of the symmetry created by the fins. Then, the physical problem is modeled using the finite volume method implemented in the commercial code Ansys Fluent.

Regarding the materials, the fins and tubes are made of aluminum, while the PCM is assumed to be RT70HC manufactured by Rubitherm GmbH [20]. A summary of the PCM properties is reported in Table 1. The enthalpy–porosity technique is used for modeling the melting and solidification process [21]. Combining this method with the hypothesis of negligible natural convection, the only governing equation to be solved is the energy equation, as expressed in Equation (1).

$$\frac{\partial}{\partial t}(\rho H) - \nabla \cdot (k \nabla T) = 0 \tag{1}$$

where *H* is the sum of the sensible enthalpy and latent heat of the material.

Adiabatic boundary conditions are applied on the external face (Γ<sup>3</sup> in Figure 1), while symmetry is considered on Γ2. The heat transfer with the water in the pipes is modeled with an implicit Robin boundary condition on Γ1, as shown by Equations (2)–(4). *Tre f* represents the average temperature between the water inlet (*Tin*) and outlet (*Tout*) temperatures. However, *Tout* depends on the efficacy of the heat transfer between the flowing HTF and the rest of the LHTS (i.e., the ensemble of fins and PCM). Therefore, if *Tout* is expressed as a function of the heat flux at the wall ( . qwall), Fluent solves Equation (2)

iteratively. The heat transfer coefficient *hconv* is computed using the Dittus–Boelter correlation (5) and the Nusselt number (6). In (5), n = 0.4 for discharge or n = 0.3 for charge. This correlation was preferred over the implicit one proposed by Sieder and Tate due to its simplicity and acceptable uncertainty. Substantially, the simulated domain represents an average cross-section of the whole LHTS unit.

$$\mathbf{q}\_{\text{wall}}^{\dagger} = h\_{\text{conv}} \left( T - T\_{ref} \right) \tag{2}$$

$$T\_{ref} = \frac{\left(T\_{in} + T\_{out}\right)}{2} \tag{3}$$

$$T\_{\rm out} = T\_{\rm in} + \frac{L\_{\rm tubc} \pi d\_{\rm i}}{\rm Gr\_{p\_w}} \mathbf{q}\_{\rm wall} \tag{4}$$

$$Nu = 0.023 Re^{\frac{4}{3}} Pr^{\text{n}} \tag{5}$$

$$h\_{conv} = \frac{Nu \ast k\_w}{d\_i} \tag{6}$$

**Figure 1.** (**a**) Three-dimensional (3D) Latent Heat Thermal Storage (LHTS) shell-and-tube concept; (**b**) 2D computational domain.


**Table 1.** Phase Change Materials (PCM) physical properties.

The numerical model follows the same approach adopted by Sciacovelli et al. [6], who validated their 2D cross-section model for PCM solidification against experimental data. For further details, the interested reader is referred to [6]. The computational grid consists of a non-structured mesh made of 8278 cells with an average cell size of 3.5 <sup>×</sup> <sup>10</sup>−<sup>4</sup> m. The selected mesh proved to be sufficiently fine not to influence the results. A Second Order Upwind scheme is used for the spatial discretization of the energy equation together with the Least Squares Cell Based method for gradient calculation. The convergence is reached when residuals are lower than 10−<sup>8</sup> for the energy equation. On the other hand, the transient nature of the problem is approached with a Second Order Implicit Euler method, with a time-step of 1 s. The selected value proved to be sufficiently fine not to influence the results.

#### *2.2. D Compact Model*

Considering the LHTSs dynamics, the thermal power released or absorbed by LHTSs varies over time, as the thermal conductivity of the liquid and solid material is different. Indeed, as the solidification or liquefaction of the PCM proceeds, the power exchanged decreases, since the melting front is further from the fins. For this reason, it was assumed that the most relevant parameter influencing the thermal power input/output of an LHTS is represented by its state of charge (*SOC*). The *SOC* can be defined as the ratio between the actual energy stored inside the LHTS (*Estored*) and its reference value when fully charged (7).

$$\text{SOC} = \frac{E\_{\text{stored}}(t) - E\_{\text{min}}}{E\_{\text{max}} - E\_{\text{min}}} \tag{7}$$

where *Emax* is the energy contained in the PCM when the LHTS is fully charged. Similarly, *Emin* is the energy contained in the PCM when fully discharged.

A mathematical correlation between the current state of charge of the LHTS (*SOC*) and its thermal power input/output was identified using the results of the 2D model as shown by Equations (8)–(10). These were obtained observing the shape of the 2D model results and searching for their best fit through the sum of exponential functions. Each expression refers to a different operational phase (i.e., charge, discharge, or idle phase).

$$\dot{Q}\_{\text{dis}} = A e^{\text{B\*SOC}\_{\text{v}}} + C e^{\text{D\*SOC}\_{\text{v}}} + \mathcal{K}(1 - \text{SOC}\_{0}) e^{-\left(\frac{\text{SOC}\_{\text{v}} - E}{\overline{I}}\right)^{2}} \tag{8}$$

$$\dot{Q\_{\text{clr}\text{g}}} = A e^{B \star S \text{OC}\_{\text{u}}} + \mathbb{C} e^{D \star S \text{OC}\_{\text{u}}} + K \ast S \text{OC}\_{0} e^{-\left(\frac{\text{SOC}\_{\text{u}} - E}{F}\right)^{2}} \tag{9}$$

$$\dot{Q}\_{\text{idle}} = 0\tag{10}$$

*SOCn* is the normalized state of charge, and it is equal to -*SOC SOC*0 during discharge and - *SOC*−*SOC*<sup>0</sup> 1−*SOC*<sup>0</sup> during charge, *SOC*<sup>0</sup> is the initial state of charge, . *Qdis* is the thermal power released by each tube during discharge, . *Qchrg* is the thermal power requested by each tube during the charging phase, and . *Qidle* is the thermal power when the LHTS is not operated. Therefore, depending on the operational phase of the LHTS, the value of the thermal power released or requested by each LHTS tube can be expressed as follows (11): .

$$
\dot{Q}\_{LHTS} = \mathbf{x}\_1 \dot{Q}\_{dis} - \mathbf{x}\_2 \dot{Q}\_{dis} + \mathbf{x}\_3 \dot{Q}\_{idk} \tag{11}
$$

where *x*1, *x*2, *x*<sup>3</sup> [0; 1] and *x*<sup>1</sup> + *x*<sup>2</sup> + *x*<sup>3</sup> = 1.

Equations (8) and (9) were obtained fitting with exponential and Gaussian functions the LHTS thermal power and states of charge resulting from the simulations of the 2D model. The initial conditions of these simulations were, respectively, a fully charged LHTS and a fully discharged LHTS. Moreover, assuming that the thermal losses from the PCM to the HTF are negligible, the thermal power is always zero during the idle phase. Although these expressions were retrieved considering a specific LHTS design with fixed operating parameters (such as water mass flow rate and temperatures), the proposed methodology can be easily tailored to different LHTS concepts and operating conditions. Indeed, the evolution of the LHTS heat rate highlighted by several studies [12,15,22] is comparable with the one expressed by Equations (8) and (9).

#### **3. Case Study**

#### *3.1. LHTS System Description*

The LHTS system considered in this study is a vertical shell-and-tube type, with a height of 1.5 m. Each tube has an inner diameter of 19.05 mm (6/8) and is equipped with 16 longitudinal fins whose extension is 30 mm (thickness: 1 mm). Instead, the interaxial distance between each tube is assumed to be 91 mm. The PCM is characterized by an average phase change temperature of 70 ◦C, as reported in Table 1. The HTF inlet temperature is assumed to be 75 ◦C during the charging phase and 48 ◦C during the discharge. These values were selected analyzing the working temperatures of the heating system of the building described in Section 3.2. The HTF mass flow rate is 0.168 kg/s per tube. This value was obtained fixing an average HTF velocity of 0.6 m/s. Finally, the size of the LHTS system (in terms of stored energy) is not fixed a priori. As a matter of fact, the LHTS can be conceived as a modular system composed by the aggregation of several tubes (*ntubes*) and their associated PCM depending on the considered application. Each unit (tube + PCM) is able to store 2637.2 kJ, considering both the sensible and latent heat content of the PCM (at the selected temperatures).

#### *3.2. Model Application to Thermal Energy Networks*

Among the strategies to improve the operational management of DH networks, the adoption of distributed thermal storage systems at the building level is gaining interest as a measure to shave peak demands [23]. Therefore, the first case study considered in this article is the reduction of the peak thermal request from a substation of the DH network of Torino (Italy), which is characterized by a climate zone E. This substation belongs to a building known as the Energy Center, which is a research center owned by Politecnico di Torino. Thanks to its monitoring system, the building energy data are available with a frequency of 15 min. For this case study, a representative day of January 2020 is taken as a reference. As can be seen in Figure 2, the building heat demand has a common shape characterized by a morning peak (from 6:30 to 9 am) and a rather steady thermal request during the remainder of the day. During the peak hours, the amount of heat requested to the DH network is 3100 MJ, which represents the 27.1% of the overall heat demand within the same day.

**Figure 2.** Building heat demand (15-min average thermal power).

Overall, an LHTS system composed of 300 units (791 MJ) is needed in order to achieve an energy demand reduction of 25.5% during the peak. However, due to the high power ceased or absorbed by each LHTS unit at the beginning of its operation, the overall storage system is subdivided into six bundles made of 50 units (tube + PCM) as depicted in Figure 3. These bundles are sequentially activated with a delay of few minutes both during charge (in the night) and discharge (in the morning peak). Moreover, each bundle follows the subsequent duty cycle: it is charged up to *SOC* = 0.97 and discharged down to *SOC* = 0.02. Further partial charge and discharge cycles are not considered for this case study. This is done since the building heat demand is rather constant for the rest of the day.

**Figure 3.** LHTS charge (**a**) and discharge (**b**) schedule for District Heating (DH) case study.

#### *3.3. Model Application to Multi-Energy Systems*

This section shows the application of the 0D model for the achievement of the best operation in a multi-energy system context, as schematized in Figure 4. The multi-energy system is conceived to supply heat, electricity, and cold to a block of flats with 15 dwellings. The production system is characterized by a group of several technologies. In this context, the adoption of a thermal storage provides significant benefits to improve the system flexibility. In particular, an LHTS can be extremely advantageous in buildings, since it requires small installation volumes. In order to easily model the time-dependent behavior of the LHTS, the 0D model is adopted.

The most convenient technologies to be operated are selected by a Mixed Integer Non-Linear Programming optimization tool that includes the 0D LHTS model. The optimization algorithm allows achieving the best operation for the multi-energy system. Therefore, its output is the optimal time evolution of (a) the power from the production/conversion technologies; (b) the power purchase from the grid; and (c) the released/absorbed power in the storages.

The optimization model includes constraints such as (a) the overall produced and purchased energy vectors must be equal to the loads and (b) the energy production of each technology must not exceed its maximum capacity (the same for the conversion technologies and the storages). Moreover, additional constraints are set to model the LHTS. The maximum thermal power that can be absorbed/released by the LHTS is imposed at each time step in relation to the outcome of the LHTS 0D model.

Concerning the energy production/conversion, both traditional and innovative technologies are considered to be installed. Among the traditional technologies that can be adopted are a heat-only boiler (HOB) for heat production, an electric chiller (EHP) for cold production, and the grid connection for the electricity/heat supply. Among the innovative technologies are a micro-cogeneration system (for combined heat and power production) fed by natural gas, an absorption chiller, fuel cells, photovoltaics (PV), thermal solar, and wind turbine. Furthermore, a series of storages (hot, cold, and electric) provide a significant level of flexibility for the multi-energy system.

**Figure 4.** Schematics of the multi-energy system optimization with a 0D LHTS model.

#### **4. Results and Discussion**

In order to assess the validity of the 0D model, this section compares the results with the ones obtained from the 2D model. Afterwards, the impact of an LHTS system on the operational strategies of two relevant case studies is quantified through the application of the compact 0D model. In the former case, it shows how the LHTS allows achieving a significant reduction of the thermal power peak in a District Heating (DH) substation. In the latter case, it is used to include the latent thermal storage in an optimizer that finds the optimal performance of a multi-energy system.

#### *4.1. Comparison between 0D and 2D Models*

The coefficients resulting from the fitting procedure described in Section 2 are reported in Table 2. Equations (8) and (9) approximate very well the 2D power discharge and charge curves when the initial conditions are, respectively, *SOC*<sup>0</sup> = 1 and *SOC*<sup>0</sup> = 0 (Figures 5 and 6). The former fitting curve is characterized by a correlation coefficient R<sup>2</sup> = 0.995, while the latter has R<sup>2</sup> = 0.992. Overall, the minimum value is R<sup>2</sup> = 0.939, while the maximum standard deviation is 0.138 kW.


**Table 2.** Fitting coefficients in Equations (8) and (9).

**Figure 5.** Latent Heat Thermal Storage (LHTS) discharge starting from different initial states of charge (after a partial/complete charge).

**Figure 6.** LHTS charge starting from different initial states of charge (after a partial/complete discharge).

The correspondence between the LHTS thermal power and its state of charge is not univocal. As a matter of fact, the LHTS thermal power depends both on the current state of charge (*SOC*) and on the initial state of charge (*SOC*0). However, *SOC*<sup>0</sup> affects the results only in case a new charging or discharging phase starts. For instance, if there are two consecutive charging or discharging phases (only interrupted by an idle period), the value of *SOC*<sup>0</sup> should not be updated at the beginning of the second charge/discharge simulation. Instead, the previous value must be retained. This issue is due to the fact that when a partial charge is performed, only the PCM close to the fins becomes liquid, while the PCM far from the fins remains solid. If the subsequent operational phase is a discharge, a peak of thermal power occurs at the beginning of this latter process because the first layer of PCM encountered in the heat propagation is liquid (Figure 5).

A similar consideration is valid when a charge is performed, starting from a partial discharge of the system (Figure 6). The 0D model for the LHTS discharge thermal power (8) is slightly less accurate when *SOC*<sup>0</sup> is much smaller than 1 (i.e., when a discharge is preceded by a very short charging phase). However, it is rather unlikely to discharge a thermal storage starting from such a low content of energy. Similarly, the 0D model for the LHTS charge thermal power (9) is less accurate when the *SOC*<sup>0</sup> is much larger than 0 (i.e., when a charge is preceded by a very short discharge phase). As a matter of fact, the 2D curves in Figure 5 are both shifted and contracted with respect to one other, while the 2D curves in Figure 6 are mainly shifted. This dissimilarity might be due to the fact that the HTF inlet temperature during the charging phase (75 ◦C) is much closer to the PCM average phase change temperature (70 ◦C) compared to the HTF inlet temperature during discharge (48 ◦C). In fact, as indicated by [24], the HTF inlet temperature affects the LHTS thermal power.

#### *4.2. Distributed LHTS in DH Networks*

The 0D LHTS model is here applied to the reduction of the thermal power peak that is observed in DH networks at the beginning of their daily operation in some areas of the world (such as the Mediterranean area). In order to achieve this goal, the LHTS system follows the subsequent duty cycle, as detailed in Figure 7:


**Figure 7.** LHTS charge and discharge duty cycle to reduce DH morning peak.

Although all the energy absorbed by the LHTS during the charging phase is released in the discharging phase, the duration of these two processes is remarkably different. Since the HTF mass flow rate is constant, the thermal power is affected only by the temperature difference between the incoming HTF (48 ◦C or 75 ◦C) and the PCM average phase change temperature (70 ◦C) in both the operational phases. Consequently, the discharging phase is faster because the temperature difference is four times higher compared to the charging phase. Moreover, the spikes registered in the LHTS thermal power does not constitute a desirable feature from the point of view of the building heating system management. However, this highly unsteady behavior could be easily changed by a finer regulation of the HTF mass flow rate [16,17] (that is out of the goals of the present analysis). Furthermore, as detailed in Section 2.2, a perfect LHTS insulation is assumed. Consequently, heat losses are negligible when the storage is fully charged.

As shown in Figure 8, the time delay between the activation of each LHTS unit strongly affects the maximum power requested by the building to the DH network, whose reference value is 385 kW for this case study. Considering the LHTS system configuration analyzed, the peak thermal power of the DH supply is minimized when the activation of each LHTS unit is delayed by 16 min with respect to its preceding unit during the discharging phase. On the contrary, the minimization of the peak DH supply during the charge of the whole LHTS system would yield the obvious solution of charging each LHTS unit separately. However, this situation would not be compatible with the available time in the night. For this reason, an activation delay of 16 min was retained also for the charging phase.

**Figure 8.** Minimization of DH peak thermal power through LHTS discharge delay.

Finally, Figure 9 shows the DH supply curve when the LHTS system presented in this study is deployed. For consistency with the data of the building energy demand, the DH supply curve describes the average thermal power over 15-min intervals. The morning peak reduction is evident from the figure. The thermal energy supplied by the DH between 6.30 and 9 am is reduced by 25.5%, while the maximum value of the DH thermal power is decreased by 18%. This result is achieved by simply shifting the supplied thermal energy in the night hours. This shifting could also be extensively advantageous for DH users in case of different tariffs during the daytime and nighttime. As an example, in Turin, the heat during the night is sold at a price that is up to 50% lower than the daily price [25]. Furthermore, the overall volume occupied by the LHTS system proposed in this study amounts to 3.4 m3, while a sensible water storage tank would require 7 m3 for the same energy content (assuming also the same temperature difference 75–48 ◦C). Therefore, this feature constitutes a great advantage when the thermal energy storage is located in the small technical rooms available in the buildings.

**Figure 9.** Building heat demand (15-min average thermal power) with and without LHTS.

#### *4.3. LHTS in Multi-Energy Systems*

This section reports the results achieved in the analysis of the Multi-Energy System. The technologies activated for supplying the electricity load are the micro-cogeneration (CHP) system and the photovoltaics (PV). The extra demand is supplied purchasing the electricity from the grid; the adoption of the batteries makes the production more flexible. The electricity consumption is due to both the building load and the electricity absorbed by the electric heat pump and the electric storage. Concerning the cooling load, this is supplied through an electric heat pump (EHP). As far as heating is concerned, Figure 10 includes the power evolution of the technologies selected by the optimizer. On the left of the figure, the production evolutions are reported; on the right, consumptions are reported. The production consists in the sum of the power due to the production/conversion technologies operated, the energy purchased from the grid, and the storage discharging phase. The consumption consists in the sum of the building load, the energy sold to the grid, the storage charging phase, and the consumption of the other technologies such as the heat consumed by the absorption heat pump. The sum of the production and consumption must be equal. As a result, the heat load is mainly provided by the micro-cogeneration system, the heat-only boiler, and the thermal solar. The LHTS allows storing energy in the valley to supply it during the peak loads.

**Figure 10.** Production/consumption evolution for the energy vectors: (**a**) hot production; (**b**) hot consumption.

#### **5. Conclusions**

In this study, a simple 0D model for a modular shell-and-tube LHTS is presented. The model was obtained considering the main quantities affecting the charging/discharging evolution and through a fitting of the LHTS thermal power achieved with a more detailed CFD model. The compact model is obtained parameterizing the 2D model results as a function of the LHTS initial states of charge in order to account for partial charge and discharge. Results show that the 0D model is extremely accurate when compared with the 2D model, as the minimum value of the correlation coefficient (R2) is 0.939, and the maximum standard deviation is 0.138 kW.

Overall, thanks to its simplicity and low computational cost, the 0D model can be easily used for simulations at the system level. It allows investigating the performance of numerous control strategies in energy systems equipped with an LHTS, even in case the LHTS is only partially charged or discharged. In this work, two specific relevant applications of the compact model are considered. The former regards the thermal power peak reduction in DH substations, while the latter concerns the optimal operational strategy of a set of production/conversion/storage technologies (included LHTS) in a multi-energy system. However, this study describes the functioning of the LHTS at the system level only on the basis of numerical results. Therefore, future work should test the LHTS also in a real case application to better quantify the model uncertainties.

**Author Contributions:** Conceptualization, A.C., E.G. and V.V.; methodology, A.C.; software, A.C. and G.M.; validation, A.C.; formal analysis, A.C.; investigation, A.C. and G.M.; resources, A.C.; data curation, A.C.; writing—original draft preparation, A.C. and G.M.; writing—review and editing, E.G. and A.L.; visualization, A.C. and G.M.; supervision, V.V., E.G. and A.L.; project administration, V.V.; funding acquisition, V.V. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by EUROPEAN COMMISSION, grant number 815301. The APC was funded by EUROPEAN COMMISSION, grant number 815301.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
