opensubscriber
   Find in this group all groups
 
The Haskell Cafe more information…

h : haskell-cafe@haskell.org 16 October 2006 • 1:14PM -0400

Re: [Haskell-cafe] Debugging Newton's method for square roots
by ajb

REPLY TO AUTHOR
 
REPLY TO GROUP



G'day all.

Quoting Tamas K Papp <tpapp@Prin...>:

> 2. Newton's method is not guaranteed to converge.

For computing square roots, it is.  The square root function is
exceedingly well-behaved.

But you can make things better by choosing a smart initial estimate.
The initial estimate in this version is so good that you only need a
fixed number of iterations, so there are no loops here.

mysqrt :: Double -> Double
mysqrt x
    | x <= 0    = if x == 0 then 0 else error "mysqrt domain error"
    | r == 1    = sqrt2 * sqrtx
    | otherwise = sqrtx
    where
        (q,r) = exponent x `divMod` 2

        sqrtx = scaleFloat q (sqrtSignificand (significand x))

        -- Feel free to add more digits if this is insufficient
        sqrt2 = 1.41421356237309504880168872420969807856967


-- Compute the square root of a number in the range [0.5,1)
sqrtSignificand :: Double -> Double
sqrtSignificand x
    = iter . iter . iter $ d0
    where
        -- d0 is a third-order Chebyshev approximation to sqrt x
        x' = x * 4.0 - 3.0
        d2 = -6.177921038278929e-3
        d1 = 2 * x' * d2 + 0.14589875298243973
        d0 = x' * d1 - d2 + 0.5 * 1.7196949654923188

        iter sx = 0.5 * (sx + x / sx)

Cheers,
Andrew Bromage
_______________________________________________
Haskell-Cafe mailing list
Haskell-Cafe@hask...
http://www.haskell.org/mailman/listinfo/haskell-cafe

Bookmark with:

Delicious   Digg   reddit   Facebook   StumbleUpon

Related Messages

opensubscriber is not affiliated with the authors of this message nor responsible for its content.