Fixed the accuracy of log1pf(float x) when `x` is small.

[ZERICO2005]

The accuracy of log1pf(float x) has been fixed. It will return x when x is small for 12.0bits of precision, or it can return x - 0.5f * (x * x) when x is small for 16.0bits of precision. It falls-back to the naive return logf(x + 1.0f) when x is large.

[parisseb]

This fix (https://github.com/CE-Programming/toolchain/commit/862e2608f03bad92a43761cd9ca098a61364f921) should be improved, because it generates a relative error that is smaller than 2^(-13), and that's far from 2^-24. It could be improved like this

    if (fabsf(x) < 2.44140625e-4f /* 0x1.0p-12f */) {
        return x*(1-ld_exp(x,-1));
    }

but for values just above 2^-12, using log(1+x) will still lead to low relative precision (about 2^-12), because computing 1+x will loose the last 12 bits of x.

Here is an improved log1p version, using Pade approximant, 4 multiplications and one division, relative error is bounded by about 2^-21=5e-7

float log1pf(float x){
  if (fabs(x)<=0.125){
    // pade(series(ln(1+x),x=0,6,polynom),x,5,3)
    // (-57*x**2-90*x)/(x**3-21*x**2-102*x-90)
    // relative error less than 1e-7
       return x*(57*x+90)/(((21-x)*x+102)*x+90);
  }
  // relative error about 2^-21 if abs(x) is just above 0.125
  return logf(1 + x);
}

FYI, I have also re-implemented exp/ln/sin/cos/tan/atan in order to speed up a little bit plotting inside KhiCAS (between 10% and 20%), because the toolchain source (from Zilog) seems to be adapted for 14 or 15 digits relative precision, about twice the 24 bits relative precision of floats. The new versions are available in the archive containing all changes I had to make on the toolchain for optimal compilation of KhiCAS: https://www-fourier.univ-grenoble-alpes.fr/~parisse/ti/celibc.tgz (files should be copied in the toolchain src/libc directory, beware that it also contains my modified version of malloc to access 2 heaps, one of them being the upper LCD part).

[adriweb]

Indeed, very nice! I'll commit that version, then.

[adriweb]

I'm closing this PR since it's now obsolete.

[ZERICO2005]

because the toolchain source (from Zilog) seems to be adapted for 14 or 15 digits relative precision

Yeah I noticed that too. We could keep those routines around should we decide on implementing Float64 (long double)