$$\rho c_p \frac{\partial T}{\partial t} = \nabla \cdot (k \nabla T) - \rho L \frac{\partial f_l}{\partial t}$$
Heat equation including phase change (enthalpy method): \(L\) latent heat of fusion [J/kg], \(f_l\) liquid fraction
$$\mathrm{Ste} = \frac{c_p (T_{wall} - T_m)}{L}$$
Stefan number: ratio of sensible heat to latent heat. When Ste < 1, latent heat dominates
$$s(t) = 2\lambda\sqrt{\alpha t}$$
Neumann solution (semi-infinite body): \(s\) melting interface position [m], \(\alpha = k/(\rho c_p)\) thermal diffusivity
What is Phase Change Material (PCM) Analysis?
Physical Model & Key Equations
The core of the simulation is the 1D heat conduction equation, reformulated using enthalpy (H) to naturally incorporate the phase change. The governing equation solved at each node in the finite difference grid is:
$$\rho \frac{\partial H}{\partial t}= \frac{\partial}{\partial x}\left( k(T) \frac{\partial T}{\partial x}\right)$$Here, $\rho$ is the density (which you can set), $H$ is the total enthalpy per unit mass, $t$ is time, $k(T)$ is the temperature-dependent thermal conductivity, and $T$ is temperature. The key is that enthalpy $H$ includes both sensible and latent heat.
The relationship between enthalpy and temperature defines the phase change. It's a piecewise function that creates the "mushy zone":
$$ H(T) = \begin{cases}c_{p,s}\cdot T & T < T_m \\ c_{p,s}\cdot T_m + L_f \cdot \phi + c_{p,l}\cdot (T - T_m) & T = T_m \\ c_{p,s}\cdot T_m + L_f + c_{p,l}\cdot (T - T_m) & T > T_m \end{cases}$$Where $T_m$ is the melting point, $L_f$ is the latent heat of fusion, $c_{p,s}$ and $c_{p,l}$ are the specific heats for solid and liquid, and $\phi$ is the liquid fraction (between 0 and 1). This is why you see a jump in the enthalpy plot—it's the latent heat $L_f$ being added at $T_m$.
Frequently Asked Questions
Real-World Applications
Building Thermal Management: Paraffin-based PCMs are integrated into wallboards or ceiling tiles. They melt during the day, absorbing excess heat to keep a room cool, and solidify at night, releasing the stored heat. This simulator helps engineers size the PCM layer thickness and select the optimal melting temperature based on local climate.
Electronic Device Cooling: High-performance electronics generate intense, intermittent heat bursts. A PCM package attached to a CPU can absorb this heat during a spike, preventing thermal throttling, while giving fans and heat sinks more time to dissipate the energy. The conductivity values (k_s, k_l) are critical design parameters here.
Concentrated Solar Power (CSP): Molten salts (like the preset in the simulator) are used as both heat transfer fluid and thermal storage medium in CSP plants. They store solar energy as latent heat during the day to generate steam and electricity at night. Accurate simulation of their melting/freezing cycles is essential for plant efficiency and reliability.
Food Processing & Cold Chain: In shipping containers for temperature-sensitive pharmaceuticals or food, ice or other PCMs are used as thermal buffers. The simulation helps predict how long the payload will stay within the safe temperature range based on the PCM's properties and the external ambient temperature (set by the "Surface Temperature T_s" parameter).
Common Misunderstandings and Points to Note
When you start using this simulator, there are a few points beginners often stumble on. A major misconception is the idea that simply increasing the latent heat will improve performance. While it's true that a larger melting latent heat L_f increases thermal storage capacity, if the thermal conductivity k is low, heat struggles to penetrate deep into the material, preventing you from utilizing that capacity. For example, paraffin has high latent heat but low thermal conductivity (about 0.2 W/mK). Try comparing "Thermal Conductivity k" values of 0.1 and 0.5 in the tool. You'll see that even with the same latent heat, the speed of heat transfer is completely different, significantly altering the half-domain melting time. In practice, you often need to incorporate fins or metal matrices to enhance thermal conductivity.
Next, consider the realism of boundary condition settings. This tool uses the simple condition of a "constant surface temperature." However, in real applications like building walls, external air temperature fluctuates with solar radiation, and convective heat transfer cannot be ignored. Simulating with a fixed "Surface Temperature T_s" in the tool is a good first step, but for actual design, you need to consider boundary conditions of the third kind (convective boundary conditions). Try settings that rapidly raise and lower the surface temperature to get a feel for the difference from reality.
Finally, understand the limitations of the 1D model. This calculation assumes shapes where one-dimensional heat flow is dominant, like slabs or infinitely long cylinders. Real thermal storage units experience three-dimensional heat spreading and the effects of natural convection. This is why when you select "Ice" to melt in the tool, it melts differently than an actual block of ice. Remember, this tool is for understanding basic principles and checking parameter sensitivity; more advanced simulation is needed for detailed design.
How to Use
- Select a PCM material: ice (Tm=0°C, Lf=334 kJ/kg), paraffin (Tm=50°C, Lf=200 kJ/kg), or molten salt (Tm=290°C, Lf=350 kJ/kg). System auto-populates melting temperature (Tm) and latent heat (Lf).
- Input material density (rho in kg/m³) and specific heat capacity (cps in kJ/kg·K). For paraffin wax: use rho=900 kg/m³, cps=2.5 kJ/kg·K.
- Run the 1D finite difference simulation. Solver applies the enthalpy method across 50 nodes over 1000 seconds. Observe temperature profile evolution and phase front propagation in the visualization panel.
Worked Example
Simulate ice melting at domain boundaries (T_boundary=5°C). Set Tm=0°C, Lf=334 kJ/kg, rho=917 kg/m³, cps=2.09 kJ/kg·K. After 600 seconds, the phase front penetrates approximately 35mm into a 100mm domain. Enthalpy H increases from 0 (solid ice) to 334 kJ/kg at the melting interface, then continues rising in liquid phase. Temperature remains pinned at 0°C across the mushy zone until complete fusion.
Practical Notes
- Latent heat dominates energy transport during phase change—cps has minimal effect near Tm. A 10% increase in boundary temperature accelerates melting by ~40% due to Stefan number effects.
- Paraffin wax in building thermal storage: Tm=47°C, Lf=190 kJ/kg. Density variation (solid 900 vs. liquid 770 kg/m³) is ignored in this 1D fixed-grid solver; account separately for expansion loads.
- Molten salt (solar thermal systems): very high latent heat requires extended simulation time. Thermal diffusivity α=k/(rho·cps)≈1.2e-7 m²/s limits front speed to ~0.3 mm/s in realistic geometries.
