[SOLVED] Identical pure expressions evaluate to two different outputs

I was reading through Iago Leal’s post on lazily solving differential equations, so I have this Stream a type which I think represents the coefficients of the Taylor series of f(x) about x=a.

{-# LANGUAGE DeriveTraversable, NoMonomorphismRestriction #-}
import Data.Foldable (toList)

data Stream a = a :> Stream a
  deriving (Functor, Foldable)

infixr 2 :>

ex = 1.0 :> ex

where ex represents the coefficients of the series of f(x)=e^x.
Then I define a helper function to return the overall coefficients of the Taylor series by dividing the first term by 0!, the second by 1!, the third by 2! etc.

toTaylor :: Fractional a => Stream a -> [a]
toTaylor f = zipWith (/) (toList f) (fromIntegral <$> facs)
  where
    facs :: [Int]
    facs = (:) 1 $ scanl1 (*) [1..]

My output from running mapM_ print $ take 25 $ zipWith (/) (toList ex) (fromIntegral <$> facs) in ghci-9.2.8 is:

1.0
1.0
0.5
0.16666666666666666
4.1666666666666664e-2
8.333333333333333e-3
1.388888888888889e-3
1.984126984126984e-4
2.48015873015873e-5
2.7557319223985893e-6
2.755731922398589e-7
2.505210838544172e-8
2.08767569878681e-9
1.6059043836821613e-10
1.1470745597729725e-11
7.647163731819816e-13
4.779477332387385e-14
2.8114572543455206e-15
1.5619206968586225e-16
8.22063524662433e-18
4.110317623312165e-19
1.9572941063391263e-20
8.896791392450574e-22
3.8681701706306835e-23
1.6117375710961184e-24

However my output from running mapM_ print $ take 25 $ toTaylor ex is:

1.0
1.0
0.5
0.16666666666666666
4.1666666666666664e-2
8.333333333333333e-3
1.388888888888889e-3
1.984126984126984e-4
2.48015873015873e-5
2.7557319223985893e-6
2.755731922398589e-7
2.505210838544172e-8
2.08767569878681e-9
1.6059043836821613e-10
1.1470745597729725e-11
7.647163731819816e-13
4.779477332387385e-14
2.8114572543455206e-15
1.5619206968586225e-16
8.22063524662433e-18
4.110317623312165e-19
-2.353334294364486e-19
-7.995773634431128e-19
1.2302708207447328e-19
-1.2762938906401676e-19

which is completely different as the first list converges to 0 whereas the above oscillates between positive and negative numbers, despite the fact that toTaylor ex = zipWith (/) (toList ex) (fromIntegral <$> facs) so they should be the same list. What’s more, both lists give the same output up to the number 4.110317623312165e-19 beyond which they diverge from each other.

I can’t think of anything that can explain this behaviour. Any help would be appreciated.

I’m 99% sure that in mapM_ print $ take 25 $ zipWith (/) (toList ex) (fromIntegral <$> facs) you either have given the type signature facs :: [Integer], or haven’t explicitly given facs one, in which case it’ll default to [Integer].

Integer is arbitrary precision, but Int is not and will overflow to a negative compare. Try running take 25 $ 1:(scanl1 (*) [1::Integer ..]) and take 25 $ 1:(scanl1 (*) [1::Int ..]) in GHCi and comparing the results.

3 Likes

I swapped [Int] for [Integer] in facs and it works, I get the first list. Thanks for the help and I’ll bear in mind that Int is fixed precision so it overflows at high precision.

I hadn’t seen that post, neat!

I’d like to add that there’s some related (subsequent) work I did ages ago that makes such (approximate numerical) solvers more practical by moving from a single taylor series to a lazy list of them and providing compositional “solver” strategies: http://gbaz.github.io/slides/ode-draft-2009.pdf

3 Likes