1

I have a question regarding numerical operations in Haskell.

I have a basic function:

derv::(Num a, Fractional a) => (a -> a) -> a -> a -> a
derv f a deltax = ((f (a+deltax))-(f a))/deltax

When I test it, this is the output I get: 

    *Main> derv (\x->x*x) 2 0.000000000000000001
    0.0
    *Main> derv (\x->x*x) 2 0.00000000000001
    4.085620730620576
    *Main> derv (\x->x*x) 2 0.0000000001
    4.000000330961484
    *Main> derv (\x->x*x) 2 0.0001
    4.0001000000078335
    *Main>

When the divisor gets smaller, it makes the answer automatically go towards zero, rather than a more refined convergence toward 4. I'm curious as to why this is happening, especially given my type definition.

Mike S
  • 41,895
  • 11
  • 89
  • 84
Eigenvalue
  • 1,093
  • 1
  • 14
  • 35
  • 7
    This has nothing whatsoever to do with haskell. Read [what every computer scientist should know about floating point](http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html). (Note that you don't get 0 in the first test if you add `:: Rational`, because ghci uses exact rationals rather than floating point for the calculations then.) — BTW, using a fixed `δx` is a pretty horrible way of differentiating – very often, [automatic differentiation](http://hackage.haskell.org/package/ad) is a much better choice. – leftaroundabout Aug 31 '14 at 20:32
  • `2 + epsilon == 2` in floating point arithmetic, for sufficiently small `epsilon`. Rounding does matter. – chi Aug 31 '14 at 20:33
  • 1
    Side note: a central difference quotient `(f (a+deltax) - f (a-deltax)) / (2*deltax)` converges much faster. – luqui Aug 31 '14 at 22:01

2 Answers2

3

In your code the 0.000000000000000001 is probably being defaulted to Double, causing a loss of precision after 2 is added because of rounding.

Using an exact representation such as Rational does not exhibit the same issue:

> import Data.Ratio
> derv (\x->x*x) 2 0.000000000000000001 :: Rational
4000000000000000001 % 1000000000000000000
> fromRational (derv (\x->x*x) 2 0.000000000000000001) :: Double
4.0

In the last line the loss of precision happens after the incremental ratio is computed, so the result is close to the exact fraction shown above.

chi
  • 111,837
  • 3
  • 133
  • 218
0

This is probably due to floating point rounding errors. You can use fixed point numbers, as is shown here.

Community
  • 1
  • 1
zmbq
  • 38,013
  • 14
  • 101
  • 171