Free42 with IEEE 754-2008 decimal floating-point -- interested in a sneak preview?

[Thomas Okken]

(03-17-2014 01:29 PM)Thomas Okken Wrote:  

(03-17-2014 11:08 AM)J-F Garnier Wrote:  I noted that there is also a similar issue in the reverse situation:
1E5 LOG 5 - --> 1E-33 (also for 1E40 and other)
This didn't happen in versions 1.4.7x

OK, I'll add code to handle powers of 10 in LOG. (I'm amazed that code isn't there already, since testing for exact powers of 10 is so easy in decimal. Continuity issues maybe?)

In order to make sure LOG doesn't just return exact results for exact powers of 10, but also satisfies

  10ⁿ ≤ x < 10ⁿ⁺¹ → n ≤ LOG(x) < n+1

which is necessary for code like FWIW to work reliably, I'm thinking of this argument reduction code:

Code:

function log10(x) {

  // These cases aren't relevant in Free42, since it never calls log10 with x ≤ 0, but still.
  if (x < 0)
    return NaN;
  if (x == 0)
    return -Infinity;

  // Extract normalized decimal mantissa and exponent
  // (make sure 1 ≤ mantissa < 10)
  m = decimal mantissa of x;
  e = decimal exponent of x;
  while (m ≥ 10) {
    m = m / 10;
    e = e + 1;
  }
  while (m < 1) {
    m = m * 10;
    e = e - 1;
  }

  // Handle the trivial case
  if (m == 1)
    return e;

  // Handle the general case, while enforcing
  // 10^n ≤ x < 10^(n+1) → n ≤ LOG(x) < n+1
  r = e + bid128_log10(m);
  if (r < e)
    return e;
  else if (r ≥ e + 1)
    return highest representable number less than e + 1;
  else
    return r;
}