```{-
Copyright (C) 2011 Dr. Alistair Ward

This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.

This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.

You should have received a copy of the GNU General Public License
along with this program.  If not, see <http://www.gnu.org/licenses/>.
-}
{- |
[@AUTHOR@]	Dr. Alistair Ward

[@DESCRIPTION@]	Defines the unit with which precision is measured, and operations on it.
-}
module Factory.Math.Precision(
-- * Types
-- ** Type-synonyms
ConvergenceOrder,
ConvergenceRate,
DecimalDigits,
-- * Constants
linearConvergence,
cubicConvergence,
quarticConvergence,
-- * Functions
getIterationsRequired,
getTermsRequired,
roundTo,
promote,
simplify
) where

import qualified	Data.Ratio

-- | The /order of convergence/; <http://en.wikipedia.org/wiki/Rate_of_convergence>.
type ConvergenceOrder	= Int

-- | The /rate of convergence/; <http://en.wikipedia.org/wiki/Rate_of_convergence>.
type ConvergenceRate	= Double

-- | A number of decimal digits; presumably positive.
type DecimalDigits	= Int

-- | /Linear/ convergence-rate; which may be qualified by the /rate of convergence/.
linearConvergence :: ConvergenceOrder
linearConvergence	= 1

-- | /Quadratic/ convergence-rate.

-- | /Cubic/ convergence-rate.
cubicConvergence :: ConvergenceOrder
cubicConvergence	= 3

-- | /Quartic/ convergence-rate.
quarticConvergence :: ConvergenceOrder
quarticConvergence	= 4

-- | The predicted number of iterations, required to achieve a specific accuracy, at a given /order of convergence/.
getIterationsRequired :: Integral i
=> ConvergenceOrder
-> DecimalDigits	-- ^ The precision of the initial estimate.
-> DecimalDigits	-- ^ The required precision.
-> i
getIterationsRequired convergenceOrder initialDecimalDigits requiredDecimalDigits
| initialDecimalDigits <= 0	= error \$ "Factory.Math.Precision.getIterationsRequired:\tinsufficient 'initialDecimalDigits'; " ++ show initialDecimalDigits
| precisionRatio <= 1		= 0
| otherwise			= ceiling \$ fromIntegral convergenceOrder `logBase` precisionRatio
where
precisionRatio :: Double
precisionRatio	= fromIntegral requiredDecimalDigits / fromIntegral initialDecimalDigits

{- |
* The predicted number of terms which must be extracted from a series,
if it is to converge to the required accuracy,
at the specified linear /convergence-rate/.

* The /convergence-rate/ of a series, is the error in the series after summation of @(n+1)th@ terms,
divided by the error after only @n@ terms, as the latter tends to infinity.
As such, for a /convergent/ series (in which the error get smaller with successive terms), it's value lies in the range @0 .. 1@.

* <http://en.wikipedia.org/wiki/Rate_of_convergence>.
-}
getTermsRequired :: Integral i
=> ConvergenceRate
-> DecimalDigits	-- ^ The additional number of correct decimal digits.
-> i
getTermsRequired _ 0		= 0
getTermsRequired convergenceRate requiredDecimalDigits
| convergenceRate <= 0 || convergenceRate >= 1	= error \$ "Factory.Math.Precision.getTermsRequired:\t(0 < convergence-rate < 1); " ++ show convergenceRate
| requiredDecimalDigits < 0			= error \$ "Factory.Math.Precision.getTermsRequired:\t'requiredDecimalDigits' must be positive; " ++ show requiredDecimalDigits
| otherwise					= ceiling \$ fromIntegral requiredDecimalDigits / negate (logBase 10 convergenceRate)

-- | Rounds the specified number, to a positive number of 'DecimalDigits'.
roundTo :: (RealFrac a, Fractional f) => DecimalDigits -> a -> f
roundTo decimals = (/ fromInteger promotionFactor) . fromInteger . round . (* fromInteger promotionFactor)	where
promotionFactor :: Integer
promotionFactor	= 10 ^ decimals

-- | Promotes the specified number, by a positive number of 'DecimalDigits'.
promote :: Num n => n -> DecimalDigits -> n
promote x	= (* x) . (10 ^)

{- |
* Reduces a 'Rational' to the minimal form required for the specified number of /fractional/ decimal places;
irrespective of the number of integral decimal places.

* A 'Rational' approximation to an irrational number, may be very long, and provide an unknown excess precision.
Whilst this doesn't sound harmful, it costs in performance and memory-requirement, and being unpredictable isn't actually useful.
-}
simplify :: RealFrac operand
=> DecimalDigits	-- ^ The number of places after the decimal point, which are required.
-> operand
-> Rational
simplify decimalDigits operand	= Data.Ratio.approxRational operand . recip \$ 4 * 10 ^ succ decimalDigits	-- Tolerate any error less than half the least significant digit required.

```