{-
   The Taylor series for the exponential function e^x:
    e^z = lim n|->inf \Sum_{n=1}^inf z^n / n!
   See: http://en.wikipedia.org/wiki/Power_series

   See the corresponding C+GMP version for comparison

   Compile:  ghc -O2 -o p_circ power_circ.hs
   Run:      ./p_circ 4000 10 
   Measure:  time ./p_circ 4000 10 1>/dev/null
-}

import System.Environment
import Control.Monad
import Data.Ratio

fact 0 = 1
fact n = n*(fact (n-1))

-- compute e^z up to d digits 
power_e :: Integer -> Integer -> Rational
power_e z d = sum (takeWhile (>(1 % (10^d))) taylor)
              where -- infinite list of entire Taylor series
            	    taylor =  [ (pow_ints!!n % 1) / (factorials!!n % 1) | n<-[0,1..] ]
	            -- 2 circular lists, to memoise earlier computations 
	            factorials = 1:1:[ (toInteger i)*(factorials!!(i-1)) | i<-[2..] ] 
                    pow_ints = 1:[ z*(pow_ints!!(i-1)) | i<-[1..] ]

main = do
         args <- getArgs
         let (z:d:_) = map read args :: [Integer]
         let res = power_e z d
         putStrLn $ "e^"++(show z)++" with a precision of "++(show d)++" digits is "++(show ((fromRational res)::Double))  
         putStrLn $ " or as Rational " ++(show res)
         putStrLn $ "Expected result: " ++(show (exp (fromInteger z)))