Thermophysical Characterization of a Thermoregulating Interior Coating Containing a Bio-Sourced Phase Change Material

Abstract

This paper presents the work carried out as part of a study of a proactive interior coating based on both plaster and a phase change material (PCM), intended to improve the energy efficiency of buildings. This bio-based PCM is composed of a mixture of vegetable oils, methyl stearate, and methyl palmitate micro-encapsulated into polymer capsules. These components with distinct thermal properties constitute a mixture that displays supercooling and proves difficult to characterize using methods known in the literature. This article focuses on the thermophysical characterization (i.e., thermal conductivities, thermal capacities, latent heat, melting temperatures) and numerical modeling of a sample of this coating tested in the laboratory. This characterization is derived from experimental measurements carried out on a fluxmeter bench and by inverse methods. A new model of PCM composite characterization is presented and simulated using Python; the output shows a high degree of accuracy in describing the thermal behavior of the coating, regardless of the thermal stress applied, even making it possible to represent the phenomenon of supercooling or partial melting/solidification.

1. Introduction

The building sector is critical in the effort to move towards carbon neutrality in France. In 2019, it was responsible for 20% of greenhouse gas emissions and accounted for 46% of all energy consumption. This consumption level grew continuously between 1990 and 2001, the year when it reached its peak. In 2018, 82% of CO₂ emissions in the residential sector were generated by the heating of buildings [1]. In the current context of energy transition and the fight against global warming, the building sector is indeed a priority with great potential for improvement. Major renovations are necessary to reduce the energy and environmental impact of built structures and comply with the requirements for an energy transition. It has thus become necessary to design tools and materials to enhance not only energy performance but also the building environment and life cycle, while preserving the health and well-being of its occupants. This study is focused on the energy storage in a coating, by using phase change materials (PCMs) to reduce energy consumption in buildings. These PCMs are materials that can transition between solid and liquid phases, reversibly, and exchange energy with their environment. They melt when reaching the corresponding melting temperature, in absorbing heat and releasing it later once the temperature has dropped and then solidifying. These materials have a high latent energy per unit volume. If this amount of energy is high, the heat storage/release properties will be more attractive. PCM stores more energy than a conventional building material, such as brick, with an equal mass. Thus, latent heat and sensible heat together constitute the total amount of energy the PCM can store.

Several types of PCM exist, including hydrated salts, paraffin, fatty acids, and eutectics of organic and non-organic compounds. A. Abhat [2] established a classification (Figure 1) of a large number of PCMs potentially used for energy storage, grouped into four families: paraffin, organic, inorganic, and eutectic. This author also detailed the main selection criteria for these materials (

Table 1. PCM selection criteria [2].

Chemical	Thermodynamic	Economic
Chemically stable	Melting point	Inexpensive
Non-flammable Non-explosive Non-corrosive	High latent heat of melting High conductivity and density Limited supercooling Volume variation during phase transition	Easy to produce Available in large quantities Durability

). However, PCMs must satisfy several technical requirements, such as chemical stability and toxicity [3,4], to ensure safe use.

As detailed in Figure 1, the range of PCMs is quite wide, featuring various thermophysical properties and paving the way for many possible uses. Hence, energy storage technologies, thanks to the use of phase change materials (PCMs), are used in many fields. M.M. Farid et al. [5] described the primary uses for these types of materials. Nasa [6] introduced certain paraffin-based PCMs (dodecane, hexadecane) in energy storage modules that served to control the temperature of the rover batteries used on Mars, by means of damping extreme temperature variations. A. Hasan et al. [7] also dealt with the application of PCMs (paraffin, salt, milk fat) for thermal control, but especially for cooling electronic packaging to keep the temperature below a safety limit (60 °C). J.D. Renteria et al. [8] used a PCM combined with graphene to reduce overheating in electric batteries; in this case, adding graphene improved the thermal conductivity of the PCM in order to enhance the thermal exchanges. Other authors, e.g., F. Frusteri et al. [9], used an inorganic salt-based PCM and combined it with a conventional material, namely carbon fiber, to increase thermal conductivity and evaluate the thermal storage properties. It is also possible to produce solar energy storage with PCMs, which constitutes one of the main advantages of this type of material. In the literature, X. Xiao et al. [10] studied mixtures with nitrates/expanded graphite and PCMs to obtain a composite material for solar energy storage; they showed that certain additives improve the thermal conductivity of the PCM and thus upgrade its storage efficiency.

For the past few decades, several research approaches have been pursued on using PCMs in the building field, primarily to improve thermal comfort and/or reduce energy consumption. M.M. Farid et al. [5] and A.M. Khudhair et al. [11] analyzed the various types of PCM used in the civil engineering sector. Their work discussed the integration of PCMs into walls, floor heating systems, and ceilings. I.O. Salyer et al. [12] even used PCMs to fireproof building materials and delay the outbreak of fire.

This study combines PCM with a plaster-based interior coating in order to store energy. The proactive coating (ETMV 8) presented here incorporates a micro-encapsulated vegetable oil-based phase change material (named Inertek 23), which is mainly methyl palmitate as well as methyl stearate (naturally present in many plant species). Similar studies were conducted in Florida by M. Shapiro [13], who demonstrated that methyl palmitate and methyl stearate mixtures are perfectly suitable for use in gypsum boards given their thermal storage and high latent heat. Other authors have employed mixtures based on methyl palmitate and stearate in combination with other construction materials such as concrete. The research conducted by L.C. Liston [14] showed that using such a fatty acid mixture as a PCM can improve the thermal performance of concrete. For their part, R. Nikolić et al. [15] performed a very detailed analysis of these two fatty acids and concluded that combining the two yields the most promising PCM for energy storage in the construction industry.

The phase change material used in this study (Intertek 23) offers the advantage of being compatible with other conventional building materials. Combined with a coating intended for an interior overlay application, it could be useful in reducing temperature fluctuations and smoothing out heat spikes by absorbing excess during summer or, if it proves to be efficient, the off-season as well.

This type of material is also attractive during winter, thanks to the storage and return of free energy obtained from solar radiation. This feature could have the effect of delaying the start of the heating system while reducing the sensation of a cold wall, which constitutes a source of discomfort in the building. To evaluate the efficiency and usefulness of such a product, it is essential to characterize and test it in situ. Generally speaking, the main thermal characteristics of PCMs, namely latent heat, phase-change temperatures, and enthalpy, are obtained by means of differential scanning calorimetry (DSC), which is one of the most widely used techniques to study reactions related to the transformation of a material subjected to temperature constraints. The underlying principle is to measure the temperature difference between an analyzed sample and an inert material, which does not undergo a phase change and serves as a reference [16]. When the sample is subjected to temperature stress, it will absorb or release heat, which is then quantified by measuring the change in heat flux. This method has some limitations [17] due, among other things, to the fact that only a very small quantity (a few milligrams) of material can be tested. Indeed, small mass samples are not representative of the thermal behavior of a material on a large scale [18,19]. In addition, the heating rate affects DSC results and can induce errors in the determination of thermal properties [20]. Non-DSC methods rely on enthalpy evolution to study the thermodynamic behavior of PCMs. According to M. Thonon et al. [21], most of these other methods, which generally consider the PCM as a binary solution, have limitations and do not always correctly represent the PCM thermal behavior. Moreover, the fatty acid mixture studied here is not a binary solution. Indeed, the two components of Inertek 23, i.e., methyl palmitate and methyl stearate, display completely different thermal properties (Figure 1b). Let’s also note that the transitions between solid/liquid phases during the state changes of the two materials do not occur over the same temperature range. The supercooling phenomena that may occur in mixtures are difficult to analyze with the classical methods reported in the literature. Supercooling reflects the ability of the PCM to lower its temperature below the crystallization point while remaining in liquid form. This supercooling effect can be a drawback if the degree of supercooling is too high by virtue of preventing energy discharge of the PCM [22]. Under these conditions, the stored heat is released at a temperature below the PCM melting temperature.

This article presents the thermophysical characterization of a coating (ETMV 8) sample tested in the laboratory. The sample’s properties are determined through experimental measurements (conducted in LGCgE laboratory on a fluxmeter bench) combined with identification techniques using inverse methods, by relying on a new numerical model. This numerical model is necessary to reliably represent the thermal behavior of the coating during solid/liquid phase transitions. This numerical model is also intended to integrate a thermal simulation software at the building scale.

The thermophysical characterization is performed during the melting and solidification phases of the PCM. The particularities of this PCM are that it is encapsulated and issued from a mixture of two PCM. This produces a specific thermal behavior of the material which is little dealt with in the literature. The melting of the PCM mixture shows a progressive melting of the two PCM in a simultaneous way. On the other hand, on cooling the material exposes the phenomenon of supercooling and shows a dissociated solidification of the two mixed PCM. To establish the numerical model of the composite material (gypsum + microencapsulated PCM), temperature ramps of varying lengths were imposed on the material, corresponding to the thermal dynamics observed in buildings, which is the aimed application for this material. The very different heating and cooling behavior of the PCM associated with the supercooling phenomenon is an issue for the choice of the analytical model. Many researchers [23] consider the assumption of a hysteresis, which implies different enthalpy curves for melting and solidification of the material. Other researchers believe that this is not thermodynamically correct, that the enthalpy curve is unique in both heating and cooling [24], and this is what we found in previous studies [25] of another micro-encapsulated MCP (BASF Micronal 5001DX) mixed with cement mortar. In the present state of our research, we have not found any other solution than to consider this hysteresis. Having two enthalpy curves for the same material introduces a problem during incomplete melting or solidification since it is necessary to “switch” from one curve to the other when the heating or cooling process is reversed. This “switch” is not to be underestimated because it can be frequently encountered when the energy gains in the building are not sufficient to completely melt the PCM in the thickness of the coating or when the night temperatures in the building are not low enough to allow the release of all the energy stored during the day by the coating. The fact of “jumping” from one curve to another has consequences from a numerical point of view by introducing either errors in the energy balance or temperature steps in the PCM. Methods (hysteresis models “curve track”, “curve switch” and “curve scale”) to deal with this problem have been provided by different authors [26,27,28]. We will specify which method we have chosen. One point that should be highlighted for this paper is the introduction of Gumbel’s law to deal with phase transitions. Long temperature ramps in cooling show supercooling of PCM and two distinct peaks of energy release at different temperatures. Gumbel’s law has the advantage of being an asymmetric law which allows by its form to model the progressive melting process and the more abrupt release of energy at the time of solidification of the PCM. This law has been associated with a binary solution formulation in order to model the thermal behavior of the two mixed PCMs.

The identification of thermophysical properties and the development of a reliable model are essential to model the thermal behavior of the coating at the building scale.

2. Thermal Characterization

2.1. ETMV 8 Sample

Let’s note that the ETMV 8 is a composite plaster-based material that incorporates a bio-based phase change material, called Inertek 23. The feasibility of such an environmentally friendly mixture has been demonstrated by many authors [29,30,31]. This PCM is a mixture of two vegetable oils encapsulated in polymer microcapsules. The sealed microcapsules, with a diameter ranging from 5 to 25 µm, are produced by polymerization. The encapsulation process allows the mixing of PCM with other building materials [32,33] and moreover provides optimal protection for PCM, thus ensuring a thermal storage capacity of 180 kJ/kg (Figure 2), as measured by the manufacturer [34].

ISO 9001-certified, Inertek 23 has a long lifespan and maintains stable thermal properties over time even after several hundred heating-cooling cycles [34]. According to its manufacturer, this PCM is perfectly resistant to temperature constraints ranging from 6° to 33 °C without affecting its performance. The melting range lies between 23° and 27 °C (Figure 2) with a peak at 25 °C, which turns out to be an attractive property for the specific objectives of saving energy and improving building comfort. On the other hand, solidification starts at 23 °C and ends below 19 °C. The ETMV 8 composition is given in

Table 2. Composition of the ETMV 8 coating.

Material	Composition [%]	Density [kg·m−3]
Coating Inertek 23 ETMV 8	30 70 100	650 350 610

.

The ETMV 8 coating (plaster + PCM (Table 2)) was poured to a thickness of 1.05 cm between two gypsum boards of dimensions 25 × 25 × 1.25 cm³ each. This coating was supplied by the manufacturer ready for use. The complete sample obtained (Figure 3), weighing 1.68 kg (dry weight, at 60 days) and sized 25 × 25 × 3.55 cm³, was dried at room temperature (approx. 20 °C, relative humidity around 50%).

The gypsum board has a density of 814 kg·m⁻³, while the ETMV 8 density equals 610 kg·m⁻³. It was decided to place the plaster between 2 gypsum boards in order to obtain homogeneous plaster compactness and constant thickness over the entire sample. The thermophysical properties of the gypsum board were determined based on heat flux and temperature measurements performed on a fluxmeter bench coupled with an inverse method, as in previous works [25].

2.2. Experimental Measurement Method

The LGCgE fluxmeter measurement bench is schematically presented in Figure 4; its various uses are detailed in [21,25,35,36,37].

The underlying principle consists of placing the sample between two exchange plates in order to impose a homogeneous surface temperature on both sides [25]. To limit lateral heat loss and ensure a 1D transfer through the coating and gypsum boards, a guard ring made of insulating material was positioned on the lateral faces of the sample. Each exchanger plate was connected to a thermostatic bath to control the temperature stresses applied to both sides of the material. The heat flux and temperature on each side of the sample were measured throughout the test using tangential gradient fluxmeters [38] integrating a thermocouple at the center.

For our case study, the temperature loadings varied between 10° and 30 °C in order to force the PCM to undergo melting and solidification cycles. From measurements carried out on the sample surface and through the use of the inverse method, the thermophysical properties of the coating, such as thermal conductivity, latent heat, specific heat capacity, and melting temperature, can be identified [37]. The temperatures imposed and measured on the sample faces constitute the input data (boundary conditions) for the numerical model, while the experimental heat fluxes serve for optimization by minimizing the difference between the numerically calculated and experimental heat fluxes. The objective is to identify, by means of inverse methods [39,40], the parameters of the numerical model that minimize an objective function. This function represents the quadratic difference between the numerically simulated heat fluxes and the experimentally measured fluxes, to the left and right of the sample; it is calculated as the difference between the experimental and theoretical heat fluxes on both sides of the sample according to the following formula:

$$f_{obj}={{\displaystyle \sum }}^{ }\left({\mathsf{\varphi }}_{lef{t}_{ exp}}-{\mathsf{\varphi }}_{lef{t}_{ cal}}\right)2+\left({\mathsf{\varphi }}_{right{ }_{exp}} – {\mathsf{\varphi }}_{righ{t}_{ cal}}\right)2$$

(1)

As such, during the first step of the characterization process, this method determines the thermophysical properties of the gypsum boards alone. At the present time, the thermal properties of gypsum boards are considered to be known (

Table 3. Properties of gypsum boards.

${\lambda }_S$ = 0.30 W·m−1·°C−1

$c_{s }$= 1071 J·kg−1·°C−1

ρ = 815 kg·m−3

).

In a second step, the system (sample) under consideration is the sandwich composed of coating surrounded by gypsum boards.

3. Numerical Modeling on Python

3.1. Binary Solution Model and Its Limitations

To analyze the coating presented in this study, it was first necessary to identify an existing model from the literature capable of describing the thermal behavior of the Inertek 23 PCM-based composite on two bio-based PCMs. Figure 5 below shows the evolution of the heat flux measured experimentally on the sandwich (gypsum board/ETMV 8 coating/gypsum board). The progressive melting, as well as the sudden flux peak during cooling, seems to be similar to the thermal behavior of a known binary solution studied previously [35] (composite material containing encapsulated-PCM).

At the beginning of melting ( ), a nearly identical curve of thermal flux measured on both the ETMV 8 sample and the binary solution can be observed. Concerning the solidification phase, two common points between the two bodies essentially stand out, namely the levels at the beginning ( ) and end of solidification ( ). It should also be noted that during the solidification phase, except for the second peak of thermal flux identified for the ETMV 8 coating, the curves are on the whole similar. All these elements serve to reinforce the choice of a binary solution model to analyze the coating behavior.

One of the most accurate existing models to describe the thermal behavior of PCMs with a binary solution behavior is detailed in [35,36,37,41]. In particular, this model has been used to accurately characterize a composite material consisting of cement mortar and paraffin-based PCM. P. Tittelein et al. [35] compared this model with the other main existing models, in deducing that it is the most precise to represent the thermal behavior of binary PCMs.

In general, the thermal behavior of PCMs can be described by expressing the evolution of enthalpy (h) as a function of temperature [21]. According to the first principle of thermodynamics, the heat equation associated with a phase change problem can be described in terms of the derivative of enthalpy (Equation (2)). For a binary solution, this derivative can be expressed by Equation (3), and the liquid fraction f_{liq_bin} (T) used to model the change of state is given by Equation (4) [42].

$$\rho .\frac{\partial h}{\partial T}.\frac{\partial T}{\partial t}=\lambda .\frac{\partial 2T}{\partial x2}$$

(2)

$$\frac{\partial h}{\partial T}=\frac{d{f}_{liq\_bin}\left(T\right)}{dT}.L_A+C_S.\left(1-f_{liq\_bin}\left(T\right)\right)+C_L.f_{liq\_bin}\left(T\right)$$

(3)

$$where f_{liq\_bin}·(T)=\{\begin{array}{c}\frac{\left(T_A-T_M\right)}{(T_A-T)} T\le T_M\\ 1 T>T_M \end{array}$$

(4)

Other existing models describe the thermal behavior of binary PCMs, although numerous authors [21,43,44,45,46] have shown that many of them do not satisfactorily model the behavior of a PCM composed of several materials with distinct thermal behavior. Some models fail to simulate a continuous evolution of the heat capacity over the entire temperature range, while others consider the heat capacity to be identical in the solid and liquid states or simply do not allow for asymmetric PCM behavior. However, the main problem with the model presented above (Equation (2) through to Equation (4)), which considers the PCM as a binary solution, is its inability to describe the phenomena caused by the presence of two oils in the PCM, resulting in the appearance of several heat flux peaks during solidification.

The use of this binary solution model, even if a different equation is considered for solidification and melting (Figure 6), is not satisfactory to accurately reproduce the thermal behavior of the ETMV 8 coating. Figure 7 shows the comparison between the experimental heat flux and the heat flux obtained numerically by the given model, for a thermal load varying between 10 °C and 30 °C over 10 h. Notable differences are visible between the two heat fluxes, such as the asymmetry between the melting/solidification phases and a 2nd flux peak during solidification. This observation is quite consistent since the PCM applied in the coating studied here is not a binary solution due to its composition (i.e., a mixture of two fatty acids). The methyl palmitate and methyl stearate included in the PCM display distinct thermal properties, e.g., melting temperatures, latent heat, and capacities.

3.2. New Hybrid Model

Since the binary model is not satisfactory to model this material, which presents a very distinct behavior in melting and solidification. The hybrid model considers that the equations characterizing the thermal behavior of PCM differ in melting and in solidification. In melting, for instance in solidification, the heat equation of a binary solution (Equation (3)) without eutectic [46] is used while being coupled to Gumbel’s law [47] in order to reproduce the effects related to dissymmetry and supercooling, as generated by the presence of several materials in the PCM. Natural convection in PCM microcapsules is neglected due to their small diameter, varying from 5 to 25 µm. The boundary conditions are of the DIRICHLET type, the imposed temperatures are T_LP on the left and T_RP on the right. The numerical model depends mainly on C_S, C_L, λ_S, λ_L, L_A, T_A, T_M, R_c-right, and R_c-left. R_c-right and R_c-left are the thermal contact resistances between the heat fluxmeter and the sample. Other easily measurable parameters, such as density and thickness, are also needed to complete the input dataset. The operation of the model implemented (simplified to a cell), i.e., homogeneous, tri-layer unidirectional conduction (2 gypsum boards, 1 ETMV 8), is presented in Figure 8, where the depiction has been simplified yet each layer features several meshes, each containing two thermal resistances and one thermal capacity.

Thus, the equations derived from this model are as follows:

$$\begin{array}{cc}m=\rho .\frac{e}{n}.S& n: \mathrm{number} \mathrm{of} \mathrm{meshes}, \mathrm{and} S: \mathrm{surface} \mathrm{area} \mathrm{(}{m}^2\mathrm{)}\end{array}$$

(5)

$$\begin{array}{cc}Q=m.C.\frac{dT}{dt}& dT=T_i-T_{i-1}, \mathrm{with} dt \mathrm{the} \mathrm{interval} \mathrm{of} \mathrm{the} \mathrm{calculation} \mathrm{time} \mathrm{step} \mathrm{(}\mathrm{s}\mathrm{)}\end{array}$$

(6)

$$\begin{array}{cc}{\mathsf{\varphi }}_{i, i-1 }=\frac{1}{R}.(T_i-T_{i-1})& \mathrm{where} \mathrm{R}=\frac{\left(e/n\right)}{\lambda }\end{array}$$

(7)

Contact resistances are integrated into the heat flux calculation at the material boundaries (left and right):

$${\mathsf{\varphi }}_{left}=\frac{1}{R+R_{c-left}}.\left(T_{LP}-T_i\right)$$

(8)

$${\mathsf{\varphi }}_{right}=\frac{1}{R+R_{c-right}}.\left(T_i-T_{RP}\right)$$

(9)

The derivative of enthalpy is calculated by integrating the Gumbel distribution (F_G) and is represented in Figure 9. This distribution is a continuous probability function generally used to fit extreme data. It can be found in several fields, e.g., finance, hydrology, geology, and meteorology [48]. In this paper, we use the Gumbel probability density function (f_G), defined as the derivative of the Gumbel function:

$$F_G\left(x,a,b\right)=\exp \left(-\exp \left(-\frac{x-a}{b}\right)\right)$$

(10)

Its derivative is given by:

$$f_G\left(x,a,b\right)=\frac{1}{b}·\exp \left(-\exp \left(-\frac{x-a}{b}\right)\right)·\exp \left(-\frac{x-a}{b}\right)$$

(11)

with x ϵ ℝ, a ∈ ℝ is a position parameter, and b > 0 is a spread parameter.

The advantage of this Gumbel derivative lies in its ability to easily obtain an asymmetric and adjustable form thanks to its position and spread parameters. Thus, the derivative of enthalpy is calculated according to the following equations.

During solidification (s):

$$\frac{d{h}_S}{dT}=\{\begin{array}{cc}{c}_S·\left(1-\frac{T_{A_s}-T_{M_s}}{T_{A_s}-T}\right)+c_L.\left(\frac{T_{A_s}-T_{M_s}}{T_{A_s}-T}\right)+\left(\frac{T_{A_s}-T_{M_s}}{T_{A_s}-T}\right).\left(\frac{{\mathsf{\psi }}_s.L_A}{T_{A_s}-T}\right) +f_{G_s}\hfill & T\le T_{M_s}\hfill \\ c_L\hfill & \hfill T\ge T_{M_s}\end{array}$$

(12)

$$f_{G_s}\left(T\right)=\left(\frac{1}{b_s}·\exp \left(-\exp \left(-\frac{T-a_s}{b_s}\right)\right)·\exp \left(-\frac{T-a_s}{b_s}\right)\right).(1-{\mathsf{\psi }}_s)·L_A$$

(13)

During melting (m):

$$\frac{d{h}_S}{dT}=\{\begin{array}{cc}{c}_S·\left(1-\frac{T_{A_m}-T_{M_m}}{T_{A_m}-T}\right)+ c_L.\left(\frac{T_{A_m}-T_{M_m}}{T_{A_m}-T}\right)+ \left(\frac{T_{A_m}-T_{M_m}}{T_{A_m}-T}\right).\left(\frac{{\psi }_m.L_A}{T_{A_m}-T}\right) + f_{G_m}\hfill & T\le T_{M_m}\hfill \\ c_L\hfill & \hfill T\ge T_{M_m}\end{array}$$

(14)

$$f_{G_m}\left(T\right)=\left(\frac{1}{b_m}·\exp \left(-\exp \left(\frac{T-a_m}{b_m}\right)\right)·\exp \left(\frac{T-a_m}{b_m}\right)\right)·(1-{\mathsf{\psi }}_m)L_A$$

(15)

${\mathsf{\psi }}_s et {\mathsf{\psi }}_m$ denote the latent heat distribution coefficients between the binary and Gumbel functions during the solidification and melting phase of PCM, respectively.

Figure 9 and Figure 10 present the evolution of the full enthalpy derivative, the enthalpy derivative without the addition of the Gumbel function, as well as with the Gumbel function. Let’s observe that the Gumbel function plays a more significant role during the melting phase in reproducing the effects of progressive melting of the material. During the solidification phase, the Gumbel function is mainly used to reproduce the effects related to supercooling, which occurs at approx. 17 °C.

This undercooling stems from the fact that the PCM used in this study is a mixture of two fatty acid-based PCMs [49], each with distinct thermal behavior. A difference of about 4% exists between the liquid fraction variation at the beginning of heating and at the end of cooling (Figure 11), as observed by K. Merlin [49] on several PCMs as well. The author assumed that this phenomenon was due to either incomplete melting of the material or measurement errors. An error analysis related to experimental uncertainties is carried out and detailed in Section 5.2.

Moreover, the choice of using an inverse Gumbel function in the melting phase is justified by the fact that this function reflects a slow and progressive increase, followed by an abrupt decrease. It presented the best results in comparison with other mathematical laws, e.g., Gauss’s law. The opposite behavior was observed during the solidification phase, thus explaining the applicability of the Gumbel function, which shows a rapid rise with a steep slope followed by a gradual decrease.

To optimize the parameters of this model and reduce the objective function to characterize the ETMV 8 sample, a differential evolution algorithm [50] available in the optimization module of Python, called Scipy, has been used. Due to its stochastic nature (i.e., the gradient method is not employed), this algorithm can search within a large candidate space area for the minimum of the objective function; its operation is detailed in [51].

To reduce the number of parameters to be optimized, the thermal conductivities and capacities, as well as the contact resistances, are determined beforehand during experiments where the PCM is in either the solid or liquid state. Yet despite this, the optimization step still requires the identification of 11 parameters (L_A, ${\mathrm{T}}_{{\mathrm{A} }_{\mathrm{m}}}$, ${\mathrm{T}}_{{\mathrm{M}}_{\mathrm{m}}}$, ${\mathrm{T}}_{{\mathrm{A}}_{\mathrm{s}}}$, ${\mathrm{T}}_{{\mathrm{M}}_s}$, ${\mathsf{\psi }}_s, { \mathsf{\psi }}_m$, ${\mathrm{a}}_{\mathrm{m}}$, b_m, a_s and ${\mathrm{b}}_{\mathrm{s}}$).

4. Results and Discussion

4.1. Optimization of the Numerical Model Parameters

The experimental set-up described above (Figure 4) has been used to characterize the ETMV 8 coating. The sample was stressed by temperature ramps ranging from 3 h to 10 h. The fluxes measured for each case are shown in Figure 12; it can be observed that the measured flux has a shape that varies greatly with the duration of the adopted temperature ramp. Similar behavior has already been observed and detailed by M. Thonon et al. [21]. The melting of the PCM is progressive during the heating ramps, the material stores sensible and latent heat. The maximum flux density is reached at the end of the temperature ramp and then returns to a zero value until the material returns to an isothermal thermal state. Let us notice in particular that the use of a temperature ramp of short duration during cooling, e.g., 3 h, does not expose the second flux peak, which is caused by the presence of several materials in the PCM. This second peak is only visible from the 5-h temperature ramp because for the faster ramp (3-h), the temperatures vary too quickly, and the two energy release peaks occur at very close times. However, when installed in a building, the material rises in temperature over periods that can be on the order of one day. It is therefore important to take these longer periods into account.

It is therefore essential to use a large temperature ramp to more accurately determine the thermal characteristics of the coating. In this study, the thermal characterization of the sample was performed with the experimental measurements of the test performed with the T-10H thermal loading ramp.

From these heat flux and temperature measurements, inverse methods allow for determining: thermal conductivity, heat capacity in the solid and liquid states, thermal contact resistances between the plates and the sample, and latent heat. This latent heat (126 kJ/kg) is most consistent since it closely corresponds to 70% of the latent heat of PCM (180 kJ/kg), as obtained by the DSC method [52,53]. The melting temperatures TA (melting temperature of the pure solvent) and TM (temperature of the liquidus) are also evaluated. The resulting parameters are as follows (

Table 4. Thermal characteristics of ETMV 8.

$T_{A{ }_m}$= 29.02	$T_{A_{ s}}$ = 23.94	${\lambda }_S$ = 0.135	$c_{s }$= 1213	LA = 126,300
$T_{M_m}{}_{ }$= 26.89	$T_{M_s}=22.71$	${\lambda }_L$ = 0.128	CL = 1124

and

Table 5. Parameters of the Gumbel function.

$\mathrm{Parameters} \mathrm{of} \mathrm{the} \mathrm{function} f_{Gumbel \mathrm{melting}}$
$(1-{\mathsf{\psi }}_{\mathrm{m}}).{\mathrm{L}}_{\mathrm{A}}$ = 81,942	am = 24.18	bm = 0.98
$\mathrm{Parameters} \mathrm{of} \mathrm{the} \mathrm{function} f_{Gumbel solidification}$
$(1-{\mathsf{\psi }}_{\mathrm{s}}).{\mathrm{L}}_{\mathrm{A}}$ = 12,300	as = 17.74	bs = 0.4

):

These parameters serve to numerically simulate the thermal behavior of the coating. Figure 13 shows a comparison between the thermal flux measured experimentally and that calculated with the numerical model and optimized parameters. This figure reveals that the model presented herein, i.e., the derivative of enthalpy coupled with a Gumbel function, yields strong results (as opposed to the model of a simple binary solution) and moreover accurately reproduces the behavior of the ETMV 8 coating, whether during melting or solidification; in particular, it allows reproducing the effect of supercooling during the solidification phase.

It has been demonstrated that the dynamics of the melting and solidification processes are well represented by the numerical model developed above. It remains to be verified that this model is efficient for various loadings and in the case of partial melting or solidification since these cases are difficult to model and can lead to errors in the calculation of energy balances.

4.2. Test of the Numerical Model with Various Thermal Sollicitations

To verify the effectiveness of the model just built, the same sample was loaded in several different ways.

Figure 14 and Figure 15 show the results obtained for the case of temperature stress in the form of a 3-h and 6-h ramp. Figure 16 displays a comparison between the experimental fluxes and those obtained numerically, for the case of an asymmetric loading. Various temperatures were imposed on the left and right sides of the sample at different heating/cooling rates. A good level of agreement is always found between the experimental results and numerical model results, regardless of the heating/cooling rate.

It was also necessary to test the model under constraints leading to partial state changes of the PCM. The problem arising with partial phase changes is mainly caused by the transition between the melting and solidification processes. In the literature, two methods are basically employed to model these types of PCM state changes. Several authors, including R. Chandrasekharan et al. [54], have suggested an absence of hysteresis until the PCM is completely melted or solidified. This assumption implies no difference between the state changes of the material as it alternates between melting and solidification. The approach proposed by R. Chandrasekharan et al. [54] thus maintains the same liquid fraction curve as a function of the PCM state (Figure 17a). When inducing changes in the partial melting states, only the liquid fraction curve is therefore used in the melting phase. On the other hand, when carrying out a partial solidification, only the liquid fraction curve in the solidification phase gets used, along the opposite direction of the curve when the temperature variation is reversed. Another method that takes PCM hysteresis into account also exists, leading J. Bony et al. [55] and E. Moreles et al. [56] to propose alternating between the melting and solidification liquid fraction curves, depending on the temperature stresses imposed on the material (Figure 17b). R. Chandrasekharan et al. [54] compared the accuracy of these two methods and found relative similarity.

Other authors, namely [57,58], have come up with different opinions; they have shown, through studies on walls containing PCM, that the hysteresis of the PCM plays an important role and has a significant impact on both the energy balance and thermal performance of the material. The studies in [26,58] compared the two methods and issued a divided opinion on the subject. They concluded that the efficiency of both methods depended on the PCM phase.

Notably, B. Delcroix et al. [26] observed that although the method taking PCM hysteresis into account is more accurate, it nonetheless does not allow for a satisfactory representation of the partial melting. The hybrid model proposed in this paper, incorporating supercooling, was developed with the first method, which does not take hysteresis into account when analyzing melting and partial solidification. Further studies will draw a comparison of these two methods. In our model, the liquid fraction (Figure 17) during solidification and melting, both complete and partial, is calculated based on Equations (12) and (14), yielding the enthalpy variations Δh_m and Δh_s (Figure 18) by means of integrating the derivative of enthalpy during melting dh_m/dT and solidification dh_s/dT.

We have:

$$\Delta h_m\left(T\right)={{\displaystyle \int }}_{T_i}^T \frac{d{h}_m \left(T\right)}{dT} dT-C_L\left(T-T_i\right)$$

(16)

where T_i denotes the initial temperature at which the PCM is assumed to be completely liquid. It is also assumed that the liquid fraction of the hybrid model, during either the melting or partial melting phase, is as follows:

$$f_{liq\_m}(\mathit{T})=\frac{\Delta h_m\left(T\right)}{\Delta h_m\left(T_{M{ }_m}\right)}$$

(17)

Similarly, as for the solidification phase, we can determine:

$$\Delta h_s\left(T\right)={{\displaystyle \int }}_{T_i}^T \frac{d{h}_s \left(T\right)}{dT} \mathit{dT}-C_L\left(T-T_i\right)$$

(18)

The liquid fraction of the hybrid model, during solidification or partial solidification, is then:

$$f_{liq\_s}(\mathit{T})=\frac{\Delta h_s\left(T\right)}{\Delta h_s\left(T_{M{ }_s}\right)}$$

(19)

In Equations (16) and (18), it was assumed that only C_L is used given its proximity to C_s.

The partial melting tests (Figure 19) were performed over a temperature range extending from T₁ (10 °C) to T₂ (24 °C) through an intermediate temperature T₃ (15 °C). The sample was also subjected to temperature variations ranging from T₄ (30 °C) to T₅ (18 °C) through T₆ (23 °C), causing partial solidification of the PCM (Figure 20). The evolution of the liquid fraction during these partial state changes is obtained from Equations (17) and (19). Figure 19 and Figure 20 present the results obtained in the case of temperature stress leading to partial melting and solidification, respectively. It can be noted that the melting temperatures obtained for the coating containing the Inertek 23 PCM (Figure 18) are very close to those obtained by DSC on PCM (Figure 2).

A good level of agreement could be observed between the experimental flux and the numerical flux obtained with the hybrid model integrating a Gumbel function. These experiments confirm that the model developed herein can correctly simulate the behavior of the ETMV 8 coating even in the presence of various thermal loadings, which generate partial changes to the state of the PCM contained in the mixture.

5. Energy Balance and Error Analysis

5.1. Energy Balance

To complete the model verification step, it was necessary to determine the amount of energy involved during the PCM melting and solidification phases, for various loadings, in order to ensure energy conversation during these storage and retrieval cycles. This determination entailed integrating the heat flux (both experimental and numerical) with respect to time.

These results indicate a very slight difference (less than 1%) between the energy calculated from the experimental flux measurements during the melting phases and that obtained during the solidification phases, for each ramp (

Table 6. Energy balance, based on the experimental heat flux.

Experimental
Loading	Melting [J]	Solidification [J]	Difference [%]
3-h ramp	136,233	135,727	0.37
6-h ramp	134,579	135,850	0.94
10-h ramp	135,770	135,937	0.12
Max difference [%]	1.2	0.2

). This difference is probably due either to measurement errors during the tests or to the supercooling, which stochastically influences the energy storage/release, especially during the solidification phase. Several studies have reported that such supercooling appears unpredictably and depends on various factors. According to J. C. Delabrouille et al. [59], supercooling essentially depends on the purity of the PCM. Indeed, impurities in the material will act as primers to accelerate the crystallization phenomenon. In addition, other authors, including J.-P. Dumas [60] and J.-P. Bedecarrats [61], have shown that all heating/cooling cycles imposed on the material during its lifetime (nature, repetitions) constitute the thermal history of the PCM, which in turn has a direct influence on the degree of supercooling. All these key factors act on the PCM supercooling and thus cause differences between the energy balances in melting and solidification for each thermal stress cycle.

Let’s also note a difference of roughly 2.4% when performing the energy balance from the heat fluxes calculated numerically by the hybrid binary-Gumbel model (

Table 7. Energy balance, from the numerical heat flux.

Numerical
Loading	Melting [J]	Solidification [J]	Difference [%]
3-h ramp	137,232	133,968	2.38
6-h ramp	137,235	133,978	2.37
10-h ramp	137,225	133,922	2.41
Difference [%]	1.0 10−2	4.0 10−2

). This difference is likely due to inaccuracies related to the optimization of the numerical model parameters. Indeed, errors in experimental measurements of heat flux or temperature lead to inaccuracies when characterizing the thermophysical parameters of the PCM. Therefore, an error analysis is performed in order to evaluate the impact of these measurement errors.

5.2. Evaluation of Errors Related to Uncertainties

In the context of this study, uncertainties may arise in the measurement of heat flux and temperature of the sample during testing, as previously mentioned. These error sources can lead to misinterpretation of the results and thus reduce the accuracy of the characterization of thermophysical properties of the material. According to the law of propagation of uncertainties, the compound uncertainty (U_c) includes several components [62]. In the present case, it is essential to take into account errors related to the T-type thermocouples (U_T), the tangential gradient fluxmeters FGT (U_φ), and the data acquisition unit (U_mod). For heat flux measurements, this compound uncertainty is given by the following equation:

$$U_{c,\mathsf{\varphi }}={\left[U_{\mathsf{\varphi }}{}^2+U_{b,\mathsf{\varphi }}{}^2\right]}^{\frac{1}{2}}$$

(20)

Similarly, the compound uncertainty related to temperature measurements is expressed as:

$$U_{\mathrm{c}, \mathrm{T}}={\left[U_T{}^2+U_{b,T}{}^2\right]}^{\frac{1}{2}}$$

(21)

The error on fluxmeter sensitivity (114 µV·W⁻¹·m⁻²) has been estimated at 2% [63], with the uncertainty values being calculated for a heat flux of approx. 100 W·m⁻² (

Table 8. Uncertainties related to the heat flux and temperature measurements.

Heat flux	Component	Designation	Uncertainty [W/m2]	Comment
	Captec fluxmeter	Uɸ	2.04	For 100 W·m−2
	Data acquisition system	Umod,ɸ	0.0009	[64]
	Compound uncertainty	Uc,ɸ	2.04	(Equation (20))
Temperature	Component	Designation	Uncertainty [°C]	Comment
	T-type thermocouple	UT	0.5	[65]
	Data acquisition system	Umod, T	0.001	[64]
	Compound uncertainty	Uc,T	0.5	(Equation (21))

). The standard uncertainty of the acquisition unit is provided by the manufacturer. Thus, the compound uncertainty for the heat flux measurement is evaluated here at 2.04 W·m⁻². The temperature measurement errors are mainly correlated with both the accuracy of thermocouples provided by the manufacturer and the accuracy of the central acquisition unit. Temperature is therefore measured with a precision of plus or minus 0.5 °C. The error range related to uncertainties of heat flux measurements during tests conducted to characterize the ETMV 8 coating is shown in Figure 21.

The actual heat flux through the coating may fall within the range of uncertainty shown in Figure 21. The accumulation of these uncertainties in heat flux and temperature measurements can affect the thermophysical characterization of the coating using inverse methods since such a characterization depends mainly on experimentally collected data. This aspect can also have an impact on the calculated energy balance during the melting and solidification phases of the material, as discussed in Section 5.1. However, these errors, equal to around 2% on the heat flux, as well as the inaccuracy of the temperature measurement (±0.5 °C), which generated about 4% error on the determination of the PCM liquid fraction, do remain quite acceptable, which allows validating the results proposed in this article.

6. Conclusions

This paper has presented a new numerical model to analyze the thermal behavior of a composite material integrating a bio-based micro-encapsulated PCM. This model has been developed by using the heat flux and temperature from measurements performed with the LGCgE experimental fluxmeter bench. Coupled with an inverse method, these data have allowed characterizing the parameters of the proposed numerical model.

Through this model, the thermodynamic analysis of the PCM is based on the enthalpy equations of a binary solution, associated with the equations of Gumbel’s Law. The latter serves to numerically describe the complex phenomena generated not only by supercooling during the material’s solidification phase but also by the presence of several components with distinct thermal behavior within the same PCM. The hybrid model developed during this work was coupled to inverse methods in order to characterize the thermophysical properties of a bio-sourced gypsum plaster incorporating a vegetable oil-based PCM, intended to improve the energy efficiency of buildings. An experimental error analysis was performed to evaluate the accuracy of this characterization, and verifications were carried out on model reliability according to various thermal dynamics and in the case of partial melting or solidification. Given the results obtained, the following conclusions could be drawn:

• This new hybrid model, associating the enthalpy of a binary solution with Gumbel’s Law equations, succeeds in describing the thermal behavior of the studied composite material (gypsum board + micro-encapsulated PCM) with high accuracy even when subjected to undercooling.
• The various dynamic thermal loads (ramps) between 10° and 30 °C depict the different thermal behaviors well represented by the hybrid model. The energy balances performed between these two temperature levels show a good level of repeatability and measurement accuracy.
• This model also takes into account the partial melting/solidification of the PCM, which can occur during its use in a building subjected to meteorological constraints.
• The model has successfully characterized the thermophysical properties of the coating incorporating PCM and, in particular, the latent heat contained in the material. By knowing that the coating is loaded with 70% PCM and in comparison, with the DSC measurements conducted by the manufacturer, the latent heat identified by characterization could be verified.
• An error calculation was used to determine uncertainties of around 4%, which appears acceptable for our application.

A subsequent work program would call for integrating this new hybrid model into a building simulation model in order to predict the energy performance of this type of coating.