ionisation schemes + equation of state

[danieljprice]

currently we have a variety of different modules designed to capture ionisation and recombination, including the following that assume ionisation/recombination equilibrium:

1. the "gasradrec" equation of state

   • assumes ionisation equilibrium
   • works with any composition
   • does H and He ionisation but does not currently do H2 formation/dissociation
   • adapts mu so that u is conserved but does not update gamma (since it would only change with H2 formation)

2. the tabulated MESA equation of state

   • assumes ionisation equilibrium
   • adapts mu/gamma so that u is conserved
   • requires fixed composition (X,Y,Z)
   • probably does not include H2

3. the non-tabulated Helmholtz equation of state

   • assumes ionisation equilibrium
   • the reverse solver for T from u is not implemented
   • composition could be local in principle
   • assumes full ionisation so valid only for high temperatures

4. NICIL library for the non-ideal MHD coefficients

   • includes cosmic ray and thermal ionisation of H and He (and maybe also H2 formation/destruction?)
   • assumes chemical equilibrium
   • currently not used to change the equation of state, though could be in principle?

5. ieos=5 which

   • does H2 formation/destruction but not ionisation
   • adapts mu/gamma so that u is conserved

6. the H2 chemistry module (h2chem.f90)

   • does not assume ionisation or chemical equilibrium
   • tracks formation/destruction of H2 and ionisation of H but not He
   • updates mu and could be used to update gamma in principle

7. cv_table for the flux limited diffusion radiation hydro

   • assumes fixed composition (X,Y,Z)
   • includes H/He ionisation and also H2 formation/destruction
   • tabulated eos, but we produce the table so can be extended in principle

The working plan is to try to consolidate these efforts into a few "master" prescriptions.

Suggested steps are:

• extend gasradrec to deal with H2 formation/destruction and variable gamma
• get rid of ieos=5 and merge into gasradrec
• keep MESA as reference but implement unit tests that check that results are the same in the overlap region
• hopefully this means we can remove the cv_table (item 8), this needs to be checked for speed in the implicit radiation solver
• output and store n_e, mu and gamma as part of the equation of state calculation, for use in MCFOST atomic line transfer
• use the inferred abundances of H/H2/He/dust etc. in the cooling rates (e.g. Bate & Keto scheme and cooling_ism module)