RCR 340: Approximate comparison of floats

Add a new operator to the Float class that supports approximate comparison of floating point values. This operator will determine if two floating point values are equal within some small epsilon.

Floating point arithmatic is not exact. This is best illustrated with the following line of ruby code:

  (0.1 + 0.2) == 0.3
  => false

The comparison is obviously true, but round-off error in the processor causes the addition operation to produce an inexact number.

The solution is to add two new methods to the Float class that perform approximate comparisons in much the same way that Float#== and Float#eql? perform exact comparisons.

The first new method would be Float#=~ It would perform approximate comparison with another object, and it would follow the same semantics of Float#== i.e. perform type conversion. The Float::EPSILON constant would be used internally by this method.

The second new method would be Float#epsilon_eql? It would perform approximate comparison with another Float object, and it would allow the user to specify the epsilon. This method would follow the same semantics of Float#eql?

Syntax:

  1.0 =~ 1
  => true
  1.0000000000000001 =~ 1.0
  => true
  1.000000000000001 =~ 1.0
  => false
  (0.1 + 0.2) =~ 0.3
  => true
  1.0.epsilon_eql? 1       # using Float::EPSILON
  => false
  1.0.epsilon_eql? 1.0     # using Float::EPSILON
  => true
  1.000001.epsilon_eql?(1.0,  0.000001)
  => true
  1.000001.epsilon_eql?(1.0,  0.0000001)
  => false

This addition would not break any existing functionality and it is fairly trivial to implement. The =~ operator in mathematics is defined to mean "approximately equal to", and so it fits in nicely with Float.

  /* in numeric.c */
  static VALUE
  flo_epsilon_eq(x, y)
      VALUE x, y;
  {
      volatile double a, b;
      switch (TYPE(y)) {
        case T_FIXNUM:
          b = FIX2LONG(y);
          break;
        case T_BIGNUM:
          b = rb_big2dbl(y);
          break;
        case T_FLOAT:
          b = RFLOAT(y)->value;
          if (isnan(b)) return Qfalse;
          break;
        default:
          /* defaults to normal equality check */
          return num_equal(x, y);
      }
      a = RFLOAT(x)->value;
      if (isnan(a)) return Qfalse;
      if (b == 0.0) return (a <= DBL_EPSILON)?Qtrue:Qfalse;
      return (fabs((b-a)/b) <= DBL_EPSILON)?Qtrue:Qfalse;                         
  }
  static VALUE
  flo_epsilon_eql(argc, argv, x)
      int argc;
      VALUE* argv;
      VALUE x;
  {
      VALUE y, e;
      rb_scan_args( argc, argv, "11", &y, &e );
      if (TYPE(y) == TFLOAT) {
          double epsilon = DBL_EPSILON;
          if (!NIL_P(e)) {
              if (TYPE(e) != TFLOAT) {
                  rb_raise(rb_eTypeError, "epsilon must be a Float");
              }
              epsilon = RFLOAT(e)->value;
              if (isnan(epsilon)) return Qfalse;
              epsilon = fabs(epsilon);
          }
          double a = RFLOAT(x)->value;
          double b = RFLOAT(y)->value;
          if (isnan(a) || isnan(b)) return Qfalse;
          if (b == 0.0) return (a <= epsilon)?Qtrue:Qfalse;
          return (fabs((b-a)/b) <= epsilon)?Qtrue:Qfalse;
      }
      return Qfalse;
  }
  rb_define_method(rb_cFloat, "=~", flo_epsilon_eq, 1);
  rb_define_method(rb_cFloat, "epsilon_eql?", flo_epsilon_eql, -1);