ITP releases energy as u relaxes to ground state. That energy is fed as a local source into a heat equation for T(r,τ):
  ∂u/∂τ = ½∇²u + (K − 2P)u  (ITP relaxation)
  ∂T/∂τ = D ∇²T + Q(r,τ)  where Q(r,τ) = −∂e_wave(r,τ)/∂τ
With e_wave(r) = ½(u′)²·r² + (K(r) − 2P(r))·u²·r² the local energy density. Total bookkeeping: E_wave + ∫T·r²·dr·4π = const (energy conservation, no statistics).

Z N_e D (thermal diffusivity)  Run Stop Reset click Run

step: 0  E_wave: — Ha  ∫T: — Ha  E_total: — Ha  drift from start: — Ha