Thermal Evolution Simulations

Updated 7 December 2025

Thermal evolution simulations are computational models that solve time-dependent heat redistribution using coupled PDEs for various physical systems.
They employ diverse numerical methods, such as finite-volume schemes, adaptive mesh refinement, and PINNs, to capture conduction, convection, and multiphysics interactions.
These simulations provide critical insights into planetary differentiation, neutron star cooling, and industrial processes while addressing challenges like scale and microphysical uncertainties.

Thermal evolution simulations encompass a wide class of computational models that solve for the time-dependent redistribution of thermal energy in complex systems, typically governed by coupled partial differential equations (PDEs) that express energy conservation, heat transport, and, where applicable, the influence of non-thermal physics (e.g., magnetic, chemical, mechanical, or compositional evolution). These simulations are central to the paper of planetary and stellar interiors, neutron stars, interstellar media, engineered materials, and more, providing predictive insight into time-dependent structural, observational, and dynamical properties.

1. Principles and Governing Equations

The fundamental mathematical foundation of thermal evolution simulations is the energy conservation equation, appropriately generalized for the composition and geometry of the physical system. For planetary and stellar interiors, the governing equation usually takes the form: $\rho\,c_p\,\frac{\partial T}{\partial t} = \nabla \cdot (\kappa \nabla T) + Q$ where $\rho$ is density, $c_p$ is specific heat capacity, $\kappa$ the (tensorial, possibly anisotropic) thermal conductivity, and $Q$ the sum of source/sink terms (e.g., radiogenic heating, chemical energy release, neutrino losses, Ohmic/Joule dissipation). Extensions incorporate additional coupled PDEs—magnetic induction for magnetized objects (Dehman, 30 Apr 2024), stress-strain relations for coupled thermo-mechanical problems (Sharma et al., 25 Dec 2024), and chemical or compositional rate equations (Gray et al., 2013).

Layered or multiphase systems (e.g., magma ocean, mush, solid mantle, and core) are discretized with dynamically evolving internal boundaries. Phase changes introduce latent heat and necessitate tracking melt fractions as in the lunar and planetary simulations (Sahijpal et al., 2020, Herath et al., 17 Sep 2024, Zhang et al., 2022, Henke et al., 2011).

Convection models (mixing-length, Nusselt–Rayleigh scaling, or full mantle convection) augment Fourier conduction by parameterized or fully resolved turbulent heat transport as a function of local gradients and thermophysical properties (Zhang et al., 2022, Plesa et al., 2022, Herath et al., 17 Sep 2024).

For non-solid systems (e.g., protoplanetary disks (Vorobyov et al., 2020)), thin-disk approximations and local energy equations incorporate radiative cooling/heating, viscous dissipation, irradiation, and, where relevant, distinct gas and dust temperatures.

2. Computational Methodologies

Thermal evolution simulations are implemented using a diverse array of computational techniques tailored to the system of interest:

Finite-volume and finite-difference schemes: Standard for spherically symmetric or multidimensional solid bodies (planetary interiors, neutron stars), with implicit time integration schemes for parabolic PDEs to handle stiff source terms and enable large evolutionary time-steps (Dehman, 30 Apr 2024, Ascenzi et al., 28 Jan 2024, Henke et al., 2011).
Moving-mesh or Lagrangian approaches: Applied when interior structure (density, pressure, phase boundaries) evolves significantly due to contraction/expansion work and phase separation; the Henyey/NR solver framework is common in planetary models (Zhang et al., 2022).
Spectral and pseudo-spectral decompositions: Used in global 3D geodynamics and magneto-thermal simulations to represent angular or spatial structures in spherical or periodic domains (e.g., neutron-star crusts (Grandis et al., 2020), Martian mantle (Plesa et al., 2022)).
Adaptive mesh refinement (AMR): Critical in resolving sharp gradients, thin boundary layers, or localized condensations in high-resolution studies of radiative-MHD systems (coronal rain, thermal instabilities) (Kohutova et al., 2020, Sen et al., 2023).
Physics-informed neural networks (PINNs): Recently advanced for rapid surrogate modeling of thermal-mechanical evolution in industrial processes, where PINNs couple neural network architectures to the underlying PDE system with automatic differentiation and data assimilation (Sharma et al., 25 Dec 2024).
Tensor network methods: For quantum many-body thermal evolution, MPO/MPS representations and tangent space time-dependent variational principle (TDVP) evolution accommodate tractable imaginary-time cooling and fixed-filling quantum statistical mechanics (Li et al., 10 Nov 2025).

3. Physical Ingredients and Parameterizations

Accurate thermal evolution simulations require detailed microphysical (and sometimes macrophysical) models:

Thermal conductivity $\kappa$ and specific heat $c_p$ : Parameterized as functions of temperature, composition, and phase (e.g., solid, partial melt, mush)—including tensorial anisotropy in magnetized or layered media (Ascenzi et al., 28 Jan 2024, Dehman, 30 Apr 2024).
Viscosity $\eta$ : Strongly temperature and phase dependent, critical in convection, creep, and sintering models (Arrhenius or piecewise-constant parameterizations) (Henke et al., 2011, Herath et al., 17 Sep 2024, Plesa et al., 2022).
Radiative and radiogenic source terms: Incorporation of all relevant isotopic heating channels ( $^{26}$ Al, $^{60}$ Fe, etc.), stellar irradiation, and background heating/cooling terms (Henke et al., 2011, Zhang et al., 2022, Vorobyov et al., 2020).
Phase transitions and latent heat: Explicit treatment of melting/freezing fronts, tracking solidus/liquidus curves and integrating latent heat into the enthalpy or energy balance (Sahijpal et al., 2020, Herath et al., 17 Sep 2024, Zhang et al., 2022, Henke et al., 2011).
Chemistry and non-equilibrium cooling: Coupled chemical networks (e.g., for collapsing filaments, protoplanetary disks) (Gray et al., 2013, Vorobyov et al., 2020) feed back into the thermal balance via molecular/atomic cooling rates.
Ohmic/Joule heating, neutrino and photon emission: Important in magnetized compact objects (neutron stars) and laboratory plasmas; microphysics drawn from first-principles electron, ion, and nucleon transport (Dehman, 30 Apr 2024, Ascenzi et al., 28 Jan 2024, Viganò et al., 2021, Potekhin et al., 2019).
Tidal and mechanical dissipation: For planets and moons subject to external forcing, tidal heating is modeled either by empirical scaling or via dynamical Love-number approaches (Herath et al., 17 Sep 2024).

4. Validation, Benchmarking, and Parameter Studies

Thermal evolution codes are rigorously validated through analytic benchmarks, laboratory data, or comparison to observed diagnostic constraints. For example:

Porous and sintering models are calibrated against laboratory compaction and thermal conductivity measurements for granular materials (Henke et al., 2011).
Planetary case studies: Simulations for chondritic parent bodies, the early Moon, Mars, and rocky exoplanets are compared to meteoritic closure ages, crustal thickness, present-day seismic or heat-flow diagnostics, and mantle processing rates (Henke et al., 2011, Sahijpal et al., 2020, Zhang et al., 2022, Plesa et al., 2022).
Magneto-thermal neutron star codes are cross-benchmarked via prototypical field configurations and compared to analytic heat-diffusion/surface-temperature solutions, and specific observational phenomena (X-ray luminosities, pulse profiles) (Dehman, 30 Apr 2024, Ascenzi et al., 28 Jan 2024, Grandis et al., 2020, Viganò et al., 2021, Potekhin et al., 2019).
Protoplanetary disk models: Testing thermal schemes (β-cooling, viscous/radiative equilibrium, two-temperature) under controlled opacities, irradiation, and initial disk conditions, and evaluating fragmentation thresholds using the Jeans criterion (Vorobyov et al., 2020).
Sensitivity studies: Systematic scans over body size, initial porosity, irradiation/insolation, viscosity contrasts, and heating rates to map out the regimes of solidification, core formation, dynamo activity, and atmospheric escape (Herath et al., 17 Sep 2024, Zhang et al., 2022, Kubyshkina et al., 2020).

5. Representative Applications and Physical Insights

Thermal evolution simulations yield insight across astrophysics, geodynamics, and materials science:

Planetesimal and planetary differentiation: Establishing timescales for core formation, magma ocean solidification, and subsequent cooling, as in the H-chondrite parent body and the early Moon (Henke et al., 2011, Sahijpal et al., 2020). For tidally locked lava planets, hemispheric thermal evolution is critically sensitive to melt viscosity, day-night coupling, and tidal dissipation (Herath et al., 17 Sep 2024). Earth's dynamo longevity and that of super-Earths depends on core-mantle boundary heat flow, mantle blanketing, and radiogenic inventory (Zhang et al., 2022).
Neutron star magneto-thermal evolution: Hall-driven cascades, Ohmic decay, and feedback from crustal heat deposition govern the observable X-ray luminosity, temperature distribution, persistence of hot spots, and timing properties over $10^3$ – $10^6$ yr, with full 3D evolution revealing non-axisymmetric patterns inaccessible to prior 2D theory (Dehman, 30 Apr 2024, Ascenzi et al., 28 Jan 2024, Grandis et al., 2020, Viganò et al., 2021).
Accreting compact objects: Evolving crustal composition in neutron stars affects the thermal quiescent luminosity and the diagnostic power of quiescent X-ray emission, discriminating between modified Urca and direct Urca cooling, and probing superfluid gap magnitudes (Potekhin et al., 2019).
Protoplanetary and circumstellar disks: The interplay of radiative, viscous, and irradiation heating, and the decoupling of dust and gas temperatures, determines disk structure, fragmentation, and the site of chemical/dust evolution (Vorobyov et al., 2020). Self-consistent coupling of thermal contraction with atmospheric escape is essential to infer true planetary radius–envelope mass tracks over Gyr timescales (Kubyshkina et al., 2020).
Plasma and astrophysical flows: Simulations of 3D radiative MHD in solar and astrophysical contexts elucidate how reconnection-driven thermal instabilities lead to condensation formation (coronal rain, prominences) and link the statistics of energy deposition with observable phenomena (Sen et al., 2023, Kohutova et al., 2020).
Quantum statistical mechanics: Advanced tensor network algorithms enable fully controlled thermal evolution at fixed filling, yielding accurate thermodynamic and correlation data for lattice quantum systems (e.g., the Hubbard model) and mapping temperature scales for emergent phases (Li et al., 10 Nov 2025).
Engineering and materials processing: Physics-informed ML surrogates (PINN frameworks) provide mesh-free, transferable predictions for thermal stress and temperature fields in metal additive manufacturing, enabling rapid “soft-sensing” within operational loops (Sharma et al., 25 Dec 2024).

6. Limitations, Challenges, and Extensions

Thermal evolution simulations are limited by:

Microphysical uncertainties: Conductivities, viscosities, yield criteria, superfluid gap magnitudes, and EOS parameters critically affect outcomes, especially under extrapolation to extreme P–T regimes (Dehman, 30 Apr 2024, Potekhin et al., 2019).
Multiphysics and scale coupling: Strongly nonlinear coupling between thermal, magnetic, mechanical, and chemical degrees of freedom can generate numerical stiffness, require implicit integration, and challenge parallel scaling (e.g., crust-core coupled magneto-thermal evolution) (Ascenzi et al., 28 Jan 2024, Viganò et al., 2021).
Dimensionality and resolution: Full 3D simulations are computationally demanding (e.g., MATINS/Parody codes for neutron stars, mantle convection for Mars), often necessitating adaptive, scalable algorithms and simplified physics; local small-scale structures (current sheets, filaments, phase fronts) demand high dynamic range (Ascenzi et al., 28 Jan 2024, Grandis et al., 2020, Sen et al., 2023).
Model assumptions: Assumptions such as instantaneous formation, zero atmosphere, constant composition, and axisymmetry may limit realism in certain regimes (Henke et al., 2011, Herath et al., 17 Sep 2024, Zhang et al., 2022).

Extensions under active research include:

Adaptive neural architectures and hybrid data–physics models for industrial and geoscience applications (Sharma et al., 25 Dec 2024).
Multi-phase, multi-component, anelastic or compressible convection for planetary/stellar interiors (Plesa et al., 2022).
Coupling to observational inference for thermophysical property inversion (surface emission, lightcurves, phase curves) (Herath et al., 17 Sep 2024, Ascenzi et al., 28 Jan 2024, Kohutova et al., 2020).
Higher fidelity chemistry/radiation transport in ISM/cloud collapse models (Gray et al., 2013).

7. Broader Implications and Observational Connection

Thermal evolution simulations provide essential quantitative predictions for the time-dependent state of astrophysical and planetary bodies, directly connecting microphysics and dynamics to observables:

Surface–interior coupling: Planetary and stellar surface observables (luminosity, crustal structure, tectonics, dynamo-generated magnetic fields) are diagnostic of deep interior thermal transport.
Population synthesis: Atmospheric retention and radius distributions in sub-Neptunes, fragmentation statistics in protoplanetary disks, and neutron star luminosity populations can be modeled only with self-consistent thermal-evolution frameworks (Kubyshkina et al., 2020, Vorobyov et al., 2020, Potekhin et al., 2019).
Laboratory and engineering translation: Transferring methodologies from astrophysical to manufacturing domains (e.g., PINNs for additive manufacturing), improving speed and transferability of high-fidelity simulations (Sharma et al., 25 Dec 2024).

The domain of thermal evolution simulations is thus central to the interpretation of present and future mission data (e.g., InSight, JWST, NICER), bridging advanced computational physics with multi-disciplinary observational programs.