A simple equation of state solver for explosive astrophysics.
"main" is commented out. If you wish to compile as a stand-alone code, uncomment main and compile with:
gcc supercrude.c
It's that easy. But then comment main back out when you want to incorporate this into your code.
MIT License. Not sticky! Feel free to add to your code and even make your code proprietary if you like. But since we're talking about it, everyone should make their code public. Take out the swear words from the comments and post it already. Coding is a language and that language dies if it is not spoken. Bill Paxton taught me that and he tried to tell Xerox that back in the day but they didn't listen.
Code paper is on the way! In the meantime, I will tell you roughly what this solver is intended for:
SuperCrude is designed for incorporation into hydrodynamics codes. If used in hydrodynamics, i.e. only the closure relation P(E) and perhaps the sound speed are needed, I can guarantee SuperCrude will be as accurate as any other solver to within a few percent. It deviates up to that few percent ONLY when the temperature is around 10^9 K and the density is about 10^8 g/cm^3. This is the point where all forms of pressure are comparable (degeneracy, gas, radiation, pairs) and also the fluid is trans-relativistic.
Everywhere else, the errors are much lower, sometimes even smaller than HelmEoS. And it's just one page of C code that can very easily be incorporated into any hydro code.
Degeneracy and Radiation are computed exactly, for any density and temperature. Gas pressure from ions is also exact. Electron gas and pair plasma are handled "crudely" using smooth analytic fitting functions.
Thermodynamic Consistency is EXACT. These are analytic functions that I designed to be self-consistent.
Everything is infinitely differentiable and the code applies over all densities and temperatures from zero to infinity! Obviously I don't know what happens at 10^30 Kelvin; maybe I need to take into account dark matter pair production at that point. So obviously we don't do all physics in those regimes. But it will return a reasonable value for pressure/energy/entropy/etc for any positive temperature and density you give it.
Temperature root find is fast and you don't even have to provide a guess value for the temperature! How nice am I?
If you want to do anything with the temperature, you should know I don't promise better than 10% precision on the temperature found. Biggest errors are typically when degeneracy pressure is comparable to gas pressure. But generally hydro doesn't require the temperature. Unless you have some cooling prescription or something.
I provide many functions that you can use to get pressure, energy density, entropy, heat capacity, etc. But I'm sure there are always functions you might like to have that I never got around to coding. Please email me if there's something specific you'd like! I am very nice but also very busy. So if I don't reply right away don't take it the wrong way.
If you use the code, please cite me! I'll have a code paper out sometime soon!
Good Luck!