ITADN

ResampledCooling does not subtract magnetic energy when computing eint

#1885OpenBenWibking 创建于 2026-05-11
bug: wrong answer/failure/crashpriority:highMHDblocking-release
B
BenWibkingcommented
### Discussed in https://github.com/quokka-astro/quokka/discussions/1883 <div type='discussions-op-text'> <sup>Originally posted by **YaoguangPei** May 11, 2026</sup> I found a cooling-related issue when testing an MHD run with resampled cooling in QUOKKA. My setup is: - is_hydro_enabled = true - is_mhd_enabled = true - is_radiation_enabled = false - periodic box - uniform initial magnetic field ``` [ \mathbf{B}0 = (5,\mu{\rm G},,0,,0) ] ``` The main issue seems to be in the cooling path: In cooling/ResampledCooling.hpp, cooling updates use ```cpp const Real Eint = RadSystem<problem_t>::ComputeEintFromEgas(rho, x1Mom, x2Mom, x3Mom, Egas); ``` However, in src/radiation/radiation_system.hpp, the implementation is currently ```cpp const double p_sq = X1GasMom * X1GasMom + X2GasMom * X2GasMom + X3GasMom * X3GasMom; const double Ekin = p_sq / (2.0 * density); const double Eint = Etot - Ekin; ``` and the inverse function is ``` cpp const double Etot = Eint + Ekin; ``` This is valid for pure hydro, but not for MHD, because magnetic energy is not included. For MHD, the gas internal energy should instead be recovered as ``` [ E{\rm int} = E_{\rm tot} - E_{\rm kin} - E_{\rm mag}, \qquad E_{\rm mag} = \frac{B^2}{8\pi}, ] and similarly the inverse relation should include magnetic energy: [ E_{\rm tot} = E_{\rm int} + E_{\rm kin} + E_{\rm mag}. ] ``` Because ResampledCooling.hpp uses ComputeEintFromEgas() to obtain the thermal energy before cooling integration, this means that in MHD runs the cooling module may be evolving an incorrect internal energy. This likely causes inconsistent thermodynamics in MHD + cooling runs. In my outputs, I observed: - magnetic field initialization looks correct - but thermal quantities become problematic after evolution - inferred gas pressure becomes negative - plasma beta becomes negative - temperature can collapse to extremely small values So the cooling-related concern is: ResampledCooling.hpp currently relies on RadSystem<problem_t>::ComputeEintFromEgas(), but that function appears to assume pure hydro and does not subtract magnetic energy. Therefore, the cooling update may be using the wrong gas internal energy in MHD runs.` This suggests that the current resampled cooling path is not MHD-safe unless internal energy recovery is modified to include magnetic energy.</div>
1 条评论