round half-integers to even instead of away from zero

This changes the behaviour to match r7rs (round x) instead of C round(x).

An answer to https://stackoverflow.com/questions/32746523/ieee-754-compliant-round-half-to-even
suggests using remainder(). The following will work if FE_TONEAREST is defined, but C11
requires FE_TONEAREST to be defined if and only if the implemenetation supports it in
fegetround() and fesetround() [Draft N1570]. On the other hand, remainder() must be defined.
C23 will have roundeven(), but this is not yet available on all platforms.
The behaviour of remainder is described in Draft N1570, page 254, footnote 239.

Alternative implementation:

  double round_to_nearest_even(double x)
  {
  #pragma STDC FENV_ACCESS ON
    int mode;
    double nearest;
    mode = fegetround();
    fesetround(FE_TONEAREST);
    nearest = nearbyint(x);
    fesetround(mode);
  #pragma STDC FENV_ACCESS OFF
    return nearest;
  }
This commit is contained in:
Yorick Hardy 2024-02-01 22:25:47 +02:00
parent 7d8f70fb07
commit 4bbceeb4d6
3 changed files with 8 additions and 2 deletions

View file

@ -504,6 +504,7 @@ int Cyc_have_mstreams();
} \
return_closcall1(data, cont, &d)
double round_to_nearest_even(double);
void Cyc_exact(void *data, object cont, object z);
object Cyc_exact_no_cps(void *data, object ptr, object z);

View file

@ -8765,6 +8765,11 @@ int num2ratio(double x, double *numerator, double *denominator)
return 0;
}
double round_to_nearest_even(double x)
{
return x-remainder(x,1.0);
}
/**
* Receive a Scheme number and pass requested portion of a rational number to
* the continuation `cont`. Pass numerator if `numerator` is true, else the

View file

@ -1372,9 +1372,9 @@
" return_double_op_no_cps(data, ptr, trunc, z);")
(define-c round
"(void *data, int argc, closure _, object k, object z)"
" return_double_op(data, k, round, z); "
" return_double_op(data, k, round_to_nearest_even, z); "
"(void *data, object ptr, object z)"
" return_double_op_no_cps(data, ptr, round, z);")
" return_double_op_no_cps(data, ptr, round_to_nearest_even, z);")
(define-c exact
"(void *data, int argc, closure _, object k, object z)"
" Cyc_exact(data, k, z); "