```{-|
Module      :  Data.Number.ER.Real.Approx.Sequence
Description :  exact reals via convergent sequences
Copyright   :  (c) Michal Konecny

Maintainer  :  mik@konecny.aow.cz
Stability   :  experimental
Portability :  portable

Types and methods related to explicit
convergent sequences of real number approximations.
-}
module Data.Number.ER.Real.Approx.Sequence
(
ConvergRealSeq,
makeFastConvergRealSeq,
convertFuncRA2Seq,
convertBinFuncRA2Seq,
convergRealSeqElem,
showConvergRealSeq,
showConvergRealSeqAuto
)
where

import qualified Data.Number.ER.Real.Approx as RA
import Data.Number.ER.BasicTypes

import Data.Maybe
import Data.Ratio

{-|
A converging sequence of real number approximations.

* Every finite subsequence has a non-empty intersection.

* The limit should be a singleton.
-}
data ConvergRealSeq ra =
ConvergRealSeq (EffortIndex -> ra)

convergRealSeqElem :: (ConvergRealSeq ra) -> EffortIndex -> ra
convergRealSeqElem (ConvergRealSeq sq) ix = sq ix

{-|
Using this operator, a unary funtion working over
approximations can be converted to one that works
over exact numbers represented through a sequence
of approximations.
-}
convertFuncRA2Seq ::
(EffortIndex -> ra -> ra) ->
(ConvergRealSeq ra) ->
(ConvergRealSeq ra)
convertFuncRA2Seq f (ConvergRealSeq argSeq) =
ConvergRealSeq resultSeq
where
resultSeq ix =
f ix (argSeq ix)

{-|
The same as above, where f is binary
-}
convertBinFuncRA2Seq ::
(EffortIndex -> ra -> ra -> ra) ->
(ConvergRealSeq ra) ->
(ConvergRealSeq ra) ->
(ConvergRealSeq ra)

convertBinFuncRA2Seq f (ConvergRealSeq arg1) (ConvergRealSeq arg2) =
ConvergRealSeq resultSeq
where
resultSeq ix =
f ix (arg1 ix) (arg2 ix)

{-|
Turn an arbitrary convergent sequence into one with
a guaranteed convergence rate - the precision (as defined
by 'RA.ERApprox.RA.getPrecision') of x_ix is at least ix.
-}
makeFastConvergRealSeq ::
(RA.ERApprox ra) =>
(ConvergRealSeq ra) ->
(ConvergRealSeq ra)
makeFastConvergRealSeq (ConvergRealSeq argSeq) =
ConvergRealSeq fastSeq
where
fastSeq ix =
head \$ catMaybes \$ map (precisionOK . argSeq) indexSeries
where
indexSeries =
--        take 5 \$ -- upper bound on iteration - for testing
binGeomSeries (max 1 ix)
precisionOK ra
| RA.getPrecision ra >= (effIx2prec ix) = Just ra
| otherwise = Nothing

{-|
binGeomSeries n is the geometric series
[ n, 2n, 4n, 8n, ...]
-}
binGeomSeries
:: (Num a)
=> a
-> [a]
binGeomSeries n =
n : (binGeomSeries (2 * n))

instance (RA.ERApprox ra) => Show (ConvergRealSeq ra)
where
show = showConvergRealSeq 6 True True 10 -- cheating here, should throw an error

{-|
Show function for ConvergRealSeq's with full arguments.
-}
showConvergRealSeq
:: (RA.ERApprox ra)
=> Int
-> Bool
-> Bool
-> Precision
-> (ConvergRealSeq ra)
-> String

showConvergRealSeq numDigits showGran showComponents prec r =
RA.showApprox numDigits showGran showComponents \$
convergRealSeqElem (makeFastConvergRealSeq r) (prec2effIx prec)

{-|
Show function for ConvergRealSeq's with all parameters fixed
except for number of digits
-}
showConvergRealSeqAuto
:: (RA.ERApprox ra)
=> Int
-> (ConvergRealSeq ra)
-> String
showConvergRealSeqAuto numDigits argSeq =
showConvergRealSeq numDigits True True prec argSeq
where
prec = effIx2prec \$ ceiling \$ (fromInteger \$ toInteger numDigits) * 3.3219280948873626

instance
(RA.ERApprox ra)
=> Eq (ConvergRealSeq ra)
where
r1 == r2 =
iterateRA_A raEq 2 [r1, r2]
where
raEq _ ([a1,a2]) = RA.equalReals a1 a2

instance
(RA.ERApprox ra)
=> Ord (ConvergRealSeq ra)
where
compare r1 r2 =
iterateRA_A eraComp 2 [r1, r2]
where
eraComp _ ([a1,a2]) = RA.compareReals a1 a2

pointwiseConvergRealSeq1 f (ConvergRealSeq sq) =
ConvergRealSeq (f . sq)
pointwiseConvergRealSeq2 f (ConvergRealSeq sq1) (ConvergRealSeq sq2) =
ConvergRealSeq (\ix -> f (sq1 ix) (sq2 ix))

instance
(RA.ERApprox ra)
=> Num (ConvergRealSeq ra)
where
fromInteger n = ConvergRealSeq sq
where
sq ix =
RA.setMinGranularity (effIx2gran ix) \$ fromInteger n
abs = pointwiseConvergRealSeq1 \$ abs
signum = pointwiseConvergRealSeq1 \$ signum
negate = pointwiseConvergRealSeq1 \$ negate
(+) = pointwiseConvergRealSeq2 \$ (+)
(*) = pointwiseConvergRealSeq2 \$ (*)

instance
(RA.ERApprox ra)
=> Fractional (ConvergRealSeq ra)
where
fromRational q = ConvergRealSeq sq
where
sq ix =
(RA.setMinGranularity (effIx2gran ix) num) / denom
num = fromInteger \$ numerator q
denom = fromInteger \$ denominator q
recip = pointwiseConvergRealSeq1 \$ recip

{-|
Take a converging sequence of partial functions F_i that operate on
real approximations and turn it into a function F that operates on converging sequences.
F looks for some members of the real approximation sequences
and an i so that F_i is defined for the chosen approximations
and returns its result.
-}
iterateRA_A
:: (EffortIndex -> [ra] -> Maybe a)
-- ^ a sequence of partial functions based on approximations
-> EffortIndex -- ^ a starting index to use when searching sequences
-> ([ConvergRealSeq ra] -> a)
-- ^ a total function based on sequences

iterateRA_A fn_RA startIx args =
head \$ catMaybes \$ map ((uncurry fn_RA) . args_Prec) indexSeries
where
indexSeries =
--        take 5 \$ -- upper bound on iteration - for testing
binGeomSeries \$ max 1 startIx
-- [(max 1 startIx)..]
args_Prec currentIndex =
(currentIndex, map (\ arg -> convergRealSeqElem arg currentIndex) args)

```