# The code of simplified heat transfer model for temperature prediction in an insulated box equipped with phase change material

Tanathep Leungtongkum, Onrawee Laguerre, Denis Flick

► To cite this version:

Tanathep Leungtongkum, Onrawee Laguerre, Denis Flick. The code of simplified heat transfer model for temperature prediction in an insulated box equipped with phase change material. Software Impacts, 2023, 17, pp.100538. 10.1016/j.simpa.2023.100538. hal-04157413

HAL Id: hal-04157413

https://hal.inrae.fr/hal-04157413

Submitted on 10 Jul 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Original software publication

# The code of simplified heat transfer model for temperature prediction in an insulated box equipped with phase change material

<figure><img src="image_3.png" /></figure>

Tanathep Leungtongkum<sup>a,b,*</sup>, Onrawee Laguerre<sup>a</sup>, Denis Flick<sup>b</sup>

## ARTICLE INFO

### Keywords:

Insulated box

Phase change material

Simplified model

Heat transfer

Temperature

## ABSTRACT

This code presents a simplified heat transfer model based on zonal approach. It allows the temperature prediction (air and product) inside an insulated box with PCM (Phase Change Material as cold source) to maintain product at low temperature during food transport to assure its quality. The model involves airflow, conduction, convection and radiation. The studied geometry is a rectangular insulated box, with PCM on a sidewall and loaded with test product. The code requires a short calculation time and can be useful to researchers and stakeholders to study the effect of box design, PCM/product properties and mass and ambient temperature.

### Code metadata

Current code version

Permanent link to code/repository used for this code version

Permanent link to reproducible capsule

Legal code license

Code versioning system used

Software code languages, tools and services used

Compilation requirements, operating environments and dependencies

If available, link to developer documentation/manual

Support email for questions

V1

https://github.com/SoftwareImpacts/SIMPAC-2023-179

https://codeocean.com/capsule/1172946/tree/v1

MIT License

None

Python

Package manager needed: conda

Package needed: numpy and matplotlib.pyplot

tanathep.leungtongkum@inrae.fr

## 1. Motivation and significance

Insulated boxes equipped with a phase change material (PCM) are commonly used to transport temperature-sensitive products. They are often used in the last-mile delivery of small product quantities when refrigeration equipment is unavailable. However, temperature abuse (spatial and temporal temperature variation) in the product loaded in an insulated box equipped with PCM is often observed mainly due to poor insulation, insufficient PCM mass and inappropriate PCM positioning.

An empirical approach is still mainly used in practice. To understand the real phenomena in the system, Computational Fluid Dynamic

models can be used to predict the temperature field and its evolution with time. However, this requires significant calculation time and expertise in fluid mechanics. As an alternative, the simplified thermal model relating to this code can provide spatial and temporal temperature variation in an insulated box equipped with PCM with a short calculation time (less than 10 s using a computer with 64 GB of RAM). This short calculation time is vital to detect temperature abuse when it occurs in the cold chain. The stakeholders can immediately determine the impact on product temperature; thus, quality and safety and the appropriate actions can be implemented to reduce the damage. This code paper is linked with the research paper [1].

DOI of original article: https://doi.org/10.1016/j.ijrefrig.2023.02.009.

The code (and data) in this article has been certified as Reproducible by Code Ocean: (https://codeocean.com/). More information on the Reproducibility Badge Initiative is available at https://www.elsevier.com/physical-sciences-and-engineering/computer-science/journals.

* Corresponding author at: Université Paris-Saclay, INRAE, FRISE, 92761 Antony, France.

E-mail address: tanathep.leungtongkum@inrae.fr (T. Leungtongkum).

https://doi.org/10.1016/j.simpa.2023.100538

Received 28 April 2023; Received in revised form 12 June 2023; Accepted 21 June 2023

2665-9638/© 2023 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).

<footnote><sup>a</sup> Université Paris-Saclay, INRAE, FRISE, 92761 Antony, France</footnote>

<footnote><sup>b</sup> Université Paris-Saclay, INRAE, AgroParisTech, UMR SayFood, 91120 Palaiseau, France</footnote>
## 2. Simplified thermal model description

This model can be applied to an insulated box loaded by the product with PCM located at a side wall. During PCM melting, it releases cold and this allows to maintain low temperatures inside the box and assure sanitary quality throughout the supply chain. The code allows the prediction of air and product temperature evolution with time at different positions inside the box. The code can be used to study the effect of box dimension, insulation, product mass and its physical properties, PCM mass and its properties and ambient temperature in the supply chain.

The model considers a temperature heterogeneity inside the box, i.e., the difference between air, product surface and core, top and bottom, cold side (near the PCM), and warm side (near opposite wall). The product is considered as four blocks of length C with shell (C/4 length) and core (3C/4 length), Fig. 1. Airflow and heat transfer (conduction, convection, and radiation) in two dimensions (y and z directions) were considered. The model was validated by comparing the predicted air and product temperatures with the ones measured in a rectangular box loaded with test product (Tylose).

The model is based on the following assumptions:

- Heat transfer coefficient between internal air and solid surfaces (box walls, PCM container wall, product surface) is constant and is applied for the whole simulation period

- Air flows downward along the PCM container wall and then circulates in a clockwise loop along the other box walls. A mass flow rate of air is assumed as approximately proportional to the heat transfer coefficient.

- Thermal inertia of air (mass multiplied by the heat capacity) is assumed as negligible.

There are 12 solid zones and the following state variables describe the system at a given time:

- Average temperature in product shell: $T_{sh,n}$ when $n \in [1,4]$

- Average temperature in product core: $T_{c,n}$ when $n \in [1,4]$

- Wall temperature: $T_{w,n}$ when $n \in [1,4]$ where $T_{w,1}$ is the PCM surface temperature and $T_{w,2}$ to $T_{w,4}$ are of the internal walls

PCM is characterized by its temperature ($T_{pcm}$) and its ice fraction ($\phi$). Air temperature evolves along the air circulation loop from $T_{a,n}$ to $T'_{a,n}$ when $n \in [1,4]$ by exchanging with the product shells and from $T'_{a,n}$ to $T_{a,n+1}$ when $n \in [1,3]$ and from $T'_{a,4}$ to $T_{a,1}$ by exchanging with the internal walls.

Some of the equations representing air or product heat balances and radiative heat exchanges are presented here.

### Air temperature estimation from product shell and wall temperatures

The air at the top of the box (position I in Fig. 1a) exchanges heat with wall 1 (PCM surface) and its temperature changes from $T'_{a,1}$ to $T_{a,2}$. The equation governing the heat balance between the adjacent air and wall 1 can be written in the same manner as that in a heat exchanger.

$$ \dot{m}_a C_{p,a} dT_a = h_w (T_{w,1} - T_a) dA $$

$$ \ln \left( \frac{T_{a,2} - T_{w,1}}{T'_{a,1} - T_{w,1}} \right) = - \frac{\dot{h}_w A_{w,1}}{\dot{m}_a C_{p,a}} $$

Finally,

$$ (T_{a,2} - T_{w,1}) = \alpha_{w,1} (T'_{a,1} - T_{w,1}) \quad \text{with} \quad \alpha_{w,1} = \exp \left( - \frac{h_w A_{w,1}}{\dot{m}_a C_{p,a}} \right) \quad (1) $$

where $\dot{m}_a$ is the mass flow rate of air [kg s⁻¹]

$C_{p,a}$ is the specific heat capacity of air [J kg⁻¹ K⁻¹]

$h_w$ is the heat transfer coefficient between air and internal wall [W m⁻² K⁻¹]

$A_{w,1}$ is the area of wall 1 [m²]

This approach is also used when air exchanges heat with the product shell, for example from $T_{a,2}$ to $T'_{a,2}$. This leads to eight equations involving the eight air temperatures. They can be expressed in a matrix form.

$$ A \cdot T_a = B \quad \Rightarrow \quad T_a = A^{-1} \cdot B \quad (2) $$

where

$$ A = \begin{bmatrix} -\alpha_p & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & -\alpha_{w,1} & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & -\alpha_p & 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & -\alpha_{w,2} & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & -\alpha_p & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & -\alpha_{w,3} & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & -\alpha_p & 1 & 0 \\ 1 & 0 & 0 & 0 & 0 & 0 & 0 & -\alpha_{w,4} & 0 \end{bmatrix} $$

$$ T_a = \begin{bmatrix} T'_{a,1} \\ T'_{a,2} \\ T'_{a,3} \\ T'_{a,4} \end{bmatrix} \quad \text{and} \quad B = \begin{bmatrix} (1 - \alpha_p) T_{sh,1} \\ (1 - \alpha_{w,1}) T_{w,1} \\ (1 - \alpha_p) T_{sh,2} \\ (1 - \alpha_{w,2}) T_{w,2} \\ (1 - \alpha_p) T_{sh,3} \\ (1 - \alpha_{w,3}) T_{w,3} \\ (1 - \alpha_p) T_{sh,4} \\ (1 - \alpha_{w,4}) T_{w,4} \end{bmatrix} $$

58

59

where $\alpha_p$ is the dimensionless convective heat transfer coefficient between internal air and product shell [-]

This allows to estimate the air temperatures from the product shell and box's wall temperatures.

Convective heat transfer coefficient can be estimated from free convection correlations or local temperature measurements. Air mass flow rate can be determined from a developed relation with heat transfer coefficient. More detail can be found in Leungtongkum et al. [1].

### Radiative heat exchange

The radiative heat exchange between the lateral surface of product block 1 ($T_{s,1}$ in K) and wall 1 (PCM surface, $T_{w,1}$ in K) is shown in Eq. (3).

$$ q_{r,s1,w1} = \varepsilon_{w,1} \sigma (T_{s,1}^4 - T_{w,1}^4) CL \quad (3) $$

Where $q_{r,s1,w1}$ is the radiative heat exchange between surface of product block 1 and wall 1 [W]

$\varepsilon_{w,1}$ is the surface emissivity of wall 1 [-]

$\sigma$ is the Stefan-Boltzmann constant = $5.67 \times 10^8 \text{ W m}^{-2} \text{ K}^{-4}$

C is the length of the cross-section of the product block [m]

L is the length of the product block [m]

More details can be found in Leungtongkum et al. [1].

### Temperature evolution of product, walls, and PCM, and PCM ice fraction evolution

Eq. (4) is the unsteady heat balance equations for the shell of product block 1 ($T_{sh,1}$).

$$ MC_{p,sh} \frac{dT_{sh,1}}{dt} = \dot{m}_a C_{p,a} (T_{a,1} - T'_{a,1}) + \frac{T_{c,1} - T_{sh,1}}{R_{sh,c}} - q_{r,s1,w1} - q_{r,s1,w4} \quad (4) $$

where $MC_{p,sh}$ is the thermal inertia of product shell [J K⁻¹]

$R_{sh,c}$ is the heat transfer resistance from shell to core of product [K W⁻¹]

$q_{r,s1,w1}$ is the radiative heat exchange between surface of product block 1 and wall 1 [W]

$q_{r,s1,w4}$ is the radiative heat exchange between surface of product block 1 and wall 4 [W]

More details can be found in Leungtongkum et al. [1].

The same approach was applied to the 12 solid zones: shell of product blocks, core of product blocks and box internal walls.

Energy balance of PCM can be written for three different PCM states starting from solid state when its initial temperature is lower than melting temperature to liquid state after PCM is completely melted as follows.
<figure><img src="image_4.png" /><figcaption>Fig. 1. Simplified 2D heat transfer and airflow diagram in an insulated box with PCM on a side wall (a) Front view and (b) Perspective view of a quarter of product.</figcaption></figure>

– if PCM is completely frozen ($T_{PCM} < T_m$ and $\varphi = 1$) or completely melted ($T_{PCM} > T_m$ and $\varphi = 0$);

$$MC_{p.pcm} \frac{dT_{pcm}}{dt} = \frac{T_{w1} - T_{pcm}}{R_{pcm.w1}} + \frac{T_{ext} - T_{pcm}}{R_{w.ext.1}} \quad (5)$$

– if PCM is partially melted ($T_{pcm} = T_m$ and $0 < \varphi < 1$);

$$M_{pcm} \frac{d\varphi}{dt} \Delta H_{fus} = \frac{T_{w1} - T_{pcm}}{R_{pcm.w1}} + \frac{T_{ext} - T_{pcm}}{R_{w.ext.1}} \quad (6)$$

Where $MC_{p.pcm}$ [J K⁻¹] is the thermal inertia of PCM (different for frozen or melted state)

$R_{pcm.w1}$ is the heat transfer resistance from PCM to its container wall [K W⁻¹]

$R_{w.ext.1}$ is the heat transfer resistance from internal wall to the external [K W⁻¹]

$\varphi$ is the ice fraction of PCM [–]

$\Delta H_{fus}$ is the Latent heat of fusion [J kg⁻¹]

## Numerical solving

At a given time t, firstly, the eight air temperatures were deduced from the wall and product shell temperatures (Eq. (2)). To solve the eight linear algebraic equations, the (8 × 8) matrix A was inverted with the numpy.linalg library. Then, the eight radiative heat exchanges (example in Eq. (3)) were calculated. Finally, the thirteen heat balance equations, i.e., 4 wall temperatures, 4 product shell temperatures, 4 product core temperatures (example in Eq. (4)) and 1 PCM (temperature or ice fraction, Eqs. (5) and (6), respectively) were applied with an explicit scheme with 5 s time intervals (shorter time intervals led to the same results).

## 3. Code application

This code serves two purposes, first, for a simulation under a transient state, i.e., temperature evolution with time which is the case from the beginning of delivery (PCM partially melted) to its end (PCM completely melted). Second, for a simulation under steady state, i.e., when temperatures in the box reach equilibrium values while PCM is not entirely melted.
The code allows the users to adjust the input parameters as follows:

- Box dimension and thermal conductivity

- PCM dimension, mass and fusion temperature

- Product dimension and physical properties

- Supply chain conditions, i.e., ambient temperature and duration

This code can be applied to a rectangular box loaded with the product of the same form and a PCM slab on a side wall. It can be applied to other configurations, but some modifications of the code are necessary.

The developed code was validated by comparing with the air and product temperatures measured (1-minute intervals) in an insulated box loaded with a test product made of methylcellulose. The experiment was undertaken in a test room with several controlled ambient temperatures. More detail on constitutive equations and experimental validation can be found in Leungtongkum et al. [1].

## 4. Impact and conclusions

Food safety and security has become an urgent issue for several years. A way to maintain food quality is to store it under low temperature along the cold chain. Transportation using insulated boxes with Phase Change Material (PCM) can play an important role, particularly when cooling devices are not available [2]. Although food transport in an insulated box is practical and cost-effective, spatial and temporal temperature variations were observed [3,4] and may cause temperature abuse. This is due to insufficient PCM mass or inappropriate PCM positioning while an empirical approach is still mainly used in practice.

This software paper presents a code of simplified thermal model to predict air and load temperatures in an insulated box equipped with PCM and exposed to different time-temperature profiles as in a supply chain. This simplified model based on a zonal approach and coded in Python language, is original comparing to the finite elements or finite volumes CFD models. This model gives good precision as it involves all three heat transfer modes, i.e., conduction, convection, radiation, and airflow inside the box. The simulation takes less than 10 s using a computer with 64 GB of RAM comparing with more than 3 days by CFD for the box with the same configuration. This model can be easily used as a real-time prediction tool to estimate the food temperature inside an insulated box in a real supply chain. This code also allows the researchers and the stakeholders to investigate the effect of box configurations and operating conditions on spatial and temporal temperature variations. It is possible to couple this model with microbial growth models to quantify the contamination load in a supply chain until consumption. It can help stakeholders to limit food loss and waste and assure food safety.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Uncited references

[5], [6]

## Acknowledgments

The authors would like to thank Royal Thai Government Scholarship and Chulalongkorn University, Bangkok, Thailand for T. Leungtongkum's PhD scholarship. This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

## Appendix A. Supplementary data

Supplementary material related to this article can be found online at https://doi.org/10.1016/j.simpa.2023.100538.

## References

[1] T. Leungtongkum, O. Laguerre, D. Flick, Simplified heat transfer model for real-time temperature prediction in insulated boxes equipped with a phase change material, Int. J. Refrig. (2023) http://dx.doi.org/10.1016/j.ijrefrig.2023.02.009.

[2] J. Robertson, L. Franzel, D. Maire, Innovations in cold chain equipment for immunization supply chains, Vaccine 35 (17) (2017) 2252-2259, http://dx.doi.org/10.1016/j.vaccine.2016.11.094.

[3] O. Laguerre, H.M. Hoang, D. Flick, Experimental investigation and modelling in the food cold chain: Thermal and quality evolution, Trends Food Sci. Technol. 29 (2) (2013) 87-97, http://dx.doi.org/10.1016/j.tifs.2012.08.001.

[4] S. Mercier, S. Villeneuve, M. Mondor, I. Uysal, Time-temperature management along the food cold chain: A review of recent developments: Food preservation along the cold chain..., Compr. Rev. Food Sci. Food Safety 16 (2017) http://dx.doi.org/10.1111/1541-4337.12269.

[5] T. Leungtongkum, O. Laguerre, D. Flick, Tanathep/ temp-predict-box-pcm-side (Version 1.0) [Software], 2023, Available from https://github.com/Tanathep/temp-predict-box-pcm-side.

[6] T. Leungtongkum, O. Laguerre, D. Flick, Simplified-model-for real-time-temperature-prediction-in-insulated-boxes-equipped-with-a-PCM [Source Code], 2023, http://dx.doi.org/10.24433/CO.1221916.v1.