{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE ExistentialQuantification #-}
module Numeric.Series (
Sequence(..)
, enumSequenceFrom
, enumSequenceFromStep
, scanSequence
, sumSeries
, sumPowerSeries
, sequenceToList
, evalContFractionB
) where
import Control.Applicative
import Data.List (unfoldr)
import Numeric.MathFunctions.Constants (m_epsilon)
data Sequence a = forall s. Sequence s (s -> (a,s))
instance Functor Sequence where
fmap :: forall a b. (a -> b) -> Sequence a -> Sequence b
fmap a -> b
f (Sequence s
s0 s -> (a, s)
step) = forall a s. s -> (s -> (a, s)) -> Sequence a
Sequence s
s0 (\s
s -> let (a
a,s
s') = s -> (a, s)
step s
s in (a -> b
f a
a, s
s'))
{-# INLINE fmap #-}
instance Applicative Sequence where
pure :: forall a. a -> Sequence a
pure a
a = forall a s. s -> (s -> (a, s)) -> Sequence a
Sequence () (\() -> (a
a,()))
Sequence s
sA s -> (a -> b, s)
fA <*> :: forall a b. Sequence (a -> b) -> Sequence a -> Sequence b
<*> Sequence s
sB s -> (a, s)
fB = forall a s. s -> (s -> (a, s)) -> Sequence a
Sequence (s
sA,s
sB) forall a b. (a -> b) -> a -> b
$ \(!s
sa,!s
sb) ->
let (a -> b
a,s
sa') = s -> (a -> b, s)
fA s
sa
(a
b,s
sb') = s -> (a, s)
fB s
sb
in (a -> b
a a
b, (s
sa',s
sb'))
{-# INLINE pure #-}
{-# INLINE (<*>) #-}
instance Num a => Num (Sequence a) where
+ :: Sequence a -> Sequence a -> Sequence a
(+) = forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2 forall a. Num a => a -> a -> a
(+)
* :: Sequence a -> Sequence a -> Sequence a
(*) = forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2 forall a. Num a => a -> a -> a
(*)
(-) = forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2 (-)
{-# INLINE (+) #-}
{-# INLINE (*) #-}
{-# INLINE (-) #-}
abs :: Sequence a -> Sequence a
abs = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap forall a. Num a => a -> a
abs
signum :: Sequence a -> Sequence a
signum = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap forall a. Num a => a -> a
signum
fromInteger :: Integer -> Sequence a
fromInteger = forall (f :: * -> *) a. Applicative f => a -> f a
pure forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. Num a => Integer -> a
fromInteger
{-# INLINE abs #-}
{-# INLINE signum #-}
{-# INLINE fromInteger #-}
instance Fractional a => Fractional (Sequence a) where
/ :: Sequence a -> Sequence a -> Sequence a
(/) = forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2 forall a. Fractional a => a -> a -> a
(/)
recip :: Sequence a -> Sequence a
recip = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap forall a. Fractional a => a -> a
recip
fromRational :: Rational -> Sequence a
fromRational = forall (f :: * -> *) a. Applicative f => a -> f a
pure forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. Fractional a => Rational -> a
fromRational
{-# INLINE (/) #-}
{-# INLINE recip #-}
{-# INLINE fromRational #-}
enumSequenceFrom :: Num a => a -> Sequence a
enumSequenceFrom :: forall a. Num a => a -> Sequence a
enumSequenceFrom a
i = forall a s. s -> (s -> (a, s)) -> Sequence a
Sequence a
i (\a
n -> (a
n,a
nforall a. Num a => a -> a -> a
+a
1))
{-# INLINE enumSequenceFrom #-}
enumSequenceFromStep :: Num a => a -> a -> Sequence a
enumSequenceFromStep :: forall a. Num a => a -> a -> Sequence a
enumSequenceFromStep a
n a
d = forall a s. s -> (s -> (a, s)) -> Sequence a
Sequence a
n (\a
i -> (a
i,a
iforall a. Num a => a -> a -> a
+a
d))
{-# INLINE enumSequenceFromStep #-}
scanSequence :: (b -> a -> b) -> b -> Sequence a -> Sequence b
{-# INLINE scanSequence #-}
scanSequence :: forall b a. (b -> a -> b) -> b -> Sequence a -> Sequence b
scanSequence b -> a -> b
f b
b0 (Sequence s
s0 s -> (a, s)
step) = forall a s. s -> (s -> (a, s)) -> Sequence a
Sequence (b
b0,s
s0) forall a b. (a -> b) -> a -> b
$ \(b
b,s
s) ->
let (a
a,s
s') = s -> (a, s)
step s
s
b' :: b
b' = b -> a -> b
f b
b a
a
in (b
b,(b
b',s
s'))
sumSeries :: Sequence Double -> Double
{-# INLINE sumSeries #-}
sumSeries :: Sequence Double -> Double
sumSeries (Sequence s
sInit s -> (Double, s)
step)
= Double -> s -> Double
go Double
x0 s
s0
where
(Double
x0,s
s0) = s -> (Double, s)
step s
sInit
go :: Double -> s -> Double
go Double
x s
s | forall a. Num a => a -> a
abs (Double
dforall a. Fractional a => a -> a -> a
/Double
x) forall a. Ord a => a -> a -> Bool
>= Double
m_epsilon = Double -> s -> Double
go Double
x' s
s'
| Bool
otherwise = Double
x'
where (Double
d,s
s') = s -> (Double, s)
step s
s
x' :: Double
x' = Double
x forall a. Num a => a -> a -> a
+ Double
d
sumPowerSeries :: Double -> Sequence Double -> Double
sumPowerSeries :: Double -> Sequence Double -> Double
sumPowerSeries Double
x Sequence Double
ser = Sequence Double -> Double
sumSeries forall a b. (a -> b) -> a -> b
$ forall (f :: * -> *) a b c.
Applicative f =>
(a -> b -> c) -> f a -> f b -> f c
liftA2 forall a. Num a => a -> a -> a
(*) (forall b a. (b -> a -> b) -> b -> Sequence a -> Sequence b
scanSequence forall a. Num a => a -> a -> a
(*) Double
1 (forall (f :: * -> *) a. Applicative f => a -> f a
pure Double
x)) Sequence Double
ser
{-# INLINE sumPowerSeries #-}
sequenceToList :: Sequence a -> [a]
sequenceToList :: forall a. Sequence a -> [a]
sequenceToList (Sequence s
s s -> (a, s)
f) = forall b a. (b -> Maybe (a, b)) -> b -> [a]
unfoldr (forall a. a -> Maybe a
Just forall b c a. (b -> c) -> (a -> b) -> a -> c
. s -> (a, s)
f) s
s
evalContFractionB :: Sequence (Double,Double) -> Double
{-# INLINE evalContFractionB #-}
evalContFractionB :: Sequence (Double, Double) -> Double
evalContFractionB (Sequence s
sInit s -> ((Double, Double), s)
step)
= let ((Double
_,Double
b0),s
s0) = s -> ((Double, Double), s)
step s
sInit
f0 :: Double
f0 = Double -> Double
maskZero Double
b0
in Double -> Double -> Double -> s -> Double
go Double
f0 Double
f0 Double
0 s
s0
where
tiny :: Double
tiny = Double
1e-60
maskZero :: Double -> Double
maskZero Double
0 = Double
tiny
maskZero Double
x = Double
x
go :: Double -> Double -> Double -> s -> Double
go Double
f Double
c Double
d s
s
| forall a. Num a => a -> a
abs (Double
delta forall a. Num a => a -> a -> a
- Double
1) forall a. Ord a => a -> a -> Bool
>= Double
m_epsilon = Double -> Double -> Double -> s -> Double
go Double
f' Double
c' Double
d' s
s'
| Bool
otherwise = Double
f'
where
((Double
a,Double
b),s
s') = s -> ((Double, Double), s)
step s
s
d' :: Double
d' = forall a. Fractional a => a -> a
recip forall a b. (a -> b) -> a -> b
$ Double -> Double
maskZero forall a b. (a -> b) -> a -> b
$ Double
b forall a. Num a => a -> a -> a
+ Double
aforall a. Num a => a -> a -> a
*Double
d
c' :: Double
c' = Double -> Double
maskZero forall a b. (a -> b) -> a -> b
$ Double
b forall a. Num a => a -> a -> a
+ Double
aforall a. Fractional a => a -> a -> a
/Double
c
delta :: Double
delta = Double
c'forall a. Num a => a -> a -> a
*Double
d'
f' :: Double
f' = Double
fforall a. Num a => a -> a -> a
*Double
delta