A regularized approach for derivative-based numerical solution of non-linearities in phase change static hysteresis modeling

Phase change materials (PCMs) represent a promising solution for thermal energy storage (TES) since they can store and release energy in the form of latent heat during solid ↔ liquid transitions. Nevertheless, accurately simulating the thermal behavior of PCMs remains challenging due to the non-linearities concerning latent heat effects and enthalpy hysteresis. This work introduces a stable and robust procedure based on the finite element method (FEM) under a mixed enthalpy–temperature formulation to address such non-linearities, which enables the numerical solutions using derivative-based algorithms such as the Newton–Raphson (NR) method. The static hysteresis model (SHM) is implemented in the FEM-based formulation via a regularization of the liquid fraction function in response to the sign of the temperature rate. This novel approach ensures a continuous and smooth heating ↔ cooling transition while retaining the SHM energy-conservative features to properly solve its non-linearities. The method is validated through a one-dimensional benchmark problem, demonstrating high performance and physical fidelity for both complete and partial phase changes. It achieves second-order convergence rates, ensures numerical stability even for large time steps, and maintains accuracy under diverse thermal boundary conditions. Finally, the method is extended to two-dimensional problems, highlighting its robustness and scalability for practical applications in TES systems.