module Data.Number.ER.RnToRm.UnitDom.ChebyshevBase.Polynom.Tests.Generate
where
import Data.Number.ER.RnToRm.UnitDom.ChebyshevBase.Polynom.Basic
import Data.Number.ER.RnToRm.UnitDom.ChebyshevBase.Polynom.Ring
import Data.Number.ER.RnToRm.UnitDom.ChebyshevBase.Polynom.Reduce
import Data.Number.ER.RnToRm.UnitDom.ChebyshevBase.Polynom.Eval
import Data.Number.ER.RnToRm.UnitDom.ChebyshevBase.Polynom.Bounds
import Data.Number.ER.RnToRm.UnitDom.ChebyshevBase.Polynom.Enclosure
import qualified Data.Number.ER.Real.Base as B
import qualified Data.Number.ER.Real.DomainBox as DBox
import Data.Number.ER.Real.DomainBox (VariableID(..), DomainBox, DomainBoxMappable, DomainIntBox)
import Data.Number.ER.Misc
import Data.Number.ER.BasicTypes
import Data.Number.ER.Real.DefaultRepr
import Data.Number.ER.Real.DomainBox.IntMap
import Data.Number.ER.Real.Approx.Interval
import qualified Data.Number.ER.Real.Approx as RA
import Test.QuickCheck hiding (two, three)
import qualified Data.Map as Map
type P = ERChebPoly (Box Int) BM
newtype PNoLimits = PNoLimits (Int, P) deriving (Show)
newtype PSize10Degree3 = PSize10Degree3 (Int, P) deriving (Show)
newtype PSize10Degree10 = PSize10Degree10 (Int, P) deriving (Show)
newtype PSize10 = PSize10 (Int, P) deriving (Show)
newtype PSize30 = PSize30 ((Int, Int), P) deriving (Show)
instance (Arbitrary PNoLimits)
where
arbitrary =
elements $ map PNoLimits $ zip [0..] $
polynomials1200ish id
coarbitrary p =
error "ERChebPoly: Generate: Arbitrary: coarbitrary not implemented for polynomials"
instance (Arbitrary PSize10Degree3)
where
arbitrary =
elements $ map PSize10Degree3 $ zip [0..] $ polynomials1200ishSize10Degree3
coarbitrary (PSize10Degree3 p) =
error "ERChebPoly: Generate: Arbitrary: coarbitrary not implemented for polynomials"
polynomials1200ishSize10Degree3 =
polynomials1200ish $ chplReduceTermCountUp 10 . chplReduceDegreeUp 3
instance (Arbitrary PSize10Degree10)
where
arbitrary =
elements $ map PSize10Degree10 $ zip [0..] $
polynomials1200ishSize10Degree10
coarbitrary (PSize10Degree10 p) =
error "ERChebPoly: Generate: Arbitrary: coarbitrary not implemented for polynomials"
polynomials1200ishSize10Degree10 =
polynomials1200ish $ chplReduceTermCountUp 10 . chplReduceDegreeUp 10
instance (Arbitrary PSize10)
where
arbitrary =
elements $ map PSize10 $ zip [0..] $ polynomials1200ishSize10
coarbitrary (PSize10 p) =
error "ERChebPoly: Generate: Arbitrary: coarbitrary not implemented for polynomials"
polynomials1200ishSize10 =
polynomials1200ish $ chplReduceTermCountUp 10
instance (Arbitrary PSize30)
where
arbitrary =
sized arbitrarySized
where
arbitrarySized n
| n <= 28 =
elements $ map PSize30 $
zip (map (\n -> (0,n)) [0..]) $
polynomials200ishSize30
| otherwise =
elements $ map PSize30 $
zip (map (\n -> (1,n)) [0..]) $
polynomials1200ishSize30
coarbitrary (PSize30 p) =
error "ERChebPoly: Generate: Arbitrary: coarbitrary not implemented for polynomials"
polynomials1200ishSize30 =
polynomials1200ish $ chplReduceTermCountUp 30
polynomials200ishSize30 =
polynomials200ishSmall $ chplReduceTermCountUp 30
data Int20 = Int20 Int deriving (Show)
instance (Arbitrary Int20)
where
arbitrary =
do
n <- choose (2,20)
return $ Int20 n
coarbitrary (Int20 n) =
error "ERChebPoly: Generate: Arbitrary: coarbitrary not implemented for EffIx20"
data Deg20Size20 = Deg20Size20 Int Int deriving (Show)
instance (Arbitrary Deg20Size20)
where
arbitrary =
do
maxDegree <- choose (2,20)
maxSize <- choose (10,20)
return $ Deg20Size20 maxDegree maxSize
coarbitrary (Deg20Size20 maxDegree maxSize) =
error "ERChebPoly: Generate: Arbitrary: coarbitrary not implemented for Deg20Size20"
data SmallRatio = SmallRatio Int Int deriving (Show)
instance (Arbitrary SmallRatio)
where
arbitrary =
do
num <- choose (1000000,1000000)
denom <- choose (1,1000000)
return $ SmallRatio num denom
coarbitrary (SmallRatio num denom) =
error "ERChebPoly: Generate: Arbitrary: coarbitrary not implemented for SmallRatio"
chplAtKeyPointsCanBeLeq ::
(B.ERRealBase b, DomainBox box varid Int, Ord box,
DomainBoxMappable boxb boxbb varid b [ERInterval b], Show boxb) =>
ERChebPoly box b ->
ERChebPoly box b ->
Bool
chplAtKeyPointsCanBeLeq p1 p2 =
and $ map testPoint points
where
points = getKeyPoints (p1 +^ p2)
testPoint point
| lower1 <= upper2 =
True
| otherwise =
unsafePrint
(
"Failure at point = " ++ (show point)
) $
False
where
lower1 = chplEvalDown p1 point
upper2 = chplEvalUp p2 point
getKeyPoints p =
getKeyPointsForVars $ chplGetVars p
getKeyPointsForVars vars =
points
where
points = map DBox.fromList $ allCombinations $ map getVarPoints varDoms
varDoms = map (\v -> (v,unitInterval)) vars
unitInterval = ERInterval (1) 1
getVarPoints (var, dom) =
(var, [domLB, domMB, domRB])
where
ERInterval domLB domRB = dom
domMB = (domLB + domRB)/2
chplAtKeyPointsPointwiseBinaryDownUpConsistent ::
(B.ERRealBase b, RealFrac b, DomainBox box varid Int, Ord box,
DomainBoxMappable boxb boxbb varid b [ERInterval b], Show boxb) =>
((ERInterval b) -> (ERInterval b) -> (ERInterval b)) ->
ERChebPoly box b ->
ERChebPoly box b ->
(ERChebPoly box b, ERChebPoly box b) ->
Bool
chplAtKeyPointsPointwiseBinaryDownUpConsistent raOp p1 p2 (resLow, resHigh) =
and $ map testPoint points
where
points = getKeyPoints (p1 +^ p2)
testPoint point
| ok = ok
| otherwise =
unsafePrint
(
"chplAtKeyPointsPointwiseBinaryDownUpConsistent failed:"
++ "\n point = " ++ show point
++ "\n raOpAtPointHigh = " ++ show raOpAtPointHigh
++ "\n raOpAtPointLow = " ++ show raOpAtPointLow
++ "\n resAtPointHigh = " ++ show resAtPointHigh
++ "\n resAtPointLow = " ++ show resAtPointLow
)
ok
where
ok =
raOpAtPointLow <= resAtPointHigh
&&
raOpAtPointHigh >= resAtPointLow
resAtPointLow = chplEvalDown resLow point
resAtPointHigh = chplEvalUp resHigh point
raOpAtPoint@(ERInterval raOpAtPointLow raOpAtPointHigh) =
raOp p1AtPoint p2AtPoint
p1AtPoint = ERInterval p1AtPointLow p1AtPointHigh
(p1AtPointLow, p1AtPointHigh) = chplEval p1 point
p2AtPoint = ERInterval p2AtPointLow p2AtPointHigh
(p2AtPointLow, p2AtPointHigh) = chplEval p2 point
makeThickEncl maxDegree maxSize p1 p2 =
(pMax q1Neg q2Neg, pMax q1 q2)
where
q1Neg = chplNeg q1
q2Neg = chplNeg q2
q1 = p1 +^ p2Mp1ScaledDown
q2 = p1 -^ p2Mp1ScaledDown
p2Mp1ScaledDown =
chplScaleUp (10/sizeB) p2Mp1
where
sizeB = max (abs upperB) (abs lowerB)
(lowerB, upperB) = chplBounds 10 p2Mp1
p2Mp1 = p2 -^ p1
pMax = chplMaxUp maxDegree maxSize
makeParalEncl p num denom =
(pNeg, p +^ cP)
where
pNeg = chplNeg p
cP = chplConst cB
cB = abs $ numB / (1000 * denomB)
numB = fromInteger $ toInteger num
denomB = fromInteger $ toInteger denom
enclRestrictRange ix (Nothing, Nothing) pEncl = (True, pEncl)
enclRestrictRange ix (maybeLower, maybeUpper) preEncl =
(succeeded, pEncl)
where
succeeded = lowerSucceeded && upperSucceeded
lowerSucceeded =
case maybeLower of
Nothing -> True
Just lower -> pLowerBound > lower
upperSucceeded =
case maybeUpper of
Nothing -> True
Just upper -> pUpperBound < upper
(pLowerBound, pUpperBound) = enclBounds ix pEncl
pEncl =
case (maybeLower, maybeUpper) of
(Just lowerB, Nothing) ->
case lowerB <= preLowerBoundB of
True -> preEncl
False ->
enclAddConst (lowerB preLowerBoundB + sepB) preEncl
(Nothing, Just upperB) ->
case preUpperBoundB <= upperB of
True -> preEncl
False ->
enclAddConst (upperB preUpperBoundB sepB) preEncl
(Just lowerB, Just upperB) ->
case lowerB <= preLowerBoundB && preUpperBoundB <= upperB of
True -> preEncl
_ ->
case preWidthB + sepB <= widthB of
True ->
enclAddConst
(lowerB preLowerBoundB + (preWidthB widthB)/2)
preEncl
_ ->
enclAddConst
(lowerB + sepB) $
enclScaleNonneg
(widthB / saferPreWidthB) $
enclAddConst
(sepB preLowerBoundB)
preEncl
where
widthB = upperB lowerB
saferPreWidthB = preWidthB + 2 * sepB
sepB = preWidthB / 1000000
preWidthB = preUpperBoundB preLowerBoundB
(preLowerBoundB, preUpperBoundB) = enclBounds ix preEncl
enclAtKeyPointsPointwiseBinaryDownUpConsistent ::
(B.ERRealBase b, RealFrac b, DomainBox box varid Int, Ord box,
DomainBoxMappable boxb boxbb varid b [ERInterval b], Show boxb, Show testId) =>
String ->
testId ->
((ERInterval b) -> (ERInterval b) -> (ERInterval b))
->
(ERChebPoly box b, ERChebPoly box b) ->
(ERChebPoly box b, ERChebPoly box b) ->
(ERChebPoly box b, ERChebPoly box b) ->
Bool
enclAtKeyPointsPointwiseBinaryDownUpConsistent
reportFileName testId
raOpInner
p1Encl@(p1LowNeg, p1High) p2Encl@(p2LowNeg, p2High) resEncl =
and $ map testPoint points
where
points = getKeyPoints (p1High +^ p2High +^ p1LowNeg +^ p2LowNeg)
testPoint point
| result =
unsafeReport reportFileName
(
show $
(testId, point, p1OpInnerP2AtPoint, resAtPoint)
)
result
| otherwise =
unsafePrint
(
"enclAtKeyPointsPointwiseBinaryDownUpConsistent failed"
++ "\n point = " ++ show point
++ "\n p1AtPoint = " ++ show p1AtPoint
++ "\n p2AtPoint = " ++ show p2AtPoint
++ "\n p1OpInnerP2AtPoint = " ++ show p1OpInnerP2AtPoint
++ "\n resAtPoint = " ++ show resAtPoint
) $
result
where
result = p1OpInnerP2AtPoint `RA.refines` resAtPoint
p1OpInnerP2AtPoint = p1AtPoint `raOpInner` p2AtPoint
resAtPoint = enclEval resEncl point
p1AtPoint = normaliseERInterval $ enclEvalInner p1Encl point
p2AtPoint = normaliseERInterval $ enclEvalInner p2Encl point
enclAtKeyPointsPointwiseUnaryDownUpConsistent ::
(B.ERRealBase b, RealFrac b, DomainBox box varid Int, Ord box,
DomainBoxMappable boxb boxbb varid b [ERInterval b], Show boxb, Show testId) =>
String ->
testId ->
((ERInterval b) -> (ERInterval b))
->
(ERChebPoly box b, ERChebPoly box b) ->
(ERChebPoly box b, ERChebPoly box b) ->
Bool
enclAtKeyPointsPointwiseUnaryDownUpConsistent
reportFileName testId
raOpInner
pEncl@(pLowNeg, pHigh) resEncl =
and $ map testPoint points
where
points = getKeyPoints (pHigh +^ pLowNeg)
testPoint point
| result =
unsafeReport reportFileName
(
show $
(testId, point, opInnerPAtPoint, resAtPoint)
)
result
| otherwise =
unsafePrint
(
"enclAtKeyPointsPointwiseUnaryDownUpConsistent failed"
++ "\n point = " ++ show point
++ "\n pAtPoint = " ++ show pAtPoint
++ "\n opInnerPAtPoint = " ++ show opInnerPAtPoint
++ "\n resAtPoint = " ++ show resAtPoint
) $
result
where
result = opInnerPAtPoint `RA.refines` resAtPoint
opInnerPAtPoint = raOpInner pAtPoint
resAtPoint = enclEval resEncl point
pAtPoint =
enclEvalInner pEncl point
enclAtKeyPointsConsistent ::
(B.ERRealBase b, RealFrac b, DomainBox box varid Int, Ord box,
DomainBoxMappable boxb boxbb varid b [ERInterval b], Show boxb, Show testId) =>
String ->
testId ->
(boxb -> (ERInterval b))
->
[varid] ->
(ERChebPoly box b, ERChebPoly box b) ->
Bool
enclAtKeyPointsConsistent
reportFileName testId
opInner allVars resEncl@(resLowNeg, resHigh) =
and $ map testPoint points
where
points = getKeyPointsForVars allVars
testPoint point
| result =
unsafeReport reportFileName
(
show $
(testId, point, opInnerAtPoint, resAtPoint)
)
result
| otherwise =
unsafePrint
(
"enclAtKeyPointsConsistent failed"
++ "\n point = " ++ show point
++ "\n opInnerAtPoint = " ++ show opInnerAtPoint
++ "\n resAtPoint = " ++ show resAtPoint
) $
result
where
result = opInnerAtPoint `RA.refines` resAtPoint
opInnerAtPoint = opInner point
resAtPoint = enclEval resEncl point
type E = (P,P)
vars :: [P]
vars = map chplVar [0..7]
varsE :: [E]
varsE = map (\p -> (chplNeg p, p)) vars
x0 = vars !! 0
x1 = vars !! 1
x2 = vars !! 2
x3 = vars !! 3
x4 = vars !! 4
x0E = varsE !! 0
x1E = varsE !! 1
x2E = varsE !! 2
x3E = varsE !! 3
x4E = varsE !! 4
one :: P
[mone, one, two, three, seven, thousand, million, tiny, huge] =
map chplConst
[1,1,2,3,7,1000,1000000,10^^(200),10^^200]
oneE :: E
[moneE, oneE, twoE, threeE, sevenE, thousandE, millionE, tinyE, hugeE] =
map (\ c -> (chplConst (c), chplConst c))
[1,1,2,3,7,1000,1000000,10^^(200),10^^200]
polynomials1200ish rdc =
concat $ map (powers10 rdc) $
concat $ map addConsts3 $
concat $ map multConsts3 $
polyBase13
polynomials200ish rdc =
concat $ map (powers4 rdc) $
concat $ map addConsts3 $
concat $ map multConsts3 $
polyBase5
polynomials40ish rdc =
concat $ map (powers2 rdc) $
concat $ map addConsts2 $
concat $ map multConsts2 $
polyBase5
polynomials200ishSmall rdc =
concat $ map (powers4Small rdc) $
concat $ map addConsts3 $
concat $ map multConsts3 $
polyBase5
polynomials40ishSmall rdc =
concat $ map (powers2Small rdc) $
concat $ map addConsts2 $
concat $ map multConsts2 $
polyBase5
polyBase5 =
[
(two *^ x0) +^ x1
,
(seven *^ x0) -^ x1
,
(tiny *^ x0) +^ x1
,
x0 -^ x1 *^ x2
,
x0 -^ x1 +^ x2 -^ x3 +^ x4
]
polyBase13 =
[
x0
,
x0 +^ x1
,
x0 -^ x1
,
(two *^ x0) +^ x1
,
(two *^ x0) -^ x1
,
(seven *^ x0) +^ x1
,
(seven *^ x0) -^ x1
,
(tiny *^ x0) +^ x1
,
(tiny *^ x0) -^ x1
,
x0 -^ x1 +^ x2
,
x0 -^ x1 *^ x2
,
x0 +^ x1 +^ x2 +^ x3 +^ x4
,
x0 -^ x1 +^ x2 -^ x3 +^ x4
]
powersAll rdc p =
powersAux [p, rdc $ p *^ p]
where
powersAux (pNHalfM1 : pNHalf : rest) =
pNHalfM1 : (powersAux $ (pNHalf : rest) ++ [pNM1, pN])
where
pNM1 = rdc $ pNHalf *^ pNHalfM1
pN = rdc $ pNHalf *^ pNHalf
powersForExps rdc p exponents =
map pw exponents
where
pw n = pws !! (n 1)
pws = powersAll rdc p
powers10 rdc p =
powersForExps rdc p [1..10]
powers4 rdc p =
powersForExps rdc p [1,3,5,7]
powers4Small rdc p =
powersForExps rdc p [1,2,3,5]
powers2 rdc p =
powersForExps rdc p [1,7]
powers2Small rdc p =
powersForExps rdc p [1,3]
addConsts3 p =
[p +^ one, p +^ three, p +^ seven]
multConsts3 p =
[p *^ two, p *^ three, p *^ seven]
addConsts2 p =
[p +^ one, p +^ three]
multConsts2 p =
[p *^ two, p *^ seven]