module Data.CDAR.Dyadic (Dyadic(..),shiftD,sqrtD,sqrtD',sqrtRecD,sqrtRecD',initSqrtRecD,initSqrtRecDoubleD,piMachinD,piBorweinD,ln2D,logD,atanhD,divD,divD') where
import Data.Ratio
import Data.Bits
import Data.CDAR.Classes
import Data.CDAR.IntegerLog
infix 8 :^
data Dyadic = Integer :^ Int deriving (ReadPrec [Dyadic]
ReadPrec Dyadic
Int -> ReadS Dyadic
ReadS [Dyadic]
(Int -> ReadS Dyadic)
-> ReadS [Dyadic]
-> ReadPrec Dyadic
-> ReadPrec [Dyadic]
-> Read Dyadic
forall a.
(Int -> ReadS a)
-> ReadS [a] -> ReadPrec a -> ReadPrec [a] -> Read a
readListPrec :: ReadPrec [Dyadic]
$creadListPrec :: ReadPrec [Dyadic]
readPrec :: ReadPrec Dyadic
$creadPrec :: ReadPrec Dyadic
readList :: ReadS [Dyadic]
$creadList :: ReadS [Dyadic]
readsPrec :: Int -> ReadS Dyadic
$creadsPrec :: Int -> ReadS Dyadic
Read,Int -> Dyadic -> ShowS
[Dyadic] -> ShowS
Dyadic -> String
(Int -> Dyadic -> ShowS)
-> (Dyadic -> String) -> ([Dyadic] -> ShowS) -> Show Dyadic
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [Dyadic] -> ShowS
$cshowList :: [Dyadic] -> ShowS
show :: Dyadic -> String
$cshow :: Dyadic -> String
showsPrec :: Int -> Dyadic -> ShowS
$cshowsPrec :: Int -> Dyadic -> ShowS
Show)
instance Eq Dyadic where
Integer
a :^ Int
s == :: Dyadic -> Dyadic -> Bool
== Integer
b :^ Int
t | Int
s Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
t = Integer
a Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer -> Int -> Integer
forall a. Bits a => a -> Int -> a
unsafeShiftL Integer
b (Int
tInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
s)
| Bool
otherwise = Integer -> Int -> Integer
forall a. Bits a => a -> Int -> a
unsafeShiftL Integer
a (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
t) Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
b
instance Ord Dyadic where
compare :: Dyadic -> Dyadic -> Ordering
compare (Integer
a :^ Int
s) (Integer
b :^ Int
t) | Int
s Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
t = Integer -> Integer -> Ordering
forall a. Ord a => a -> a -> Ordering
compare Integer
a (Integer -> Int -> Integer
forall a. Bits a => a -> Int -> a
unsafeShiftL Integer
b (Int
tInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
s))
| Bool
otherwise = Integer -> Integer -> Ordering
forall a. Ord a => a -> a -> Ordering
compare (Integer -> Int -> Integer
forall a. Bits a => a -> Int -> a
unsafeShiftL Integer
a (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
t)) Integer
b
instance Scalable Dyadic where
scale :: Dyadic -> Int -> Dyadic
scale (Integer
a :^ Int
s) Int
n = Integer
a Integer -> Int -> Dyadic
:^ (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
n)
instance Num Dyadic where
Integer
a :^ Int
s + :: Dyadic -> Dyadic -> Dyadic
+ Integer
b :^ Int
t
| Int
s Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
t = (Integer
a Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
+ Integer -> Int -> Integer
forall a. Scalable a => a -> Int -> a
scale Integer
b (Int
tInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
s)) Integer -> Int -> Dyadic
:^ Int
s
| Bool
otherwise = (Integer -> Int -> Integer
forall a. Scalable a => a -> Int -> a
scale Integer
a (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
t) Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
+ Integer
b) Integer -> Int -> Dyadic
:^ Int
t
Integer
a :^ Int
s * :: Dyadic -> Dyadic -> Dyadic
* Integer
b :^ Int
t = (Integer
a Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
* Integer
b) Integer -> Int -> Dyadic
:^ (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
t)
negate :: Dyadic -> Dyadic
negate (Integer
a :^ Int
s) = (Integer -> Integer
forall a. Num a => a -> a
negate Integer
a) Integer -> Int -> Dyadic
:^ Int
s
abs :: Dyadic -> Dyadic
abs (Integer
a :^ Int
s) = (Integer -> Integer
forall a. Num a => a -> a
abs Integer
a) Integer -> Int -> Dyadic
:^ Int
s
signum :: Dyadic -> Dyadic
signum (Integer
a :^ Int
_) = (Integer -> Integer
forall a. Num a => a -> a
signum Integer
a)Integer -> Int -> Dyadic
:^ Int
0
fromInteger :: Integer -> Dyadic
fromInteger Integer
i = Integer
i Integer -> Int -> Dyadic
:^ Int
0
instance Real Dyadic where
toRational :: Dyadic -> Rational
toRational (Integer
a :^ Int
s) = (Integer -> Rational
forall a. Real a => a -> Rational
toRational Integer
a)Rational -> Rational -> Rational
forall a. Num a => a -> a -> a
*Rational
2Rational -> Int -> Rational
forall a b. (Fractional a, Integral b) => a -> b -> a
^^Int
s
shiftD :: Int -> Dyadic -> Dyadic
shiftD :: Int -> Dyadic -> Dyadic
shiftD Int
t (Integer
m:^Int
s) =
if Int
t Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
s
then Integer -> Int -> Integer
forall a. Bits a => a -> Int -> a
unsafeShiftL Integer
m (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
t) Integer -> Int -> Dyadic
:^ Int
t
else Integer -> Int -> Integer
forall a. Bits a => a -> Int -> a
unsafeShiftR (Integer
m Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
+ Int -> Integer
forall a. Bits a => Int -> a
bit (Int
tInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1)) (Int
tInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
s) Integer -> Int -> Dyadic
:^ Int
t
divD :: Dyadic -> Dyadic -> Dyadic
divD :: Dyadic -> Dyadic -> Dyadic
divD Dyadic
a (Integer
n:^Int
t) = let (Integer
m:^Int
_) = Int -> Dyadic -> Dyadic
shiftD (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
t) Dyadic
a
in Rational -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round (Integer
m Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
% Integer
n) Integer -> Int -> Dyadic
:^ Int
t
sqrtD :: Int -> Dyadic -> Dyadic
sqrtD :: Int -> Dyadic -> Dyadic
sqrtD Int
t Dyadic
x = Int -> Dyadic -> Dyadic -> Dyadic
sqrtD' Int
t Dyadic
x (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ Dyadic -> Dyadic
initSqrtD Dyadic
x
where
initSqrtD :: Dyadic -> Dyadic
initSqrtD :: Dyadic -> Dyadic
initSqrtD (Integer
0:^Int
_) = (Integer
0Integer -> Int -> Dyadic
:^Int
0)
initSqrtD (Integer
m:^Int
s) = let i :: Int
i = Integer -> Int
integerLog2 Integer
m
n :: Integer
n = Integer -> Int -> Integer
forall a. Bits a => a -> Int -> a
shift Integer
m (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
i)
s' :: Int
s' = (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
i) Int -> Int -> Int
forall a. Integral a => a -> a -> a
`div` Int
2 Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
3
in if Int -> Bool
forall a. Integral a => a -> Bool
odd (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
i)
then (Integer
nInteger -> Integer -> Integer
forall a. Num a => a -> a -> a
+Integer
8)Integer -> Int -> Dyadic
:^Int
s'
else (Integer
nInteger -> Integer -> Integer
forall a. Num a => a -> a -> a
+Integer
4)Integer -> Int -> Dyadic
:^Int
s'
sqrtD' :: Int -> Dyadic -> Dyadic -> Dyadic
sqrtD' :: Int -> Dyadic -> Dyadic -> Dyadic
sqrtD' Int
t x :: Dyadic
x@(Integer
m:^Int
_) Dyadic
y
| Integer
m Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
0 = Integer
0Integer -> Int -> Dyadic
:^Int
0
| Integer
m Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
> Integer
0 = [Dyadic] -> Dyadic
converge ([Dyadic] -> Dyadic) -> (Dyadic -> [Dyadic]) -> Dyadic -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Dyadic -> Dyadic) -> Dyadic -> [Dyadic]
forall a. (a -> a) -> a -> [a]
iterate (Dyadic -> Int -> Dyadic -> Dyadic
newtonStep Dyadic
x Int
t) (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ Int -> Dyadic -> Dyadic
shiftD Int
t Dyadic
y
| Bool
otherwise = String -> Dyadic
forall a. HasCallStack => String -> a
error String
"Attempting sqrt of negative dyadic number."
where
newtonStep :: Dyadic -> Int -> Dyadic -> Dyadic
newtonStep :: Dyadic -> Int -> Dyadic -> Dyadic
newtonStep Dyadic
d Int
t' Dyadic
a = Int -> Dyadic -> Dyadic
shiftD Int
t' (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ (Integer
1Integer -> Int -> Dyadic
:^(-Int
1)) Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
* (Dyadic
a Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
+ Dyadic -> Dyadic -> Dyadic
divD Dyadic
d Dyadic
a)
sqrtRecD :: Int -> Dyadic -> Dyadic
sqrtRecD :: Int -> Dyadic -> Dyadic
sqrtRecD Int
t Dyadic
a = Int -> Dyadic -> Dyadic -> Dyadic
sqrtRecD' Int
t Dyadic
a (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ Dyadic -> Dyadic
initSqrtRecDoubleD Dyadic
a
initSqrtRecD :: Dyadic -> Dyadic
initSqrtRecD :: Dyadic -> Dyadic
initSqrtRecD (Integer
0 :^ Int
_) = String -> Dyadic
forall a. HasCallStack => String -> a
error String
"Divide by zero in initSqrtRecD"
initSqrtRecD (Integer
m :^ Int
s) =
let i :: Int
i = Integer -> Int
integerLog2 Integer
m
n :: Integer
n = Integer -> Int -> Integer
forall a. Bits a => a -> Int -> a
shift Integer
m (Int
3Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
i)
s' :: Int
s' = (-Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1) Int -> Int -> Int
forall a. Integral a => a -> a -> a
`div` Int
2 Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
5
in if Int -> Bool
forall a. Integral a => a -> Bool
even (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
i)
then ([Integer
62,Integer
59,Integer
56,Integer
53,Integer
51,Integer
49,Integer
48,Integer
46][Integer] -> Int -> Integer
forall a. [a] -> Int -> a
!!(Integer -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral Integer
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
8)) Integer -> Int -> Dyadic
:^ Int
s'
else ([Integer
44,Integer
42,Integer
40,Integer
38,Integer
36,Integer
35,Integer
34,Integer
33][Integer] -> Int -> Integer
forall a. [a] -> Int -> a
!!(Integer -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral Integer
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
8)) Integer -> Int -> Dyadic
:^ Int
s'
initSqrtRecDoubleD :: Dyadic -> Dyadic
initSqrtRecDoubleD :: Dyadic -> Dyadic
initSqrtRecDoubleD (Integer
0 :^ Int
_) = String -> Dyadic
forall a. HasCallStack => String -> a
error String
"Divide by zero in initSqrtRecDoubleD"
initSqrtRecDoubleD (Integer
m :^ Int
s) =
let i :: Int
i = Integer -> Int
integerLog2 Integer
m
n :: Integer
n = Integer -> Int -> Integer
forall a. Bits a => a -> Int -> a
shift Integer
m (Int
52Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
i)
s' :: Int
s' = (-Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
s) Int -> Int -> Int
forall a. Integral a => a -> a -> a
`div` Int
2 Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
52
(Integer
m',Int
_) = Double -> (Integer, Int)
forall a. RealFloat a => a -> (Integer, Int)
decodeFloat (Double -> (Integer, Int))
-> (Double -> Double) -> Double -> (Integer, Int)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Double
1Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/) (Double -> Double) -> (Double -> Double) -> Double -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Double -> Double
forall a. Floating a => a -> a
sqrt (Double -> (Integer, Int)) -> Double -> (Integer, Int)
forall a b. (a -> b) -> a -> b
$ Integer -> Int -> Double
forall a. RealFloat a => Integer -> Int -> a
encodeFloat Integer
n (if Int -> Bool
forall a. Integral a => a -> Bool
even (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
i) then (-Int
52) else (-Int
51))
t :: Int
t = if Integer
m' Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
/= Int -> Integer
forall a. Bits a => Int -> a
bit Int
52 Bool -> Bool -> Bool
&& Int -> Bool
forall a. Integral a => a -> Bool
even (Int
sInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
i) then Int
s'Int -> Int -> Int
forall a. Num a => a -> a -> a
-Int
1 else Int
s'
in Integer
m' Integer -> Int -> Dyadic
:^ Int
t
sqrtRecD' :: Int -> Dyadic -> Dyadic -> Dyadic
sqrtRecD' :: Int -> Dyadic -> Dyadic -> Dyadic
sqrtRecD' Int
t Dyadic
a Dyadic
x0 =
let step :: Dyadic -> Dyadic
step Dyadic
x = Int -> Dyadic -> Dyadic
shiftD Int
t (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ Dyadic
x Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
+ Int -> Dyadic -> Dyadic
shiftD Int
t (Dyadic
x Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
* (Dyadic
1 Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
- Int -> Dyadic -> Dyadic
shiftD Int
t (Dyadic
x Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
* Dyadic
x Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
* Dyadic
a)) Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
* (Integer
1 Integer -> Int -> Dyadic
:^ (-Int
1)))
xs :: [Dyadic]
xs = (Dyadic -> Dyadic) -> Dyadic -> [Dyadic]
forall a. (a -> a) -> a -> [a]
iterate Dyadic -> Dyadic
step Dyadic
x0
in [Dyadic] -> Dyadic
converge [Dyadic]
xs
converge :: [Dyadic] -> Dyadic
converge :: [Dyadic] -> Dyadic
converge ((Integer
0:^Int
s):[Dyadic]
_) = (Integer
0Integer -> Int -> Dyadic
:^Int
s)
converge ((Integer
m:^Int
_):ds :: [Dyadic]
ds@(d :: Dyadic
d@(Integer
n:^Int
_):[Dyadic]
_)) =
let a :: Int
a = Integer -> Int
integerLog2 Integer
m
b :: Int
b = Integer -> Int
integerLog2 (Integer -> Integer
forall a. Num a => a -> a
abs (Integer
m Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
- Integer
n))
in if Int
a Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
6
then if Integer -> Integer
forall a. Num a => a -> a
abs (Integer
mInteger -> Integer -> Integer
forall a. Num a => a -> a -> a
-Integer
n) Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
<= Integer
1 then Dyadic
d else [Dyadic] -> Dyadic
converge [Dyadic]
ds
else if Int
2 Int -> Int -> Int
forall a. Num a => a -> a -> a
* Int
b Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
a Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
3 then Dyadic
d else [Dyadic] -> Dyadic
converge [Dyadic]
ds
converge [Dyadic]
_ = String -> Dyadic
forall a. HasCallStack => String -> a
error String
"List terminating in converge."
divD' :: Int -> Dyadic -> Dyadic -> Dyadic
divD' :: Int -> Dyadic -> Dyadic -> Dyadic
divD' Int
p Dyadic
a Dyadic
b = let (Integer
m:^Int
_) = Int -> Dyadic -> Dyadic
shiftD (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
p) Dyadic
a
(Integer
n:^Int
_) = Int -> Dyadic -> Dyadic
shiftD Int
p Dyadic
b
in Rational -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round (Integer
mInteger -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
%Integer
n) Integer -> Int -> Dyadic
:^ Int
p
piMachinD :: Int -> Dyadic
piMachinD :: Int -> Dyadic
piMachinD Int
t = let t' :: Int
t' = Int
tInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
10Int -> Int -> Int
forall a. Num a => a -> a -> a
-Integer -> Int
integerLog2 (Int -> Integer
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Int
forall a. Num a => a -> a
abs Int
t))
a :: [Dyadic]
a = (Integer -> Dyadic) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map ((Integer -> Int -> Dyadic
:^Int
t') (Integer -> Dyadic) -> (Integer -> Integer) -> Integer -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round (Rational -> Integer)
-> (Integer -> Rational) -> Integer -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Int -> Integer
forall a. Bits a => Int -> a
bit (-Int
t') Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
%)) ([Integer] -> [Dyadic]) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> a -> b
$ (Integer -> Integer) -> Integer -> [Integer]
forall a. (a -> a) -> a -> [a]
iterate ((-Integer
25)Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
*) (Integer
5 :: Integer)
b :: [Dyadic]
b = (Integer -> Dyadic) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map ((Integer -> Int -> Dyadic
:^Int
t') (Integer -> Dyadic) -> (Integer -> Integer) -> Integer -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round (Rational -> Integer)
-> (Integer -> Rational) -> Integer -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Int -> Integer
forall a. Bits a => Int -> a
bit (-Int
t') Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
%)) ([Integer] -> [Dyadic]) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> a -> b
$ (Integer -> Integer) -> Integer -> [Integer]
forall a. (a -> a) -> a -> [a]
iterate ((-Integer
57121)Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
*) (Integer
239 :: Integer)
c :: [Dyadic]
c = (Integer -> Dyadic) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map ((Integer -> Int -> Dyadic
:^Int
t') (Integer -> Dyadic) -> (Integer -> Integer) -> Integer -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round (Rational -> Integer)
-> (Integer -> Rational) -> Integer -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Int -> Integer
forall a. Bits a => Int -> a
bit (-Int
t') Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
%)) ([Integer
1,Integer
3..] :: [Integer])
d :: [Dyadic]
d = (Dyadic -> Bool) -> [Dyadic] -> [Dyadic]
forall a. (a -> Bool) -> [a] -> [a]
takeWhile (Dyadic -> Dyadic -> Bool
forall a. Eq a => a -> a -> Bool
/= Dyadic
0) ([Dyadic] -> [Dyadic])
-> ([Dyadic] -> [Dyadic]) -> [Dyadic] -> [Dyadic]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Dyadic -> Dyadic) -> [Dyadic] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map (Int -> Dyadic -> Dyadic
shiftD Int
t') ([Dyadic] -> [Dyadic]) -> [Dyadic] -> [Dyadic]
forall a b. (a -> b) -> a -> b
$ (Dyadic -> Dyadic -> Dyadic) -> [Dyadic] -> [Dyadic] -> [Dyadic]
forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
(*) [Dyadic]
a [Dyadic]
c
e :: [Dyadic]
e = (Dyadic -> Bool) -> [Dyadic] -> [Dyadic]
forall a. (a -> Bool) -> [a] -> [a]
takeWhile (Dyadic -> Dyadic -> Bool
forall a. Eq a => a -> a -> Bool
/= Dyadic
0) ([Dyadic] -> [Dyadic])
-> ([Dyadic] -> [Dyadic]) -> [Dyadic] -> [Dyadic]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Dyadic -> Dyadic) -> [Dyadic] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map (Int -> Dyadic -> Dyadic
shiftD Int
t') ([Dyadic] -> [Dyadic]) -> [Dyadic] -> [Dyadic]
forall a b. (a -> b) -> a -> b
$ (Dyadic -> Dyadic -> Dyadic) -> [Dyadic] -> [Dyadic] -> [Dyadic]
forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
(*) [Dyadic]
b [Dyadic]
c
in Int -> Dyadic -> Dyadic
shiftD Int
t (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ Dyadic
4 Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
* (Dyadic
4 Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
* [Dyadic] -> Dyadic
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [Dyadic]
d Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
- [Dyadic] -> Dyadic
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [Dyadic]
e)
piBorweinD :: Int -> Dyadic
piBorweinD :: Int -> Dyadic
piBorweinD Int
t = let t' :: Int
t' = Int
tInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
10Int -> Int -> Int
forall a. Num a => a -> a -> a
-Integer -> Int
integerLog2 (Int -> Integer
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Int
forall a. Num a => a -> a
abs Int
t))
s :: Dyadic
s = Int -> Dyadic -> Dyadic
sqrtD Int
t' Dyadic
2
d :: (Dyadic, Dyadic, Int)
d = [(Dyadic, Dyadic, Int)] -> (Dyadic, Dyadic, Int)
forall a. [a] -> a
head ([(Dyadic, Dyadic, Int)] -> (Dyadic, Dyadic, Int))
-> ([(Dyadic, Dyadic, Int)] -> [(Dyadic, Dyadic, Int)])
-> [(Dyadic, Dyadic, Int)]
-> (Dyadic, Dyadic, Int)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ((Dyadic, Dyadic, Int) -> Bool)
-> [(Dyadic, Dyadic, Int)] -> [(Dyadic, Dyadic, Int)]
forall a. (a -> Bool) -> [a] -> [a]
dropWhile (\(Dyadic
_,Dyadic
y,Int
_) -> Dyadic
y Dyadic -> Dyadic -> Bool
forall a. Eq a => a -> a -> Bool
/= Dyadic
0) ([(Dyadic, Dyadic, Int)] -> (Dyadic, Dyadic, Int))
-> [(Dyadic, Dyadic, Int)] -> (Dyadic, Dyadic, Int)
forall a b. (a -> b) -> a -> b
$ ((Dyadic, Dyadic, Int) -> (Dyadic, Dyadic, Int))
-> (Dyadic, Dyadic, Int) -> [(Dyadic, Dyadic, Int)]
forall a. (a -> a) -> a -> [a]
iterate (Int -> (Dyadic, Dyadic, Int) -> (Dyadic, Dyadic, Int)
f Int
t') (Dyadic
6Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
-Dyadic
4Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
*Dyadic
s,Dyadic
sDyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
-Dyadic
1,Int
1 :: Int)
in Int -> Dyadic -> Dyadic
shiftD Int
t (Dyadic -> Dyadic)
-> ((Dyadic, Dyadic, Int) -> Dyadic)
-> (Dyadic, Dyadic, Int)
-> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (\(Dyadic
a,Dyadic
_,Int
_) -> Int -> Dyadic -> Dyadic -> Dyadic
divD' Int
t' Dyadic
1 Dyadic
a) ((Dyadic, Dyadic, Int) -> Dyadic)
-> (Dyadic, Dyadic, Int) -> Dyadic
forall a b. (a -> b) -> a -> b
$ (Dyadic, Dyadic, Int)
d
where
f :: Int -> (Dyadic, Dyadic, Int) -> (Dyadic, Dyadic, Int)
f :: Int -> (Dyadic, Dyadic, Int) -> (Dyadic, Dyadic, Int)
f Int
l (Dyadic
a,Dyadic
y,Int
k) = let u :: Dyadic
u = Int -> Dyadic -> Dyadic
sqrtD Int
l (Dyadic -> Dyadic) -> (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Dyadic -> Dyadic
sqrtD Int
l (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ Dyadic
1Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
-Dyadic
yDyadic -> Int -> Dyadic
forall a b. (Num a, Integral b) => a -> b -> a
^(Int
4::Int)
y' :: Dyadic
y' = Int -> Dyadic -> Dyadic
shiftD Int
l (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ Dyadic -> Dyadic -> Dyadic
divD (Dyadic
1Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
-Dyadic
u) (Dyadic
1Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
+Dyadic
u)
a' :: Dyadic
a' = Int -> Dyadic -> Dyadic
shiftD Int
l (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ Dyadic
aDyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
*(Dyadic
1Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
+Dyadic
y')Dyadic -> Int -> Dyadic
forall a b. (Num a, Integral b) => a -> b -> a
^(Int
4::Int)Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
-(Int -> Integer
forall a. Bits a => Int -> a
bit (Int
2Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
kInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)Integer -> Int -> Dyadic
:^Int
0)Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
*Dyadic
y'Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
*(Dyadic
1Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
+Dyadic
y'Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
+Dyadic
y'Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
*Dyadic
y')
in (Dyadic
a',Dyadic
y',Int
kInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)
ln2D :: Int -> Dyadic
ln2D :: Int -> Dyadic
ln2D Int
t = let t' :: Int
t' = Int
t Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
10 Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
2 Int -> Int -> Int
forall a. Num a => a -> a -> a
* Integer -> Int
integerLog2 (Int -> Integer
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> Int
forall a. Num a => a -> a
abs Int
t))
a :: [Dyadic]
a = (Integer -> Dyadic) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map ((Integer -> Int -> Dyadic
:^Int
t') (Integer -> Dyadic) -> (Integer -> Integer) -> Integer -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round (Rational -> Integer)
-> (Integer -> Rational) -> Integer -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Int -> Integer
forall a. Bits a => Int -> a
bit (-Int
t') Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
%)) ([Integer] -> [Dyadic]) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> a -> b
$ (Integer -> Integer) -> Integer -> [Integer]
forall a. (a -> a) -> a -> [a]
iterate (Integer
3Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
*) (Integer
3 :: Integer)
b :: [Dyadic]
b = (Integer -> Dyadic) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map ((Integer -> Int -> Dyadic
:^Int
t') (Integer -> Dyadic) -> (Integer -> Integer) -> Integer -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round (Rational -> Integer)
-> (Integer -> Rational) -> Integer -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Int -> Integer
forall a. Bits a => Int -> a
bit (-Int
t') Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
%)) ([Integer] -> [Dyadic]) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> a -> b
$ (Integer -> Integer) -> Integer -> [Integer]
forall a. (a -> a) -> a -> [a]
iterate (Integer
4Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
*) (Integer
4 :: Integer)
c :: [Dyadic]
c = (Integer -> Dyadic) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map ((Integer -> Int -> Dyadic
:^Int
t') (Integer -> Dyadic) -> (Integer -> Integer) -> Integer -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Rational -> Integer
forall a b. (RealFrac a, Integral b) => a -> b
round (Rational -> Integer)
-> (Integer -> Rational) -> Integer -> Integer
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Int -> Integer
forall a. Bits a => Int -> a
bit (-Int
t') Integer -> Integer -> Rational
forall a. Integral a => a -> a -> Ratio a
%)) ([Integer
1..] :: [Integer])
d :: [Dyadic]
d = (Dyadic -> Dyadic -> Dyadic) -> [Dyadic] -> [Dyadic] -> [Dyadic]
forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
(+) [Dyadic]
a [Dyadic]
b
e :: [Dyadic]
e = (Dyadic -> Bool) -> [Dyadic] -> [Dyadic]
forall a. (a -> Bool) -> [a] -> [a]
takeWhile (Dyadic -> Dyadic -> Bool
forall a. Eq a => a -> a -> Bool
/= Dyadic
0) ([Dyadic] -> [Dyadic])
-> ([Dyadic] -> [Dyadic]) -> [Dyadic] -> [Dyadic]
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Dyadic -> Dyadic) -> [Dyadic] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map (Int -> Dyadic -> Dyadic
shiftD Int
t') ([Dyadic] -> [Dyadic]) -> [Dyadic] -> [Dyadic]
forall a b. (a -> b) -> a -> b
$ (Dyadic -> Dyadic -> Dyadic) -> [Dyadic] -> [Dyadic] -> [Dyadic]
forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
(*) [Dyadic]
d [Dyadic]
c
in Int -> Dyadic -> Dyadic
shiftD Int
t (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ [Dyadic] -> Dyadic
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [Dyadic]
e
logD :: Int -> Dyadic -> Dyadic
logD :: Int -> Dyadic -> Dyadic
logD Int
t x :: Dyadic
x@(Integer
a :^ Int
s) =
if Integer
a Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
<= Integer
0 then String -> Dyadic
forall a. HasCallStack => String -> a
error String
"logD: Non-positive argument"
else let t' :: Int
t' = Int
t Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
10 Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int -> Int -> Int
forall a. Ord a => a -> a -> a
max Int
0 (Integer -> Int
integerLog2 Integer
a Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
s)
r :: Int
r = Int
s Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Integer -> Int
integerLog2 (Integer
3Integer -> Integer -> Integer
forall a. Num a => a -> a -> a
*Integer
a) Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1
x' :: Dyadic
x' = Dyadic -> Int -> Dyadic
forall a. Scalable a => a -> Int -> a
scale Dyadic
x (-Int
r)
y :: Dyadic
y = Int -> Dyadic -> Dyadic -> Dyadic
divD' Int
t' (Dyadic
x' Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
- Dyadic
1) (Dyadic
x' Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
+ Dyadic
1)
z :: Dyadic
z = Int -> Dyadic
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
r Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
* Int -> Dyadic
ln2D Int
t'
in Int -> Dyadic -> Dyadic
shiftD Int
t (Dyadic -> Dyadic) -> (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
+Dyadic
z) (Dyadic -> Dyadic) -> (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Dyadic -> Int -> Dyadic) -> Int -> Dyadic -> Dyadic
forall a b c. (a -> b -> c) -> b -> a -> c
flip Dyadic -> Int -> Dyadic
forall a. Scalable a => a -> Int -> a
scale Int
1 (Dyadic -> Dyadic) -> Dyadic -> Dyadic
forall a b. (a -> b) -> a -> b
$ Int -> Dyadic -> Dyadic
atanhD Int
t' Dyadic
y
atanhD :: Int -> Dyadic -> Dyadic
atanhD :: Int -> Dyadic -> Dyadic
atanhD Int
t x :: Dyadic
x@(Integer
a :^ Int
s) =
if Integer -> Int
integerLog2 (Integer -> Integer
forall a. Num a => a -> a
abs Integer
a) Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
s Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
>= Int
0 then String -> Dyadic
forall a. HasCallStack => String -> a
error String
"atanhD: Argument outside domain, (-1,1)"
else let
t' :: Int
t' = Int
t Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
5 Int -> Int -> Int
forall a. Num a => a -> a -> a
- Integer -> Int
integerLog2 (Integer -> Integer
forall a. Num a => a -> a
abs (Integer -> Integer) -> Integer -> Integer
forall a b. (a -> b) -> a -> b
$ Int -> Integer
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
t) Int -> Int -> Int
forall a. Integral a => a -> a -> a
`div` Int
2
g :: Dyadic -> Dyadic -> Dyadic
g Dyadic
_x Dyadic
_y = Int -> Dyadic -> Dyadic
shiftD Int
t' (Dyadic
_x Dyadic -> Dyadic -> Dyadic
forall a. Num a => a -> a -> a
* Dyadic
_y)
x2 :: Dyadic
x2 = Dyadic -> Dyadic -> Dyadic
g Dyadic
x Dyadic
x
b :: [Dyadic]
b = (Dyadic -> Dyadic) -> Dyadic -> [Dyadic]
forall a. (a -> a) -> a -> [a]
iterate (Dyadic -> Dyadic -> Dyadic
g Dyadic
x2) Dyadic
1
c :: [Dyadic]
c = ((Integer -> Dyadic) -> [Integer] -> [Dyadic]
forall a b. (a -> b) -> [a] -> [b]
map (Int -> Dyadic -> Dyadic -> Dyadic
divD' Int
t' Dyadic
1 (Dyadic -> Dyadic) -> (Integer -> Dyadic) -> Integer -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Integer -> Dyadic
forall a b. (Integral a, Num b) => a -> b
fromIntegral) [Integer
1,Integer
3..])
in Int -> Dyadic -> Dyadic
shiftD Int
t (Dyadic -> Dyadic) -> ([Dyadic] -> Dyadic) -> [Dyadic] -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Dyadic -> Dyadic -> Dyadic
g Dyadic
x (Dyadic -> Dyadic) -> ([Dyadic] -> Dyadic) -> [Dyadic] -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Dyadic] -> Dyadic
forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum ([Dyadic] -> Dyadic)
-> ([Dyadic] -> [Dyadic]) -> [Dyadic] -> Dyadic
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Dyadic -> Bool) -> [Dyadic] -> [Dyadic]
forall a. (a -> Bool) -> [a] -> [a]
takeWhile (Dyadic -> Dyadic -> Bool
forall a. Eq a => a -> a -> Bool
/= Dyadic
0) ([Dyadic] -> Dyadic) -> [Dyadic] -> Dyadic
forall a b. (a -> b) -> a -> b
$ (Dyadic -> Dyadic -> Dyadic) -> [Dyadic] -> [Dyadic] -> [Dyadic]
forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith Dyadic -> Dyadic -> Dyadic
g [Dyadic]
b [Dyadic]
c