UDK 669.017.3:536.65 ISSN 1580-2949 Original scientific article/Izvirni znanstveni članek MTAEC9, 46(4)335(2012) CHALLENGES IN THE COMPUTER MODELING OF PHASE CHANGE MATERIALS IZZIVI V RAČUNALNIŠKEM MODELIRANJU MATERIALOV S FAZNO SPREMEMBO Lubomir Klimes1, Pavel Charvat1, Milan Ostry2 1Brno University of Technology, Faculty of Mechanical Engineering, Technicka 2896/2, 616 69 Brno, Czech Republic 2Brno University of Technology, Faculty of Civil Engineering, Vevefi 331/95, 602 00 Brno, Czech Republic yklime17@stud.fme.vutbr.cz Prejem rokopisa - received: 2011-10-20; sprejem za objavo - accepted for publication: 2012-03-01 Phase-change materials (PCMs) are a well-established category of materials with many possible applications, ranging from the stabilization of temperature to heat or cold storage. The main principle of PCMs is the utilization of the latent heat of a phase change for energy storage. Though many pure chemical elements can be used as PCMs, a PCM very often consists of a number of substances. The main reason for creating a PCM as a mixture of various substances is to achieve a desirable melting temperature for a particular application. However, these mixed PCMs require accurate and reliable methods for determining their physical properties, since for numerical modeling the thermal properties of materials and their proper determination represent a significant issue that considerably affects the accuracy and credibility of numerical simulations and their outcomes. The thermal properties of PCMs are usually obtained by differential scanning calorimetry (DSC), based on temperature and heat measurements of the investigated and reference materials. There are several approaches to modeling materials with phase changes. In this paper, the focus is on the enthalpy approach and the effective heat capacity method. Both techniques, which utilize the results from a particular DSC measurement, allow a treatment with the latent heat, and as a result the desired heat or cold storage may be effectively simulated. The presented methods were implemented using the results from DSC measurements in order to simulate a solar air collector with a PCM. The obtained results are presented and discussed. The paper also concerns the problems associated with the uncertainty of material properties and their influence on numerical simulations. Keywords: phase change material (PCM), enthalpy method, effective heat capacity method, solar air collector Materiali s faznimi premenami (PCM) so uveljavljena kategorija materialov z mnogimi področji uporabe: od stabilizacije temperature do shranjevanja toplote ali hladu. Glavni princip PMS je uporaba latentne toplote pri fazni premeni za shranjevanje energije. Čeprav je mogoče uporabiti kot PCM mnoge čiste kemijske elemente, se pogosto uporablja mešanica različnih elementov. Glavni razlog za uporabo mešanice za PCM je zagotavljanje želene temperature tališča za določen namen. Vendar pa mešanica PCM-ov zahteva natančne in izvedljive metode za določanje fizikalnih lastnosti. Numerično modeliranje termičnih lastnosti materialov zahteva zanesljivo določene termične lastnosti, kar bistveno vpliva na natančnost in zanesljivost simulacij in njihovih rezultatov. Termične lastnosti PCM-materialov navadno določamo z DSC-metodo, ki temelji na merjenju temperature in toplote preiskovanega in referenčnega materiala. Obstaja več načinov, kako modelirati materiale s faznimi premenami. V tem članku je glavna pozornost posvečena entalpijskemu približku in metodi efektivne toplotne kapacitete. Obe tehniki, ki uporabljata rezultate DSC-meritev, omogočata uporabo latentne toplote, kar omogoča učinkovito simuliranje shranjevanja toplote ali hladu. Predstavljeni metodi sta bili uporabljeni z uporabo rezultatov iz DSC-meritev za simulacijo sončnega kolektorja s PCM-materiali. Predstavljeni so rezultati z diskusijo. Članek obravnava tudi problem nezanesljivosti podatkov o lastnostih materialov in njihov vpliv na numerično simulacijo. Ključne besede: materiali s faznimi premenami (PCM), entalpijska metoda, metoda učinkovite toplotne kapacitete, sončni kolektor PCMs are usually composed of several substances 1 INTRODUCTION rather than being a pure chemical element. The reason is that a particular application requires a specific range of Due to the increasing consumption of energy, the phase change according to the operating conditions and a demands for thermal comfort and the need to decrease desired melting temperature can be achieved with a temperature fluctuations and energy losses, new mixture of several components with a certain chemical construction and thermal storage materials have been composition. developed over the past I0 years. Phase change materials The optimal layout, usage, integration, and design of (PCMs) are substances that are used to absorb, store and PCMs in buildings, constructions and other equipment release energy. This is accomplished by utilizing the cannot be carried out without using numerical and com-latent heat of the phase change from the liquid to the puter simulations, optimization methods and experi-solid state and vice versa, since the amount of energy mental measurements. In order to obtain accurate and included in a phase change is usually considerably larger credible results, it is essential to properly determine all than the sensible heat1. Nowadays, PCMs are being used the physical properties of the used materials and to in a wide range of applications in architecture and civil precisely describe the operating conditions2. For this engineering, from heat or cold storage to temperature purpose differential scanning calorimetry (DSC) is very stabilization and balancing. often used. There are several methods and approaches as to how the materials embracing phase change and related phenomena can be modeled. In this paper the attention is focused on the enthalpy approach and the effective heat capacity method that enable the simulation of the heat transfer in materials with phase changes. As the inputs and parameters (e.g., the physical properties of a material) of the simulation the results from a DSC measurement may be utilized. The enthalpy method4,5 introduces the volume enthalpy H/(J m-3), which includes the phase change phenomenon, to the equation since the enthalpy represents the sum of the sensible heat and the latent heat of the phase change. The function of the volume enthalpy is defined as follows: H( T) = ^ , d/s PC -pLf dT de (2) 2 DIFFERENTIAL SCANNING CALORIMETRY Differential scanning calorimetry is a thermo-analyti-cal method and is widely used to assess and construct phase diagrams and to determine the thermophysical properties of materials, e.g., the crystallization and melting temperatures, heat capacity and the enthalpy3. The principle of differential scanning calorimetry is a measurement of the difference in the amount of heat needed to raise the temperature of the sample and the reference materials. The aim is to keep the temperature of both the sample and the reference materials at the same level during the experiment, while the temperature increases linearly. In order to precisely determine the properties and parameters of a sample material, the thermophysical properties of the reference material have to be known. The particular outcome of a DSC measurement for a PCM is shown in Figure 1. 3 ENTHALPY METHOD The numerical modeling of PCMs requires solving the heat transfer problem in order to obtain the transient temperature distribution for the examined material. The heat transfer is described by the heat equation that is a partial differential equation of the second order: dT —-aAT = 0 (1) dt where T/K is the temperature, t/s is the time and a/(m2 s-1) is the thermal diffusivity of the material. where p/(kg m-3) is the density, c/(J kg-1 K-1) is the specific heat, Lf/J is the latent heat of the phase transformation and /s is the solid fraction. By substituting the enthalpy relation (2) into equation (1) we can obtain the reformulated equation describing the unsteady heat transfer in a material with phase changes: dH dt = V( AVT) (3) where A/(W m-1 K-1) is the thermal conductivity. As can be seen, in equation (3) there are two unknown quantities - the enthalpy and the temperature. However, they are coupled together through the enthalpy-temperature (H-T) relationship that can be obtained from a DSC measurement. An example of the H-T relationship for a particular material with a phase change range from 25 °C to 27 °C is shown in Figure 2. In order to solve equation (3) for a particular application the corresponding initial and boundary conditions have to be provided. To solve the system numerically, an appropriate numerical method is needed to assemble a system of discretized equations. For this purpose many numerical methods can be utilized, e.g., the finite-difference method, the finite-volume method, the finite-element method or the mesh-free methods6. For simplicity, let us assume the finite-difference method. Thus, the partial derivatives in equation (3) are replaced by their finite-difference approximations and the original system is reformulated into a system of linear algebraic equations: Hn+1 _ Hn I /(T" T^ ) ■"i, j, k ^ i, j, ^ W V-^ i, j, ^ ^ i±Ax, j±&x , k ±A7' (4) Figure 1: Results of a DSC measurement of a PCM Slika 1: Reznltatr merrtev DSC s PCM However, due to the presence of unknown values of both the enthalpy and the temperature, the explicit scheme of the finite-difference method for the discretization in time has to be used. The temperature field is calculated in two successive steps: in the first step the enthalpy H(n + 1) in time t + 1 is calculated from the known enthalpy and temperature values Ht, Tt in time t according to equation (4). After that the temperature T(n +1) in time t + 1 is calculated from the enthalpy H(n +1) according to the enthalpy-temperature relationship obtained, for instance, from a DSC measurement. Repeating the described procedure over all the time nodes leads to the targeted transient temperature field. 4 EFFECTIVE CAPACITY METHOD The effective capacity method is based on the idea of incorporating the phase change phenomenon into the heat capacity calculations4,5. This physical property represents the amount of heat needed for the matter to increase its temperature. It means that the effective heat capacity Ceff/(J kg-1 K-1) has to include the latent heat of the phase change. Therefore, the function of the effective heat capacity increases and decreases sharply with an apparent peak when the material undergoes the phase change. An example of the effective heat capacity function corresponding to the enthalpy curve is shown in the same Figure (Figure 2). In order to calculate congruous and trustworthy results it is crucial to properly determine the effective heat capacity. For this purpose differential scanning calorimetry is one of the tools and techniques available. By introducing the effective heat capacity, equation (1) becomes: dT PCeff = V(AVT) (5) In equation (5) there is only one unknown variable, the temperature T, and therefore a wider class of discretization approaches can be used. One of the most significant features is that the implicit scheme can be utilized. Hence, using the temperatures T n from the foregoing time step for the determination of the effective heat capacity in a successive time n + 1, the discretized equation for each node i, j, k can be written in the form: -ffpn +1 rp n +1 \ _ rp n J (T i, j, k , T i ±Ax, j ±Ax, k ±Az , ...) ~ T i, j, j, k (6) The system of algebraic equations can be assembled by using equation (6) applied to all the nodes of the discretized domain and an appropriate numerical method has to be used to solve it. X 10 xlO' ^ 6 1 a. 1 ^ I 3 "o > 2 1 -Effective heat capacity 1 1 1 1 1 1 1 1 1 1 1 1 1 3 I S ■S 2 o ■B 1 S 21 22 23 24 25 26 27 28 29 30 31 Temperature fCJ Figure 2: The enthalpy temperature relationship and the effective heat capacity function Slika 2: Sprememba entalpije v odvisnosti od temperature in efektivne toplotne kapacitete 5 COMPARISON OF METHODS The enthalpy method is based on the introduction of an additional quantity to the governing equation in order to include the phase change. From the computing point of view the volume enthalpy function introduces difficulties. Namely, (a) the desired quantity, i.e., the temperature, cannot be computed directly and the computation is split into two steps with an intermediate calculation, (b) the temperature search from a known enthalpy is demanding and requires a powerful search algorithm and (c) without additional techniques the explicit scheme, which is conditionally stable, must be used. The advantages of the enthalpy method are a better numerical stability when an explicit discretization scheme in time is used for both methods and its better accuracy in comparison with the effective heat capacity method4,5. The effective heat capacity method calculates the temperature field directly without intermediate steps and it allows the use of an implicit discretization scheme that is unconditionally stable6. However, without additional techniques an assumption is required that the effective heat capacity is calculated using the temperature from a previous step. In the case that the effective heat capacity is a function of an unknown time, this makes the system of equations considerably more difficult to solve, and therefore the possibility of using an arbitrary time step (owing to the implicit scheme) is strictly confined. Moreover, when a problem with a large number of nodes (a 3D case with a large and intricate geometry) the computational demands of solving the system of equations becomes more arduous and time consuming. 6 RESULTS AND DISCUSSION Both discussed methods were implemented in order to simulate and analyze the solar air collector, which is a solar thermal device that is used for heating in buildings. An issue with solar thermal systems is the general need for thermal storage in order to balance the heat demands over a certain time period. Therefore, the utilization of Figure 3: The results of the numerical simulation based on the effective heat capacity for the solar air collector Slika 3: Rezultati numericne simulacije na osnovi efektivne kapacitete za soncni kolektor zraka phase change materials and their feature of energy storage can overcome the mentioned problem. The TRNSYS software, which is a simulation tool designed for the analysis of transient systems, was used for modeling the solar air collector. Owing to the phase change material, in the case of the enthalpy method, the numerical model of heat transfer in a PCM has been implemented in MATLAB and both software packages, TRNSYS and MATLAB, were coupled together. In the case of the effective heat capacity method the numerical model of the PCM was implemented in the form of a stand-alone module (a DLL library) for TRNSYS in the programming language C++. Both the software implementations utilized a result obtained from a DSC measurement of the PCM with the phase change temperature ranging from 30 °C to 38 °C and the latent heat of the phase change being 70 kJ kg-1. The simulations performed with the described methods produced comparable results. Figure 3 shows the outcomes of the simulation (meteorological data for 24 h in June, using the effective capacity method, the finite-difference discretization and the explicit scheme, the time step 2 seconds) for the so-called front-pass solar air collector that was split longitudinally into five separate sections. As seen in Figure 3 the courses of the PCM temperatures show the phase change in the given temperature ranges, the absorption and the release of the latent heat. However, even though the material properties are properly determined by the DSC there are still quantities under uncertainty that may considerably influence the result of the simulations. A typical example is the heat transfer coefficient (HTC), since it can randomly fluctuate according to various factors. The HTC belongs to the class of operating conditions that can significantly affect the behavior of the system through the initial and boundary conditions. In order to properly solve the systems with parameters under uncertainty, additional techniques and methods, e.g., stochastic modeling and stochastic optimization, must be considered. Moreover, the presented methods and their implementations can be utilized in a wide range of applications, e.g., for steel solidification and optimal control in the steel industry7. 7 CONCLUSION The presented paper introduces two possible approaches to modeling materials undergoing phase changes. Both techniques assume some particular thermophysical properties obtained, e.g., from the differential scanning calorimetry method. The concepts of the enthalpy approach and the effective heat capacity are described and the advantages, drawbacks and the usage are discussed. The implemented methods were successfully tested on the simulation of the solar air collector with a PCM absorber. The obtained results suggest that the described approaches of numerical modeling and their implementations are valuable tools for solving problems in optimizing the design, usage and the integration of materials with phase changes. Acknowledgement The presented research was supported by the project OC10051 of the Czech Ministry of Education and by the BUT project BD13102003 for young researchers. The main author, the holder of the Brno PhD Talent Financial Aid sponsored by Brno City Municipality, also gratefully acknowledges the financial support. 8 REFERENCES 1F. Kuznik, D. David, K. Johannes, J. J. Roux, A review on phase change materials integrated in building walls, Renewable and Sustainable Energy Reviews, 15 (2011), 379-391 2F. Kuznik, J. Virgone, J. J. Roux, Energetic efficiency of room wall containing PCM wallboard: A full-scale experimental investigation, Energy and Buildings, 40 (2008), 148-156 3 G. H. Zhang, C. Y. Zhao, Thermal and rheological properties of microencapsulated phase change materials, Renewable Energy, 36 (2011), 2959-2966 4 D. M. Stefanescu, Science and Engineering of Casting Solidification, 2nd ed., Springer, 2009 5 M. Muhieddine, E. Canot, R. March, Various approaches for solving problems in heat conduction with phase change, International Journal On Finite Volumes, 6 (2009) 1 6 Handbook of numerical heat transfer. Editors W. J. Minkowycz, E. M. Sparrow, J. Y. Murthy, 2°d ed., Wiley, 2006 7 J. Stetina, F. Kavicka, The influence of the chemical composition of steel on the numerical simulation of a continuously cast slab, Mater. Tehnol., 45 (2011) 4, 363-367