aern2-mp 0.1.4 → 0.2.0.0
raw patch · 33 files changed
+1597/−1650 lines, 33 filesdep +cdar-mBounddep +collect-errorsdep +deepseqdep −cdardep −convertibledep −lensdep ~QuickCheckdep ~hspecdep ~mixed-types-num
Dependencies added: cdar-mBound, collect-errors, deepseq, reflection
Dependencies removed: cdar, convertible, lens, rounded
Dependency ranges changed: QuickCheck, hspec, mixed-types-num
Files
- LICENSE +1/−1
- README.md +87/−0
- aern2-mp.cabal +89/−111
- changelog.md +6/−0
- src-cdar/AERN2/MP/Float/Arithmetic.hs +0/−104
- src-cdar/AERN2/MP/Float/Conversions.hs +0/−152
- src-cdar/AERN2/MP/Float/Type.hs +0/−79
- src-rounded/AERN2/MP/Float/Arithmetic.hs +0/−117
- src-rounded/AERN2/MP/Float/Conversions.hs +0/−178
- src-rounded/AERN2/MP/Float/RoundedAdaptor.hs +0/−84
- src-rounded/AERN2/MP/Float/Type.hs +0/−76
- src/AERN2/Limit.hs +26/−0
- src/AERN2/MP.hs +2/−0
- src/AERN2/MP/Accuracy.hs +43/−62
- src/AERN2/MP/Ball.hs +9/−2
- src/AERN2/MP/Ball/Comparisons.hs +133/−132
- src/AERN2/MP/Ball/Conversions.hs +39/−5
- src/AERN2/MP/Ball/Elementary.hs +33/−21
- src/AERN2/MP/Ball/Field.hs +63/−128
- src/AERN2/MP/Ball/Limit.hs +155/−0
- src/AERN2/MP/Ball/PreludeOps.hs +13/−11
- src/AERN2/MP/Ball/Tests.hs +4/−2
- src/AERN2/MP/Ball/Type.hs +75/−72
- src/AERN2/MP/Dyadic.hs +86/−136
- src/AERN2/MP/Enclosure.hs +158/−131
- src/AERN2/MP/ErrorBound.hs +20/−8
- src/AERN2/MP/Float/Arithmetic.hs +104/−0
- src/AERN2/MP/Float/Conversions.hs +152/−0
- src/AERN2/MP/Float/Tests.hs +25/−25
- src/AERN2/MP/Float/Type.hs +79/−0
- src/AERN2/MP/Precision.hs +13/−12
- src/AERN2/MP/WithCurrentPrec.hs +181/−0
- src/AERN2/Norm.hs +1/−1
LICENSE view
@@ -1,4 +1,4 @@-Copyright (c) 2015-2017 Michal Konecny+Copyright (c) 2015-2021 Michal Konecny All rights reserved. Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
+ README.md view
@@ -0,0 +1,87 @@+# aern2-mp++Variable-precision interval arithmetic++## Numeric data types++This package provides the following two data types:++* `Dyadic`: variable-precision floats with exact ring operations+ +* `MPBall`: variable-precision interval (float centre ± error bound) with field & elementary interval operations+ +The type `MPBall` has instances of both [mixed-types-num](https://hackage.haskell.org/package/mixed-types-num) type classes such as `CanAdd`, `CanSqrt` as well as with traditional Prelude type classes such as `Ord`, `Num` and `Floating`.+The type `Dyadic` also has an appropriate subset of such instances.++Package [aern2-real](../aern2-real/README.md) provides an arithmetic of exact real numbers as converging lazy sequences of `MPBall`s of increasing precision. Exact real numbers offer additional convenience and readability to validated numeric programming.++### Examples++First, let us test interval arithmetic with **Prelude** operations:++ $ stack ghci aern2-mp:lib --no-load --ghci-options AERN2.MP+ *AERN2.MP> import Prelude+ *AERN2.MP Prelude>++ ...> pi100 = piBallP (prec 100)+ ...> pi10000 = piBallP (prec 10000)++ ...> pi100 ^ 2+ [9.8696044010893586188344909998725639610631902560... ± ~8.1120e-30 ~2^(-96)]++ ...> pi100 ^ pi100+ <interactive>:18:1: error:+ • No instance for (Integral MPBall) arising from a use of ‘^’++ ...> sin pi100+ [0.0000000000000000000000000000001694818351060767... ± ~7.8925e-31 ~2^(-99)]++ ...> sin pi10000+ [0.0000000000000000000000000000000000000000000000... ± ~0.0000 ~2^(-9999)]+ (0.07 secs, 64,466,432 bytes)++ ...> pi100 > 0+ True++ ...> pi100 == pi100+ *** Exception: Failed to decide equality of MPBalls. If you switch to MixedTypesNumPrelude instead of Prelude, comparison of MPBalls returns Kleenean instead of Bool.++Some things do not work with Prelude. Let us try using **MixedTypesNumPrelude** operations:++ $ stack ghci aern2-mp:lib --no-load --ghci-options AERN2.MP+ *AERN2.MP> import MixedTypesNumPrelude+ *AERN2.MP MixedTypesNumPrelude>++ ...> pi100 = piBallP (prec 100)+ ...> pi10000 = piBallP (prec 10000)++ ...> pi100 ^ pi100+ [36.4621596072079117709908260226198218149834948802... ± ~1.8696e-28 ~2^(-92)]++ ...> pi10000 ^ pi10000+ [36.4621596072079117709908260226921236663655084022... ± ~0.0000 ~2^(-9992)]+ (0.27 secs, 204,657,248 bytes)++ ...> pi100 > 0+ CertainTrue++ ...> pi100 == pi100+ TrueOrFalse++### Internal types and backends++The type `MPBall` internally uses the type:++* `MPFloat`: arbitrary-precision floats with both upwards and downwards-rounded arithmetic operations such as `*^` and `*.`++The package uses [cdar-mBound](https://hackage.haskell.org/package/cdar-mBound), a fork of [cdar](https://github.com/jensblanck/cdar) as its backend for `Dyadic` and `MPFloat`.++In previous versions, there was an MPFR backend via [rounded](https://hackage.haskell.org/package/rounded). This may be added again in future.++## Specifications and tests++This package also provides a fairly complete hspec/QuickCheck specification of algebraic properties for the above types. ++For `MPFloat`, the properties are given mostly as approximate versions of algebraic equalities with a small rounding error tolerance. ++For `MPBall`, the properties are given mostly as (interval) set over-approximations of the usual algebraic equalities.
aern2-mp.cabal view
@@ -1,128 +1,106 @@+cabal-version: 1.12++-- This file has been generated from package.yaml by hpack version 0.33.0.+--+-- see: https://github.com/sol/hpack+--+-- hash: 210a1027cd95f3ead272fe18229ed64036a8b1408d9ed7558c82f9655f6b8248+ name: aern2-mp-version: 0.1.4-cabal-version: >= 1.9.2-build-type: Simple-homepage: https://github.com/michalkonecny/aern2+version: 0.2.0.0+synopsis: Multi-precision ball (interval) arithmetic+description: Please see the README on GitHub at <https://github.com/michalkonecny/aern2/#readme>+category: Math+homepage: https://github.com/michalkonecny/aern2#readme+bug-reports: https://github.com/michalkonecny/aern2/issues author: Michal Konecny-maintainer: Michal Konecny <mikkonecny@gmail.com>-copyright: (c) 2015-2019 Michal Konecny+maintainer: mikkonecny@gmail.com+copyright: 2015-2021 Michal Konecny license: BSD3 license-file: LICENSE-extra-source-files: changelog.md-stability: experimental-category: Math-synopsis: Multi-precision ball (interval) arithmetic-Description:- This package provides the following types:- .- * Dyadic: variable-precision floats with exact ring operations- .- * MPBall: float ± error bound with field & elementary interval-like operations- .- The types have instances of both <https://hackage.haskell.org/package/mixed-types-num MixedTypeNumPrelude> - type classes as well as with traditional Prelude type classes.- .- There is a plan to add an Integer-only backend so that aern2-mp can- be used without MPFR.+build-type: Simple+extra-source-files:+ README.md+ changelog.md source-repository head- type: git- location: https://github.com/mikkonecny/aern2.git- subdir: aern2-mp--flag UseCDAR- Description: Use CDAR (mBound branch) as an Integer-only backend instead of MPFR- Default: False+ type: git+ location: https://github.com/michalkonecny/aern2 library- hs-source-dirs: src+ exposed-modules:+ AERN2.Limit+ AERN2.MP+ AERN2.MP.Accuracy+ AERN2.MP.Ball+ AERN2.MP.Ball.Comparisons+ AERN2.MP.Ball.Conversions+ AERN2.MP.Ball.Elementary+ AERN2.MP.Ball.Field+ AERN2.MP.Ball.Limit+ AERN2.MP.Ball.PreludeOps+ AERN2.MP.Ball.Tests+ AERN2.MP.Ball.Type+ AERN2.MP.Dyadic+ AERN2.MP.Enclosure+ AERN2.MP.ErrorBound+ AERN2.MP.Float+ AERN2.MP.Float.Arithmetic+ AERN2.MP.Float.Auxi+ AERN2.MP.Float.Conversions+ AERN2.MP.Float.Operators+ AERN2.MP.Float.Tests+ AERN2.MP.Float.Type+ AERN2.MP.Precision+ AERN2.MP.WithCurrentPrec+ AERN2.Norm+ AERN2.Normalize+ AERN2.Utils.Bench+ other-modules:+ Paths_aern2_mp+ hs-source-dirs:+ src+ default-extensions: RebindableSyntax, ScopedTypeVariables, DeriveGeneric, TypeFamilies, TypeOperators, ConstraintKinds, DefaultSignatures, MultiParamTypeClasses, FlexibleContexts, FlexibleInstances, UndecidableInstances+ other-extensions: TemplateHaskell+ ghc-options: -Wall build-depends:- base == 4.*+ QuickCheck+ , base ==4.*+ , cdar-mBound >=0.1.0.0+ , collect-errors >=0.1.3+ , deepseq+ , hspec , integer-logarithms- , convertible+ , mixed-types-num >=0.5.1+ , reflection , regex-tdfa- , hspec- -- , hspec-smallcheck >= 0.3- , QuickCheck- , lens , template-haskell- , mixed-types-num >= 0.3.2- if flag(UseCDAR)- hs-source-dirs: src-cdar- build-depends:- cdar- else- hs-source-dirs: src-rounded- build-depends:- rounded == 0.1.*- ghc-options: -Wall -fno-warn-orphans- extensions:- RebindableSyntax,- PostfixOperators,- ScopedTypeVariables,- DeriveGeneric,- GeneralizedNewtypeDeriving,- TypeFamilies,- TypeOperators,- ConstraintKinds,- DefaultSignatures,- MultiParamTypeClasses,- FlexibleContexts,- FlexibleInstances,- UndecidableInstances- if !flag(UseCDAR)- exposed-modules:- AERN2.MP.Float.RoundedAdaptor- exposed-modules:- -- modules that depend on backend choice:- AERN2.MP.Float.Type- AERN2.MP.Float.Arithmetic- AERN2.MP.Float.Conversions- -- modules common to all backends:- AERN2.Utils.Bench- AERN2.Normalize- AERN2.Norm- AERN2.MP.Precision- AERN2.MP.Accuracy- AERN2.MP.Enclosure- AERN2.MP.ErrorBound- AERN2.MP.Float.Auxi- AERN2.MP.Float.Operators- AERN2.MP.Float.Tests- AERN2.MP.Float- AERN2.MP.Dyadic- AERN2.MP.Ball.Type- AERN2.MP.Ball.Conversions- AERN2.MP.Ball.Comparisons- AERN2.MP.Ball.Field- AERN2.MP.Ball.Elementary- AERN2.MP.Ball.PreludeOps- AERN2.MP.Ball- AERN2.MP.Ball.Tests- AERN2.MP+ default-language: Haskell2010 -test-suite spec- type:- exitcode-stdio-1.0- ghc-options:- -Wall- extensions:- RebindableSyntax,- PostfixOperators,- ScopedTypeVariables,- FlexibleContexts+test-suite aern2-mp-test+ type: exitcode-stdio-1.0+ main-is: Spec.hs+ other-modules:+ AERN2.MP.BallSpec+ AERN2.MP.DyadicSpec+ AERN2.MP.FloatSpec+ Paths_aern2_mp hs-source-dirs: test- main-is:- Spec.hs- other-modules:- AERN2.MP.BallSpec- AERN2.MP.DyadicSpec- AERN2.MP.FloatSpec+ default-extensions: RebindableSyntax, ScopedTypeVariables, DeriveGeneric, TypeFamilies, TypeOperators, ConstraintKinds, DefaultSignatures, MultiParamTypeClasses, FlexibleContexts, FlexibleInstances, UndecidableInstances+ other-extensions: TemplateHaskell+ ghc-options: -threaded -rtsopts -with-rtsopts=-N -Wall build-depends:- base == 4.*- -- , mixed-types-num >= 0.3.1 && < 0.4+ QuickCheck , aern2-mp- , hspec >= 2.1- -- , hspec-smallcheck >= 0.3- , QuickCheck >= 2.7+ , base ==4.*+ , cdar-mBound >=0.1.0.0+ , collect-errors >=0.1.3+ , deepseq+ , hspec+ , integer-logarithms+ , mixed-types-num >=0.5.1+ , reflection+ , regex-tdfa+ , template-haskell+ default-language: Haskell2010
changelog.md view
@@ -1,5 +1,11 @@ # Change log for aern2-mp +* v 0.2.0 2021-05-17+ * switch to new simplified collect-errors, mixed-types-num 0.5.0+ * got rid of EnsureCE etc.+ * not introducing CN wrapper unless at least one parameter is already CN+ * using CDAR backend only, no MPFR for now+ * WithCurrentPrec for specifying default precision via types * v 0.1.4 2019-03-19 * CDAR-based Integer-only backend * needs the mBound branch of CDAR
− src-cdar/AERN2/MP/Float/Arithmetic.hs
@@ -1,104 +0,0 @@-{-|- Module : AERN2.MP.Float.Arithmetic- Description : Arbitrary precision floating point numbers- Copyright : (c) Michal Konecny- License : BSD3-- Maintainer : mikkonecny@gmail.com- Stability : experimental- Portability : portable-- Arbitrary precision floating-point numbers with up/down-rounded operations.--}--module AERN2.MP.Float.Arithmetic- (- -- * MPFloat basic arithmetic- addCEDU, subCEDU- , mulCEDU, divCEDU, recipCEDU- -- * MPFloat selected constants and operations- , piCEDU- , cosCEDU, sinCEDU- , sqrtCEDU, expCEDU, logCEDU- )-where--import MixedTypesNumPrelude-import qualified Prelude as P--import AERN2.MP.Precision--import qualified Data.CDAR as MPLow--import AERN2.MP.Float.Auxi-import AERN2.MP.Float.Type--{- common functions -}--instance CanNeg MPFloat where- negate = P.negate--instance CanAbs MPFloat where- abs = P.abs--addCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat-addCEDU = binaryCEDU (P.+)--subCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat-subCEDU = binaryCEDU (P.-)--mulCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat-mulCEDU = binaryCEDU (P.*)--divCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat-divCEDU x y - | y P.== (P.fromInteger 0) = getBoundsCEDU MPLow.Bottom- | otherwise = binaryCEDU (P./) x y--recipCEDU :: MPFloat -> BoundsCEDU MPFloat-recipCEDU = unaryCEDU P.recip--{- special constants and functions -}--piCEDU :: Precision -> BoundsCEDU MPFloat-piCEDU pp = - getBoundsCEDU $ MPLow.piA (p2cdarPrec pp)--cosCEDU :: MPFloat -> BoundsCEDU MPFloat-cosCEDU = unaryPrecCEDU 0 MPLow.cosA--sinCEDU :: MPFloat -> BoundsCEDU MPFloat-sinCEDU = unaryPrecCEDU 0 MPLow.sinA- -sqrtCEDU :: MPFloat -> BoundsCEDU MPFloat-sqrtCEDU = unaryCEDU MPLow.sqrtA- -expCEDU :: MPFloat -> BoundsCEDU MPFloat-expCEDU = unaryCEDU MPLow.expA--logCEDU :: MPFloat -> BoundsCEDU MPFloat-logCEDU = unaryCEDU MPLow.logA--{- auxiliary functions to automatically determine result precision from operand precisions -}--binaryCEDU ::- (MPFloat -> MPFloat -> MPFloat) ->- (MPFloat -> MPFloat -> BoundsCEDU MPFloat)-binaryCEDU op x y =- getBoundsCEDU $ op x y--unaryCEDU ::- (MPFloat -> MPFloat) ->- (MPFloat -> BoundsCEDU MPFloat)-unaryCEDU op x =- getBoundsCEDU $ op x--unaryPrecCEDU ::- Integer ->- (MPLow.Precision -> MPFloat -> MPFloat) ->- (MPFloat -> BoundsCEDU MPFloat)-unaryPrecCEDU addPrec op x@(MPLow.Approx mb _ _ s) =- getBoundsCEDU $ op ((-s P.+ mb) P.+ (int addPrec)) x-unaryPrecCEDU addPrec op MPLow.Bottom =- getBoundsCEDU $ op ((int $ integer defaultPrecision) P.+ (int addPrec)) MPLow.Bottom-
− src-cdar/AERN2/MP/Float/Conversions.hs
@@ -1,152 +0,0 @@-{-|- Module : AERN2.MP.Float.Conversions- Description : Conversions and comparisons of arbitrary precision floats- Copyright : (c) Michal Konecny- License : BSD3-- Maintainer : mikkonecny@gmail.com- Stability : experimental- Portability : portable-- Conversions and comparisons of arbitrary precision floating point numbers--}--module AERN2.MP.Float.Conversions- (- -- * MPFloat to other types (see also instances)- toDouble- -- * MPFloat constructors (see also instances)- , CanBeMPFloat, mpFloat- , fromIntegerCEDU- , fromRationalCEDU- -- * comparisons and constants (see also instances)- , zero, one, two- , nan, infinity- )-where--import MixedTypesNumPrelude-import qualified Prelude as P--import Data.Ratio-import Data.Convertible---- import AERN2.Norm-import AERN2.MP.Precision--import qualified Data.CDAR as MPLow--import AERN2.MP.Float.Auxi-import AERN2.MP.Float.Type-import AERN2.MP.Float.Arithmetic---{- conversions to MPFloat -}--type CanBeMPFloat t = ConvertibleExactly t MPFloat-mpFloat :: (CanBeMPFloat t) => t -> MPFloat-mpFloat = convertExactly--instance ConvertibleExactly Integer MPFloat where- safeConvertExactly =- Right . P.fromInteger--instance ConvertibleExactly Int MPFloat where- safeConvertExactly = safeConvertExactly . integer--fromIntegerCEDU :: Precision -> Integer -> BoundsCEDU MPFloat-fromIntegerCEDU pp =- setPrecisionCEDU pp . P.fromInteger--fromRationalCEDU :: Precision -> Rational -> BoundsCEDU MPFloat-fromRationalCEDU pp =- setPrecisionCEDU pp . (MPLow.toApprox (p2cdarPrec pp))--{- conversions from MPFloat -}--instance ConvertibleExactly MPFloat Rational where- safeConvertExactly = Right . P.toRational- -toDouble :: MPFloat -> Double-toDouble = P.fromRational . rational--instance Convertible MPFloat Double where- safeConvert x- | isFinite dbl = Right dbl- | otherwise = convError "conversion to double: out of bounds" x- where- dbl = toDouble x---instance CanRound MPFloat where- properFraction x = (n,f)- where- r = rational x- n = (numerator r) `P.quot` (denominator r)- f = ceduCentre $ x `subCEDU` (P.fromInteger n)- -{- comparisons -}--instance HasEqAsymmetric MPFloat MPFloat-instance HasEqAsymmetric MPFloat Integer where- equalTo = convertSecond equalTo-instance HasEqAsymmetric Integer MPFloat where- equalTo = convertFirst equalTo-instance HasEqAsymmetric MPFloat Int where- equalTo = convertSecond equalTo-instance HasEqAsymmetric Int MPFloat where- equalTo = convertFirst equalTo-instance HasEqAsymmetric MPFloat Rational where- equalTo = convertFirst equalTo-instance HasEqAsymmetric Rational MPFloat where- equalTo = convertSecond equalTo--instance CanTestZero MPFloat--instance HasOrderAsymmetric MPFloat MPFloat-instance HasOrderAsymmetric MPFloat Integer where- lessThan = convertSecond lessThan- leq = convertSecond leq-instance HasOrderAsymmetric Integer MPFloat where- lessThan = convertFirst lessThan- leq = convertFirst leq-instance HasOrderAsymmetric MPFloat Int where- lessThan = convertSecond lessThan- leq = convertSecond leq-instance HasOrderAsymmetric Int MPFloat where- lessThan = convertFirst lessThan- leq = convertFirst leq-instance HasOrderAsymmetric Rational MPFloat where- lessThan = convertSecond lessThan- leq = convertSecond leq-instance HasOrderAsymmetric MPFloat Rational where- lessThan = convertFirst lessThan- leq = convertFirst leq--instance CanTestPosNeg MPFloat--{- min, max -}--instance CanMinMaxAsymmetric MPFloat MPFloat--{- constants -}--zero, one, two :: MPFloat-zero = mpFloat 0-one = mpFloat 1-two = mpFloat 2--nan, infinity :: MPFloat-nan = MPLow.Bottom-infinity = nan--itisNaN :: MPFloat -> Bool-itisNaN MPLow.Bottom = True-itisNaN _ = False--instance CanTestFinite MPFloat where- isInfinite = itisNaN- isFinite = not . itisNaN--instance CanTestNaN MPFloat where- isNaN = itisNaN
− src-cdar/AERN2/MP/Float/Type.hs
@@ -1,79 +0,0 @@-{-# LANGUAGE DeriveGeneric, DeriveDataTypeable, StandaloneDeriving #-}-{-|- Module : AERN2.MP.Float.Type- Description : Arbitrary precision floating point numbers (via cdar)- Copyright : (c) Michal Konecny- License : BSD3-- Maintainer : mikkonecny@gmail.com- Stability : experimental- Portability : portable-- Arbitrary precision floating-point numbers, re-using CDAR Approx type.--}--module AERN2.MP.Float.Type- (- -- * MPFloat numbers and their basic operations- MPFloat- , showMPFloat- , setPrecisionCEDU- , p2cdarPrec- , getBoundsCEDU- )-where--import MixedTypesNumPrelude-import qualified Prelude as P---- import Data.Bits (unsafeShiftL)-import Data.Typeable--import AERN2.Norm-import AERN2.MP.Precision-import AERN2.MP.Float.Auxi--import qualified Data.CDAR as MPLow--{-| Multiple-precision floating-point type based on CDAR.Approx with 0 radius. -}-type MPFloat = MPLow.Approx--showMPFloat :: MPFloat -> String-showMPFloat x = MPLow.showA x--deriving instance (Typeable MPFloat)--p2cdarPrec :: Precision -> MPLow.Precision-p2cdarPrec = P.fromInteger . integer--getBoundsCEDU :: MPFloat -> BoundsCEDU MPFloat-getBoundsCEDU (MPLow.Approx mb m e s) = - BoundsCEDU - (MPLow.Approx mb m 0 s) (MPLow.approxMB eb_mb e 0 s)- (MPLow.Approx mb (m-e) 0 s) (MPLow.Approx mb (m+e) 0 s)-getBoundsCEDU MPLow.Bottom =- BoundsCEDU- MPLow.Bottom MPLow.Bottom MPLow.Bottom MPLow.Bottom--{-| The bit-size bound for the error bound in CEDU -}-eb_prec :: Precision-eb_prec = prec 63--{-| The bit-size bound for the error bound in CEDU -}-eb_mb :: Int-eb_mb = int $ integer eb_prec--instance HasPrecision MPFloat where- getPrecision (MPLow.Approx mb _ _ _) = prec (P.toInteger $ mb)- getPrecision MPLow.Bottom = error "illegal MPFloat (Bottom)"- --instance CanSetPrecision MPFloat where- setPrecision p = ceduCentre . setPrecisionCEDU p--setPrecisionCEDU :: Precision -> MPFloat -> BoundsCEDU MPFloat-setPrecisionCEDU pp = getBoundsCEDU . MPLow.enforceMB . MPLow.setMB (p2cdarPrec pp)--instance HasNorm MPFloat where- getNormLog (MPLow.Approx _ m _ s) = (getNormLog m) + (integer s)- getNormLog MPLow.Bottom = error "getNormLog undefined for Bottom"
− src-rounded/AERN2/MP/Float/Arithmetic.hs
@@ -1,117 +0,0 @@-{-|- Module : AERN2.MP.Float.Arithmetic- Description : Arbitrary precision floating point numbers- Copyright : (c) Michal Konecny- License : BSD3-- Maintainer : mikkonecny@gmail.com- Stability : experimental- Portability : portable-- Arbitrary precision floating-point numbers with up/down-rounded operations.--}--module AERN2.MP.Float.Arithmetic- (- -- * MPFloat basic arithmetic- addCEDU, subCEDU- , mulCEDU, divCEDU, recipCEDU- -- * MPFloat selected constants and operations- , piCEDU- , cosCEDU, sinCEDU- , sqrtCEDU, expCEDU, logCEDU- -- * auxiliary functions- , constCEDU, unaryCEDU, binaryCEDU- )-where--import MixedTypesNumPrelude-import qualified Prelude as P--import AERN2.MP.Precision--import qualified AERN2.MP.Float.RoundedAdaptor as MPLow--import AERN2.MP.Float.Auxi-import AERN2.MP.Float.Type--one :: MPFloat-one = MPLow.one--{- common functions -}--instance CanNeg MPFloat where- negate = ceduUp . unaryCEDU MPLow.neg--instance CanAbs MPFloat where- abs x- | x P.< MPLow.zero = negate x- | otherwise = x---addCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat-addCEDU = binaryCEDU MPLow.add--subCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat-subCEDU = binaryCEDU MPLow.sub--mulCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat-mulCEDU = binaryCEDU MPLow.mul--divCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat-divCEDU = binaryCEDU MPLow.div--recipCEDU :: MPFloat -> BoundsCEDU MPFloat-recipCEDU x = divCEDU one x--{- special constants and functions -}--piCEDU :: Precision -> BoundsCEDU MPFloat-piCEDU pp = - constCEDU MPLow.pi (p2mpfrPrec pp)--cosCEDU :: MPFloat -> BoundsCEDU MPFloat-cosCEDU = unaryCEDU MPLow.cos--sinCEDU :: MPFloat -> BoundsCEDU MPFloat-sinCEDU = unaryCEDU MPLow.sin- -sqrtCEDU :: MPFloat -> BoundsCEDU MPFloat-sqrtCEDU = unaryCEDU MPLow.sqrt- -expCEDU :: MPFloat -> BoundsCEDU MPFloat-expCEDU = unaryCEDU MPLow.exp- -logCEDU :: MPFloat -> BoundsCEDU MPFloat-logCEDU = unaryCEDU MPLow.log--{- auxiliary functions to automatically determine result precision from operand precisions -}--binaryCEDU :: - (MPLow.RoundMode -> MPLow.Precision -> MPFloat -> MPFloat -> MPFloat) -> - MPFloat -> MPFloat -> BoundsCEDU MPFloat-binaryCEDU op x y =- getCEDU d u- where- d = op MPLow.Down p x y- u = op MPLow.Up p x y- p = p2mpfrPrec $ (getPrecision x) `max` (getPrecision y)--unaryCEDU :: - (MPLow.RoundMode -> MPLow.Precision -> MPFloat -> MPFloat) -> - MPFloat -> BoundsCEDU MPFloat-unaryCEDU op x =- getCEDU d u- where- d = op MPLow.Down p x- u = op MPLow.Up p x- p = p2mpfrPrec $ getPrecision x--constCEDU :: - (MPLow.RoundMode -> MPLow.Precision -> MPFloat) -> - MPLow.Precision -> BoundsCEDU MPFloat-constCEDU op p =- getCEDU d u- where- d = op MPLow.Down p- u = op MPLow.Up p
− src-rounded/AERN2/MP/Float/Conversions.hs
@@ -1,178 +0,0 @@-{-|- Module : AERN2.MP.Float.Conversions- Description : Conversions and comparisons of arbitrary precision floats- Copyright : (c) Michal Konecny- License : BSD3-- Maintainer : mikkonecny@gmail.com- Stability : experimental- Portability : portable-- Conversions and comparisons of arbitrary precision floating point numbers--}--module AERN2.MP.Float.Conversions- (- -- * MPFloat to other types (see also instances)- toDouble- -- * MPFloat constructors (see also instances)- , CanBeMPFloat, mpFloat- , fromIntegerCEDU- , fromRationalCEDU- -- * comparisons and constants (see also instances)- , zero, one, two- , nan, infinity- )-where--import MixedTypesNumPrelude-import qualified Prelude as P--import Data.Ratio-import Data.Convertible--import AERN2.Norm-import AERN2.MP.Precision--import AERN2.MP.Float.Auxi-import AERN2.MP.Float.Type-import AERN2.MP.Float.Arithmetic--import qualified AERN2.MP.Float.RoundedAdaptor as MPLow--mpToDouble :: MPLow.RoundMode -> MPFloat -> Double-mpToDouble = MPLow.toDoubleA--mpToRational :: MPFloat -> Rational-mpToRational x- | x == 0 = 0.0- | otherwise = MPLow.toRationalA x--mpFromRationalA :: MPLow.RoundMode -> MPLow.Precision -> Rational -> MPFloat-mpFromRationalA = MPLow.fromRationalA--{- conversions to MPFloat -}--type CanBeMPFloat t = ConvertibleExactly t MPFloat-mpFloat :: (CanBeMPFloat t) => t -> MPFloat-mpFloat = convertExactly--instance ConvertibleExactly Integer MPFloat where- safeConvertExactly n =- findExact $ map (flip fromIntegerCEDU n) $ standardPrecisions initPrec- where- initPrec =- case getNormLog n of- NormBits b -> prec (b + 8)- _ -> prec 8- findExact [] =- convError "integer too high to represent exactly" n- findExact (cedu : rest)- | ceduErr cedu P.> zero = findExact rest- | otherwise = Right (ceduCentre cedu)--instance ConvertibleExactly Int MPFloat where- safeConvertExactly = safeConvertExactly . integer--fromIntegerCEDU :: Precision -> Integer -> BoundsCEDU MPFloat-fromIntegerCEDU pp n =- constCEDU (\r p -> MPLow.fromIntegerA r p n) (p2mpfrPrec pp)--fromRationalCEDU :: Precision -> Rational -> BoundsCEDU MPFloat-fromRationalCEDU pp q =- constCEDU (\r p -> mpFromRationalA r p q) (p2mpfrPrec pp)--{- conversions from MPFloat -}--instance ConvertibleExactly MPFloat Rational where- safeConvertExactly = Right . mpToRational--toDouble :: MPFloat -> Double-toDouble = mpToDouble MPLow.Up--instance Convertible MPFloat Double where- safeConvert x- | isFinite dbl = Right dbl- | otherwise = convError "conversion to double: out of bounds" x- where- dbl = toDouble x--instance CanRound MPFloat where- properFraction x = (n,f)- where- r = rational x- n = (numerator r) `P.quot` (denominator r)- f = ceduCentre $ x `subCEDU` (mpFloat n)--{- comparisons -}--instance HasEqAsymmetric MPFloat MPFloat-instance HasEqAsymmetric MPFloat Integer where- equalTo = convertSecond equalTo-instance HasEqAsymmetric Integer MPFloat where- equalTo = convertFirst equalTo-instance HasEqAsymmetric MPFloat Int where- equalTo = convertSecond equalTo-instance HasEqAsymmetric Int MPFloat where- equalTo = convertFirst equalTo-instance HasEqAsymmetric MPFloat Rational where- equalTo = convertFirst equalTo-instance HasEqAsymmetric Rational MPFloat where- equalTo = convertSecond equalTo--instance CanTestZero MPFloat--instance HasOrderAsymmetric MPFloat MPFloat-instance HasOrderAsymmetric MPFloat Integer where- lessThan = convertSecond lessThan- leq = convertSecond leq-instance HasOrderAsymmetric Integer MPFloat where- lessThan = convertFirst lessThan- leq = convertFirst leq-instance HasOrderAsymmetric MPFloat Int where- lessThan = convertSecond lessThan- leq = convertSecond leq-instance HasOrderAsymmetric Int MPFloat where- lessThan = convertFirst lessThan- leq = convertFirst leq-instance HasOrderAsymmetric Rational MPFloat where- lessThan = convertSecond lessThan- leq = convertSecond leq-instance HasOrderAsymmetric MPFloat Rational where- lessThan = convertFirst lessThan- leq = convertFirst leq--instance CanTestPosNeg MPFloat--{- min, max -}--instance CanMinMaxAsymmetric MPFloat MPFloat--{- constants -}--zero, one, two :: MPFloat-zero = MPLow.zero-one = MPLow.one-two = MPLow.add MPLow.Up (MPLow.getPrec one) one one--nan, infinity :: MPFloat-nan = ceduCentre $ divCEDU zero zero -infinity = ceduCentre $ divCEDU one zero --itisNaN :: MPFloat -> Bool-itisNaN x = not $ x P.== x--itisInfinite :: MPFloat -> Bool-itisInfinite x =- ceduCentre (mulCEDU x two) P.== x- &&- x P./= zero--instance CanTestFinite MPFloat where- isInfinite = itisInfinite- isFinite x = not (itisInfinite x || itisNaN x)--instance CanTestNaN MPFloat where- isNaN = itisNaN--
− src-rounded/AERN2/MP/Float/RoundedAdaptor.hs
@@ -1,84 +0,0 @@-{-# LANGUAGE DataKinds, ExistentialQuantification, RankNTypes #-}--- {-# LANGUAGE DeriveGeneric, DeriveDataTypeable, StandaloneDeriving #-}-{-|- Module : AERN2.MP.Float.UseRounded.RoundedAdaptor- Description : Numeric.Rounded + variable precision- Copyright : (c) Michal Konecny- License : BSD3-- Maintainer : mikkonecny@gmail.com- Stability : experimental- Portability : portable-- Numeric.Rounded + variable precision--}-module AERN2.MP.Float.RoundedAdaptor-(- module AERN2.MP.Float.RoundedAdaptor-, module Numeric.Rounded.Simple-)-where--import Prelude hiding (div, pi)--- import qualified Prelude as P--import Numeric.Rounded.Simple--- import qualified Numeric.RoundedSimple as R--instance Show Rounded where- show = show'--getPrec :: Rounded -> Int-getPrec = precision--getExp :: Rounded -> Int-getExp = exponent'--data RoundMode = Up | Down--withRoundMode :: (RoundingMode -> t) -> (RoundMode -> t)-withRoundMode op Up = op TowardInf-withRoundMode op Down = op TowardNegInf-{-# INLINE withRoundMode #-}--set :: RoundMode -> Precision -> Rounded -> Rounded-set = withRoundMode precRound--defaultPrecision :: Precision-defaultPrecision = 53--pi :: RoundMode -> Precision -> Rounded-pi = withRoundMode kPi--fromIntegerA :: RoundMode -> Precision -> Integer -> Rounded-fromIntegerA = withRoundMode fromInteger'--zero, one :: Rounded-zero = fromIntegerA Up defaultPrecision 0-one = fromIntegerA Up defaultPrecision 1--toDoubleA :: RoundMode -> Rounded -> Double-toDoubleA = withRoundMode toDouble--fromRationalA :: RoundMode -> Precision -> Rational -> Rounded-fromRationalA = withRoundMode fromRational'--toRationalA :: Rounded -> Rational-toRationalA = toRational' TowardNearest--add, sub, mul, div, atan2 :: RoundMode -> Precision -> Rounded -> Rounded -> Rounded-add = withRoundMode add_-sub = withRoundMode sub_-mul = withRoundMode mul_-div = withRoundMode div_-atan2 = withRoundMode atan2_--neg, abs, sqrt, exp, log, sin, cos :: RoundMode -> Precision -> Rounded -> Rounded-neg = withRoundMode negate_-abs = withRoundMode abs_-sqrt = withRoundMode sqrt_-exp = withRoundMode exp_-log = withRoundMode log_-sin = withRoundMode sin_-cos = withRoundMode cos_--- TODO: add more ops
− src-rounded/AERN2/MP/Float/Type.hs
@@ -1,76 +0,0 @@-{-# LANGUAGE DeriveGeneric, DeriveDataTypeable, StandaloneDeriving #-}-{-|- Module : AERN2.MP.Float.Type- Description : Arbitrary precision floating point numbers (MPFR)- Copyright : (c) Michal Konecny- License : BSD3-- Maintainer : mikkonecny@gmail.com- Stability : experimental- Portability : portable-- Arbitrary precision floating-point numbers using MPFR via package rounded.--}--module AERN2.MP.Float.Type- (- -- * MPFloat numbers and their basic operations- MPFloat- , showMPFloat- , setPrecisionCEDU- , getCEDU- , p2mpfrPrec- )-where--import MixedTypesNumPrelude-import qualified Prelude as P--import AERN2.Norm-import AERN2.MP.Precision-import AERN2.MP.Float.Auxi--import qualified AERN2.MP.Float.RoundedAdaptor as MPLow-import Data.Typeable--{-| Multiple-precision floating-point type based on MPFR via rounded. -}-type MPFloat = MPLow.Rounded--showMPFloat :: MPFloat -> String-showMPFloat = show--deriving instance (Typeable MPFloat)--p2mpfrPrec :: Precision -> MPLow.Precision-p2mpfrPrec = P.fromInteger . integer--getCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat-getCEDU d u = BoundsCEDU c e d u- where- c = u- e = MPLow.sub MPLow.Up (MPLow.getPrec c) u d--instance HasPrecision MPFloat where- getPrecision x = prec (P.toInteger $ MPLow.getPrec x)--instance CanSetPrecision MPFloat where- setPrecision = setPrecisionUp--setPrecisionUp :: Precision -> MPFloat -> MPFloat-setPrecisionUp p = MPLow.set MPLow.Up (p2mpfrPrec p)--setPrecisionDown :: Precision -> MPFloat -> MPFloat-setPrecisionDown p = MPLow.set MPLow.Down (p2mpfrPrec p)--setPrecisionCEDU :: Precision -> MPFloat -> BoundsCEDU MPFloat-setPrecisionCEDU p x =- getCEDU d u - where- d = setPrecisionDown p x- u = setPrecisionUp p x--instance HasNorm MPFloat where- getNormLog x- | x P.== MPLow.zero = NormZero- | otherwise = NormBits (P.toInteger $ MPLow.getExp x)-
+ src/AERN2/Limit.hs view
@@ -0,0 +1,26 @@+{-|+ Module : AERN2.Limit+ Description : limit operation+ Copyright : (c) Michal Konecny+ License : BSD3++ Maintainer : mikkonecny@gmail.com+ Stability : experimental+ Portability : portable++ Generic limit operation+-}+module AERN2.Limit where++-- import MixedTypesNumPrelude++-- import qualified Numeric.CollectErrors as CN++---------+-- limit+---------++class HasLimits ix s where+ type LimitType ix s+ limit :: (ix -> s) -> LimitType ix s+
src/AERN2/MP.hs view
@@ -18,6 +18,8 @@ , module AERN2.MP.ErrorBound , module AERN2.MP.Enclosure , MPBall(..), CanBeMPBall, mpBall, CanBeMPBallP, mpBallP+, reducePrecionIfInaccurate+, giveUpIfVeryInaccurate ) where
src/AERN2/MP/Accuracy.hs view
@@ -8,20 +8,19 @@ Stability : experimental Portability : portable - A type for roughly measuring the accuracy of an enclosure.+ A type for measuring the accuracy of an enclosing set,+ roughly corresponding to the maximum absolute error in some distance metric+ approximately measured in bits. -} module AERN2.MP.Accuracy (Accuracy(NoInformation, Exact), bits, fromAccuracy, HasAccuracy(..),- HasAccuracyGuide(..), CanSetAccuracyGuide(..), adjustAccuracyGuide,- getFiniteAccuracy, ac2prec,- CanReduceSizeUsingAccuracyGuide(..),- specCanReduceSizeUsingAccuracyGuide, iterateUntilAccurate,- convergentList2CauchySeq,- seqByPrecision2CauchySeq,+ convergentList2seqByAccuracy,+ seqByPrecision2seqByAccuracy, setPrecisionAtLeastAccuracy,+ ShowWithAccuracy(..), HasApproximate(..)) where @@ -32,8 +31,8 @@ import Data.Complex -import Test.Hspec-import Test.QuickCheck+-- import Test.Hspec+import Test.QuickCheck ( Arbitrary(arbitrary), frequency ) import AERN2.Norm import AERN2.MP.Precision@@ -174,52 +173,48 @@ class HasAccuracy a where getAccuracy :: a -> Accuracy+ {-| Return accuracy, except when the element is Exact, return its nominal Precision dressed as Accuracy.+ This function is useful when we have a convergent sequence where all elements happen to be+ actually equal to the limit and we need the property that the sequence elements keep improving.+ -}+ getFiniteAccuracy :: a -> Accuracy+ default getFiniteAccuracy :: (HasPrecision a) => a -> Accuracy+ getFiniteAccuracy b =+ case getAccuracy b of+ Exact -> bits $ getPrecision b+ a -> a -instance (HasAccuracy a, SuitableForCE es) => HasAccuracy (CollectErrors es a) where+instance (HasAccuracy a, CanBeErrors es) => HasAccuracy (CollectErrors es a) where getAccuracy (CollectErrors ma es) = case ma of Just a | not (hasCertainError es) -> getAccuracy a _ -> NoInformation+ getFiniteAccuracy (CollectErrors ma es) = + case ma of+ Just a | not (hasCertainError es) -> getFiniteAccuracy a+ _ -> NoInformation -instance HasAccuracy Int where getAccuracy _ = Exact-instance HasAccuracy Integer where getAccuracy _ = Exact-instance HasAccuracy Rational where getAccuracy _ = Exact-instance HasAccuracy Bool where getAccuracy _ = Exact+instance HasAccuracy Int where getAccuracy _ = Exact; getFiniteAccuracy _ = NoInformation+instance HasAccuracy Integer where getAccuracy _ = Exact; getFiniteAccuracy _ = NoInformation+instance HasAccuracy Rational where getAccuracy _ = Exact; getFiniteAccuracy _ = NoInformation+instance HasAccuracy Bool where getAccuracy _ = Exact; getFiniteAccuracy _ = NoInformation+instance HasAccuracy Kleenean where getAccuracy _ = Exact; getFiniteAccuracy _ = NoInformation instance HasAccuracy t => HasAccuracy (Complex t) where getAccuracy (a :+ i) = (getAccuracy a) `min` (getAccuracy i)-+ getFiniteAccuracy (a :+ i) =+ (getFiniteAccuracy a) `min` (getFiniteAccuracy i)+ instance HasAccuracy t => HasAccuracy [t] where getAccuracy xs = foldl min Exact $ map getAccuracy xs+ getFiniteAccuracy xs = foldl min Exact $ map getFiniteAccuracy xs instance HasAccuracy t => HasAccuracy (Maybe t) where getAccuracy (Just x) = getAccuracy x getAccuracy _ = NoInformation--class HasAccuracyGuide a where- getAccuracyGuide :: a -> Accuracy--class HasAccuracyGuide a => CanSetAccuracyGuide a where- setAccuracyGuide :: Accuracy -> a -> a--adjustAccuracyGuide ::- (CanSetAccuracyGuide a) =>- (Accuracy -> Accuracy) -> a -> a-adjustAccuracyGuide adj_acG a =- setAccuracyGuide (adj_acG (getAccuracyGuide a)) a--{-| Return accuracy, except when the element is Exact, return its nominal Precision dressed as Accuracy.- This function is useful when we have a convergent sequence where all elements happen to be- actually equal to the limit and we need the property that the sequence elements keep improving.--}-getFiniteAccuracy ::- (HasAccuracy t, HasPrecision t) =>- t -> Accuracy-getFiniteAccuracy b =- case getAccuracy b of- Exact -> bits $ getPrecision b- a -> a+ getFiniteAccuracy (Just x) = getFiniteAccuracy x+ getFiniteAccuracy _ = NoInformation iterateUntilAccurate :: (HasAccuracy t) =>@@ -238,22 +233,22 @@ Bits b -> prec (max 2 $ b + 50) _ -> prec 100 -seqByPrecision2CauchySeq ::+seqByPrecision2seqByAccuracy :: (HasAccuracy t) => (Precision -> t) -> (Accuracy -> t)-seqByPrecision2CauchySeq seqByPrecision ac =- convergentList2CauchySeq list ac+seqByPrecision2seqByAccuracy seqByPrecision ac =+ convergentList2seqByAccuracy list ac where list = map seqByPrecision $ dropWhile (lowPrec ac) (standardPrecisions (ac2prec ac)) lowPrec Exact _ = False lowPrec _ p = bits p < ac -convergentList2CauchySeq :: (HasAccuracy t) => [t] -> (Accuracy -> t)-convergentList2CauchySeq list ac = findAccurate list+convergentList2seqByAccuracy :: (HasAccuracy t) => [t] -> (Accuracy -> t)+convergentList2seqByAccuracy list ac = findAccurate list where findAccurate [] =- error "convergentList2CauchySeq: the sequence either converges too slowly or it does not converge"+ error "convergentList2seqByAccuracy: the sequence either converges too slowly or it does not converge" findAccurate (b : rest) | getAccuracy b >= ac = b | otherwise = findAccurate rest@@ -263,7 +258,7 @@ it is at least as high as the supplied accuracy (assuming the accuracy is finite). -}-setPrecisionAtLeastAccuracy :: (CanSetPrecision t) => Accuracy -> t -> t+setPrecisionAtLeastAccuracy :: (HasPrecision t, CanSetPrecision t) => Accuracy -> t -> t setPrecisionAtLeastAccuracy acc b | p_b < p_acc = setPrecision p_acc b | otherwise = b@@ -275,22 +270,8 @@ _ -> prec $ max 2 (fromAccuracy acc) p_b = getPrecision b --class CanReduceSizeUsingAccuracyGuide t where- reduceSizeUsingAccuracyGuide :: Accuracy -> t -> t--specCanReduceSizeUsingAccuracyGuide ::- ( CanReduceSizeUsingAccuracyGuide t- , HasEqCertainly t t- , Arbitrary t, Show t)- =>- (T t) -> Spec-specCanReduceSizeUsingAccuracyGuide (T tName :: T t) =- describe ("CanReduceSizeUsingAccuracyGuide " ++ tName) $ do- it "is safe" $- property $- \ (t :: t) (ac :: Accuracy) ->- reduceSizeUsingAccuracyGuide ac t ?==? t+class ShowWithAccuracy t where+ showWithAccuracy :: Accuracy -> t -> String {-| An unsafe approximation of an enclosure or exact value, useful mainly for showing something brief and readable to humans.
src/AERN2/MP/Ball.hs view
@@ -1,3 +1,4 @@+{-# OPTIONS_GHC -Wno-orphans #-} {-| Module : AERN2.MP.Ball Description : Arbitrary precision ball arithmetic@@ -19,12 +20,17 @@ , module AERN2.MP.ErrorBound , module AERN2.MP.Enclosure -- * The Ball type- , MPBall(..), CanBeMPBall, mpBall, CanBeMPBallP, mpBallP+ , MPBall(..), CanBeMPBall, mpBall, cnMPBall, CanBeMPBallP, mpBallP, cnMPBallP , reducePrecionIfInaccurate+ , giveUpIfVeryInaccurate+ {- -- * Ball construction/extraction functions- -- , endpointsMP, fromEndpointsMP+ , endpointsMP, fromEndpointsMP+ -} -- * Ball operations (see also instances) , piBallP+ , intersectCNMPBall+ , hullMPBall -- * Helpers for constructing ball functions , byEndpointsMP , fromApproxWithLipschitz@@ -45,6 +51,7 @@ import AERN2.MP.Ball.Conversions () import AERN2.MP.Ball.Comparisons import AERN2.MP.Ball.Field ()+import AERN2.MP.Ball.Limit () import AERN2.MP.Ball.Elementary import AERN2.MP.Ball.PreludeOps ()
src/AERN2/MP/Ball/Comparisons.hs view
@@ -1,3 +1,4 @@+{-# OPTIONS_GHC -Wno-orphans #-} {-| Module : AERN2.MP.Ball.Comparisons Description : Comparisons of arbitrary precision dyadic balls@@ -18,13 +19,19 @@ , reducePrecionIfInaccurate -- * Helpers for constructing ball functions , byEndpointsMP+ -- * intersection and hull+ , intersectCNMPBall+ , hullMPBall ) where import MixedTypesNumPrelude -- import qualified Prelude as P +import qualified Control.CollectErrors as CE import Control.CollectErrors+ ( CollectErrors(getMaybeValue), CanBeErrors )+import qualified Numeric.CollectErrors as CN import AERN2.Norm import AERN2.MP.Dyadic (Dyadic)@@ -38,176 +45,168 @@ {- comparisons -} instance HasEqAsymmetric MPBall MPBall where- type EqCompareType MPBall MPBall = Maybe Bool+ type EqCompareType MPBall MPBall = Kleenean b1 `equalTo` b2 = b1 >= b2 && b1 <= b2 instance HasEqAsymmetric MPBall Integer where- type EqCompareType MPBall Integer = Maybe Bool+ type EqCompareType MPBall Integer = Kleenean b1 `equalTo` b2 = b1 >= b2 && b1 <= b2 instance HasEqAsymmetric Integer MPBall where- type EqCompareType Integer MPBall = Maybe Bool+ type EqCompareType Integer MPBall = Kleenean b1 `equalTo` b2 = b1 >= b2 && b1 <= b2 instance HasEqAsymmetric MPBall Int where- type EqCompareType MPBall Int = Maybe Bool+ type EqCompareType MPBall Int = Kleenean b1 `equalTo` b2 = b1 >= b2 && b1 <= b2 instance HasEqAsymmetric Int MPBall where- type EqCompareType Int MPBall = Maybe Bool+ type EqCompareType Int MPBall = Kleenean b1 `equalTo` b2 = b1 >= b2 && b1 <= b2 instance HasEqAsymmetric MPBall Rational where- type EqCompareType MPBall Rational = Maybe Bool+ type EqCompareType MPBall Rational = Kleenean b1 `equalTo` b2 = b1 >= b2 && b1 <= b2 instance HasEqAsymmetric Rational MPBall where- type EqCompareType Rational MPBall = Maybe Bool+ type EqCompareType Rational MPBall = Kleenean b1 `equalTo` b2 = b1 >= b2 && b1 <= b2 instance HasEqAsymmetric MPBall Dyadic where- type EqCompareType MPBall Dyadic = Maybe Bool+ type EqCompareType MPBall Dyadic = Kleenean b1 `equalTo` b2 = b1 >= b2 && b1 <= b2 instance HasEqAsymmetric Dyadic MPBall where- type EqCompareType Dyadic MPBall = Maybe Bool+ type EqCompareType Dyadic MPBall = Kleenean b1 `equalTo` b2 = b1 >= b2 && b1 <= b2 instance (HasEqAsymmetric MPBall b- , CanEnsureCE es b- , CanEnsureCE es (EqCompareType MPBall b)- , IsBool (EnsureCE es (EqCompareType MPBall b))- , SuitableForCE es)+ , IsBool (EqCompareType MPBall b)+ , CanBeErrors es) =>- HasEqAsymmetric MPBall (CollectErrors es b)+ HasEqAsymmetric MPBall (CollectErrors es b) where- type EqCompareType MPBall (CollectErrors es b) =- EnsureCE es (EqCompareType MPBall b)- equalTo = lift2TLCE equalTo+ type EqCompareType MPBall (CollectErrors es b) =+ CollectErrors es (EqCompareType MPBall b)+ equalTo = CE.liftT1 equalTo instance (HasEqAsymmetric a MPBall- , CanEnsureCE es a- , CanEnsureCE es (EqCompareType a MPBall)- , IsBool (EnsureCE es (EqCompareType a MPBall))- , SuitableForCE es)+ , IsBool (EqCompareType a MPBall)+ , CanBeErrors es) => HasEqAsymmetric (CollectErrors es a) MPBall where- type EqCompareType (CollectErrors es a) MPBall =- EnsureCE es (EqCompareType a MPBall)- equalTo = lift2TCE equalTo+ type EqCompareType (CollectErrors es a) MPBall =+ CollectErrors es (EqCompareType a MPBall)+ equalTo = CE.lift1T equalTo instance HasOrderAsymmetric MPBall MPBall where- type OrderCompareType MPBall MPBall = Maybe Bool+ type OrderCompareType MPBall MPBall = Kleenean lessThan b1 b2- | r1 < l2 = Just True- | r2 <= l1 = Just False- | otherwise = Nothing+ | r1 < l2 = CertainTrue+ | r2 <= l1 = CertainFalse+ | otherwise = TrueOrFalse where- (l1, r1) = endpointsMP b1- (l2, r2) = endpointsMP b2+ (l1, r1) = endpoints b1+ (l2, r2) = endpoints b2 leq b1 b2- | r1 <= l2 = Just True- | r2 < l1 = Just False- | otherwise = Nothing+ | r1 <= l2 = CertainTrue+ | r2 < l1 = CertainFalse+ | otherwise = TrueOrFalse where- (l1, r1) = endpointsMP b1- (l2, r2) = endpointsMP b2+ (l1, r1) = endpoints b1+ (l2, r2) = endpoints b2 instance HasOrderAsymmetric Integer MPBall where- type OrderCompareType Integer MPBall = Maybe Bool+ type OrderCompareType Integer MPBall = Kleenean lessThan = convertFirst lessThan leq = convertFirst leq instance HasOrderAsymmetric MPBall Integer where- type OrderCompareType MPBall Integer = Maybe Bool+ type OrderCompareType MPBall Integer = Kleenean lessThan = convertSecond lessThan leq = convertSecond leq instance HasOrderAsymmetric Int MPBall where- type OrderCompareType Int MPBall = Maybe Bool+ type OrderCompareType Int MPBall = Kleenean lessThan = convertFirst lessThan leq = convertFirst leq instance HasOrderAsymmetric MPBall Int where- type OrderCompareType MPBall Int = Maybe Bool+ type OrderCompareType MPBall Int = Kleenean lessThan = convertSecond lessThan leq = convertSecond leq instance HasOrderAsymmetric Dyadic MPBall where- type OrderCompareType Dyadic MPBall = Maybe Bool+ type OrderCompareType Dyadic MPBall = Kleenean lessThan = convertFirst lessThan leq = convertFirst leq instance HasOrderAsymmetric MPBall Dyadic where- type OrderCompareType MPBall Dyadic = Maybe Bool+ type OrderCompareType MPBall Dyadic = Kleenean lessThan = convertSecond lessThan leq = convertSecond leq instance HasOrderAsymmetric MPBall Rational where- type OrderCompareType MPBall Rational = Maybe Bool+ type OrderCompareType MPBall Rational = Kleenean lessThan b1 q2- | r1 < l2 = Just True- | r2 <= l1 = Just False- | otherwise = Nothing+ | r1 < l2 = CertainTrue+ | r2 <= l1 = CertainFalse+ | otherwise = TrueOrFalse where- (l1, r1) = endpointsMP b1+ (l1, r1) = endpoints b1 l2 = q2 r2 = q2 leq b1 q2- | r1 <= l2 = Just True- | r2 < l1 = Just False- | otherwise = Nothing+ | r1 <= l2 = CertainTrue+ | r2 < l1 = CertainFalse+ | otherwise = TrueOrFalse where- (l1, r1) = endpointsMP b1+ (l1, r1) = endpoints b1 l2 = q2 r2 = q2 instance HasOrderAsymmetric Rational MPBall where- type OrderCompareType Rational MPBall = Maybe Bool+ type OrderCompareType Rational MPBall = Kleenean lessThan q1 b2- | r1 < l2 = Just True- | r2 <= l1 = Just False- | otherwise = Nothing+ | r1 < l2 = CertainTrue+ | r2 <= l1 = CertainFalse+ | otherwise = TrueOrFalse where- (l2, r2) = endpointsMP b2+ (l2, r2) = endpoints b2 l1 = q1 r1 = q1 leq q1 b2- | r1 <= l2 = Just True- | r2 < l1 = Just False- | otherwise = Nothing+ | r1 <= l2 = CertainTrue+ | r2 < l1 = CertainFalse+ | otherwise = TrueOrFalse where- (l2, r2) = endpointsMP b2+ (l2, r2) = endpoints b2 l1 = q1 r1 = q1 instance (HasOrderAsymmetric MPBall b- , CanEnsureCE es b- , CanEnsureCE es (OrderCompareType MPBall b)- , IsBool (EnsureCE es (OrderCompareType MPBall b))- , SuitableForCE es)+ , IsBool (OrderCompareType MPBall b)+ , CanBeErrors es) => HasOrderAsymmetric MPBall (CollectErrors es b) where type OrderCompareType MPBall (CollectErrors es b) =- EnsureCE es (OrderCompareType MPBall b)- lessThan = lift2TLCE lessThan- leq = lift2TLCE leq- greaterThan = lift2TLCE greaterThan- geq = lift2TLCE geq+ CollectErrors es (OrderCompareType MPBall b)+ lessThan = CE.liftT1 lessThan+ leq = CE.liftT1 leq+ greaterThan = CE.liftT1 greaterThan+ geq = CE.liftT1 geq instance (HasOrderAsymmetric a MPBall- , CanEnsureCE es a- , CanEnsureCE es (OrderCompareType a MPBall)- , IsBool (EnsureCE es (OrderCompareType a MPBall))- , SuitableForCE es)+ , IsBool (OrderCompareType a MPBall)+ , CanBeErrors es) => HasOrderAsymmetric (CollectErrors es a) MPBall where type OrderCompareType (CollectErrors es a) MPBall =- EnsureCE es (OrderCompareType a MPBall)- lessThan = lift2TCE lessThan- leq = lift2TCE leq- greaterThan = lift2TCE greaterThan- geq = lift2TCE geq+ CollectErrors es (OrderCompareType a MPBall)+ lessThan = CE.lift1T lessThan+ leq = CE.lift1T leq+ greaterThan = CE.lift1T greaterThan+ geq = CE.lift1T geq instance CanTestZero MPBall instance CanTestPosNeg MPBall@@ -265,117 +264,119 @@ instance (CanMinMaxAsymmetric MPBall b- , CanEnsureCE es b- , CanEnsureCE es (MinMaxType MPBall b)- , SuitableForCE es)+ , CanBeErrors es) => CanMinMaxAsymmetric MPBall (CollectErrors es b) where type MinMaxType MPBall (CollectErrors es b) =- EnsureCE es (MinMaxType MPBall b)- min = lift2TLCE min- max = lift2TLCE max+ CollectErrors es (MinMaxType MPBall b)+ min = CE.liftT1 min+ max = CE.liftT1 max instance (CanMinMaxAsymmetric a MPBall- , CanEnsureCE es a- , CanEnsureCE es (MinMaxType a MPBall)- , SuitableForCE es)+ , CanBeErrors es) => CanMinMaxAsymmetric (CollectErrors es a) MPBall where type MinMaxType (CollectErrors es a) MPBall =- EnsureCE es (MinMaxType a MPBall)- min = lift2TCE min- max = lift2TCE max+ CollectErrors es (MinMaxType a MPBall)+ min = CE.lift1T min+ max = CE.lift1T max {- intersection -} instance CanIntersectAsymmetric MPBall MPBall where intersect a b- | rL > rR =- noValueNumErrorCertainCN $ NumError $ "intersect: empty intersection: " ++ show a ++ "; " ++ show b- | otherwise = cn $ fromEndpointsMP rL rR+ | l > r =+ CN.noValueNumErrorCertain $ CN.NumError $ "intersect: empty intersection: " ++ show a ++ "; " ++ show b+ | otherwise = cn $ setPrecision p $ fromMPFloatEndpoints l r where- rL = max aL bL- rR = min aR bR- (aL,aR) = endpointsMP a- (bL,bR) = endpointsMP b+ p = getPrecision a+ l = max aL bL+ r = min aR bR+ (aL,aR) = endpoints a+ (bL,bR) = endpoints b +intersectCNMPBall :: CN MPBall -> CN MPBall -> CN MPBall+intersectCNMPBall = intersect+ -- case (fst $ ensureNoCN x, fst $ ensureNoCN y) of + -- (Nothing, Nothing) -> x+ -- (Just _ , Nothing) -> x+ -- (Nothing, Just _ ) -> y+ -- (Just _ , Just _ ) -> lift2CE intersect x y+ instance (CanIntersectAsymmetric MPBall b- , CanEnsureCE es b- , CanEnsureCE es (IntersectionType MPBall b)- , SuitableForCE es)+ , CanBeErrors es) => CanIntersectAsymmetric MPBall (CollectErrors es b) where type IntersectionType MPBall (CollectErrors es b) =- EnsureCE es (IntersectionType MPBall b)- intersect = lift2TLCE intersect+ CollectErrors es (IntersectionType MPBall b)+ intersect = CE.liftT1 intersect instance (CanIntersectAsymmetric a MPBall- , CanEnsureCE es a- , CanEnsureCE es (IntersectionType a MPBall)- , SuitableForCE es)+ , CanBeErrors es) => CanIntersectAsymmetric (CollectErrors es a) MPBall where type IntersectionType (CollectErrors es a) MPBall =- EnsureCE es (IntersectionType a MPBall)- intersect = lift2TCE intersect+ CollectErrors es (IntersectionType a MPBall)+ intersect = CE.lift1T intersect +{- hull -}++hullMPBall :: MPBall -> MPBall -> MPBall+hullMPBall a b = + fromEndpoints rL rR+ where+ rL = min aL bL+ rR = max aR bR+ (aL,aR) = endpoints a+ (bL,bR) = endpoints b+ {- union -} instance CanUnionAsymmetric MPBall MPBall where union a b =- case getMaybeValueCN (a `intersect` b) of- Just _ -> prependErrorsCN [(ErrorPotential, err)] r- _ -> prependErrorsCN [(ErrorCertain, err)] r+ case getMaybeValue (a `intersect` b) of+ Just _ -> r+ _ -> CN.prependErrorCertain err r where- err = NumError $ "union of enclosures: not enclosing the same value"- r = cn $ fromEndpointsMP rL rR- rL = min aL bL- rR = max aR bR- (aL,aR) = endpointsMP a- (bL,bR) = endpointsMP b-+ err = CN.NumError $ "union of enclosures: not enclosing the same value"+ r = cn $ hullMPBall a b instance (CanUnionAsymmetric MPBall b- , CanEnsureCE es b- , CanEnsureCE es (UnionType MPBall b)- , SuitableForCE es)+ , CanBeErrors es) => CanUnionAsymmetric MPBall (CollectErrors es b) where type UnionType MPBall (CollectErrors es b) =- EnsureCE es (UnionType MPBall b)- union = lift2TLCE union+ CollectErrors es (UnionType MPBall b)+ union = CE.liftT1 union instance (CanUnionAsymmetric a MPBall- , CanEnsureCE es a- , CanEnsureCE es (UnionType a MPBall)- , SuitableForCE es)+ , CanBeErrors es) => CanUnionAsymmetric (CollectErrors es a) MPBall where type UnionType (CollectErrors es a) MPBall =- EnsureCE es (UnionType a MPBall)- union = lift2TCE union+ CollectErrors es (UnionType a MPBall)+ union = CE.lift1T union {-|- Computes an *increasing* ball fucntion @f@ from *exact* MPFR operations.+ Compute an MPBall function from *exact* MPFloat operations on interval endpoints.+ This works only for *non-decreasing* operations, eg addition, min, max. -} byEndpointsMP :: (MPFloat -> MPFloat -> MPFloat) -> (MPBall -> MPBall -> MPBall) byEndpointsMP op b1 b2 =- fromEndpointsMP (l1 `op` l2) (r1 `op` r2)+ fromEndpoints (l1 `op` l2) (r1 `op` r2) where- (l1,r1) = endpointsMP b1- (l2,r2) = endpointsMP b2--{- random generation -}+ (l1,r1) = endpoints b1+ (l2,r2) = endpoints b2
src/AERN2/MP/Ball/Conversions.hs view
@@ -1,3 +1,4 @@+{-# OPTIONS_GHC -Wno-orphans #-} {-| Module : AERN2.MP.Ball.Conversions Description : Conversions of arbitrary precision dyadic balls@@ -19,16 +20,19 @@ import MixedTypesNumPrelude -- import qualified Prelude as P +-- import qualified Numeric.CollectErrors as CN+ import Data.Typeable-import Data.Convertible+-- import Data.Convertible import AERN2.MP.Dyadic (Dyadic, dyadic) import qualified AERN2.MP.Float as MPFloat-import AERN2.MP.Float (mpFloat)+import AERN2.MP.Float (MPFloat, mpFloat) -- import AERN2.MP.Float.Operators import AERN2.MP.Precision -- import qualified AERN2.MP.ErrorBound as EB-import AERN2.MP.ErrorBound (ErrorBound, errorBound)+-- import qualified AERN2.MP.ErrorBound as EB+import AERN2.MP.ErrorBound (ErrorBound, errorBound, CanBeErrorBound) import AERN2.MP.Ball.Type @@ -38,13 +42,13 @@ integerBounds b = (floor l, ceiling r) where- (l,r) = endpointsMP b+ (l,r) = endpoints b instance Convertible MPBall ErrorBound where safeConvert b = Right (errorBound (max (abs l) (abs r))) where- (l,r) = endpointsMP b+ (l,r) = endpoints b {--- constructing an exact ball ---} @@ -115,3 +119,33 @@ b = MPBall xFlt (xe + eUp) -- beware, precision may be too high relative to accuracy (MPBall xFlt xe) = mpBallP p x eUp = errorBound e++{--- constructing a fat ball ---}++instance (CanBeErrorBound t) => CanPlusMinus MPBall t where+ plusMinus b e = updateRadius (+ (errorBound e)) b++instance (CanBeErrorBound t) => CanPlusMinus (CN MPBall) t where+ plusMinus b e = updateRadius (+ (errorBound e)) b++instance (CanBeErrorBound t) => CanPlusMinus MPFloat t where+ type PlusMinusType MPFloat t = MPBall+ plusMinus b e = MPBall b (errorBound e)++instance (CanBeErrorBound t) => CanPlusMinus Dyadic t where+ type PlusMinusType Dyadic t = MPBall+ plusMinus b e = MPBall (mpFloat b) (errorBound e)++instance (CanBeErrorBound t) => CanPlusMinus Integer t where+ type PlusMinusType Integer t = MPBall+ plusMinus b e = MPBall (mpFloat b) (errorBound e)++instance (CanBeErrorBound t) => CanPlusMinus Int t where+ type PlusMinusType Int t = MPBall+ plusMinus b e = MPBall (mpFloat b) (errorBound e)++instance (CanBeErrorBound t) => CanPlusMinus Rational t where+ type PlusMinusType Rational t = MPBall+ plusMinus b e = (mpBallP p b) +- e+ where+ p = prec 100
src/AERN2/MP/Ball/Elementary.hs view
@@ -1,3 +1,5 @@+{-# LANGUAGE TemplateHaskell #-}+{-# OPTIONS_GHC -Wno-orphans #-} {-| Module : AERN2.MP.Ball.Elementary Description : Elementary operations on arbitrary precision dyadic balls@@ -52,20 +54,22 @@ sinB :: Integer -> MPBall -> MPBall sinB i x = -- increasingPrecisionUntilNotImproving (fromApproxWithLipschitz MPFloat.sinDown MPFloat.sinUp lip) x+ max (-1) $ min 1 $ fromApproxWithLipschitz MPFloat.sinCEDU lip x where lip | i == 0 = mpFloat 1- | otherwise = snd $ endpointsMP $ abs $ cosB (i - 1) x+ | otherwise = endpointR $ abs $ cosB (i - 1) x cosB :: Integer -> MPBall -> MPBall cosB i x = -- increasingPrecisionUntilNotImproving (fromApproxWithLipschitz MPFloat.cosDown MPFloat.cosUp lip) x+ max (-1) $ min 1 $ fromApproxWithLipschitz MPFloat.cosCEDU lip x where lip | i == 0 = mpFloat 1- | otherwise = snd $ endpointsMP $ abs $ sinB (i - 1) x+ | otherwise = endpointR $ abs $ sinB (i - 1) x -- increasingPrecisionUntilNotImproving :: (MPBall -> MPBall) -> (MPBall -> MPBall) -- increasingPrecisionUntilNotImproving f x =@@ -89,49 +93,45 @@ exp = intervalFunctionByEndpointsUpDown MPFloat.expDown MPFloat.expUp instance CanLog MPBall where- type LogType MPBall = CN MPBall+ type LogType MPBall = MPBall log x | x_!>! 1 =- cn $ setPrecision p $ ballFunctionUsingLipschitz log_ logLip x_+ setPrecision p $ ballFunctionUsingLipschitz log_ logLip x_ | x_!>! 0 =- cn $ setPrecision p $ intervalFunctionByEndpoints log_ x_+ setPrecision p $ intervalFunctionByEndpoints log_ x_ | x !>! 0 =- cn $ setPrecision p $ intervalFunctionByEndpoints log_ x- | x !<=! 0 = noValueNumErrorCertainCN err- | otherwise = noValueNumErrorPotentialCN err+ intervalFunctionByEndpoints log_ x+ | otherwise = err where p = getPrecision x x_ = reducePrecionIfInaccurate x- err = OutOfRange $ "log: argument must be > 0: " ++ show x+ err = error $ "log: argument must be > 0: " ++ show x log_ (MPBall c e) = MPBall lc (e + (errorBound le)) where (lc, le) = ceduCentreErr $ MPFloat.logCEDU c- logLip y = errorBound $ (1/!y)+ logLip y = errorBound $ (1/y) instance CanPow MPBall MPBall where- powNoCN b e = (~!) $ pow b e- pow = powUsingExpLog (mpBall 0) (mpBall 1)+ pow = powUsingExpLog (mpBall 1) instance CanPow MPBall Dyadic where- powNoCN b e = (~!) $ pow b e- pow b e = powUsingExpLog (mpBall 0) (mpBall 1) b (mpBall e)+ pow b e = powUsingExpLog (mpBall 1) b (mpBall e) instance CanPow MPBall Rational where- powNoCN b e = (~!) $ pow b e- pow b e = powUsingExpLog (mpBall 0) (mpBall 1) b (mpBallP (getPrecision b) e)+ pow b e = powUsingExpLog (mpBall 1) b (mpBallP (getPrecision b) e) instance CanSqrt MPBall where- type SqrtType MPBall = CN MPBall+ type SqrtType MPBall = MPBall sqrt x- | x !>=! 0 = cn $ aux x- | x !<! 0 = noValueNumErrorCertainCN err- | otherwise = prependErrorsCN [(ErrorPotential, err)] $ cn $ aux (max 0 x)+ | x !>=! 0 = aux x+ | x ?>=? 0 = aux $ max 0 x+ | otherwise = err where aux = intervalFunctionByEndpointsUpDown (\ e -> MPFloat.sqrtDown (P.max (mpFloat 0) e)) (\ e -> MPFloat.sqrtUp (P.max (mpFloat 0) e))- err = OutOfRange $ "sqrt: argument must be >= 0: " ++ show x+ err = error $ "sqrt: argument must be >= 0: " ++ show x {- generic methods for computing real functions from MPFR-approximations -} @@ -151,3 +151,15 @@ setPrecision (getPrecision xc) $ -- beware, some MPFloat functions may increase precision, eg sine and cosine (MPBall fxC (errorBound fxErr)) err = (errorBound lip) * xe + fxe++$(declForTypes+ [[t| Integer |], [t| Int |], [t| Rational |]]+ (\ b -> [d|++ instance + CanPow $b MPBall + where+ type PowType $b MPBall = MPBall+ pow x e = pow (mpBallP (getPrecision e) x) e+ |]))+
src/AERN2/MP/Ball/Field.hs view
@@ -1,4 +1,7 @@ {-# LANGUAGE TemplateHaskell #-}+{-# LANGUAGE PartialTypeSignatures #-}+{-# OPTIONS_GHC -Wno-orphans #-}+{-# OPTIONS_GHC -Wno-partial-type-signatures #-} {-| Module : AERN2.MP.Ball.Field Description : Field operations on arbitrary precision dyadic balls@@ -18,7 +21,7 @@ import MixedTypesNumPrelude -- import qualified Prelude as P -import Control.CollectErrors+import qualified Numeric.CollectErrors as CN import AERN2.Normalize @@ -71,28 +74,20 @@ add = convertPFirst add instance- (CanAddAsymmetric MPBall b- , CanEnsureCE es b- , CanEnsureCE es (AddType MPBall b)- , SuitableForCE es)+ (CanAddAsymmetric MPBall b) =>- CanAddAsymmetric MPBall (CollectErrors es b)+ CanAddAsymmetric MPBall (CN b) where- type AddType MPBall (CollectErrors es b) =- EnsureCE es (AddType MPBall b)- add = lift2TLCE add+ type AddType MPBall (CN b) = CN (AddType MPBall b)+ add = CN.liftT1 add instance- (CanAddAsymmetric a MPBall- , CanEnsureCE es a- , CanEnsureCE es (AddType a MPBall)- , SuitableForCE es)+ (CanAddAsymmetric a MPBall) =>- CanAddAsymmetric (CollectErrors es a) MPBall+ CanAddAsymmetric (CN a) MPBall where- type AddType (CollectErrors es a) MPBall =- EnsureCE es (AddType a MPBall)- add = lift2TCE add+ type AddType (CN a) MPBall = CN (AddType a MPBall)+ add = CN.lift1T add {- subtraction -} @@ -111,28 +106,20 @@ instance CanSub Dyadic MPBall instance- (CanSub MPBall b- , CanEnsureCE es b- , CanEnsureCE es (SubType MPBall b)- , SuitableForCE es)+ (CanSub MPBall b) =>- CanSub MPBall (CollectErrors es b)+ CanSub MPBall (CN b) where- type SubType MPBall (CollectErrors es b) =- EnsureCE es (SubType MPBall b)- sub = lift2TLCE sub+ type SubType MPBall (CN b) = CN (SubType MPBall b)+ sub = CN.liftT1 sub instance- (CanSub a MPBall- , CanEnsureCE es a- , CanEnsureCE es (SubType a MPBall)- , SuitableForCE es)+ (CanSub a MPBall) =>- CanSub (CollectErrors es a) MPBall+ CanSub (CN a) MPBall where- type SubType (CollectErrors es a) MPBall =- EnsureCE es (SubType a MPBall)- sub = lift2TCE sub+ type SubType (CN a) MPBall = CN (SubType a MPBall)+ sub = CN.lift1T sub {- multiplication -} @@ -173,43 +160,27 @@ mul = convertPFirst mul instance- (CanMulAsymmetric MPBall b- , CanEnsureCE es b- , CanEnsureCE es (MulType MPBall b)- , SuitableForCE es)+ (CanMulAsymmetric MPBall b) =>- CanMulAsymmetric MPBall (CollectErrors es b)+ CanMulAsymmetric MPBall (CN b) where- type MulType MPBall (CollectErrors es b) =- EnsureCE es (MulType MPBall b)- mul = lift2TLCE mul+ type MulType MPBall (CN b) = CN (MulType MPBall b)+ mul = CN.liftT1 mul instance- (CanMulAsymmetric a MPBall- , CanEnsureCE es a- , CanEnsureCE es (MulType a MPBall)- , SuitableForCE es)+ (CanMulAsymmetric a MPBall) =>- CanMulAsymmetric (CollectErrors es a) MPBall+ CanMulAsymmetric (CN a) MPBall where- type MulType (CollectErrors es a) MPBall =- EnsureCE es (MulType a MPBall)- mul = lift2TCE mul+ type MulType (CN a) MPBall = CN (MulType a MPBall)+ mul = CN.lift1T mul {- division -} instance CanDiv MPBall MPBall where- type DivTypeNoCN MPBall MPBall = MPBall- divideNoCN b1 b2 = (~!) (divide b1 b2)- type DivType MPBall MPBall = CN MPBall- divide (MPBall x1 e1) b2@(MPBall x2 e2)- | isCertainlyNonZero b2 =- cn $ normalize $ MPBall x12C err- | isCertainlyZero b2 =- noValueNumErrorCertainCN DivByZero- | otherwise =- noValueNumErrorPotentialCN DivByZero+ type DivType MPBall MPBall = MPBall+ divide (MPBall x1 e1) (MPBall x2 e2) = normalize $ MPBall x12C err where (x12C, e12) = MPFloat.ceduCentreErr $ MPFloat.divCEDU x1 x2 x12AbsUp = (abs x12C) +^ e12@@ -241,113 +212,77 @@ [[t| Integer |], [t| Int |], [t| Dyadic |]] (\ t -> [d| instance CanDiv MPBall $t where- type DivType MPBall $t = CN MPBall+ type DivType MPBall $t = MPBall divide = convertSecond divide- type DivTypeNoCN MPBall $t = MPBall- divideNoCN = convertSecond divideNoCN instance CanDiv $t MPBall where- type DivType $t MPBall = CN MPBall+ type DivType $t MPBall = MPBall divide = convertFirst divide- type DivTypeNoCN $t MPBall = MPBall- divideNoCN = convertFirst divideNoCN |])) instance CanDiv Dyadic Dyadic where- type DivTypeNoCN Dyadic Dyadic = MPBall- divideNoCN a b = divideNoCN (mpBall a) (mpBall b)+ type DivType Dyadic Dyadic = MPBall divide a b = divide (mpBall a) (mpBall b) instance CanDiv MPBall Rational where- type DivTypeNoCN MPBall Rational = MPBall- divideNoCN = convertPSecond divideNoCN+ type DivType MPBall Rational = MPBall divide = convertPSecond divide instance CanDiv Rational MPBall where- type DivTypeNoCN Rational MPBall = MPBall- divideNoCN = convertPFirst divideNoCN+ type DivType Rational MPBall = MPBall divide = convertPFirst divide instance- (CanDiv MPBall b- , CanEnsureCE es b- , CanEnsureCE es (DivType MPBall b)- , CanEnsureCE es (DivTypeNoCN MPBall b)- , SuitableForCE es)+ (CanDiv MPBall b, CanTestZero b) =>- CanDiv MPBall (CollectErrors es b)+ CanDiv MPBall (CN b) where- type DivType MPBall (CollectErrors es b) =- EnsureCE es (DivType MPBall b)- divide = lift2TLCE divide- type DivTypeNoCN MPBall (CollectErrors es b) =- EnsureCE es (DivTypeNoCN MPBall b)- divideNoCN = lift2TLCE divideNoCN+ type DivType MPBall (CN b) = CN (DivType MPBall b)+ divide a b = divide (cn a) b instance- (CanDiv a MPBall- , CanEnsureCE es a- , CanEnsureCE es (DivType a MPBall)- , CanEnsureCE es (DivTypeNoCN a MPBall)- , SuitableForCE es)+ (CanDiv a MPBall) =>- CanDiv (CollectErrors es a) MPBall+ CanDiv (CN a) MPBall where- type DivType (CollectErrors es a) MPBall =- EnsureCE es (DivType a MPBall)- divide = lift2TCE divide- type DivTypeNoCN (CollectErrors es a) MPBall =- EnsureCE es (DivTypeNoCN a MPBall)- divideNoCN = lift2TCE divideNoCN+ type DivType (CN a) MPBall = CN (DivType a MPBall)+ divide a b = divide a (cn b) {- integer power -} instance CanPow MPBall Integer where- powNoCN b e = (~!) $ powUsingMulRecip (mpBall 1) b e- pow = powUsingMulRecip (mpBall 1)+ pow = powUsingMulRecipCutNeg (mpBall 1) instance CanPow MPBall Int where- powNoCN b e = (~!) $ powUsingMulRecip (mpBall 1) b e- pow = powUsingMulRecip (mpBall 1)+ pow = powUsingMulRecipCutNeg (mpBall 1) +powUsingMulRecipCutNeg :: _ => t -> t -> e -> DivType Integer t+powUsingMulRecipCutNeg one x e + | even e = + max 0 $ powUsingMulRecip one x e+ | otherwise = powUsingMulRecip one x e+ instance- (CanPow MPBall b- , CanEnsureCE es b- , CanEnsureCE es (PowType MPBall b)- , CanEnsureCE es (PowTypeNoCN MPBall b)- , SuitableForCE es)+ (CanPow MPBall b) =>- CanPow MPBall (CollectErrors es b)+ CanPow MPBall (CN b) where- type PowTypeNoCN MPBall (CollectErrors es b) =- EnsureCE es (PowTypeNoCN MPBall b)- powNoCN = lift2TLCE powNoCN- type PowType MPBall (CollectErrors es b) =- EnsureCE es (PowType MPBall b)- pow = lift2TLCE pow+ type PowType MPBall (CN b) = CN (PowType MPBall b)+ pow = CN.liftT1 pow instance- (CanPow a MPBall- , CanEnsureCE es a- , CanEnsureCE es (PowType a MPBall)- , CanEnsureCE es (PowTypeNoCN a MPBall)- , SuitableForCE es)+ (CanPow a MPBall) =>- CanPow (CollectErrors es a) MPBall+ CanPow (CN a) MPBall where- type PowTypeNoCN (CollectErrors es a) MPBall =- EnsureCE es (PowTypeNoCN a MPBall)- powNoCN = lift2TCE powNoCN- type PowType (CollectErrors es a) MPBall =- EnsureCE es (PowType a MPBall)- pow = lift2TCE pow+ type PowType (CN a) MPBall = CN (PowType a MPBall)+ pow = CN.lift1T pow instance CanDivIMod MPBall MPBall where+ type DivIType MPBall MPBall = Integer divIMod x m - | m !>! 0 = (cn d, cn xm)- | otherwise = (err (0 :: Integer), err xm)+ | m !>! 0 = (d, xm)+ | otherwise = error $ "modulus not positive: " ++ show m where- d = floor $ centre $ (centreAsBall x) /! (centreAsBall m)+ d = floor $ centre $ (centreAsBall x) / (centreAsBall m) xm = x - m*d- err :: (CanEnsureCN t) => t -> EnsureCN t- err s = noValueNumErrorCertainECN (Just s) $ OutOfRange $ "modulus not positive: " ++ show m
+ src/AERN2/MP/Ball/Limit.hs view
@@ -0,0 +1,155 @@+{-# LANGUAGE CPP #-}+-- #define DEBUG++{-# OPTIONS_GHC -Wno-orphans #-}+{-|+ Module : AERN2.MP.Ball.Limit+ Description : limits of MPBall sequences+ Copyright : (c) Michal Konecny+ License : BSD3++ Maintainer : mikkonecny@gmail.com+ Stability : experimental+ Portability : portable++ Limits of MPBall sequences.+-}+module AERN2.MP.Ball.Limit where++#ifdef DEBUG+import Debug.Trace (trace)+#define maybeTrace trace+#define maybeTraceIO putStrLn+#else+#define maybeTrace (\ (_ :: String) t -> t)+#define maybeTraceIO (\ (_ :: String) -> return ())+#endif++import MixedTypesNumPrelude++import qualified Numeric.CollectErrors as CN++-- import Math.NumberTheory.Logarithms (integerLog2)++import AERN2.Limit+import AERN2.MP.Ball.Type+import AERN2.MP.Ball.Field ()++{-+ The following strategies are inspired by+ Mueller: The iRRAM: Exact Arithmetic in C++, Section 10.1+ https://link.springer.com/chapter/10.1007/3-540-45335-0_14+-}++instance HasLimits Rational (CN MPBall -> CN MPBall) where+ type LimitType Rational (CN MPBall -> CN MPBall) = (CN MPBall -> CN MPBall)+ limit = limitMPBall (\ac -> 0.5^(fromAccuracy ac))++instance HasLimits Integer (CN MPBall -> CN MPBall) where+ type LimitType Integer (CN MPBall -> CN MPBall) = (CN MPBall -> CN MPBall)+ limit = limitMPBall (\ac -> (fromAccuracy ac))+ -- limit s = limit (s . getN)+ -- where+ -- -- q > 0 => 0 < 0.5^(getN q) <= q+ -- getN :: Rational -> Integer+ -- getN r = integer $ integerLog2 (max 1 $ (ceiling $ 1/r) - 1) + 1++instance HasLimits Int (CN MPBall -> CN MPBall) where+ type LimitType Int (CN MPBall -> CN MPBall) = (CN MPBall -> CN MPBall)+ limit s = limit (s . c)+ where+ c :: Integer -> Int+ c = int++ +limitMPBall :: (Accuracy -> ix) -> (ix -> (CN MPBall -> CN MPBall)) -> (CN MPBall -> CN MPBall)+limitMPBall ac2ix fs x =+ maybeTrace ("limit (CN MPBall -> CN MPBall): x = " ++ show x) $+ maybeTrace ("limit (CN MPBall -> CN MPBall): xPNext = " ++ show xPNext) $+ maybeTrace ("limit (CN MPBall -> CN MPBall): accuracies = " ++ show accuracies) $+ tryAccuracies accuracies+ where+ acX = getFiniteAccuracy x+ accuracies = aux (fromAccuracy acX)+ where+ aux a+ | a >= 4 = bits a : aux ((100 * a) `divI` 105)+ | otherwise = [bits a]+ xPNext = setPrecision (increaseP $ getPrecision x) x+ increaseP p =+ prec $ (integer p) + 10+ -- prec $ ((101 * (integer p)) `div` 100) + 1++ tryAccuracies [] =+ CN.noValueNumErrorPotential $ CN.NumError "limit (CN MPBall -> CN MPBall) failed"+ tryAccuracies (ac : rest) =+ let result = withAccuracy ac in+ case CN.hasError result of+ False -> result+ _ -> tryAccuracies rest++ withAccuracy ac =+ maybeTrace ("limit (CN MPBall -> CN MPBall): withAccuracy: ac = " ++ show ac) $+ (fs ix xPNext) +- e+ where+ ix = ac2ix ac+ e = 0.5^(fromAccuracy ac)++-- data WithLipschitz f = WithLipschitz f f++-- instance HasLimits Rational (WithLipschitz (MPBall -> CN MPBall)) where+-- type LimitType Rational (WithLipschitz (MPBall -> CN MPBall)) = (MPBall -> CN MPBall)+-- limit ffs xPre =+-- maybeTrace ("limit (MPBall -> CN MPBall): x = " ++ show x) $+-- maybeTrace ("limit (MPBall -> CN MPBall): xC = " ++ show xC) $+-- maybeTrace ("limit (MPBall -> CN MPBall): xE = " ++ show xE) $+-- maybeTrace ("limit (MPBall -> CN MPBall): accuracies = " ++ show accuracies) $+-- tryAccuracies accuracies+-- where+-- acX = getFiniteAccuracy x+-- accuracies = aux (fromAccuracy acX)+-- where+-- aux a+-- | a >= 4 = bits a : aux ((100 * a) `divI` 105)+-- | otherwise = [bits a]+-- x = increasePrec xPre+-- where+-- increasePrec z = setPrecision (inc $ getPrecision xPre) z+-- inc p =+-- prec $ (integer p) + 10+-- -- prec $ ((101 * (integer p)) `div` 100) + 1+-- xC = centreAsBall x+-- xE = radius x++-- tryAccuracies [] =+-- noValueNumErrorPotentialECN (Nothing :: Maybe MPBall) $+-- NumError "limit (MPBall -> CN MPBall) failed"+-- tryAccuracies (ac : rest) =+-- let result = withAccuracy ac in+-- case getErrorsCN result of+-- [] -> result+-- _ -> tryAccuracies rest++-- withAccuracy ac =+-- maybeTrace ("limit (MPBall -> CN MPBall): withAccuracy: ac = " ++ show ac) $+-- lift1CE (updateRadius (+ (errorBound e))) fx+-- where+-- e = 0.5^!(fromAccuracy ac)+-- WithLipschitz f f' = ffs e+-- fxC_CN = f xC+-- f'x_CN = f' x+-- fx :: CN MPBall+-- fx =+-- case (ensureNoCN fxC_CN, ensureNoCN f'x_CN) of+-- ((Just fxC, []), (Just f'x, [])) ->+-- cn $ updateRadius (+ (xE * (errorBound f'x))) fxC+-- _ ->+-- f x -- fallback++-- instance HasLimits Rational (WithLipschitz (CReal -> CReal)) where+-- type LimitType Rational (WithLipschitz (CReal -> CReal)) = (CReal -> CReal)+-- limit ffs x = limit fs x+-- where+-- fs e = f+-- where+-- WithLipschitz f _ = ffs e
src/AERN2/MP/Ball/PreludeOps.hs view
@@ -1,3 +1,4 @@+{-# OPTIONS_GHC -Wno-orphans #-} {-# LANGUAGE CPP #-} {-| Module : AERN2.MP.Ball.PreludeOps@@ -33,18 +34,19 @@ instance P.Eq MPBall where a == b = case a == b of- Just t -> t+ CertainTrue -> True+ CertainFalse -> False _ ->- error "Failed to decide equality of MPBalls. If you switch to MixedTypesNumPrelude instead of Prelude, comparison of MPBalls returns Maybe Bool instead of Bool."+ error "Failed to decide equality of MPBalls. If you switch to MixedTypesNumPrelude instead of Prelude, comparison of MPBalls returns Kleenean instead of Bool." instance P.Ord MPBall where compare a b = case (a < b, a == b, a > b) of- (Just True, _, _) -> P.LT- (_, Just True, _) -> P.EQ- (_, _, Just True) -> P.GT+ (CertainTrue, _, _) -> P.LT+ (_, CertainTrue, _) -> P.EQ+ (_, _, CertainTrue) -> P.GT _ ->- error "Failed to decide order of MPBalls. If you switch to MixedTypesNumPrelude instead of Prelude, comparison of MPBalls returns Maybe Bool instead of Bool."+ error "Failed to decide order of MPBalls. If you switch to MixedTypesNumPrelude instead of Prelude, comparison of MPBalls returns Kleenean instead of Bool." instance P.Num MPBall where fromInteger = convertExactly@@ -56,16 +58,16 @@ instance P.Fractional MPBall where fromRational = convertExactly . dyadic -- will work only for dyadic rationals- recip = (~!) . recip- (/) = (/!)+ recip = recip+ (/) = (/) instance P.Floating MPBall where- pi = error "MPBall: no pi :: MPBall, use pi ? (bitsS n) instead"- sqrt = (~!) . sqrt+ pi = error "There is no pi :: MPBall, use pi :: Real instead"+ sqrt = sqrt exp = exp sin = sin cos = cos- log = (~!) . log+ log = log atan = error "MPBall: atan not implemented yet" atanh = error "MPBall: atanh not implemented yet" asin = error "MPBall: asin not implemented yet"
src/AERN2/MP/Ball/Tests.hs view
@@ -1,4 +1,5 @@ {-# LANGUAGE CPP #-}+{-# OPTIONS_GHC -Wno-orphans #-} {-| Module : AERN2.MP.Ball.Tests Description : Tests for operations on arbitrary precision balls@@ -17,7 +18,6 @@ stack test aern2-mp --test-arguments "-a 1000 -m MPBall" @ -}- module AERN2.MP.Ball.Tests ( specMPBall, tMPBall@@ -69,11 +69,13 @@ tMPBall :: T MPBall tMPBall = T "MPBall" +-- tCNMPBall :: T (CN MPBall)+-- tCNMPBall = T "(CN MPBall)"+ specMPBall :: Spec specMPBall = describe ("MPBall") $ do specCanSetPrecision tMPBall (printArgsIfFails2 "`contains`" contains)- specCanReduceSizeUsingAccuracyGuide tMPBall specConversion tInteger tMPBall mpBall (fst . integerBounds) describe "order" $ do specHasEqNotMixed tMPBall
src/AERN2/MP/Ball/Type.hs view
@@ -1,4 +1,5 @@ {-# LANGUAGE TemplateHaskell #-}+{-# OPTIONS_GHC -Wno-orphans #-} {-| Module : AERN2.MP.Ball.Type Description : Arbitrary precision dyadic balls@@ -18,20 +19,25 @@ , module AERN2.MP.Accuracy , module AERN2.MP.Enclosure -- * The Ball type- , MPBall(..), CanBeMPBall, mpBall, CanBeMPBallP, mpBallP+ , MPBall(..), CanBeMPBall, mpBall, cnMPBall+ , CanBeMPBallP, mpBallP, cnMPBallP , reducePrecionIfInaccurate -- * Ball construction/extraction functions- , endpointsMP, fromEndpointsMP+ , fromMPFloatEndpoints+ , mpBallEndpoints, fromMPBallEndpoints ) where import MixedTypesNumPrelude -- import qualified Prelude as P -import Control.CollectErrors+import qualified Numeric.CollectErrors as CN import GHC.Generics (Generic)+import Control.DeepSeq +import qualified Data.List as List+ import Text.Printf import AERN2.Normalize@@ -56,20 +62,43 @@ -- } deriving (Generic) -instance Show MPBall- where- show b@(MPBall x _e) =- -- printf "[%s ± %s](prec=%s)" (show x) (showAC $ getAccuracy b) (show $ integer $ getPrecision b)- printf "[%s ± %s]" (showMPFloat x) (showAC $ getAccuracy b)- -- "[" ++ show x ++ " ± " ++ show e ++ "](prec=" ++ (show $ integer $ getPrecision x) ++ ")"- where- showAC Exact = "0"- showAC NoInformation = "oo"- showAC ac = "<2^(" ++ show (negate $ fromAccuracy ac) ++ ")"-+instance NFData MPBall -instance (SuitableForCE es) => CanEnsureCE es MPBall where+instance Show MPBall where+ show = showWithAccuracy (bits 50) +instance ShowWithAccuracy MPBall where+ showWithAccuracy displayAC b@(MPBall x e) =+ -- printf "[%s ± %s](prec=%s)" (show x) (showAC $ getAccuracy b) (show $ integer $ getPrecision b)+ printf "[%s ± %s%s]" (dropSomeDigits $ showMPFloat x) eDS (showAC $ getAccuracy b)+ -- "[" ++ show x ++ " ± " ++ show e ++ "](prec=" ++ (show $ integer $ getPrecision x) ++ ")"+ where+ eDS + | e == 0 = "0"+ | otherwise =+ case safeConvert (dyadic e) of+ Right (eD :: Double) -> printf "~%.4g" $ eD+ _ -> ""+ dropSomeDigits s =+ case List.findIndex (== '.') s of+ Nothing -> s+ Just ix -> withDotIx ix+ where+ withDotIx ix =+ let maxLength = ix + displayAC_n in+ let sTrimmed = take maxLength s in+ if length sTrimmed < maxLength+ then sTrimmed+ else (take (maxLength - 3) sTrimmed) <> "..."+ displayAC_n = + case displayAC of+ Exact -> 1000000000+ NoInformation -> 0+ _ -> integer $ ac2prec displayAC+ showAC Exact = ""+ showAC NoInformation = "(oo)"+ showAC ac = " ~2^(" ++ show (negate $ fromAccuracy ac) ++ ")"+ -- instance CanTestValid MPBall where -- isValid = isFinite @@ -108,6 +137,13 @@ p_e_nb = prec $ max 2 (10 + nb + fromAccuracy bAcc) (NormBits nb) = bNorm +instance CanGiveUpIfVeryInaccurate MPBall where+ giveUpIfVeryInaccurate = (aux =<<)+ where+ aux b@(MPBall _ e)+ | e > 1000 = CN.noValueNumErrorCertain $ numErrorVeryInaccurate "MPBall" ""+ | otherwise = cn b+ instance CanTestContains MPBall MPBall where contains (MPBall xLarge eLarge) (MPBall xSmall eSmall) = xLargeDy - eLargeDy <= xSmallDy - eSmallDy@@ -134,42 +170,29 @@ {- ball construction/extraction functions -} -instance IsInterval MPBall MPFloat where+instance IsInterval MPBall where+ type IntervalEndpoint MPBall = MPFloat fromEndpoints l u | u < l = fromEndpoints u l | otherwise =- MPBall (mpFloat cDy) (errorBound $ mpFloat eDy)+ MPBall c (errorBound e) where- lDy = dyadic l- uDy = dyadic u- cDy = (lDy + uDy) * (dyadic 0.5)- eDy = (uDy - cDy) `max` (cDy - lDy)- endpoints (MPBall x e) = (mpFloat lDy, mpFloat uDy)+ c = (l +. u) *. (mpFloat $ dyadic 0.5)+ e = (u -^ c) `max` (c -^ l)+ endpoints (MPBall x e) = (l, u) where- xDy = dyadic x- eDy = dyadic e- lDy = xDy - eDy- uDy = xDy + eDy+ eFl = mpFloat e+ l = x -. eFl+ u = x +^ eFl -fromEndpointsMP :: MPFloat -> MPFloat -> MPBall-fromEndpointsMP = fromEndpoints+fromMPFloatEndpoints :: MPFloat -> MPFloat -> MPBall+fromMPFloatEndpoints = fromEndpoints -endpointsMP :: MPBall -> (MPFloat, MPFloat)-endpointsMP = endpoints+fromMPBallEndpoints :: MPBall -> MPBall -> MPBall+fromMPBallEndpoints = fromEndpointsAsIntervals -instance IsInterval MPBall MPBall where- fromEndpoints l r = -- works as union even when r < l- fromEndpointsMP lMP uMP- where- lMP = min llMP rlMP- uMP = max luMP ruMP- (llMP, luMP) = endpointsMP l- (rlMP, ruMP) = endpointsMP r- endpoints x = (l,u)- where- l = MPBall lMP (errorBound 0)- u = MPBall uMP (errorBound 0)- (lMP, uMP) = endpointsMP x+mpBallEndpoints :: MPBall -> (MPBall, MPBall)+mpBallEndpoints = endpointsAsIntervals instance IsBall MPBall where type CentreType MPBall = Dyadic@@ -188,6 +211,8 @@ mpBallP :: (CanBeMPBallP t) => Precision -> t -> MPBall mpBallP = convertP +cnMPBallP :: (CanBeMPBallP a) => Precision -> CN a -> CN MPBall+cnMPBallP p = fmap (mpBallP p) {--- constructing an exact ball ---} @@ -196,41 +221,19 @@ mpBall :: (CanBeMPBall t) => t -> MPBall mpBall = convertExactly +cnMPBall :: (CanBeMPBall a) => CN a -> CN MPBall+cnMPBall = fmap mpBall++ {-- extracting approximate information about a ball --} instance HasAccuracy MPBall where getAccuracy = getAccuracy . ball_error -instance CanReduceSizeUsingAccuracyGuide MPBall where- reduceSizeUsingAccuracyGuide acGuide b@(MPBall x _e) =- case acGuide of- NoInformation -> lowerPrecisionIfAbove (prec 2) b- _ | getAccuracy b > acGuide -> tryPrec newPrec- _ -> b- where- tryPrec p- | getAccuracy bP >= acGuide = bP- | otherwise = tryPrec (p + 10)- where- bP = lowerPrecisionIfAbove p b- queryBits = 1 + fromAccuracy acGuide- newPrec =- case (getNormLog x) of- NormBits xNormBits ->- prec (max 2 (queryBits + xNormBits + 2))- NormZero ->- prec $ max 2 queryBits- -- bWithLowAC =- -- case acGuide of- -- Exact -> b- -- NoInformation -> b- -- _ -> normalize $- -- MPBall x (errorBound ((0.5^(fromAccuracy acGuide))⚡))- instance HasNorm MPBall where getNormLog ball = getNormLog boundMP where- (_, MPBall boundMP _) = endpoints $ absRaw ball+ (_, MPBall boundMP _) = mpBallEndpoints $ absRaw ball instance HasApproximate MPBall where type Approximate MPBall = (MPFloat, Bool)@@ -267,8 +270,8 @@ absRaw :: MPBall -> MPBall absRaw b | l < 0 && 0 < r =- fromEndpointsMP (mpFloat 0) (max (-l) r)+ fromEndpoints (mpFloat 0) (max (-l) r) | 0 <= l = b | otherwise = -b where- (l,r) = endpointsMP b+ (l,r) = endpoints b
src/AERN2/MP/Dyadic.hs view
@@ -1,6 +1,6 @@ {-# LANGUAGE CPP #-} -- #define DEBUG-{-# LANGUAGE DeriveDataTypeable #-}+{-# LANGUAGE GeneralizedNewtypeDeriving #-} {-| Module : AERN2.MP.Dyadic Description : Dyadics with exact ring operations@@ -16,7 +16,6 @@ Currently, we use hmpfr when compiling with ghc 7.10 and higher and haskell-mpfr when compiling with ghc 7.8. -}- module AERN2.MP.Dyadic ( -- * Dyadic numbers and their basic operations@@ -40,13 +39,15 @@ import MixedTypesNumPrelude import qualified Prelude as P -import Control.CollectErrors+import Control.CollectErrors (CollectErrors(..), CanBeErrors)+import qualified Control.CollectErrors as CE+-- import qualified Numeric.CollectErrors as CN import Text.Printf import Text.Regex.TDFA import Data.Typeable-import Data.Convertible+-- import Data.Convertible import Test.Hspec import Test.QuickCheck@@ -75,6 +76,7 @@ instance OrderedCertainlyRing (CN Dyadic) instance HasAccuracy Dyadic where getAccuracy _ = Exact+instance CanGiveUpIfVeryInaccurate Dyadic -- ie, never give up instance Show Dyadic where show (Dyadic x)@@ -106,7 +108,7 @@ [nS,eS] -> case (reads nS, reads eS) of ([(n,"")],[(e,"")]) ->- [((n :: Integer)*(dyadic 0.5)^!(e :: Integer), afterS)]+ [((n :: Integer)*(dyadic 0.5)^(e :: Integer), afterS)] _ -> tryNext _ -> tryNext where@@ -114,9 +116,7 @@ dyadicS =~ "\\`dyadic \\(([-0-9]*)\\*0.5\\^([0-9]*)\\)" :: (String, String, String, [String]) -instance (SuitableForCE es) => CanEnsureCE es Dyadic - {-- conversions --} type HasDyadics t = ConvertibleExactly Dyadic t@@ -151,7 +151,7 @@ | isDyadic = Right $ Dyadic (ceduCentre $ fromRationalCEDU (prec $ max 2 (dp + np + 1)) q) | otherwise = convError "this number is not dyadic" q where- isDyadic = d == 2^!dp+ isDyadic = d == 2^dp dp = integerLog2 d d = denominator q np = integerLog2 (max 1 $ abs $ numerator q)@@ -180,29 +180,25 @@ instance (HasEqAsymmetric Dyadic b- , CanEnsureCE es b- , CanEnsureCE es (EqCompareType Dyadic b)- , IsBool (EnsureCE es (EqCompareType Dyadic b))- , SuitableForCE es)+ , IsBool (CollectErrors es (EqCompareType Dyadic b))+ , CanBeErrors es) =>- HasEqAsymmetric Dyadic (CollectErrors es b)+ HasEqAsymmetric Dyadic (CollectErrors es b) where- type EqCompareType Dyadic (CollectErrors es b) =- EnsureCE es (EqCompareType Dyadic b)- equalTo = lift2TLCE equalTo+ type EqCompareType Dyadic (CollectErrors es b) =+ CollectErrors es (EqCompareType Dyadic b)+ equalTo = CE.liftT1 equalTo instance (HasEqAsymmetric a Dyadic- , CanEnsureCE es a- , CanEnsureCE es (EqCompareType a Dyadic)- , IsBool (EnsureCE es (EqCompareType a Dyadic))- , SuitableForCE es)+ , IsBool (CollectErrors es (EqCompareType a Dyadic))+ , CanBeErrors es) => HasEqAsymmetric (CollectErrors es a) Dyadic where type EqCompareType (CollectErrors es a) Dyadic =- EnsureCE es (EqCompareType a Dyadic)- equalTo = lift2TCE equalTo+ CollectErrors es (EqCompareType a Dyadic)+ equalTo = CE.lift1T equalTo instance CanTestZero Dyadic @@ -228,35 +224,31 @@ instance (HasOrderAsymmetric Dyadic b- , CanEnsureCE es b- , CanEnsureCE es (OrderCompareType Dyadic b)- , IsBool (EnsureCE es (OrderCompareType Dyadic b))- , SuitableForCE es)+ , IsBool (CollectErrors es (OrderCompareType Dyadic b))+ , CanBeErrors es) => HasOrderAsymmetric Dyadic (CollectErrors es b) where type OrderCompareType Dyadic (CollectErrors es b) =- EnsureCE es (OrderCompareType Dyadic b)- lessThan = lift2TLCE lessThan- leq = lift2TLCE leq- greaterThan = lift2TLCE greaterThan- geq = lift2TLCE geq+ CollectErrors es (OrderCompareType Dyadic b)+ lessThan = CE.liftT1 lessThan+ leq = CE.liftT1 leq+ greaterThan = CE.liftT1 greaterThan+ geq = CE.liftT1 geq instance (HasOrderAsymmetric a Dyadic- , CanEnsureCE es a- , CanEnsureCE es (OrderCompareType a Dyadic)- , IsBool (EnsureCE es (OrderCompareType a Dyadic))- , SuitableForCE es)+ , IsBool (CollectErrors es (OrderCompareType a Dyadic))+ , CanBeErrors es) => HasOrderAsymmetric (CollectErrors es a) Dyadic where type OrderCompareType (CollectErrors es a) Dyadic =- EnsureCE es (OrderCompareType a Dyadic)- lessThan = lift2TCE lessThan- leq = lift2TCE leq- greaterThan = lift2TCE greaterThan- geq = lift2TCE geq+ CollectErrors es (OrderCompareType a Dyadic)+ lessThan = CE.lift1T lessThan+ leq = CE.lift1T leq+ greaterThan = CE.lift1T greaterThan+ geq = CE.lift1T geq instance CanTestPosNeg Dyadic@@ -306,29 +298,25 @@ instance (CanMinMaxAsymmetric Dyadic b- , CanEnsureCE es b- , CanEnsureCE es (MinMaxType Dyadic b)- , SuitableForCE es)+ , CanBeErrors es) => CanMinMaxAsymmetric Dyadic (CollectErrors es b) where type MinMaxType Dyadic (CollectErrors es b) =- EnsureCE es (MinMaxType Dyadic b)- min = lift2TLCE min- max = lift2TLCE max+ CollectErrors es (MinMaxType Dyadic b)+ min = CE.liftT1 min+ max = CE.liftT1 max instance (CanMinMaxAsymmetric a Dyadic- , CanEnsureCE es a- , CanEnsureCE es (MinMaxType a Dyadic)- , SuitableForCE es)+ , CanBeErrors es) => CanMinMaxAsymmetric (CollectErrors es a) Dyadic where type MinMaxType (CollectErrors es a) Dyadic =- EnsureCE es (MinMaxType a Dyadic)- min = lift2TCE min- max = lift2TCE max+ CollectErrors es (MinMaxType a Dyadic)+ min = CE.lift1T min+ max = CE.lift1T max {- addition -} @@ -358,27 +346,23 @@ instance (CanAddAsymmetric Dyadic b- , CanEnsureCE es b- , CanEnsureCE es (AddType Dyadic b)- , SuitableForCE es)+ , CanBeErrors es) => CanAddAsymmetric Dyadic (CollectErrors es b) where type AddType Dyadic (CollectErrors es b) =- EnsureCE es (AddType Dyadic b)- add = lift2TLCE add+ CollectErrors es (AddType Dyadic b)+ add = CE.liftT1 add instance (CanAddAsymmetric a Dyadic- , CanEnsureCE es a- , CanEnsureCE es (AddType a Dyadic)- , SuitableForCE es)+ , CanBeErrors es) => CanAddAsymmetric (CollectErrors es a) Dyadic where type AddType (CollectErrors es a) Dyadic =- EnsureCE es (AddType a Dyadic)- add = lift2TCE add+ CollectErrors es (AddType a Dyadic)+ add = CE.lift1T add {- subtraction -} @@ -408,27 +392,23 @@ instance (CanSub Dyadic b- , CanEnsureCE es b- , CanEnsureCE es (SubType Dyadic b)- , SuitableForCE es)+ , CanBeErrors es) => CanSub Dyadic (CollectErrors es b) where type SubType Dyadic (CollectErrors es b) =- EnsureCE es (SubType Dyadic b)- sub = lift2TLCE sub+ CollectErrors es (SubType Dyadic b)+ sub = CE.liftT1 sub instance (CanSub a Dyadic- , CanEnsureCE es a- , CanEnsureCE es (SubType a Dyadic)- , SuitableForCE es)+ , CanBeErrors es) => CanSub (CollectErrors es a) Dyadic where type SubType (CollectErrors es a) Dyadic =- EnsureCE es (SubType a Dyadic)- sub = lift2TCE sub+ CollectErrors es (SubType a Dyadic)+ sub = CE.lift1T sub {- multiplication -}@@ -459,119 +439,89 @@ instance (CanMulAsymmetric Dyadic b- , CanEnsureCE es b- , CanEnsureCE es (MulType Dyadic b)- , SuitableForCE es)+ , CanBeErrors es) => CanMulAsymmetric Dyadic (CollectErrors es b) where type MulType Dyadic (CollectErrors es b) =- EnsureCE es (MulType Dyadic b)- mul = lift2TLCE mul+ CollectErrors es (MulType Dyadic b)+ mul = CE.liftT1 mul instance (CanMulAsymmetric a Dyadic- , CanEnsureCE es a- , CanEnsureCE es (MulType a Dyadic)- , SuitableForCE es)+ , CanBeErrors es) => CanMulAsymmetric (CollectErrors es a) Dyadic where type MulType (CollectErrors es a) Dyadic =- EnsureCE es (MulType a Dyadic)- mul = lift2TCE mul+ CollectErrors es (MulType a Dyadic)+ mul = CE.lift1T mul instance CanPow Dyadic Integer where- powNoCN = powUsingMul (dyadic 1)- pow = integerPowCN (powUsingMul (dyadic 1))+ pow = powUsingMul (dyadic 1) instance CanPow Dyadic Int where- powNoCN = powUsingMul (dyadic 1)- pow = integerPowCN (powUsingMul (dyadic 1))+ pow = powUsingMul (dyadic 1) instance (CanDiv a Dyadic- , CanEnsureCE es a- , CanEnsureCE es (DivType a Dyadic)- , CanEnsureCE es (DivTypeNoCN a Dyadic)- , SuitableForCE es)+ , CanBeErrors es) => CanDiv (CollectErrors es a) Dyadic where type DivType (CollectErrors es a) Dyadic =- EnsureCE es (DivType a Dyadic)- divide = lift2TCE divide- type DivTypeNoCN (CollectErrors es a) Dyadic =- EnsureCE es (DivTypeNoCN a Dyadic)- divideNoCN = lift2TCE divideNoCN+ CollectErrors es (DivType a Dyadic)+ divide = CE.lift1T divide instance CanDiv Integer Dyadic where- type DivTypeNoCN Integer Dyadic = Rational- divideNoCN a b = divideNoCN a (rational b)+ type DivType Integer Dyadic = Rational+ divide a b = divide a (rational b) instance CanDiv Dyadic Integer where- type DivTypeNoCN Dyadic Integer = Rational- divideNoCN a b = divideNoCN (rational a) b+ type DivType Dyadic Integer = Rational+ divide a b = divide (rational a) b instance CanDiv Int Dyadic where- type DivTypeNoCN Int Dyadic = Rational- divideNoCN a b = divideNoCN a (rational b)+ type DivType Int Dyadic = Rational+ divide a b = divide a (rational b) instance CanDiv Dyadic Int where- type DivTypeNoCN Dyadic Int = Rational- divideNoCN a b = divideNoCN (rational a) b+ type DivType Dyadic Int = Rational+ divide a b = divide (rational a) b instance CanDiv Rational Dyadic where- type DivTypeNoCN Rational Dyadic = Rational- divideNoCN = convertSecond divideNoCN+ type DivType Rational Dyadic = Rational+ divide = convertSecond divide instance CanDiv Dyadic Rational where- type DivTypeNoCN Dyadic Rational = Rational- divideNoCN = convertFirst divideNoCN+ type DivType Dyadic Rational = Rational+ divide = convertFirst divide instance (CanDiv Dyadic b- , CanEnsureCE es b- , CanEnsureCE es (DivType Dyadic b)- , CanEnsureCE es (DivTypeNoCN Dyadic b)- , SuitableForCE es)+ , CanBeErrors es) =>- CanDiv Dyadic (CollectErrors es b)+ CanDiv Dyadic (CollectErrors es b) where- type DivType Dyadic (CollectErrors es b) =- EnsureCE es (DivType Dyadic b)- divide = lift2TLCE divide- type DivTypeNoCN Dyadic (CollectErrors es b) =- EnsureCE es (DivTypeNoCN Dyadic b)- divideNoCN = lift2TLCE divideNoCN+ type DivType Dyadic (CollectErrors es b) =+ CollectErrors es (DivType Dyadic b)+ divide = CE.liftT1 divide instance (CanPow Dyadic b- , CanEnsureCE es b- , CanEnsureCE es (PowTypeNoCN Dyadic b)- , CanEnsureCE es (PowType Dyadic b)- , SuitableForCE es)+ , CanBeErrors es) => CanPow Dyadic (CollectErrors es b) where- type PowTypeNoCN Dyadic (CollectErrors es b) =- EnsureCE es (PowTypeNoCN Dyadic b)- powNoCN = lift2TLCE powNoCN type PowType Dyadic (CollectErrors es b) =- EnsureCE es (PowType Dyadic b)- pow = lift2TLCE pow+ CollectErrors es (PowType Dyadic b)+ pow = CE.liftT1 pow instance (CanPow a Dyadic- , CanEnsureCE es a- , CanEnsureCE es (PowType a Dyadic)- , CanEnsureCE es (PowTypeNoCN a Dyadic)- , SuitableForCE es)+ , CanBeErrors es) => CanPow (CollectErrors es a) Dyadic where- type PowTypeNoCN (CollectErrors es a) Dyadic =- EnsureCE es (PowTypeNoCN a Dyadic)- powNoCN = lift2TCE powNoCN type PowType (CollectErrors es a) Dyadic =- EnsureCE es (PowType a Dyadic)- pow = lift2TCE pow+ CollectErrors es (PowType a Dyadic)+ pow = CE.lift1T pow instance CanTestFinite Dyadic where isFinite = isFinite . dyadicMPFloat
src/AERN2/MP/Enclosure.hs view
@@ -1,4 +1,4 @@-{-# LANGUAGE Arrows #-}+{-# OPTIONS_GHC -Wno-orphans #-} {-| Module : AERN2.MP.Enclosure Description : Enclosure operations@@ -14,11 +14,16 @@ module AERN2.MP.Enclosure ( IsBall(..), ballFunctionUsingLipschitz- , IsInterval(..), intervalFunctionByEndpoints, intervalFunctionByEndpointsUpDown+ , IsInterval(..), endpointL, endpointR+ , fromEndpointsAsIntervals, endpointsAsIntervals, endpointLAsInterval, endpointRAsInterval+ , intervalFunctionByEndpoints, intervalFunctionByEndpointsUpDown+ , CanPlusMinus(..), (+-) , CanTestContains(..), CanMapInside(..), specCanMapInside , CanIntersectAsymmetric(..), CanIntersect+ , CanIntersectBy, CanIntersectSameType , CanIntersectCNBy, CanIntersectCNSameType , CanUnionAsymmetric(..), CanUnion+ , CanUnionBy, CanUnionSameType , CanUnionCNBy, CanUnionCNSameType ) where@@ -26,15 +31,15 @@ import MixedTypesNumPrelude -- import qualified Prelude as P -import Control.Arrow+-- import Control.Arrow import Test.Hspec import Test.QuickCheck -import Control.CollectErrors+import qualified Numeric.CollectErrors as CN import AERN2.MP.ErrorBound-import AERN2.MP.Accuracy+-- import AERN2.MP.Accuracy {- ball-specific operations -} @@ -52,12 +57,18 @@ to the explicit radius, making the centre exact. This may lose some information, but as a ball is equivalent to the original. For MPBall this function is pointless because it is equivalent to the identity. -}- makeExactCentre :: (IsBall t) => t -> t+ makeExactCentre :: t -> t makeExactCentre v = updateRadius (+r) c where (c, r) = centreAsBallAndRadius v +instance (IsBall t => IsBall (CN t)) where+ type CentreType (CN t) = CN (CentreType t)+ centre = fmap centre+ updateRadius f = fmap (updateRadius f)+ centreAsBallAndRadius = error $ "centreAsBallAndRadius not defined for CN types"+ {-| Computes a ball function @f@ on the centre and updating the error bound using a Lipschitz constant. -}@@ -75,33 +86,84 @@ {- interval-specific operations -}+class IsInterval i where+ type IntervalEndpoint i+ endpoints :: i -> (IntervalEndpoint i, IntervalEndpoint i)+ fromEndpoints :: IntervalEndpoint i -> IntervalEndpoint i -> i -class IsInterval i e where- fromEndpoints :: e -> e -> i- endpoints :: i -> (e,e)+instance (IsInterval t) => (IsInterval (CN t)) where+ type (IntervalEndpoint (CN t)) = CN (IntervalEndpoint t)+ fromEndpoints l u = CN.lift2 fromEndpoints l u+ endpoints = CN.liftPair endpoints +endpointL :: (IsInterval i) => i -> IntervalEndpoint i+endpointL = fst . endpoints++endpointR :: (IsInterval i) => i -> IntervalEndpoint i+endpointR = snd . endpoints++endpointsAsIntervals :: + (IsInterval i) => i -> (i,i)+endpointsAsIntervals x = (lI,rI)+ where+ lI = fromEndpoints l l+ rI = fromEndpoints r r+ (l,r) = endpoints x++endpointLAsInterval :: (IsInterval i) => i -> i+endpointLAsInterval = fst . endpointsAsIntervals++endpointRAsInterval :: (IsInterval i) => i -> i+endpointRAsInterval = snd . endpointsAsIntervals++fromEndpointsAsIntervals :: + (IsInterval i, CanMinMaxSameType (IntervalEndpoint i)) =>+ i -> i -> i+fromEndpointsAsIntervals l r = + fromEndpoints lMP uMP+ where+ lMP = min llMP rlMP+ uMP = max luMP ruMP+ (llMP, luMP) = endpoints l+ (rlMP, ruMP) = endpoints r++{- plusMinus (+-) operator -}++class CanPlusMinus t1 t2 where+ type PlusMinusType t1 t2+ type PlusMinusType t1 t2 = t1+ {-| Operator for constructing or enlarging enclosures such as balls or intervals -}+ plusMinus :: t1 -> t2 -> PlusMinusType t1 t2++infixl 6 +-++{-| Operator for constructing or enlarging enclosures such as balls or intervals -}+(+-) :: (CanPlusMinus t1 t2) => t1 -> t2 -> PlusMinusType t1 t2+(+-) = plusMinus++ {-| Computes a *monotone* ball function @f@ on intervals using the interval endpoints. -} intervalFunctionByEndpoints ::- (IsInterval t t, HasEqCertainly t t)+ (IsInterval t, CanMinMaxSameType (IntervalEndpoint t), HasEqCertainly t t) => (t -> t) {-^ @fThin@: a version of @f@ that works well on thin intervals -} -> (t -> t) {-^ @f@ on *large* intervals -} intervalFunctionByEndpoints fThin x | l !==! u = fThin l- | otherwise = fromEndpoints (fThin l) (fThin u)+ | otherwise = fromEndpointsAsIntervals (fThin l) (fThin u) where- (l,u) = endpoints x+ (l,u) = endpointsAsIntervals x {-| Computes a *monotone* ball function @f@ on intervals using the interval endpoints. -} intervalFunctionByEndpointsUpDown ::- (IsInterval t e)+ (IsInterval t) =>- (e -> e) {-^ @fDown@: a version of @f@ working on endpoints, rounded down -} ->- (e -> e) {-^ @fUp@: a version of @f@ working on endpoints, rounded up -} ->+ (IntervalEndpoint t -> IntervalEndpoint t) {-^ @fDown@: a version of @f@ working on endpoints, rounded down -} ->+ (IntervalEndpoint t -> IntervalEndpoint t) {-^ @fUp@: a version of @f@ working on endpoints, rounded up -} -> (t -> t) {-^ @f@ on intervals rounding *outwards* -} intervalFunctionByEndpointsUpDown fDown fUp x = fromEndpoints (fDown l) (fUp u)@@ -144,170 +206,135 @@ (CanIntersectAsymmetric e1 e2, CanIntersectAsymmetric e1 e2 , IntersectionType e1 e2 ~ IntersectionType e2 e1) +{-| A set intersection (usually partial) -} class CanIntersectAsymmetric e1 e2 where type IntersectionType e1 e2- type IntersectionType e1 e2 = EnsureCN e1+ type IntersectionType e1 e2 = CN e1 intersect :: e1 -> e2 -> IntersectionType e1 e2 +type CanIntersectBy e1 e2 =+ (CanIntersect e1 e2, IntersectionType e1 e2 ~ e1)++type CanIntersectSameType e1 =+ (CanIntersectBy e1 e1)+ type CanIntersectCNBy e1 e2 =- (CanIntersect e1 e2, IntersectionType e1 e2 ~ EnsureCN e1- , CanIntersect (EnsureCN e1) e2, IntersectionType (EnsureCN e1) e2 ~ EnsureCN e1)+ (CanIntersect e1 e2, IntersectionType e1 e2 ~ CN e1)+ type CanIntersectCNSameType e1 =- (CanIntersectCNBy e1 e1- , CanIntersect (EnsureCN e1) (EnsureCN e1), IntersectionType (EnsureCN e1) (EnsureCN e1) ~ EnsureCN e1)+ (CanIntersectCNBy e1 e1) -instance CanIntersectAsymmetric Bool Bool where+instance+ CanIntersectAsymmetric Bool Bool+ where intersect b1 b2 | b1 == b2 = cn b1 | otherwise =- noValueNumErrorCertainCN $ NumError "empty Boolean intersection"+ CN.noValueNumErrorCertain $ CN.NumError "empty Boolean intersection" instance- (CanIntersectAsymmetric Bool b- , CanEnsureCE es b- , CanEnsureCE es (IntersectionType Bool b)- , SuitableForCE es)- =>- CanIntersectAsymmetric Bool (CollectErrors es b)+ CanIntersectAsymmetric Kleenean Kleenean where- type IntersectionType Bool (CollectErrors es b) =- EnsureCE es (IntersectionType Bool b)- intersect = lift2TLCE intersect+ intersect CertainTrue CertainFalse =+ CN.noValueNumErrorCertain $ CN.NumError "empty Kleenean intersection"+ intersect CertainFalse CertainTrue =+ CN.noValueNumErrorCertain $ CN.NumError "empty Kleenean intersection"+ intersect TrueOrFalse k2 = cn k2+ intersect k1 _ = cn k1 -instance- (CanIntersectAsymmetric a Bool- , CanEnsureCE es a- , CanEnsureCE es (IntersectionType a Bool)- , SuitableForCE es)+instance + (CanIntersectAsymmetric a b, IntersectionType a b ~ CN c) =>- CanIntersectAsymmetric (CollectErrors es a) Bool+ CanIntersectAsymmetric (CN a) (CN b) where- type IntersectionType (CollectErrors es a) Bool =- EnsureCE es (IntersectionType a Bool)- intersect = lift2TCE intersect+ type IntersectionType (CN a) (CN b) = IntersectionType a b+ intersect = CN.lift2CN intersect instance- (CanIntersectAsymmetric (Maybe a) b- , CanEnsureCE es b- , CanEnsureCE es (IntersectionType (Maybe a) b)- , SuitableForCE es)+ (CanIntersectAsymmetric (CN Bool) (CN b)) =>- CanIntersectAsymmetric (Maybe a) (CollectErrors es b)+ CanIntersectAsymmetric Bool (CN b) where- type IntersectionType (Maybe a) (CollectErrors es b) =- EnsureCE es (IntersectionType (Maybe a) b)- intersect = lift2TLCE intersect+ type IntersectionType Bool (CN b) = IntersectionType (CN Bool) (CN b)+ intersect b1 = intersect (cn b1) instance- (CanIntersectAsymmetric a (Maybe b)- , CanEnsureCE es a- , CanEnsureCE es (IntersectionType a (Maybe b))- , SuitableForCE es)+ (CanIntersectAsymmetric (CN a) (CN Bool)) =>- CanIntersectAsymmetric (CollectErrors es a) (Maybe b)+ CanIntersectAsymmetric (CN a) Bool where- type IntersectionType (CollectErrors es a) (Maybe b) =- EnsureCE es (IntersectionType a (Maybe b))- intersect = lift2TCE intersect-+ type IntersectionType (CN a) Bool = IntersectionType (CN a) (CN Bool)+ intersect b1 b2 = intersect b1 (cn b2) instance- (CanIntersectAsymmetric a b- , CanEnsureCN a, IntersectionType a b ~ EnsureCN a- , CanEnsureCN (EnsureCN a)- , CanEnsureCN b, EnsureCN b ~ EnsureCN a)+ (CanIntersectAsymmetric (CN Kleenean) (CN b)) =>- CanIntersectAsymmetric (Maybe a) (Maybe b)+ CanIntersectAsymmetric Kleenean (CN b) where- type IntersectionType (Maybe a) (Maybe b) = EnsureCN (Maybe (IntersectionType a b))- intersect (ma :: Maybe a) (mb :: Maybe b) =- case (ma, mb) of- (Just a, Just b) -> justCN sample_r (intersect a b)- (Just a, Nothing) -> justCN sample_r (ensureCN a)- (Nothing, Just b) -> justCN sample_r (ensureCN b)- _ -> cn (Nothing :: Maybe a)- where- sample_r = Nothing :: EnsureCN (Maybe (IntersectionType a b))--justCN :: (CanEnsureCN a) => Maybe a -> EnsureCN a -> EnsureCN (Maybe a)-justCN (_sample_a :: Maybe a) aCN =- case deEnsureCN aCN of- Right a -> cn (Just (a :: a))- _ -> cn (Nothing :: Maybe a)-+ type IntersectionType Kleenean (CN b) = IntersectionType (CN Kleenean) (CN b)+ intersect k1 = intersect (cn k1) --- --- Version that removes inner CN:--- instance--- (CanIntersectCNSameType a, CanEnsureCN a)--- =>--- CanIntersectAsymmetric (Maybe a) (Maybe a)--- where--- type IntersectionType (Maybe a) (Maybe a) = CN (Maybe (WithoutCN (IntersectionType a a)))--- intersect (Just a) (Just b) = fmap Just (intersect a b)--- intersect (Just a) Nothing = fmap Just (ensureCN a)--- intersect Nothing (Just b) = fmap Just (ensureCN b)--- intersect Nothing Nothing = cn Nothing--- instance- (CanIntersectAsymmetric e1 e2, SuitableForCE es- , CanEnsureCE es e1, CanEnsureCE es e2- , CanEnsureCE es (IntersectionType e1 e2))+ (CanIntersectAsymmetric (CN a) (CN Kleenean)) =>- CanIntersectAsymmetric (CollectErrors es e1) (CollectErrors es e2)+ CanIntersectAsymmetric (CN a) Kleenean where- type IntersectionType (CollectErrors es e1) (CollectErrors es e2) =- EnsureCE es (IntersectionType e1 e2)- intersect = lift2CE intersect+ type IntersectionType (CN a) Kleenean = IntersectionType (CN a) (CN Kleenean)+ intersect k1 k2 = intersect k1 (cn k2) -{- union -}+{- set union -} type CanUnion e1 e2 = (CanUnionAsymmetric e1 e2, CanUnionAsymmetric e1 e2 , UnionType e1 e2 ~ UnionType e2 e1) +{-| A set union (usually partial) -} class CanUnionAsymmetric e1 e2 where type UnionType e1 e2- type UnionType e1 e2 = EnsureCN e1+ type UnionType e1 e2 = CN e1 union :: e1 -> e2 -> UnionType e1 e2 +type CanUnionBy e1 e2 =+ (CanUnion e1 e2, UnionType e1 e2 ~ e1)++type CanUnionSameType e1 =+ (CanUnionBy e1 e1)+ type CanUnionCNBy e1 e2 =- (CanUnion e1 e2, UnionType e1 e2 ~ EnsureCN e1- , CanUnion (EnsureCN e1) e2, UnionType (EnsureCN e1) e2 ~ EnsureCN e1)+ (CanUnion e1 e2, UnionType e1 e2 ~ CN e1) type CanUnionCNSameType e1 =- (CanUnionCNBy e1 e1- , CanUnion (EnsureCN e1) (EnsureCN e1), UnionType (EnsureCN e1) (EnsureCN e1) ~ EnsureCN e1)+ (CanUnionCNBy e1 e1) -instance- (CanUnionAsymmetric e1 e2- , CanEnsureCN e1, CanEnsureCN e2- , CanEnsureCN (UnionType e1 e2))+instance + (CanUnionAsymmetric a b, UnionType a b ~ CN c) =>- CanUnionAsymmetric (CN e1) (CN e2)- -- a more general `CollectErrors` instance conflicts with the Arrow instance below+ CanUnionAsymmetric (CN a) (CN b) where- type UnionType (CN e1) (CN e2) =- EnsureCN (UnionType e1 e2)- union = lift2CE union+ type UnionType (CN a) (CN b) = UnionType a b+ union = CN.lift2CN union -instance- (Arrow to, CanUnionAsymmetric e1 e2)- =>- CanUnionAsymmetric (to Accuracy e1) (to Accuracy e2)- -- this instance is important for "parallel if"- where- type UnionType (to Accuracy e1) (to Accuracy e2) =- to Accuracy (UnionType e1 e2)- union xA yA =- proc ac ->- do- x <- xA -< ac- y <- yA -< ac- returnA -< union x y+-- TODO: move to aern2-real (or aern2-net)+-- instance+-- (Arrow to, CanUnionAsymmetric e1 e2)+-- =>+-- CanUnionAsymmetric (to Accuracy e1) (to Accuracy e2)+-- -- this instance is important for "parallel if"+-- where+-- type UnionType (to Accuracy e1) (to Accuracy e2) =+-- to Accuracy (UnionType e1 e2)+-- union xA yA =+-- proc ac ->+-- do+-- x <- xA -< ac+-- y <- yA -< ac+-- returnA -< union x y -instance (CanUnionCNSameType t, CanEnsureCN t) =>- HasIfThenElse (Maybe Bool) t+instance (CanUnionSameType t, CN.CanTakeCNErrors t) =>+ HasIfThenElse Kleenean t where- type IfThenElseType (Maybe Bool) t = EnsureCN t- ifThenElse (Just b) e1 e2 = cn $ if b then e1 else e2- ifThenElse Nothing e1 e2 = e1 `union` e2+ type IfThenElseType Kleenean t = t+ ifThenElse CertainTrue e1 _ = e1+ ifThenElse CertainFalse _ e2 = e2+ ifThenElse TrueOrFalse e1 e2 = e1 `union` e2+
src/AERN2/MP/ErrorBound.hs view
@@ -22,11 +22,15 @@ import qualified Prelude as P import Data.Typeable+import GHC.Generics (Generic)+import Control.DeepSeq import Test.QuickCheck -import Data.Convertible+-- import Data.Convertible +import Data.Ratio (numerator)+ import Math.NumberTheory.Logarithms (integerLog2) import AERN2.MP.Precision@@ -45,8 +49,10 @@ {-| A non-negative Double value to serve as an error bound. Arithmetic is rounded towards +infinity. -} newtype ErrorBound = ErrorBound { er2mp :: MPFloat }- deriving (P.Eq, P.Ord, Typeable)+ deriving (P.Eq, P.Ord, Typeable, Generic) +instance NFData ErrorBound+ instance Show ErrorBound where show (ErrorBound d) = show d @@ -64,9 +70,16 @@ where eN = floor $ rational e eRecipN = ceiling $ rational $ one /. e+ getFiniteAccuracy eb@(ErrorBound e) + | e == 0 = bits errorBoundPrecision+ | otherwise = getAccuracy eb + {- conversions -} +instance ConvertibleExactly ErrorBound ErrorBound where+ safeConvertExactly = Right+ instance ConvertibleExactly ErrorBound MPFloat where safeConvertExactly = Right . er2mp @@ -197,12 +210,11 @@ | otherwise = error "trying to multiply ErrorBound by a negative integer" instance CanDiv ErrorBound Integer where- type DivTypeNoCN ErrorBound Integer = ErrorBound type DivType ErrorBound Integer = ErrorBound- divideNoCN = divide- divide (ErrorBound a) i- | i > 0 = ErrorBound $ a /^ (MPFloat.fromIntegerUp errorBoundPrecision i)- | otherwise = error "trying to multiply ErrorBound by a non-positive integer"+ divide = divide+ -- divide (ErrorBound a) i+ -- | i > 0 = ErrorBound $ a /^ (MPFloat.fromIntegerUp errorBoundPrecision i)+ -- | otherwise = error "trying to multiply ErrorBound by a non-positive integer" instance Arbitrary ErrorBound where arbitrary =@@ -216,6 +228,6 @@ | otherwise = do (s :: Integer) <- arbitrary- let resultR = ((abs s) `P.mod` (2^!35))/!(2^!32)+ let resultR = ((abs s) `P.mod` (numerator $ 2^35))/(2^32) let result = convert resultR return result
+ src/AERN2/MP/Float/Arithmetic.hs view
@@ -0,0 +1,104 @@+{-# OPTIONS_GHC -Wno-orphans #-}+{-|+ Module : AERN2.MP.Float.Arithmetic+ Description : Arbitrary precision floating point numbers+ Copyright : (c) Michal Konecny+ License : BSD3++ Maintainer : mikkonecny@gmail.com+ Stability : experimental+ Portability : portable++ Arbitrary precision floating-point numbers with up/down-rounded operations.+-}+module AERN2.MP.Float.Arithmetic+ (+ -- * MPFloat basic arithmetic+ addCEDU, subCEDU+ , mulCEDU, divCEDU, recipCEDU+ -- * MPFloat selected constants and operations+ , piCEDU+ , cosCEDU, sinCEDU+ , sqrtCEDU, expCEDU, logCEDU+ )+where++import MixedTypesNumPrelude+import qualified Prelude as P++import AERN2.MP.Precision++import qualified Data.CDAR as MPLow++import AERN2.MP.Float.Auxi+import AERN2.MP.Float.Type++{- common functions -}++instance CanNeg MPFloat where+ negate = P.negate++instance CanAbs MPFloat where+ abs = P.abs++addCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat+addCEDU = binaryCEDU (P.+)++subCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat+subCEDU = binaryCEDU (P.-)++mulCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat+mulCEDU = binaryCEDU (P.*)++divCEDU :: MPFloat -> MPFloat -> BoundsCEDU MPFloat+divCEDU x y + | y P.== (P.fromInteger 0) = getBoundsCEDU MPLow.Bottom+ | otherwise = binaryCEDU (P./) x y++recipCEDU :: MPFloat -> BoundsCEDU MPFloat+recipCEDU = unaryCEDU P.recip++{- special constants and functions -}++piCEDU :: Precision -> BoundsCEDU MPFloat+piCEDU pp = + getBoundsCEDU $ MPLow.piA (p2cdarPrec pp)++cosCEDU :: MPFloat -> BoundsCEDU MPFloat+cosCEDU = unaryCEDU MPLow.cosA++sinCEDU :: MPFloat -> BoundsCEDU MPFloat+sinCEDU = unaryCEDU MPLow.sinA+ +sqrtCEDU :: MPFloat -> BoundsCEDU MPFloat+sqrtCEDU = unaryCEDU MPLow.sqrtA+ +expCEDU :: MPFloat -> BoundsCEDU MPFloat+expCEDU = unaryCEDU MPLow.expA++logCEDU :: MPFloat -> BoundsCEDU MPFloat+logCEDU = unaryCEDU MPLow.logA++{- auxiliary functions to automatically determine result precision from operand precisions -}++binaryCEDU ::+ (MPFloat -> MPFloat -> MPFloat) ->+ (MPFloat -> MPFloat -> BoundsCEDU MPFloat)+binaryCEDU op x y =+ getBoundsCEDU $ op x y++unaryCEDU ::+ (MPFloat -> MPFloat) ->+ (MPFloat -> BoundsCEDU MPFloat)+unaryCEDU op x =+ getBoundsCEDU $ op x++-- unaryPrecCEDU ::+-- Integer ->+-- (MPLow.Precision -> MPFloat -> MPFloat) ->+-- (MPFloat -> BoundsCEDU MPFloat)+-- unaryPrecCEDU addPrec op x@(MPLow.Approx mb _ _ s) =+-- getBoundsCEDU $ op ((-s P.+ mb) P.+ (int addPrec)) x+-- unaryPrecCEDU addPrec op MPLow.Bottom =+-- getBoundsCEDU $ op ((int $ integer defaultPrecision) P.+ (int addPrec)) MPLow.Bottom+
+ src/AERN2/MP/Float/Conversions.hs view
@@ -0,0 +1,152 @@+{-# OPTIONS_GHC -Wno-orphans #-}+{-|+ Module : AERN2.MP.Float.Conversions+ Description : Conversions and comparisons of arbitrary precision floats+ Copyright : (c) Michal Konecny+ License : BSD3++ Maintainer : mikkonecny@gmail.com+ Stability : experimental+ Portability : portable++ Conversions and comparisons of arbitrary precision floating point numbers+-}+module AERN2.MP.Float.Conversions+ (+ -- * MPFloat to other types (see also instances)+ toDouble+ -- * MPFloat constructors (see also instances)+ , CanBeMPFloat, mpFloat+ , fromIntegerCEDU+ , fromRationalCEDU+ -- * comparisons and constants (see also instances)+ , zero, one, two+ , nan, infinity+ )+where++import MixedTypesNumPrelude+import qualified Prelude as P++import Data.Ratio+-- import Data.Convertible++-- import AERN2.Norm+import AERN2.MP.Precision++import qualified Data.CDAR as MPLow++import AERN2.MP.Float.Auxi+import AERN2.MP.Float.Type+import AERN2.MP.Float.Arithmetic+++{- conversions to MPFloat -}++type CanBeMPFloat t = ConvertibleExactly t MPFloat+mpFloat :: (CanBeMPFloat t) => t -> MPFloat+mpFloat = convertExactly++instance ConvertibleExactly Integer MPFloat where+ safeConvertExactly =+ Right . P.fromInteger++instance ConvertibleExactly Int MPFloat where+ safeConvertExactly = safeConvertExactly . integer++fromIntegerCEDU :: Precision -> Integer -> BoundsCEDU MPFloat+fromIntegerCEDU pp =+ setPrecisionCEDU pp . P.fromInteger++fromRationalCEDU :: Precision -> Rational -> BoundsCEDU MPFloat+fromRationalCEDU pp =+ setPrecisionCEDU pp . (MPLow.toApproxMB (p2cdarPrec pp))++{- conversions from MPFloat -}++instance ConvertibleExactly MPFloat Rational where+ safeConvertExactly = Right . P.toRational+ +toDouble :: MPFloat -> Double+toDouble = P.fromRational . rational++instance Convertible MPFloat Double where+ safeConvert x+ | isFinite dbl = Right dbl+ | otherwise = convError "conversion to double: out of bounds" x+ where+ dbl = toDouble x+++instance CanRound MPFloat where+ properFraction x = (n,f)+ where+ r = rational x+ n = (numerator r) `P.quot` (denominator r)+ f = ceduCentre $ x `subCEDU` (P.fromInteger n)+ +{- comparisons -}++instance HasEqAsymmetric MPFloat MPFloat+instance HasEqAsymmetric MPFloat Integer where+ equalTo = convertSecond equalTo+instance HasEqAsymmetric Integer MPFloat where+ equalTo = convertFirst equalTo+instance HasEqAsymmetric MPFloat Int where+ equalTo = convertSecond equalTo+instance HasEqAsymmetric Int MPFloat where+ equalTo = convertFirst equalTo+instance HasEqAsymmetric MPFloat Rational where+ equalTo = convertFirst equalTo+instance HasEqAsymmetric Rational MPFloat where+ equalTo = convertSecond equalTo++instance CanTestZero MPFloat++instance HasOrderAsymmetric MPFloat MPFloat+instance HasOrderAsymmetric MPFloat Integer where+ lessThan = convertSecond lessThan+ leq = convertSecond leq+instance HasOrderAsymmetric Integer MPFloat where+ lessThan = convertFirst lessThan+ leq = convertFirst leq+instance HasOrderAsymmetric MPFloat Int where+ lessThan = convertSecond lessThan+ leq = convertSecond leq+instance HasOrderAsymmetric Int MPFloat where+ lessThan = convertFirst lessThan+ leq = convertFirst leq+instance HasOrderAsymmetric Rational MPFloat where+ lessThan = convertSecond lessThan+ leq = convertSecond leq+instance HasOrderAsymmetric MPFloat Rational where+ lessThan = convertFirst lessThan+ leq = convertFirst leq++instance CanTestPosNeg MPFloat++{- min, max -}++instance CanMinMaxAsymmetric MPFloat MPFloat++{- constants -}++zero, one, two :: MPFloat+zero = mpFloat 0+one = mpFloat 1+two = mpFloat 2++nan, infinity :: MPFloat+nan = MPLow.Bottom+infinity = nan++itisNaN :: MPFloat -> Bool+itisNaN MPLow.Bottom = True+itisNaN _ = False++instance CanTestFinite MPFloat where+ isInfinite = itisNaN+ isFinite = not . itisNaN++instance CanTestNaN MPFloat where+ isNaN = itisNaN
src/AERN2/MP/Float/Tests.hs view
@@ -1,3 +1,4 @@+{-# OPTIONS_GHC -Wno-orphans #-} {-| Module : AERN2.MP.Float.Tests Description : Tests for operations on arbitrary precision floats@@ -16,7 +17,6 @@ stack test aern2-mp --test-arguments "-a 1000 -m MPFloat" @ -}- module AERN2.MP.Float.Tests ( specMPFloat, tMPFloat@@ -36,14 +36,14 @@ import Test.QuickCheck -- import qualified Test.Hspec.SmallCheck as SC -import Control.CollectErrors+-- import qualified Numeric.CollectErrors as CN import AERN2.Norm import AERN2.MP.Precision import AERN2.MP.Float.Auxi import AERN2.MP.Float.Type-import AERN2.MP.Float.Arithmetic+-- import AERN2.MP.Float.Arithmetic import AERN2.MP.Float.Conversions import AERN2.MP.Float.Operators@@ -62,7 +62,7 @@ (p :: Precision) <- arbitrary (s :: Integer) <- arbitrary ex <- choose (-20,10)- let resultR = s * (10.0^!ex)+ let resultR = s * (10.0^ex) let result = ceduCentre $ fromRationalCEDU p resultR return result @@ -88,7 +88,7 @@ where l = mpFloat l_ u = mpFloat u_- b = l +^ ((abs a) `modNoCN` (u-^l))+ b = l +^ ((abs a) `mod` (u-^l)) enforceRangeMP (Just l_, _) a | isInfinite a = abs a | l < a = a@@ -105,21 +105,20 @@ u = mpFloat u_ enforceRangeMP _ a = a -instance CanEnsureCE NumErrors MPFloat- instance CanDivIMod MPFloat MPFloat where- divIMod x m - | (not (isFinite m)) = (errM (d :: Integer), errM xm)- | (not (isFinite x)) = (errX (d :: Integer), errX xm)- | m > zero = (cn d, cn xm)- | otherwise = (errM (d :: Integer), errM xm)+ type DivIType MPFloat MPFloat = Integer+ divIMod x m+ | (not (isFinite m)) = ((d :: Integer), xm)+ | (not (isFinite x)) = ((d :: Integer), xm)+ | m > zero = (d, xm)+ | otherwise = ((d :: Integer), xm) where d = floor (x /^ m) xm = x -^ (mpFloat d)*^m- errM :: (CanEnsureCN t) => t -> EnsureCN t- errM s = noValueNumErrorCertainECN (Just s) $ OutOfRange $ "modulus not finite and positive: " ++ show m- errX :: (CanEnsureCN t) => t -> EnsureCN t- errX s = noValueNumErrorCertainECN (Just s) $ OutOfRange $ "modulus input not finite: " ++ show x+ -- errM :: (CanEnsureCN t) => t -> EnsureCN t+ -- errM s = CN.noValueNumErrorCertain $ CN.OutOfDomain $ "modulus not finite and positive: " ++ show m+ -- errX :: (CanEnsureCN t) => t -> EnsureCN t+ -- errX s = CN.noValueNumErrorCertain$ CN.OutOfDomain $ "modulus input not finite: " ++ show x {- approximate comparison -}@@ -145,7 +144,7 @@ | isNaN x || isNaN y = False | isInfinite x || isInfinite y = x == y | otherwise =- abs (x -. y) <= 0.5^!e+ abs (x -. y) <= 0.5^e {-| Assert equality of two MPFloat's with tolerance derived from the size and precision@@ -492,14 +491,15 @@ cosDown x <=% 1 && cosUp x >=% -1- it "cos(pi)=-1" $ do- property $ \ (p :: Precision) ->- let- piA = ceduCentre $ piCEDU p- (=~~=) = approxEqualWithArgs 1 [(piA,"pi")]- infix 4 =~~=- in- cosUp(piA) =~~= (-one)+ -- TODO: fix accuracy of CDAR mBounds cosine+ -- it "cos(pi)=-1" $ do+ -- property $ \ (p :: Precision) ->+ -- let+ -- piA = ceduCentre $ piCEDU p+ -- (=~~=) = approxEqualWithArgs 2 [(piA,"pi")]+ -- infix 4 =~~=+ -- in+ -- cosUp(piA) =~~= (-one) it "cos(x)^2 + sin(x)^2 = 1" $ do property $ \ (x_ :: MPFloat) -> let x = enforceRangeMP (Just (-1000000), Just 1000000) x_ in
+ src/AERN2/MP/Float/Type.hs view
@@ -0,0 +1,79 @@+{-# OPTIONS_GHC -Wno-orphans #-}+{-# LANGUAGE StandaloneDeriving #-}+{-|+ Module : AERN2.MP.Float.Type+ Description : Arbitrary precision floating point numbers (via cdar)+ Copyright : (c) Michal Konecny+ License : BSD3++ Maintainer : mikkonecny@gmail.com+ Stability : experimental+ Portability : portable++ Arbitrary precision floating-point numbers, re-using CDAR Approx type.+-}+module AERN2.MP.Float.Type+ (+ -- * MPFloat numbers and their basic operations+ MPFloat+ , showMPFloat+ , setPrecisionCEDU+ , p2cdarPrec+ , getBoundsCEDU+ )+where++import MixedTypesNumPrelude+import qualified Prelude as P++-- import Data.Bits (unsafeShiftL)+import Data.Typeable++import AERN2.Norm+import AERN2.MP.Precision+import AERN2.MP.Float.Auxi++import qualified Data.CDAR as MPLow++{-| Multiple-precision floating-point type based on CDAR.Approx with 0 radius. -}+type MPFloat = MPLow.Approx++showMPFloat :: MPFloat -> String+showMPFloat x = MPLow.showA x++deriving instance (Typeable MPFloat)++p2cdarPrec :: Precision -> MPLow.Precision+p2cdarPrec = P.fromInteger . integer++getBoundsCEDU :: MPFloat -> BoundsCEDU MPFloat+getBoundsCEDU (MPLow.Approx mb m e s) = + BoundsCEDU + (MPLow.Approx mb m 0 s) (MPLow.approxMB eb_mb e 0 s)+ (MPLow.Approx mb (m-e) 0 s) (MPLow.Approx mb (m+e) 0 s)+getBoundsCEDU MPLow.Bottom =+ BoundsCEDU+ MPLow.Bottom MPLow.Bottom MPLow.Bottom MPLow.Bottom++{-| The bit-size bound for the error bound in CEDU -}+eb_prec :: Precision+eb_prec = prec 63++{-| The bit-size bound for the error bound in CEDU -}+eb_mb :: Int+eb_mb = int $ integer eb_prec++instance HasPrecision MPFloat where+ getPrecision (MPLow.Approx mb _ _ _) = prec (P.toInteger $ mb)+ getPrecision MPLow.Bottom = error "illegal MPFloat (Bottom)"+ ++instance CanSetPrecision MPFloat where+ setPrecision p = ceduCentre . setPrecisionCEDU p++setPrecisionCEDU :: Precision -> MPFloat -> BoundsCEDU MPFloat+setPrecisionCEDU pp = getBoundsCEDU . MPLow.enforceMB . MPLow.setMB (p2cdarPrec pp)++instance HasNorm MPFloat where+ getNormLog (MPLow.Approx _ m _ s) = (getNormLog m) + (integer s)+ getNormLog MPLow.Bottom = error "getNormLog undefined for Bottom"
src/AERN2/MP/Precision.hs view
@@ -1,8 +1,8 @@ {-# LANGUAGE TemplateHaskell #-}-{-# LANGUAGE DeriveDataTypeable #-}+{-# LANGUAGE GeneralizedNewtypeDeriving #-} {-| Module : AERN2.MP.Precision- Description : Floating-point precision+ Description : number of significant binary digits Copyright : (c) Michal Konecny License : BSD3 @@ -10,7 +10,8 @@ Stability : experimental Portability : portable - Floating-point precision type and its operations+ Precision type and its operations. Precision expresses a limit on the size+ of an approximation, roughly corresponding to the number of significant bits for a floating-point number and its generalisations to other types. -} module AERN2.MP.Precision (@@ -27,7 +28,8 @@ import qualified Prelude as P import Text.Printf -import Control.CollectErrors+import Control.CollectErrors (CollectErrors(..), CanBeErrors)+-- import qualified Control.CollectErrors as CE import Data.Complex @@ -100,7 +102,7 @@ class HasPrecision t where getPrecision :: t -> Precision -class (HasPrecision t) => CanSetPrecision t where+class CanSetPrecision t where setPrecision :: Precision -> t -> t instance HasPrecision t => HasPrecision (Complex t) where@@ -124,24 +126,24 @@ instance HasPrecision t => HasPrecision (CollectErrors es t) where getPrecision vCE =- case getMaybeValueCE vCE of+ case getMaybeValue vCE of Just v -> getPrecision v _ -> defaultPrecision instance CanSetPrecision t => CanSetPrecision (CollectErrors es t) where setPrecision p = fmap (setPrecision p) -lowerPrecisionIfAbove :: (CanSetPrecision t) => Precision -> t -> t+lowerPrecisionIfAbove :: (HasPrecision t, CanSetPrecision t) => Precision -> t -> t lowerPrecisionIfAbove p x | getPrecision x > p = setPrecision p x | otherwise = x -raisePrecisionIfBelow :: (CanSetPrecision t) => Precision -> t -> t+raisePrecisionIfBelow :: (HasPrecision t, CanSetPrecision t) => Precision -> t -> t raisePrecisionIfBelow p x | getPrecision x < p = setPrecision p x | otherwise = x specCanSetPrecision ::- (CanSetPrecision t, CanTestFinite t, Arbitrary t, Show t, Testable prop)+ (HasPrecision t, CanSetPrecision t, CanTestFinite t, Arbitrary t, Show t, Testable prop) => (T t) -> (t -> t -> prop) -> Spec specCanSetPrecision (T typeName :: T t) check =@@ -213,7 +215,6 @@ (t1 -> t1 -> c) -> (t1 -> t2 -> c) convertPSecond = convertSecondUsing (\ b q -> convertP (getPrecision b) q) - instance Arbitrary Precision where arbitrary = sized $ \size -> choose (4*(size+1),10*(size+1)) >>= return . prec@@ -222,6 +223,6 @@ [[t| Bool |], [t| Integer |], [t| Int |], [t| Rational |], [t| Double |]] (\ t -> [d| - instance (ConvertibleWithPrecision $t t, Monoid es) => ConvertibleWithPrecision $t (CollectErrors es t) where- safeConvertP p = fmap (\v -> CollectErrors (Just v) mempty) . safeConvertP p+ instance (ConvertibleWithPrecision $t t, CanBeErrors es) => ConvertibleWithPrecision $t (CollectErrors es t) where+ safeConvertP p = fmap pure . safeConvertP p |]))
+ src/AERN2/MP/WithCurrentPrec.hs view
@@ -0,0 +1,181 @@+{-# LANGUAGE EmptyDataDecls #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE PolyKinds #-}+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE GeneralizedNewtypeDeriving #-}+{-# LANGUAGE StandaloneDeriving #-}+{-|+ Module : AERN2.MP.WithCurrentPrec+ Description : Type wrapper setting default precision+ Copyright : (c) Michal Konecny+ License : BSD3++ Maintainer : mikkonecny@gmail.com+ Stability : experimental+ Portability : portable++ Type wrapper setting default precision.++ Not suitable for use with MixedTypesNumPrelude since we need binary operators to enforce+ the same precision on both operands via the equality of their types.++ Borrowed some tricks from https://github.com/ekmett/rounded/blob/master/src/Numeric/Rounded/Precision.hs+-}+module AERN2.MP.WithCurrentPrec+(+ WithCurrentPrec(..), runWithPrec, HasCurrentPrecision(..)+ -- , _example1 , _example2 , _example3+)+where++import qualified MixedTypesNumPrelude as MxP+import Prelude+-- import Text.Printf++-- import Text.Printf+import Numeric.CollectErrors (CN, cn, NumErrors, CanTakeErrors(..))+-- import qualified Numeric.CollectErrors as CN++import Data.Proxy+import Data.Reflection+import GHC.TypeLits++-- import Data.Complex++import AERN2.Limit+import AERN2.MP.Precision+import AERN2.MP.Ball++class HasCurrentPrecision p where+ getCurrentPrecision :: proxy p -> Precision++instance KnownNat n => HasCurrentPrecision n where+ getCurrentPrecision p = max (prec 2) . min maximumPrecision $ prec (natVal p)++-- data PrecAdd10 (p :: *)++-- instance (HasCurrentPrecision p) => HasCurrentPrecision (PrecAdd10 p) where+-- isPrecision (_ :: proxy _) = 10 + isPrecision (undefined :: proxy p)++newtype WithCurrentPrec t p = WithCurrentPrec { unWithCurrentPrec :: t }+ deriving (Show)++deriving instance (CanTakeErrors NumErrors t) => (CanTakeErrors NumErrors (WithCurrentPrec t p))++runWithPrec :: Precision -> (forall n. (KnownNat n) => WithCurrentPrec t n) -> t+runWithPrec p (wfp :: (forall n. (KnownNat n) => WithCurrentPrec t n)) = + reifyNat (MxP.integer p) withNat+ where+ withNat :: KnownNat n => Proxy n -> t+ withNat (_ :: Proxy n) = + unWithCurrentPrec (wfp :: WithCurrentPrec t n)++-- -- The following does not work:+-- instance (CanAddAsymmetric t1 t2) => (CanAddAsymmetric (WithCurrentPrec t1 p) (WithCurrentPrec t2 p)) where+-- type AddType (WithCurrentPrec t1 p) (WithCurrentPrec t2 p) = WithCurrentPrec (AddType t1 t2) p+-- add (WithCurrentPrec a1) (WithCurrentPrec a2) = WithCurrentPrec $ a1 + a2++instance + (MxP.HasOrderAsymmetric t1 t2)+ =>+ MxP.HasOrderAsymmetric (WithCurrentPrec t1 p1) (WithCurrentPrec t2 p2) + where+ type OrderCompareType (WithCurrentPrec t1 p1) (WithCurrentPrec t2 p2) = MxP.OrderCompareType t1 t2+ greaterThan (WithCurrentPrec v1) (WithCurrentPrec v2) = MxP.greaterThan v1 v2+ lessThan (WithCurrentPrec v1) (WithCurrentPrec v2) = MxP.lessThan v1 v2+ geq (WithCurrentPrec v1) (WithCurrentPrec v2) = MxP.geq v1 v2+ leq (WithCurrentPrec v1) (WithCurrentPrec v2) = MxP.leq v1 v2++instance + (MxP.CanMinMaxAsymmetric t1 t2, p1 ~ p2)+ =>+ MxP.CanMinMaxAsymmetric (WithCurrentPrec t1 p1) (WithCurrentPrec t2 p2) + where+ type MinMaxType (WithCurrentPrec t1 p1) (WithCurrentPrec t2 p2) = WithCurrentPrec (MxP.MinMaxType t1 t2) p1+ min (WithCurrentPrec v1) (WithCurrentPrec v2) = WithCurrentPrec $ MxP.min v1 v2+ max (WithCurrentPrec v1) (WithCurrentPrec v2) = WithCurrentPrec $ MxP.max v1 v2++instance Eq t => Eq (WithCurrentPrec t p) where+ (==) = lift2P (==)+instance Ord t => Ord (WithCurrentPrec t p) where+ compare = lift2P compare++instance + (HasCurrentPrecision p, Num t, ConvertibleWithPrecision Integer t) + => + Num (WithCurrentPrec t p) + where+ fromInteger n = r+ where + r = WithCurrentPrec $ convertP (getCurrentPrecision r) n+ negate = lift1 negate+ abs = lift1 abs+ (+) = lift2 (+)+ (*) = lift2 (*)+ signum = lift1 signum++instance + (HasCurrentPrecision p, Fractional t+ , ConvertibleWithPrecision Integer t, ConvertibleWithPrecision Rational t) + => + Fractional (WithCurrentPrec t p) + where+ fromRational q = r+ where + r = WithCurrentPrec $ convertP (getCurrentPrecision r) q+ recip = lift1 recip+ (/) = lift2 (/)++instance (HasCurrentPrecision p) => Floating (WithCurrentPrec (CN MPBall) p) where+ pi = r + where+ r = WithCurrentPrec $ cn $ piBallP (getCurrentPrecision r)+ sqrt = lift1 sqrt+ exp = lift1 exp+ log = lift1 log+ sin = lift1 sin+ cos = lift1 cos+ asin = lift1 asin+ acos = lift1 acos+ atan = lift1 atan+ sinh = lift1 sinh+ cosh = lift1 cosh+ asinh = lift1 asinh+ acosh = lift1 acosh+ atanh = lift1 atanh++instance + (HasLimits ix (CN MPBall -> CN MPBall)+ , LimitType ix (CN MPBall -> CN MPBall) ~ (CN MPBall -> CN MPBall)+ ,HasCurrentPrecision p)+ => + HasLimits ix (WithCurrentPrec (CN MPBall) p) + where+ type LimitType ix (WithCurrentPrec (CN MPBall) p) = WithCurrentPrec (CN MPBall) p+ limit (s :: ix -> (WithCurrentPrec (CN MPBall) p)) = + WithCurrentPrec $ limit (snop) $ sample+ where+ sample :: CN MPBall+ sample = setPrecision (getCurrentPrecision sampleP) 0 + sampleP :: WithCurrentPrec MPBall p+ sampleP = error "sampleP is not defined, it is only a type proxy"+ snop :: ix -> (CN MPBall -> CN MPBall)+ snop ix _sample = unWithCurrentPrec $ s ix++lift1 :: (t1 -> t2) -> (WithCurrentPrec t1 p) -> (WithCurrentPrec t2 p)+lift1 f (WithCurrentPrec v1) = WithCurrentPrec (f v1)++lift2 :: (t1 -> t2 -> t3) -> (WithCurrentPrec t1 p) -> (WithCurrentPrec t2 p) -> (WithCurrentPrec t3 p)+lift2 f (WithCurrentPrec v1) (WithCurrentPrec v2) = WithCurrentPrec (f v1 v2)++lift2P :: (t1 -> t2 -> t3) -> (WithCurrentPrec t1 p) -> (WithCurrentPrec t2 p) -> t3+lift2P f (WithCurrentPrec v1) (WithCurrentPrec v2) = f v1 v2++_example1 :: CN MPBall+_example1 = runWithPrec (prec 1000) pi++_example2 :: CN MPBall+_example2 = runWithPrec (prec 1000) $ pi - pi++_example3 :: CN MPBall+_example3 = runWithPrec (prec 1000) $ sqrt 2
src/AERN2/Norm.hs view
@@ -56,7 +56,7 @@ getNormLog x | x == 0.0 = NormZero | abs x >= 1.0 = getNormLog $ ceiling $ abs x- | otherwise = NormBits $ negate $ integer $ integerLog2 $ floor $ (1 /! (abs x))+ | otherwise = NormBits $ negate $ integer $ integerLog2 $ floor $ (1 / (abs x)) instance (HasNorm t)