Builtin Float Epsilon? (was: Re: Does python suck or I am just stupid? )

[Dan Bishop]

Stephen Horne <intentionally at blank.co.uk> wrote in message news:<mkpf5v4ekln433ob1b0v464qf4nl0njj44 at 4ax.com>...
> On Sat, 22 Feb 2003 18:33:59 GMT, djw <dwelch91 at no.spam.attbi.com>
> wrote:
> 
> >sys.floatepsilon = 0.000000001 # or maybe some nice default value
> 
> Debatable - most people expect '==' to only return true if the values
> are precisely equal. I suspect an 'approximately equal' operator has
> been suggested in the past, though.
> 
> One problem is - just how close do two values need to be to count as
> equal. There is no single good answer to this - it depends on the
> context. If you are dealing with numbers that range from zero to
> several billion, you might consider 0 and 0.1 to be equal - if your
> numbers range from 0 to 0.00000001 you will probably take quite a
> different view.

Here's a class that solves this particular problem.  Every object of
this type stores its uncertainty, which is used in computations.  x ==
y returns true whenever (x-y).relativeUncertainty >= 1.

Also note that str (but not repr) limits the number of significant
digits according to the uncertainty.

The disadvantage of this approach is that all operations are more
computationally expensive.

------------------------------------------------------------------------

from __future__ import division

import math

_epsilon = 1.
while 1. + _epsilon != 1.:
   _epsilon /= 2

def _machineUncertainty(value):
   mantissa, exponent = math.frexp(value)
   return math.ldexp(_epsilon, exponent)

class Uncertain(object):
   def __init__(self, value, uncertainty=None):
      if uncertainty is None:
         uncertainty = _machineUncertainty(value)
      else:
         uncertainty = max(float(uncertainty),
_machineUncertainty(value))
      self.__v = float(value)
      self.__u = abs(uncertainty)
   value = property(lambda self: self.__v)
   uncertainty = property(lambda self: self.__u)
   relativeUncertainty = property(lambda self: abs(self.__u /
self.__v))
   def __float__(self):
      return self.__v
   def __repr__(self):
      return 'Uncertain(value=%s, uncertainty=%s)' % (`self.__v`,
`self.__u`)
   def __str__(self):
      return '%%.%dg' % self.significantDigits() % self.__v
   def __nonzero__(self):
      return self != 0
   def __pos__(self):
      return self
   def __neg__(self):
      return Uncertain(-self.__v, self.__u)
   def __abs__(self):
      return Uncertain(abs(self.__v), self.__u)
   def __add__(self, other):
      if not isinstance(other, Uncertain):
         other = Uncertain(other)
      return Uncertain(self.__v + other.__v, math.hypot(self.__u,
other.__u))
   __radd__ = __add__
   def __sub__(self, other):
      if not isinstance(other, Uncertain):
         other = Uncertain(other)
      return Uncertain(self.__v - other.__v, math.hypot(self.__u,
other.__u))
   def __rsub__(self, other):
      return Uncertain(other) - self
   def __mul__(self, other):
      if not isinstance(other, Uncertain):
         other = Uncertain(other)
      return Uncertain(self.__v * other.__v,
                       math.hypot(other.__v * self.__u, self.__v *
other.__u))
   __rmul__ = __mul__
   def __truediv__(self, other):
      if not isinstance(other, Uncertain):
         other = Uncertain(other)
      return Uncertain(self.__v / other.__v, math.hypot(self.__u /
other.__v,
                       other.__u * self.__v / other.__v ** 2))
   def __rtruediv__(self, other):
      return Uncertain(other) - self
   __div__ = __truediv__
   __rdiv__ = __rtruediv__
   def __pow__(self, other):
      if not isinstance(other, Uncertain):
         other = Uncertain(other)
      return Uncertain(self.__v ** other.__v, math.hypot(
                       self.__v ** (other.__v - 1) * other.__v *
self.__u,
                       self.__v ** other.__v * math.log(self.__v) *
other.__u))
   def __rpow__(self, other):
      return Uncertain(other) ** self
   def __cmp__(self, other):
      if not isinstance(other, Uncertain):
         other = Uncertain(other)
      diff = self - other
      if diff.__v <= diff.__u:
         return 0
      return cmp(self.__v, other.__v)
   def significantDigits(self, radix=10):
      return max(0, int(math.ceil(-math.log(self.relativeUncertainty)
/
                                            math.log(10))))