Need a bit-for-bit safe error function in math

`erf` (e.g. https://gcc.gnu.org/onlinedocs/gfortran/ERF.html) does not produce bit-for-bit identical results across compilers, so is not currently appropriate for general usage in MESA in the occasional circumstance where one might want to use it as part of a fitting formula.

See discussion in #399, where an error function approximation is used in the fitting formula for the diffusion coefficients. That approximation is fine for the specific diffusion coefficient application there, but is not accurate to floating point precision, so is not appropriate as a general solution to use in `math`.

Start by reading the discussion in issue #399 and the existing math implementation where the diffusion-coefficient approximation is used. Define an approach for a general error function that produces bit-for-bit identical results across compilers, then verify that it is accurate to floating-point precision and suitable for reuse in MESA fitting formulas.

Written by the indexing model from the issue text.