module Data.Number.ER.Real.Approx.Interval
(
ERInterval(..),
normaliseERInterval,
intervalTimesInner,
intervalPlusInner,
intervalDivideInner
)
where
import qualified Data.Number.ER.Real.Approx as RA
import qualified Data.Number.ER.Real.Approx.Elementary as RAEL
import qualified Data.Number.ER.Real.Base as B
import qualified Data.Number.ER.ExtendedInteger as EI
import Data.Number.ER.BasicTypes
import Data.Number.ER.Misc
import Data.Ratio
import qualified Text.Html as H
import Data.Typeable
import Data.Generics.Basics
import Data.Binary
data ERInterval base =
ERIntervalEmpty
| ERIntervalAny
| ERInterval
{
erintv_left :: base,
erintv_right :: base
}
deriving (Typeable, Data)
instance (Binary a) => Binary (ERInterval a) where
put ERIntervalEmpty = putWord8 0
put ERIntervalAny = putWord8 1
put (ERInterval a b) = putWord8 2 >> put a >> put b
get = do
tag_ <- getWord8
case tag_ of
0 -> return ERIntervalEmpty
1 -> return ERIntervalAny
2 -> get >>= \a -> get >>= \b -> return (ERInterval a b)
_ -> fail "no parse"
normaliseERInterval ::
(B.ERRealBase b) =>
ERInterval b -> ERInterval b
normaliseERInterval (ERInterval minusInfty plusInfty)
| B.isPlusInfinity plusInfty && B.isPlusInfinity ( minusInfty) =
ERIntervalAny
normaliseERInterval (ERInterval nan1 nan2)
| B.isERNaN nan1 && B.isERNaN nan2 =
ERIntervalAny
normaliseERInterval (ERInterval nan r)
| B.isERNaN nan =
ERInterval ( B.plusInfinity) r
normaliseERInterval (ERInterval l nan)
| B.isERNaN nan =
ERInterval l (B.plusInfinity)
normaliseERInterval (ERInterval l r)
| l > r = ERIntervalEmpty
normaliseERInterval i = i
erintvPrecision ::
(B.ERRealBase b) =>
ERInterval b -> EI.ExtendedInteger
erintvPrecision (ERInterval l r) =
1 (B.getApproxBinaryLog $ (r l))
erintvPrecision ERIntervalEmpty = EI.PlusInfinity
erintvPrecision ERIntervalAny = EI.MinusInfinity
erintvGranularity ::
(B.ERRealBase b) =>
ERInterval b -> Int
erintvGranularity ERIntervalAny = 0
erintvGranularity ERIntervalEmpty = 0
erintvGranularity (ERInterval l r) =
min (B.getGranularity l) (B.getGranularity r)
erintvEqualApprox ::
(B.ERRealBase b) =>
ERInterval b -> ERInterval b -> Bool
erintvEqualApprox (ERInterval l1 r1) (ERInterval l2 r2) =
l1 == l2 && r1 == r2
erintvEqualApprox ERIntervalEmpty ERIntervalEmpty = True
erintvEqualApprox ERIntervalAny ERIntervalAny = True
erintvEqualApprox _ _ = False
erintvCompareApprox ::
(B.ERRealBase b) =>
ERInterval b -> ERInterval b -> Ordering
erintvCompareApprox ERIntervalEmpty ERIntervalEmpty = EQ
erintvCompareApprox ERIntervalEmpty _ = LT
erintvCompareApprox _ ERIntervalEmpty = GT
erintvCompareApprox ERIntervalAny ERIntervalAny = EQ
erintvCompareApprox ERIntervalAny _ = LT
erintvCompareApprox _ ERIntervalAny = GT
erintvCompareApprox (ERInterval l1 r1) (ERInterval l2 r2) =
case compare l1 l2 of
EQ -> compare r1 r2
res -> res
erintvEqualReals ::
(B.ERRealBase b) =>
ERInterval b ->
ERInterval b ->
Maybe Bool
erintvEqualReals ERIntervalEmpty _ = Nothing
erintvEqualReals _ ERIntervalEmpty = Nothing
erintvEqualReals ERIntervalAny _ = Nothing
erintvEqualReals _ ERIntervalAny = Nothing
erintvEqualReals (ERInterval l1 r1) (ERInterval l2 r2)
| l1 == r1 && l2 == r2 && l1 == l2 = Just True
| r1 < l2 || l1 > r2 = Just False
| otherwise = Nothing
erintvCompareReals ::
(B.ERRealBase b) =>
ERInterval b ->
ERInterval b ->
Maybe Ordering
erintvCompareReals ERIntervalEmpty _ = Nothing
erintvCompareReals _ ERIntervalEmpty = Nothing
erintvCompareReals ERIntervalAny _ = Nothing
erintvCompareReals _ ERIntervalAny = Nothing
erintvCompareReals i1@(ERInterval l1 r1) i2@(ERInterval l2 r2)
| r1 < l2 = Just LT
| l1 > r2 = Just GT
| l1 == r1 && l2 == r2 && l1 == l2 = Just EQ
| otherwise = Nothing
erintvLeqReals ::
(B.ERRealBase b) =>
ERInterval b ->
ERInterval b ->
Maybe Bool
erintvLeqReals ERIntervalEmpty _ = Nothing
erintvLeqReals _ ERIntervalEmpty = Nothing
erintvLeqReals ERIntervalAny _ = Nothing
erintvLeqReals _ ERIntervalAny = Nothing
erintvLeqReals i1@(ERInterval l1 r1) i2@(ERInterval l2 r2)
| r1 <= l2 = Just True
| l1 > r2 = Just False
| otherwise = Nothing
erintvDefaultBisectPt ::
(B.ERRealBase b) =>
Granularity ->
(ERInterval b) ->
(ERInterval b)
erintvDefaultBisectPt gran ERIntervalAny = 0
erintvDefaultBisectPt gran ERIntervalEmpty = ERIntervalEmpty
erintvDefaultBisectPt gran (ERInterval l r) =
ERInterval m m
where
m
| B.isPlusInfinity r =
if l < 0
then 0
else 2 * (B.setMinGranularity gran l) + 1
| B.isPlusInfinity (l) =
if r > 0
then 0
else 2 * (B.setMinGranularity gran r) 1
| otherwise =
((B.setMinGranularity gran l) + r)/2
erintvBisect ::
(B.ERRealBase b, RealFrac b) =>
Granularity ->
(Maybe (ERInterval b)) ->
(ERInterval b) ->
(ERInterval b, ERInterval b)
erintvBisect gran maybePt i =
(l RA.\/ m, m RA.\/ r)
where
(l,r) = RA.bounds i
m =
case maybePt of
Just m -> m
Nothing -> erintvDefaultBisectPt gran i
instance (B.ERRealBase b) => Eq (ERInterval b) where
i1 == i2 =
case erintvEqualReals i1 i2 of
Nothing ->
error $
"ERInterval: Eq: comparing overlapping intervals:\n" ++
show i1 ++ "\n" ++
show i2
Just b -> b
instance (B.ERRealBase b) => Ord (ERInterval b) where
compare i1 i2 =
case erintvCompareReals i1 i2 of
Nothing ->
error $
"ERInterval: Ord: comparing overlapping intervals:\n" ++
show i1 ++ "\n" ++
show i2
Just r -> r
max i1@(ERInterval l1 r1) i2@(ERInterval l2 r2) =
normaliseERInterval $ ERInterval (max l1 l2) (max r1 r2)
max ERIntervalEmpty _ = ERIntervalEmpty
max _ ERIntervalEmpty = ERIntervalEmpty
max ERIntervalAny ERIntervalAny = ERIntervalAny
max ERIntervalAny (ERInterval l r) = ERInterval l B.plusInfinity
max (ERInterval l r) ERIntervalAny = ERInterval l B.plusInfinity
min i1@(ERInterval l1 r1) i2@(ERInterval l2 r2) =
normaliseERInterval $ ERInterval (min l1 l2) (min r1 r2)
min ERIntervalEmpty _ = ERIntervalEmpty
min _ ERIntervalEmpty = ERIntervalEmpty
min ERIntervalAny ERIntervalAny = ERIntervalAny
min ERIntervalAny (ERInterval l r) = ERInterval ( B.plusInfinity) r
min (ERInterval l r) ERIntervalAny = ERInterval ( B.plusInfinity) r
instance (B.ERRealBase b) => Show (ERInterval b)
where
show = erintvShow 16 True False
erintvShow numDigits showGran showComponents interval =
showERI interval
where
showERI ERIntervalEmpty = "[NONE]"
showERI ERIntervalAny = "[ANY]"
showERI (ERInterval l r)
| l == r = "<" ++ showBase l ++ ">"
| otherwise =
"[" ++ showBase l ++ "," ++ showBase r ++ "]"
showBase = B.showDiGrCmp numDigits showGran showComponents
instance (B.ERRealBase b, H.HTML b) => H.HTML (ERInterval b)
where
toHtml (ERInterval l r)
| l == r =
H.toHtml $ show l
| otherwise =
H.simpleTable [] [] [[H.toHtml l],[H.toHtml r]]
toHtml i = H.toHtml $ show i
instance (B.ERRealBase b) => Num (ERInterval b) where
fromInteger n =
normaliseERInterval $ ERInterval (fromInteger n) (fromInteger n)
abs (ERInterval l r)
| l < 0 && r > 0 = ERInterval 0 (max (l) r)
| r <= 0 = ERInterval (r) (l)
| otherwise = ERInterval l r
abs ERIntervalAny = ERInterval 0 B.plusInfinity
abs ERIntervalEmpty = ERIntervalEmpty
signum i@(ERInterval l r)
| l < 0 && r > 0 = ERInterval (1) 1
| r < 0 = ERInterval (1) (1)
| l > 0 = ERInterval 1 1
| l == 0 && r == 0 = i
| l == 0 = ERInterval 0 1
| r == 0 = ERInterval (1) 0
signum ERIntervalAny = ERInterval (1) 1
signum ERIntervalEmpty = ERIntervalEmpty
negate (ERInterval l r) = (ERInterval (r) (l))
negate ERIntervalEmpty = ERIntervalEmpty
negate ERIntervalAny = ERIntervalAny
i1 + i2 = fst $ intervalPlusOuterInner i1 i2
i1 * i2 = fst $ intervalTimesOuterInner i1 i2
intervalPlusInner ::
(B.ERRealBase b) =>
(ERInterval b) ->
(ERInterval b) ->
(ERInterval b)
intervalPlusInner i1 i2 = snd $ intervalPlusOuterInner i1 i2
intervalTimesInner ::
(B.ERRealBase b) =>
(ERInterval b) ->
(ERInterval b) ->
(ERInterval b)
intervalTimesInner i1 i2 = snd $ intervalTimesOuterInner i1 i2
intervalPlusOuterInner (ERInterval l1 r1) (ERInterval l2 r2) =
(normaliseERInterval $
ERInterval (l1 `plusDown` l2) (r1 `plusUp` r2),
ERInterval (l1 `plusUp` l2) (r1 `plusDown` r2))
intervalPlusOuterInner ERIntervalAny i2 = (ERIntervalAny, ERIntervalAny)
intervalPlusOuterInner l1 ERIntervalAny = (ERIntervalAny, ERIntervalAny)
intervalPlusOuterInner ERIntervalEmpty i2 = (ERIntervalEmpty, ERIntervalEmpty)
intervalPlusOuterInner l1 ERIntervalEmpty = (ERIntervalEmpty, ERIntervalEmpty)
intervalTimesOuterInner (ERInterval l1 r1) (ERInterval l2 r2)
| haveNan = (ERIntervalAny, ERIntervalAny)
| otherwise =
(normaliseERInterval $
ERInterval minProdOuter maxProdOuter,
ERInterval minProdInner maxProdInner)
where
haveNan = or $ map B.isERNaN (prodsUp ++ prodsDown)
minProdOuter = foldl1 min prodsDown
maxProdOuter = foldl1 max prodsUp
minProdInner = foldl1 min prodsUp
maxProdInner = foldl1 max prodsDown
prodsDown = [l1 `timesDown` l2, l1 `timesDown` r2, r1 `timesDown` l2, r1 `timesDown` r2]
prodsUp = [l1 `timesUp` l2, l1 `timesUp` r2, r1 `timesUp` l2, r1 `timesUp` r2]
intervalTimesOuterInner ERIntervalAny i2 = (ERIntervalAny, ERIntervalAny)
intervalTimesOuterInner l1 ERIntervalAny = (ERIntervalAny, ERIntervalAny)
intervalTimesOuterInner ERIntervalEmpty i2 = (ERIntervalEmpty, ERIntervalEmpty)
intervalTimesOuterInner l1 ERIntervalEmpty = (ERIntervalEmpty, ERIntervalEmpty)
instance (B.ERRealBase b) => Fractional (ERInterval b) where
fromRational rat =
(fromInteger $ numerator rat)
/ (fromInteger $ denominator rat)
i1 / i2 =
fst $ intervalDivideOuterInner i1 i2
intervalDivideInner ::
(B.ERRealBase b) =>
(ERInterval b) ->
(ERInterval b) ->
(ERInterval b)
intervalDivideInner i1 i2 = snd $ intervalDivideOuterInner i1 i2
intervalDivideOuterInner (ERInterval l1 r1) (ERInterval l2 r2)
| l2 < 0 && r2 > 0 = (ERIntervalAny, ERIntervalAny)
| haveNan =
(ERIntervalAny, ERIntervalAny)
| l2 == 0 && r2 > 0 && 1/l2 < 0 =
intervalDivideOuterInner (ERInterval l1 r1) (ERInterval (l2) r2)
| r2 == 0 && l2 < 0 && 1/r2 > 0 =
intervalDivideOuterInner (ERInterval l1 r1) (ERInterval l2 (r2))
| otherwise =
(
normaliseERInterval $
ERInterval minDivOuter maxDivOuter
,
ERInterval minDivInner maxDivInner
)
where
haveNan = or $ map B.isERNaN (divsUp ++ divsDown)
minDivOuter = foldl1 min divsDown
maxDivOuter = foldl1 max divsUp
minDivInner = foldl1 min divsUp
maxDivInner = foldl1 max divsDown
divsDown = [l1 `divideDown` l2, l1 `divideDown` r2, r1 `divideDown` l2, r1 `divideDown` r2]
divsUp = [l1 `divideUp` l2, l1 `divideUp` r2, r1 `divideUp` l2, r1 `divideUp` r2]
intervalDivideOuterInner ERIntervalAny i2 = (ERIntervalAny, ERIntervalAny)
intervalDivideOuterInner i1 ERIntervalAny = (ERIntervalAny, ERIntervalAny)
intervalDivideOuterInner ERIntervalEmpty i2 = (ERIntervalEmpty, ERIntervalEmpty)
intervalDivideOuterInner i1 ERIntervalEmpty = (ERIntervalEmpty, ERIntervalEmpty)
instance (B.ERRealBase b, RealFrac b) => RA.ERApprox (ERInterval b) where
initialiseBaseArithmetic _ =
B.initialiseBaseArithmetic (0 :: b)
getPrecision i = erintvPrecision i
getGranularity i = erintvGranularity i
setMinGranularity gr (ERInterval l r) =
normaliseERInterval $
(ERInterval ( (B.setMinGranularity gr (l))) (B.setMinGranularity gr r))
setMinGranularity _ i = i
setGranularity gr (ERInterval l r) =
normaliseERInterval $
(ERInterval ( (B.setGranularity gr (l))) (B.setGranularity gr r))
setGranularity _ i = i
bottomApprox = ERIntervalAny
emptyApprox = ERIntervalEmpty
isEmpty ERIntervalEmpty = True
isEmpty _ = False
isBottom ERIntervalAny = True
isBottom (ERInterval l r) =
B.isPlusInfinity r && B.isPlusInfinity (l)
isBottom _ = False
isExact ERIntervalEmpty = False
isExact ERIntervalAny = False
isExact (ERInterval l r) = l == r
isBounded ERIntervalEmpty = True
isBounded ERIntervalAny = False
isBounded (ERInterval l r) =
( B.plusInfinity) < l && r < B.plusInfinity
ERIntervalEmpty /\ i = ERIntervalEmpty
i /\ ERIntervalEmpty = ERIntervalEmpty
ERIntervalAny /\ i = i
i /\ ERIntervalAny = i
(ERInterval l1 r1) /\ (ERInterval l2 r2) =
normaliseERInterval $
ERInterval (max l1 l2) (min r1 r2)
intersectMeasureImprovement _ ERIntervalEmpty i = (ERIntervalEmpty, 1)
intersectMeasureImprovement _ i ERIntervalEmpty = (ERIntervalEmpty, 1)
intersectMeasureImprovement _ ERIntervalAny i = (i, 1)
intersectMeasureImprovement _ i ERIntervalAny = (i, 1)
intersectMeasureImprovement ix i1 i2 =
(isec, impr)
where
isec = i1 RA./\ i2
impr
| 0 `RA.refines` isecWidth && 0 `RA.refines` i1Width = 1
| otherwise = i1Width / isecWidth
i1Width = i1H i1L
isecWidth = isecH isecL
(isecL, isecH) = RA.bounds $ RA.setMinGranularity gran isec
(i1L, i1H) = RA.bounds $ RA.setMinGranularity gran i1
gran = effIx2gran ix
refines _ ERIntervalAny = True
refines ERIntervalEmpty _ = True
refines ERIntervalAny (ERInterval l r)
| B.isPlusInfinity r && B.isPlusInfinity (l) = True
refines ERIntervalAny _ = False
refines _ ERIntervalEmpty = False
refines (ERInterval l1 r1) (ERInterval l2 r2) =
l2 <= l1 && r1 <= r2
equalReals = erintvEqualReals
compareReals = erintvCompareReals
leqReals = erintvLeqReals
equalApprox = erintvEqualApprox
compareApprox = erintvCompareApprox
double2ra d =
ERInterval b b
where
b = B.fromDouble d
showApprox = erintvShow
instance (B.ERRealBase b, RealFrac b) => RA.ERIntApprox (ERInterval b)
where
doubleBounds ERIntervalAny = ( infinity, infinity)
where
infinity = 1/0
doubleBounds ERIntervalEmpty =
error "ERInterval: doubleBounds: empty interval"
doubleBounds (ERInterval l r) =
(B.toDouble l, B.toDouble r)
floatBounds ERIntervalAny = ( infinity, infinity)
where
infinity = 1/0
floatBounds ERIntervalEmpty =
error "ERInterval: floatBounds: empty interval"
floatBounds (ERInterval l r) =
(B.toFloat l, B.toFloat r)
integerBounds ERIntervalAny = ( infinity, infinity)
where
infinity = EI.PlusInfinity
integerBounds ERIntervalEmpty =
error "ERInterval: integerBounds: empty interval"
integerBounds (ERInterval l r) =
( (mkEI ( l)), mkEI r)
where
mkEI f
| B.isPlusInfinity f = EI.PlusInfinity
| B.isPlusInfinity (f) = EI.MinusInfinity
| otherwise = ceiling f
defaultBisectPt dom = erintvDefaultBisectPt (RA.getGranularity dom + 1) dom
bisectDomain maybePt dom =
erintvBisect (RA.getGranularity dom + 1) maybePt dom
ERIntervalEmpty \/ i = i
i \/ ERIntervalEmpty = i
ERIntervalAny \/ _ = ERIntervalAny
_ \/ ERIntervalAny = ERIntervalAny
(ERInterval l1 r1) \/ (ERInterval l2 r2) =
normaliseERInterval $
ERInterval (min l1 l2) (max r1 r2)
bounds ERIntervalAny =
(ERInterval (B.plusInfinity) (B.plusInfinity),
ERInterval B.plusInfinity B.plusInfinity)
bounds ERIntervalEmpty = (ERIntervalEmpty, ERIntervalEmpty)
bounds (ERInterval l r) =
(ERInterval l l, ERInterval r r)
fromBounds (ERInterval l1 r1, ERInterval l2 r2)
| l1 == r1 && l2 == r2 = ERInterval l1 l2
fromBounds i1i2 =
error $
"ER.Real.Approx.Interval: fromBounds: bounds not thin: "
++ show i1i2
instance (B.ERRealBase b, RealFrac b) => RAEL.ERApproxElementary (ERInterval b)