summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorMichalKonecny <>2008-12-01 23:14:34 (GMT)
committerLuite Stegeman <luite@luite.com>2008-12-01 23:14:34 (GMT)
commit207aa45f6ada4a47828546ac9580b67555382326 (patch)
tree5e1eed84bdd247d3136413f03221d4c0fcaaac22
parent41e99376b779ad918d1569a3747c34bf0e1b9b08 (diff)
version 0.4.20.4.2
-rw-r--r--AERN-RnToRm.cabal23
-rw-r--r--ChangeLog3
-rw-r--r--src/Data/Number/ER/RnToRm/Approx.hs2
-rw-r--r--src/Data/Number/ER/RnToRm/Approx/DomEdges.hs90
-rw-r--r--src/Data/Number/ER/RnToRm/Approx/DomTransl.hs98
-rw-r--r--src/Data/Number/ER/RnToRm/Approx/PieceWise.hs60
-rw-r--r--src/Data/Number/ER/RnToRm/Approx/Tuple.hs22
-rw-r--r--src/Data/Number/ER/RnToRm/BisectionTree.hs105
-rw-r--r--src/Data/Number/ER/RnToRm/UnitDom/Approx/Interval.hs141
-rw-r--r--src/Data/Number/ER/RnToRm/UnitDom/Base.hs10
-rw-r--r--src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom.hs2
-rw-r--r--src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Basic.hs45
-rw-r--r--src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Elementary.hs25
-rw-r--r--src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Field.hs6
-rw-r--r--tests/ISin3.hs2
15 files changed, 510 insertions, 124 deletions
diff --git a/AERN-RnToRm.cabal b/AERN-RnToRm.cabal
index a4c84d8..726e7d3 100644
--- a/AERN-RnToRm.cabal
+++ b/AERN-RnToRm.cabal
@@ -1,5 +1,5 @@
Name: AERN-RnToRm
-Version: 0.4.1
+Version: 0.4.2
Cabal-Version: >= 1.2
Build-Type: Simple
License: BSD3
@@ -10,7 +10,7 @@ Maintainer: mik@konecny.aow.cz
Stability: experimental
Category: Data, Math
Synopsis: polynomial function enclosures (PFEs) approximating exact real functions
-Tested-with: GHC ==6.8.2
+Tested-with: GHC ==6.8.3
Description:
AERN-RnToRm provides
datatypes and abstractions for approximating functions
@@ -34,7 +34,9 @@ Description:
.
Simple examples of usage can be found in module @Demo.hs@ in folder @tests@.
Extra-source-files:
- ChangeLog tests/Demo.hs tests/ISin3.hs
+ tests/Demo.hs tests/ISin3.hs
+Data-files:
+ ChangeLog
Flag containers-in-base
Default: False
@@ -43,10 +45,10 @@ Library
hs-source-dirs: src
if flag(containers-in-base)
Build-Depends:
- base < 3, binary >= 0.4, AERN-Real >= 0.9.7
+ base < 3, binary >= 0.4, html >= 1.0, AERN-Real >= 0.9.7
else
Build-Depends:
- base >= 3, containers, binary >= 0.4, AERN-Real >= 0.9.7
+ base >= 3, containers, binary >= 0.4, html >= 1.0, AERN-Real >= 0.9.7
Exposed-modules:
Data.Number.ER.RnToRm,
Data.Number.ER.RnToRm.BisectionTree.Path,
@@ -69,13 +71,4 @@ Library
Data.Number.ER.RnToRm.Approx.Tuple,
Data.Number.ER.RnToRm.Approx,
Data.Number.ER.RnToRm.TestingDefs
-
- Extensions:
- CPP,
- DeriveDataTypeable,
- FlexibleContexts,
- FlexibleInstances,
- FunctionalDependencies,
- MultiParamTypeClasses,
- UndecidableInstances
-
+
diff --git a/ChangeLog b/ChangeLog
index 2118449..4798a59 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,6 @@
+0.4.2: 1 December 2008
+ * substantially improved division by a constant PFE
+ * added proper handling of overflown coefficients
0.4.1: 30 September 2008
* updated to work with AERN-Real 0.9.7
0.4: 20 August 2008
diff --git a/src/Data/Number/ER/RnToRm/Approx.hs b/src/Data/Number/ER/RnToRm/Approx.hs
index 7f41bba..b5dc932 100644
--- a/src/Data/Number/ER/RnToRm/Approx.hs
+++ b/src/Data/Number/ER/RnToRm/Approx.hs
@@ -120,8 +120,8 @@ class
partialIntersect ::
EffortIndex ->
box {-^ the subdomain; defined by clipping the range of some variables -} ->
- fa {-^ function to improve by intersecting its subdomain -} ->
fa {-^ the enclosure to be used on the subdomain (but defined on the whole domain) -} ->
+ fa {-^ function to improve by intersecting its subdomain -} ->
fa
{-|
Intersect two enclosures and measure the global improvement as one number.
diff --git a/src/Data/Number/ER/RnToRm/Approx/DomEdges.hs b/src/Data/Number/ER/RnToRm/Approx/DomEdges.hs
index 88e395d..6a42104 100644
--- a/src/Data/Number/ER/RnToRm/Approx/DomEdges.hs
+++ b/src/Data/Number/ER/RnToRm/Approx/DomEdges.hs
@@ -31,6 +31,9 @@ import Data.Number.ER.BasicTypes
import Data.Number.ER.Misc
import Data.Number.ER.PlusMinus
+import Data.Number.ER.ShowHTML
+import qualified Text.Html as H
+
import Data.Typeable
import Data.Generics.Basics
import Data.Binary
@@ -152,6 +155,39 @@ instance
domB = FA.dom fa
instance
+ (FA.ERFnDomApprox box varid domra ranra fa, H.HTML fa) =>
+ H.HTML (ERFnDomEdgesApprox varid fa)
+ where
+ toHtml f@(ERFnDomEdgesApprox mainEncl edges)
+ | Map.null edges =
+ H.toHtml mainEncl
+ | otherwise =
+ H.toHtml $
+ abovesTable [H.border 2] [leftEdgesHtml, H.toHtml mainEncl, rightEdgesHtml]
+ where
+ leftEdgesHtml =
+ H.toHtml $ H.simpleTable [H.border 1] [] leftRows
+ rightEdgesHtml =
+ H.toHtml $ H.simpleTable [H.border 1] [] rightRows
+ leftRows =
+ map (makeRow fst) leftEdges
+ rightRows =
+ map (makeRow snd) rightEdges
+ (leftEdges, rightEdges) =
+ partition (\((_, pm),_) -> pm == Minus) $ Map.toList edges
+ makeRow pickEndpoint ((varId, _), faEdge) =
+ [
+ H.stringToHtml $ showVar varId ++ "=" ++ show varValueLeft
+ ,
+ H.toHtml faEdge
+ ]
+ where
+ varValueLeft =
+ pickEndpoint $
+ RA.bounds $ DBox.lookup "RnToRm.Approx.DomEdges: " varId domB
+ domB = FA.dom mainEncl
+
+instance
(FA.ERFnApprox box varid domra ranra fa) =>
Eq (ERFnDomEdgesApprox varid fa)
where
@@ -183,7 +219,7 @@ instance
instance
- (FA.ERFnDomApprox box varid domra ranra fa, VariableID varid) =>
+ (FA.ERFnDomApprox box varid domra ranra fa, VariableID varid, Show box) =>
RA.ERApprox (ERFnDomEdgesApprox varid fa)
where
initialiseBaseArithmetic _ =
@@ -211,9 +247,44 @@ instance
Map.intersectionWith (RA.intersectMeasureImprovement ix) edges1 edges2
leqReals fa1 fa2 =
RA.leqReals (erfnMainVolume fa1) (erfnMainVolume fa2)
-
+ refines
+ f1@(ERFnDomEdgesApprox mainEncl1 edges1)
+ f2@(ERFnDomEdgesApprox mainEncl2 edges2) =
+-- unsafePrint
+-- (
+-- "ERFnDomEdgesApprox: refines: "
+-- ++ "\n domB = " ++ show (FA.dom mainEncl1)
+-- ++ "\n mainRefines = " ++ show mainRefines
+-- ++ "\n mainEncl1 = " ++ show mainEncl1
+-- ++ "\n mainEncl2 = " ++ show mainEncl2
+-- )$
+ mainRefines && edgesRefine
+ where
+ mainRefines =
+ mainEncl1 `RA.refines` mainEncl2
+ edgesRefine =
+ Map.fold (&&) True $
+ Map.intersectionWith RA.refines edges1 edges2
+ compareApprox
+ f1@(ERFnDomEdgesApprox mainEncl1 edges1)
+ f2@(ERFnDomEdgesApprox mainEncl2 edges2) =
+ compareComposeMany
+ [
+ RA.compareApprox mainEncl1 mainEncl2
+ ,
+ compareListsWith compareEdges (Map.toAscList edges1) (Map.toAscList edges2)
+ ]
+ where
+ compareEdges ((v1,pm1),fa1) ((v2,pm2),fa2) =
+ compareComposeMany
+ [
+ compare v1 v2,
+ compare pm1 pm2,
+ RA.compareApprox fa1 fa2
+ ]
instance
- (FA.ERFnDomApprox box varid domra ranra fa, RA.ERIntApprox fa, VariableID varid) =>
+ (FA.ERFnDomApprox box varid domra ranra fa, RA.ERIntApprox fa,
+ VariableID varid, Show box) =>
RA.ERIntApprox (ERFnDomEdgesApprox varid fa)
where
-- doubleBounds = :: ira -> (Double, Double)
@@ -246,7 +317,8 @@ instance
f1 \/ f2 = edgesLift2 (RA.\/) f1 f2
instance
- (FA.ERFnDomApprox box varid domra ranra fa, RAEL.ERApproxElementary fa, VariableID varid) =>
+ (FA.ERFnDomApprox box varid domra ranra fa, RAEL.ERApproxElementary fa,
+ VariableID varid, Show box) =>
RAEL.ERApproxElementary (ERFnDomEdgesApprox varid fa)
where
abs ix = edgesLift1 $ RAEL.abs ix
@@ -257,7 +329,8 @@ instance
atan ix = edgesLift1 $ RAEL.atan ix
instance
- (FA.ERFnDomApprox box varid domra ranra fa, VariableID varid) =>
+ (FA.ERFnDomApprox box varid domra ranra fa,
+ VariableID varid, Show box) =>
FA.ERFnApprox box varid domra ranra (ERFnDomEdgesApprox varid fa)
where
check prgLocation (ERFnDomEdgesApprox mainEncl edges) =
@@ -312,7 +385,7 @@ instance
partialIntersectEdge (var, pm) edge1 edge2
| withinSubstitutions =
FA.partialIntersect ix substitutions edge1 edge2
- | otherwise = edge1
+ | otherwise = edge2
where
withinSubstitutions =
(varDomEndpoint pm) `RA.refines` varVal
@@ -322,7 +395,7 @@ instance
varDomEndpoint Minus = varDomLO
varDomEndpoint Plus = varDomHI
(varDomLO, varDomHI) = RA.bounds varDom
- varDom = DBox.lookup "DomEdges: partialIntersect: " var $ FA.dom mainEncl1
+ varDom = DBox.lookup "DomEdges: partialIntersect: " var $ FA.dom mainEncl2
eval ptBox (ERFnDomEdgesApprox mainEncl edges)
| null edgeVals =
mainVal
@@ -364,7 +437,8 @@ instance
substitutionsNew = DBox.delete varID substitutionsPrev
instance
- (FA.ERFnDomApprox box varid domra ranra fa, VariableID varid) =>
+ (FA.ERFnDomApprox box varid domra ranra fa,
+ VariableID varid, Show box) =>
FA.ERFnDomApprox box varid domra ranra (ERFnDomEdgesApprox varid fa)
where
dom (ERFnDomEdgesApprox mainEncl edges) = FA.dom mainEncl
diff --git a/src/Data/Number/ER/RnToRm/Approx/DomTransl.hs b/src/Data/Number/ER/RnToRm/Approx/DomTransl.hs
index a14ddf2..c01795e 100644
--- a/src/Data/Number/ER/RnToRm/Approx/DomTransl.hs
+++ b/src/Data/Number/ER/RnToRm/Approx/DomTransl.hs
@@ -35,6 +35,8 @@ import Data.Number.ER.Real.DomainBox (VariableID(..), DomainIntBox, DomainBoxMap
import Data.Number.ER.BasicTypes
import Data.Number.ER.Misc
+import qualified Text.Html as H
+
import Data.Typeable
import Data.Generics.Basics
import Data.Binary
@@ -90,6 +92,13 @@ instance
instance
(RA.ERIntApprox domra) =>
+ Ord (DomTransl domra)
+ where
+ compare (DomTransl _ _ _ _ dom1) (DomTransl _ _ _ _ dom2) =
+ RA.compareApprox dom1 dom2
+
+instance
+ (RA.ERIntApprox domra) =>
Show (DomTransl domra)
where
show (DomTransl fromA fromB toA toB dom) =
@@ -98,6 +107,15 @@ instance
" fromUnit = " ++ show fromA ++ " * x + " ++ show fromB ++ "\n" ++
" toUnit = " ++ show toA ++ " * x + " ++ show toB ++ "\n"
+instance
+ (RA.ERIntApprox domra, H.HTML domra) =>
+ H.HTML (DomTransl domra)
+ where
+ toHtml (DomTransl fromA fromB toA toB dom) =
+ error "DomTransl: toHtml not implemented yet"
+
+
+
dtrIdentity ::
(RA.ERIntApprox ira) =>
DomTransl ira
@@ -177,9 +195,9 @@ dtrUnion msg dtr1 dtr2
#endif
-dtrBShow dtrs =
+dtrBShow dtrB =
concatWith "," $
- map showOne $ DBox.toList dtrs
+ map showOne $ DBox.toList dtrB
where
showOne (var, dtr) =
showVar var ++ " in " ++ show (dtrDom dtr)
@@ -194,23 +212,41 @@ instance
"\nERFnDomTranslApprox" ++
show ufaDom ++
-- show ufa ++
- "\n dom = [" ++
- (concatWith ", " $ map showVarDom $ DBox.toList $ dtrBToDomB dtrB)
- ++ "]"
+ "\n dom = [" ++ dtrBShow dtrB ++ "]"
where
ufaDom =
- FA.composeThin ufa $
- Map.fromAscList $
- map mkToUnitUFA $
- DBox.toAscList dtrB
+ translateUfaToDom ufa dtrB
-- gr = 20 + (RA.getGranularity ufa)
- mkToUnitUFA (var, tr) =
- (var, UFA.affine [co] (Map.singleton var [sl]))
- where
- sl = FA.domra2ranra ufa $ dtrToUnitSlope tr
- co = FA.domra2ranra ufa $ dtrToUnitConst tr
- showVarDom (varID, varDom) =
- showVar varID ++ " -> " ++ show varDom
+
+translateUfaToDom ufa dtrB =
+ FA.composeThin ufa $
+ Map.fromAscList $
+ map mkToUnitUFA $
+ DBox.toAscList dtrB
+ where
+ mkToUnitUFA (var, tr) =
+ (var, UFA.affine [co] (Map.singleton var [sl]))
+ where
+ sl = FA.domra2ranra ufa $ dtrToUnitSlope tr
+ co = FA.domra2ranra ufa $ dtrToUnitConst tr
+
+
+instance
+ (UFA.ERUnitFnApprox box varid domra ranra ufa,
+ DomainBoxMappable dtrbox box varid (DomTransl domra) domra,
+ DomainBoxMappable box dtrbox varid domra (DomTransl domra),
+ H.HTML ufa) =>
+ H.HTML (ERFnDomTranslApprox dtrbox varid ufa domra)
+ where
+ toHtml (ERFnDomTranslApprox ufa dtrB) =
+ H.toHtml $ translateUfaToDom ufa dtrB
+-- H.toHtml $
+-- abovesTable []
+-- [
+-- H.toHtml $ "enclosure on dom = [" ++ dtrBShow dtrB ++ "]",
+-- H.toHtml $ translateUfaToDom ufa dtrB
+-- ]
+
instance
(UFA.ERUnitFnApprox box varid domra ranra ufa,
Eq dtrbox) =>
@@ -258,7 +294,7 @@ instance
instance
(UFA.ERUnitFnApprox box varid domra ranra ufa
- , DomainBoxMappable dtrbox box varid (DomTransl domra) domra, Eq dtrbox) =>
+ , DomainBoxMappable dtrbox box varid (DomTransl domra) domra, Eq dtrbox, Ord dtrbox) =>
RA.ERApprox (ERFnDomTranslApprox dtrbox varid ufa domra)
where
initialiseBaseArithmetic _ =
@@ -286,11 +322,18 @@ instance
vars = DBox.keys dtrB
leqReals fa1 fa2 =
RA.leqReals (erfnUnitApprox fa1) (erfnUnitApprox fa2)
-
+ refines fa1 fa2 =
+ RA.refines (erfnUnitApprox fa1) (erfnUnitApprox fa2)
+ compareApprox (ERFnDomTranslApprox ufa1 dtrB1) (ERFnDomTranslApprox ufa2 dtrB2) =
+ compareComposeMany
+ [
+ RA.compareApprox ufa1 ufa2,
+ compare dtrB1 dtrB2
+ ]
instance
(UFA.ERUnitFnApprox box varid domra ranra ufa, RA.ERIntApprox ufa
- , DomainBoxMappable dtrbox box varid (DomTransl domra) domra, Eq dtrbox) =>
+ , DomainBoxMappable dtrbox box varid (DomTransl domra) domra, Eq dtrbox, Ord dtrbox) =>
RA.ERIntApprox (ERFnDomTranslApprox dtrbox varid ufa domra)
where
-- doubleBounds = :: ira -> (Double, Double)
@@ -313,7 +356,7 @@ instance
instance
(UFA.ERUnitFnApprox box varid domra ranra ufa, RAEL.ERApproxElementary ufa
- , DomainBoxMappable dtrbox box varid (DomTransl domra) domra, Eq dtrbox) =>
+ , DomainBoxMappable dtrbox box varid (DomTransl domra) domra, Eq dtrbox, Ord dtrbox) =>
RAEL.ERApproxElementary (ERFnDomTranslApprox dtrbox varid ufa domra)
where
abs ix (ERFnDomTranslApprox ufa dtrB) =
@@ -331,10 +374,10 @@ instance
instance
(UFA.ERUnitFnApprox box varid domra ranra ufa,
- DomainIntBox box varid domra,
DomainBoxMappable dtrbox box varid (DomTransl domra) domra,
+ DomainIntBox box varid domra,
DomainBoxMappable box dtrbox varid domra (DomTransl domra),
- Eq dtrbox) =>
+ Eq dtrbox, Ord dtrbox) =>
FA.ERFnApprox box varid domra ranra (ERFnDomTranslApprox dtrbox varid ufa domra)
where
check prgLocation (ERFnDomTranslApprox ufa dtrB) =
@@ -345,6 +388,10 @@ instance
FA.ranra2domra (erfnUnitApprox fa) r
setMaxDegree maxDegree (ERFnDomTranslApprox ufa dtrB) =
ERFnDomTranslApprox (FA.setMaxDegree maxDegree ufa) dtrB
+ getMaxDegree (ERFnDomTranslApprox ufa _) =
+ FA.getMaxDegree ufa
+ getRangeApprox (ERFnDomTranslApprox ufa dtrB) =
+ FA.getRangeApprox ufa
volume (ERFnDomTranslApprox ufa dtrB) =
DBox.fold
(\tr vol -> vol * (FA.domra2ranra ufa $ dtrFromUnitSlope tr))
@@ -355,12 +402,12 @@ instance
(ERFnDomTranslApprox (FA.scale ratio ufa) dtrB)
partialIntersect ix substitutions f1 f2
| insideSubstitutions = f1 RA./\ f2
- | otherwise = f1
+ | otherwise = f2
where
insideSubstitutions =
and $ map snd $
DBox.zipWith (RA.refines) dom1 substitutions
- dom1 = FA.dom f1
+ dom1 = FA.dom f2
eval ptBox (ERFnDomTranslApprox ufa dtrB) =
FA.eval (domToUnit dtrB ptBox) ufa
partialEval substitutions (ERFnDomTranslApprox ufa dtrB) =
@@ -368,6 +415,7 @@ instance
where
dtrBNoVars =
DBox.difference dtrB substitutions
+
--instance
-- (UFA.ERUnitFnApprox box varid domra ranra ufa,
@@ -393,7 +441,7 @@ instance
DomainIntBox box varid domra,
DomainBoxMappable dtrbox box varid (DomTransl domra) domra,
DomainBoxMappable box dtrbox varid domra (DomTransl domra),
- Eq dtrbox) =>
+ Eq dtrbox, Ord dtrbox) =>
FA.ERFnDomApprox box varid domra ranra (ERFnDomTranslApprox dtrbox varid ufa domra)
where
dom (ERFnDomTranslApprox ufa dtrB) = dtrBToDomB dtrB
diff --git a/src/Data/Number/ER/RnToRm/Approx/PieceWise.hs b/src/Data/Number/ER/RnToRm/Approx/PieceWise.hs
index 6eec389..a5d34ec 100644
--- a/src/Data/Number/ER/RnToRm/Approx/PieceWise.hs
+++ b/src/Data/Number/ER/RnToRm/Approx/PieceWise.hs
@@ -41,6 +41,8 @@ import Data.Number.ER.Real.DomainBox (VariableID(..), DomainBox, DomainBoxMappab
import Data.Number.ER.BasicTypes
import Data.Number.ER.Misc
+import qualified Text.Html as H
+
import Data.Typeable
import Data.Generics.Basics
import Data.Binary
@@ -133,6 +135,13 @@ instance
show f@(ERFnPiecewise bistr) =
"\nERFnPiecewise:" ++ show bistr
+instance
+ (FA.ERFnDomApprox box varid domra ranra fa, RA.ERIntApprox fa, H.HTML fa) =>
+ H.HTML (ERFnPiecewise box varid domra fa)
+ where
+ toHtml (ERFnPiecewise bistr) =
+ H.toHtml bistr
+
instance
(FA.ERFnDomApprox box varid domra ranra fa) =>
Eq (ERFnPiecewise box varid domra fa)
@@ -202,7 +211,10 @@ instance
leqTuple [] = Just True
leqTuple _ =
error $ "ERFnTuple: leqReals not implemented"
-
+ refines f1@(ERFnPiecewise bistr1) f2@(ERFnPiecewise bistr2) =
+ and $ BISTR.collectValues $ pwbistrZipWith (RA.refines) 10 bistr1 bistr2
+ compareApprox f1@(ERFnPiecewise bistr1) f2@(ERFnPiecewise bistr2) =
+ BISTR.compare RA.compareApprox RA.compareApprox bistr1 bistr2
instance
(FA.ERFnDomApprox box varid domra ranra fa, RA.ERIntApprox fa, VariableID varid) =>
RA.ERIntApprox (ERFnPiecewise box varid domra fa)
@@ -287,17 +299,31 @@ instance
[FA.partialIntersect ix substitutions val1 val2]
updateTouch = updateInside
updateAway dom [val1, val2] =
- [val1]
+ [val2]
eval ptBox (ERFnPiecewise bistr) =
- foldl1 (zipWith (RA.\/)) $
+ foldl1 (zipWith (RA./\)) $
+-- unsafePrintReturn ("ERFnPiecewise: eval: vals = ")$
map (\fa -> FA.eval ptBox fa) $
- BISTR.collectValues $ BISTR.lookupSubtreeDom bistr ptBox
+ concat $ map BISTR.collectValues $
+ BISTR.lookupSubtreeDoms bistr ptBox
partialEval substitutions f@(ERFnPiecewise bistr) =
pwLift1 (FA.partialEval substitutions) (ERFnPiecewise bistrNoVars)
where
bistrNoVars =
BISTR.removeVars substitutions bistr
-
+ intersectMeasureImprovement ix f1@(ERFnPiecewise bistr1) f2@(ERFnPiecewise bistr2) =
+ (intersection, improvementRA)
+ where
+ (intersection, _) = RA.intersectMeasureImprovement ix f1 f2
+ improvementRA
+ | 0 `RA.refines` intersectionVolume && 0 `RA.refines` f1Volume = 1
+-- error $
+-- "ERFnInterval: intersectMeasureImprovement: inconsistent result: "
+-- ++ show intersection
+ | otherwise =
+ f1Volume / intersectionVolume
+ intersectionVolume = FA.volume intersection
+ f1Volume = FA.volume f1
instance
(FA.ERFnDomApprox box varid domra ranra fa, RA.ERIntApprox fa, Show box,
DomainBoxMappable box box varid domra domra) =>
@@ -323,9 +349,19 @@ instance
Just pt -> pt
unBisect var (ERFnPiecewise bistr1, ERFnPiecewise bistr2) =
ERFnPiecewise $
- BISTR.Node (depth1 - 1) dom var domVarMid bistr1 bistr2
+ BISTR.Node (depth - 1) dom var domVarMid bistr1Dp bistr2Dp
where
+ depth = max depth1 depth2
depth1 = BISTR.bistrDepth bistr1
+ depth2 = BISTR.bistrDepth bistr2
+ bistr1Dp
+ | depth1 == depth = bistr1
+ | otherwise =
+ BISTR.setDepth depth bistr1
+ bistr2Dp
+ | depth2 == depth = bistr2
+ | otherwise =
+ BISTR.setDepth depth bistr2
dom1 = BISTR.bistrDom bistr1
dom2 = BISTR.bistrDom bistr2
dom = DBox.unionWith (RA.\/) dom1 dom2
@@ -360,13 +396,13 @@ instance
(Just faLO, Nothing) ->
FA.partialIntersect ix
(DBox.singleton x domLO)
- bottomFn
(ERFnPiecewise (BISTR.Leaf depth domB faLO))
+ bottomFn
(Nothing, Just faHI) ->
FA.partialIntersect ix
(DBox.singleton x domHI)
- bottomFn
(ERFnPiecewise (BISTR.Leaf depth domB faHI))
+ bottomFn
bottomFn =
ERFnPiecewise $ BISTR.Leaf depth domB $ FA.bottomApprox domB (FA.getTupleSize fD)
(domLO, domHI) =
@@ -442,13 +478,13 @@ instance
(Just faLO, Nothing) ->
FA.partialIntersect ix
(DBox.singleton x domLO)
- (ERFnPiecewise bistrP)
(ERFnPiecewise (BISTR.Leaf depth domB faLO))
+ (ERFnPiecewise bistrP)
(Nothing, Just faHI) ->
FA.partialIntersect ix
(DBox.singleton x domHI)
- (ERFnPiecewise bistrP)
(ERFnPiecewise (BISTR.Leaf depth domB faHI))
+ (ERFnPiecewise bistrP)
(domLO, domHI) =
RA.bounds $
DBox.lookup "ERFnPieceWise: integrateMeasureImprovement: zipOutsideRange: " x domB
@@ -497,7 +533,7 @@ instance
(faIsect, faImpr) =
FA.integrateMeasureImprovement ix faD x integdomBox domLO faPadj
faPadj =
- FA.partialIntersect ix (DBox.singleton x domLO) faP faLO
+ FA.partialIntersect ix (DBox.singleton x domLO) faLO faP
faHI =
FA.partialEval (DBox.singleton x domHI) faIsect
(domLO, domHI) =
@@ -516,7 +552,7 @@ instance
(faIsect, faImpr) =
FA.integrateMeasureImprovement ix faD x integdomBox domHI faPadj
faPadj =
- FA.partialIntersect ix (DBox.singleton x domHI) faP faHI
+ FA.partialIntersect ix (DBox.singleton x domHI) faHI faP
faLO =
FA.partialEval (DBox.singleton x domLO) faIsect
(domLO, domHI) =
diff --git a/src/Data/Number/ER/RnToRm/Approx/Tuple.hs b/src/Data/Number/ER/RnToRm/Approx/Tuple.hs
index 9a43105..bb16ab4 100644
--- a/src/Data/Number/ER/RnToRm/Approx/Tuple.hs
+++ b/src/Data/Number/ER/RnToRm/Approx/Tuple.hs
@@ -27,6 +27,10 @@ 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.DomainBox as DBox
import Data.Number.ER.BasicTypes
+import Data.Number.ER.Misc
+
+import Data.Number.ER.ShowHTML
+import qualified Text.Html as H
import Data.Typeable
import Data.Generics.Basics
@@ -117,6 +121,18 @@ instance
showFA (fnname, fa) =
"\n>>> Function " ++ show fnname ++ ":" ++ show fa
+instance
+ (FA.ERFnDomApprox box varid domra ranra fa, H.HTML fa) =>
+ H.HTML (ERFnTuple fa)
+ where
+ toHtml (ERFnTuple []) = H.toHtml "[]"
+-- toHtml (ERFnTuple [f]) = H.toHtml f
+ toHtml (ERFnTuple fs) =
+ H.toHtml $
+ abovesTable [H.border 2] $
+-- [H.toHtml "Function tuple"] ++
+ map H.toHtml fs
+
instance
(FA.ERFnApprox box varid domra ranra fa) =>
Eq (ERFnTuple fa)
@@ -159,6 +175,8 @@ instance
setGranularity gran = tuplesLift1 (RA.setGranularity gran)
setMinGranularity gran = tuplesLift1 (RA.setMinGranularity gran)
f1 /\ f2 = tuplesLift2 "ERFnTuple: /\\: " (RA./\) f1 f2
+ refines f1@(ERFnTuple fas1) f2@(ERFnTuple fas2) =
+ and $ zipWith RA.refines fas1 fas2
intersectMeasureImprovement ix f1@(ERFnTuple fas1) f2@(ERFnTuple fas2)
| length fas1 == length fas2 =
(ERFnTuple fasIsect, ERFnTuple fasImpr)
@@ -174,7 +192,9 @@ instance
where
leqTuple [] = Just True
leqTuple _ =
- error $ "ERFnTuple: leqReals not implemented"
+ error $ "ERFnTuple: leqReals not implemented"
+ compareApprox f1@(ERFnTuple fas1) f2@(ERFnTuple fas2) =
+ compareListsWith RA.compareApprox fas1 fas2
instance
(FA.ERFnDomApprox box varid domra ranra fa, RA.ERIntApprox fa) =>
diff --git a/src/Data/Number/ER/RnToRm/BisectionTree.hs b/src/Data/Number/ER/RnToRm/BisectionTree.hs
index 093a854..c73059f 100644
--- a/src/Data/Number/ER/RnToRm/BisectionTree.hs
+++ b/src/Data/Number/ER/RnToRm/BisectionTree.hs
@@ -27,6 +27,7 @@ module Data.Number.ER.RnToRm.BisectionTree
removeVars,
sync2,
syncMany,
+ setDepth,
split,
mapWithDom,
mapLeaves,
@@ -36,7 +37,8 @@ module Data.Number.ER.RnToRm.BisectionTree
combineWith,
collectValues,
collectDomValues,
- lookupSubtreeDom
+ compare,
+ lookupSubtreeDoms
)
where
@@ -50,6 +52,9 @@ import Data.Number.ER.BasicTypes
import Data.Number.ER.Misc
+import Data.Number.ER.ShowHTML
+import qualified Text.Html as H
+
import Data.Typeable
import Data.Generics.Basics
import Data.Binary
@@ -127,6 +132,31 @@ showBisectionTree showValue =
showVD (v,d) =
showVar v ++ "->" ++ show d
+instance (Show d, H.HTML v, DomainBox box varid d) =>
+ H.HTML (BisectionTree box varid d v)
+ where
+ toHtml (Leaf depth dom val) =
+ H.toHtmlFromList $
+ [
+ H.toHtml $ concatWith "," (Prelude.map showVD $ DBox.toList dom)
+ ,
+ H.primHtml " &rarr; "
+ ,
+ H.toHtml val
+ ]
+ where
+ showVD (v,d) =
+ showVar v ++ " in " ++ show d
+ toHtml (Node depth dom dir pt lo hi) =
+ H.toHtml $
+ besidesTable [H.border 2]
+ [
+ abovesTable [] [H.toHtml $ "(" ++ showVar dir ++ ")"]
+ ,
+ abovesTable [] [H.toHtml lo, H.toHtml hi]
+ ]
+
+
isLeaf ::
BisectionTree box varid d v ->
Bool
@@ -153,6 +183,23 @@ type ValueCombiner box varid d v =
(EffortIndex -> Depth -> (BisectionTree box varid d v) -> v)
+setDepth ::
+ Depth ->
+ BisectionTree box varid d v ->
+ BisectionTree box varid d v
+setDepth depth bistr
+ | isLeaf bistr =
+ bistr { bistrDepth = depth }
+ | otherwise =
+ bistr
+ {
+ bistrLO = setDepth depthInc $ bistrLO bistr,
+ bistrHI = setDepth depthInc $ bistrHI bistr
+ }
+ where
+ depthInc = depth + 1
+
+
split ::
(RA.ERIntApprox d, DomainBox box varid d) =>
ValueSplitter box varid d v ->
@@ -559,41 +606,55 @@ collectDomValues (Node _ _ _ _ cLO cHI) =
linear ordering on bisection trees
-}
compare ::
- (Ord d, Ord varid) =>
+ (Ord varid, DomainBox box varid d) =>
+ (d -> d -> Ordering) ->
(v -> v -> Ordering) ->
(BisectionTree box varid d v) ->
(BisectionTree box varid d v) ->
Ordering
-compare compValues (Leaf _ _ _) (Node _ _ _ _ _ _) = LT
-compare compValues (Node _ _ _ _ _ _) (Leaf _ _ _) = GT
-compare compValues (Leaf _ _ val1) (Leaf _ _ val2) =
- compValues val1 val2
-compare compValues
- (Node _ _ dir1 pt1 lo1 hi1)
- (Node _ _ dir2 pt2 lo2 hi2) =
- compareComposeMany $
- [Prelude.compare dir1 dir2,
- Prelude.compare pt1 pt2,
- compare compValues lo1 lo2,
- compare compValues hi1 hi2]
+compare compDoms compValues (Leaf _ _ _) (Node _ _ _ _ _ _) = LT
+compare compDoms compValues (Node _ _ _ _ _ _) (Leaf _ _ _) = GT
+compare compDoms compValues (Leaf depth1 dom1 val1) (Leaf depth2 dom2 val2) =
+ compareComposeMany
+ [
+ Prelude.compare depth1 depth2,
+ DBox.compare compDoms dom1 dom2,
+ compValues val1 val2
+ ]
+compare compDoms compValues
+ (Node depth1 dom1 dir1 pt1 lo1 hi1)
+ (Node depth2 dom2 dir2 pt2 lo2 hi2) =
+ compareComposeMany
+ [
+ Prelude.compare dir1 dir2,
+ compDoms pt1 pt2,
+ compare compDoms compValues lo1 lo2,
+ compare compDoms compValues hi1 hi2
+ ]
{-|
- lookup the smallest subtree whose domain covers the given rectangle
+ lookup all maximal subtrees whose domain intersect the given rectangle
-}
-lookupSubtreeDom ::
+lookupSubtreeDoms ::
(RA.ERIntApprox d, DomainBox box varid d) =>
(BisectionTree box varid d v) ->
box {-^ domain to look up within the tree -} ->
- (BisectionTree box varid d v)
-lookupSubtreeDom origBistr dom =
+ [BisectionTree box varid d v]
+lookupSubtreeDoms origBistr dom =
lk origBistr
where
- lk bistr@(Leaf _ _ _) = bistr
+ lk bistr@(Leaf _ _ _) = [bistr]
lk bistr@(Node _ _ _ _ lo hi)
- | and $ Prelude.map snd $ DBox.zipWithDefault RA.bottomApprox (RA.refines) dom domHI = lk hi
- | and $ Prelude.map snd $ DBox.zipWithDefault RA.bottomApprox (RA.refines) dom domLO = lk lo
- | otherwise = bistr
+ | loDisjoint = lk hi
+ | hiDisjoint = lk lo
+ | otherwise = (lk lo) ++ (lk hi)
where
+ loDisjoint =
+ and $ Prelude.map snd $
+ DBox.zipWithDefault RA.bottomApprox (RA.isDisjoint) dom domLO
+ hiDisjoint =
+ and $ Prelude.map snd $
+ DBox.zipWithDefault RA.bottomApprox (RA.isDisjoint) dom domHI
domLO = bistrDom lo
domHI = bistrDom hi
diff --git a/src/Data/Number/ER/RnToRm/UnitDom/Approx/Interval.hs b/src/Data/Number/ER/RnToRm/UnitDom/Approx/Interval.hs
index 160b17e..48f9f69 100644
--- a/src/Data/Number/ER/RnToRm/UnitDom/Approx/Interval.hs
+++ b/src/Data/Number/ER/RnToRm/UnitDom/Approx/Interval.hs
@@ -41,6 +41,9 @@ import Data.Number.ER.BasicTypes
import Data.Number.ER.Misc
+import Data.Number.ER.ShowHTML
+import qualified Text.Html as H
+
import qualified Data.Map as Map
import Data.Typeable
@@ -120,6 +123,24 @@ instance
instance
(UFB.ERUnitFnBase boxb boxra varid b ra fb) =>
+ H.HTML (ERFnInterval fb ra)
+ where
+ toHtml (ERFnIntervalAny ctxt) =
+ H.toHtml "ERFnIntervalAny"
+ toHtml (ERFnInterval h ln ctxt gl) =
+-- H.toHtml $
+-- abovesTable
+-- [
+-- H.toHtml "ERFnInterval",
+ H.toHtml $ H.simpleTable [H.border 2] []
+ [
+ [H.toHtml "upper = ", H.toHtml $ show h],
+ [H.toHtml "lower = ", H.toHtml $ show (- ln)]
+ ]
+-- ]
+
+instance
+ (UFB.ERUnitFnBase boxb boxra varid b ra fb) =>
Eq (ERFnInterval fb ra)
where
(ERFnInterval h1 ln1 ctxt1 gl1)
@@ -135,9 +156,9 @@ instance
compare
(ERFnInterval h1 ln1 ctxt1 gl1)
(ERFnInterval h2 ln2 ctxt2 gl2) =
- error "ERFnInterval: comparison not implemented; consider leqReals from class ERApprox instead"
+ error "ERFnInterval: comparison not implemented; consider leqReals or compareApprox from class ERApprox instead"
compare _ _ =
- error "ERFnInterval: comparison not implemented; consider leqReals from class ERApprox instead"
+ error "ERFnInterval: comparison not implemented; consider leqReals or compareApprox from class ERApprox instead"
instance
@@ -149,6 +170,7 @@ instance
negate (ERFnInterval h ln ctxt gl) =
(ERFnInterval ln h ctxt (negate gl))
(ERFnInterval h1 ln1 ctxt1 gl1) + (ERFnInterval h2 ln2 ctxt2 gl2) =
+ normalise $
ERFnInterval (h1 + h2) (ln1 + ln2) ctxt (gl1 + gl2)
where
ctxt = erfnContextUnify ctxt1 ctxt2
@@ -156,6 +178,7 @@ instance
where
ctxt = erfnContextUnify (erfnContext f1) (erfnContext f2)
(ERFnInterval h1 ln1 ctxt1 gl1) * (ERFnInterval h2 ln2 ctxt2 gl2) =
+ normalise $
ERFnInterval h ln ctxt (gl1 * gl2)
where
(h, ln) =
@@ -196,6 +219,7 @@ instance
recip f@(ERFnIntervalAny _) = f
recip (ERFnInterval h ln ctxt gl)
| certainNoZero =
+ normalise $
ERFnInterval lRecipUp hnRecipUp ctxt (recip gl)
| otherwise = ERFnIntervalAny ctxt
where
@@ -237,7 +261,11 @@ instance
(ctxt { erfnCoeffGranularity = max gran (erfnCoeffGranularity ctxt) }) gl
-- getPrecision (ERFnIntervalAny _) = 0
-- getPrecision f = intLog 2 (1 + (fst $ RA.integerBounds (FA.volume f))) -- wrong!
- (ERFnInterval h1 ln1 ctxt1 gl1) /\ (ERFnInterval h2 ln2 ctxt2 gl2) =
+ f1@(ERFnInterval h1 ln1 ctxt1 gl1) /\ f2@(ERFnInterval h2 ln2 ctxt2 gl2) =
+---- #ifdef RUNTIME_CHECKS
+---- check ("ERFnInterval: /\\:\n f1:\n" ++ show f1 ++ " f2:\n" ++ show f2 ++ "\n result:\n") $
+---- #endif
+ normalise $
ERFnInterval (snd $ UFB.min maxDegr h1 h2) (snd $ UFB.min maxDegr ln1 ln2) ctxt (gl1 RA./\ gl2)
where
ctxt = erfnContextUnify ctxt1 ctxt2
@@ -254,6 +282,21 @@ instance
where
ctxt = erfnContextUnify (erfnContext f1) (erfnContext f2)
leqReals = erfnintLeq
+ refines _ (ERFnIntervalAny _) = True
+ refines (ERFnIntervalAny _) _ = False
+ refines (ERFnInterval h1 ln1 _ _) (ERFnInterval h2 ln2 _ _) =
+ (UFB.upperBound 10 (h2 - h1) >= 0)
+ &&
+ (UFB.upperBound 10 (ln2 - ln1) >= 0)
+ compareApprox (ERFnIntervalAny _) (ERFnIntervalAny _) = EQ
+ compareApprox (ERFnIntervalAny _) _ = LT
+ compareApprox _ (ERFnIntervalAny _) = GT
+ compareApprox (ERFnInterval h1 ln1 _ _) (ERFnInterval h2 ln2 _ _) =
+ compareComposeMany
+ [
+ UFB.compareApprox h1 h2,
+ UFB.compareApprox ln1 ln2
+ ]
erfnintLeq left right
| left `isClearlyBelow` right = Just True
@@ -295,7 +338,11 @@ instance
bounds (ERFnInterval u ln c g) =
(ERFnInterval (negate ln) ln c g,
ERFnInterval u (negate u) c g)
- (ERFnInterval u1 ln1 c1 g1) \/ (ERFnInterval u2 ln2 c2 g2) =
+ f1@(ERFnInterval u1 ln1 c1 g1) \/ f2@(ERFnInterval u2 ln2 c2 g2) =
+---- #ifdef RUNTIME_CHECKS
+---- check ("ERFnInterval: abs:\n f1:\n" ++ show f1 ++ " f2:\n" ++ show f2 ++ "\n result:\n") $
+---- #endif
+ normalise $
ERFnInterval u ln c (g1 RA.\/ g2)
where
u = UFB.maxUp maxDegree u1 u2
@@ -320,7 +367,11 @@ instance
where
-- default abs does not work because we do not have Prelude.abs
abs _ f@(ERFnIntervalAny _) = f
- abs _ (ERFnInterval u ln c g) =
+ abs _ f@(ERFnInterval u ln c g) =
+---- #ifdef RUNTIME_CHECKS
+---- check ("ERFnInterval: abs:\n f:\n" ++ show f ++ "\n result:\n") $
+---- #endif
+ normalise $
ERFnInterval maxulnUp maxunl0Dn c (abs g)
where
maxDegree = erfnMaxDegree c
@@ -329,17 +380,17 @@ instance
fst $ UFB.max maxDegree 0 $
fst $ UFB.max maxDegree (- u) (- ln)
exp ix f@(ERFnIntervalAny _) = f
- exp ix (ERFnInterval u ln c g) =
+ exp ix f@(ERFnInterval u ln c g) =
+ normalise $
ERFnInterval uExp lExpNeg c (RAEL.exp ix g)
where
maxDegree = erfnMaxDegree c
--- ix = int2effIx maxDegree
uExp = snd $ UFB.exp maxDegree ix u
lExpNeg =
negate $ fst $ UFB.exp maxDegree ix (negate ln)
sin ix f@(ERFnIntervalAny c) =
ERFnInterval 1 1 c ((-1) RA.\/ 1)
- sin ix (ERFnInterval u ln c g) =
+ sin ix f@(ERFnInterval u ln c g) =
-- unsafePrint
-- (
-- "ERFnInterval: RAEL.sin: "
@@ -348,13 +399,17 @@ instance
-- ++ "\n uSin = " ++ show uSin
-- ++ "\n lSinNeg = " ++ show lSinNeg
-- ) $
+---- #ifdef RUNTIME_CHECKS
+---- check ("ERFnInterval: sin:\n f:\n" ++ show f ++ "\n result:\n") $
+---- #endif
+ normalise $
ERFnInterval uSin (- lSin) c (RAEL.sin ix g)
where
(lSin, uSin) = sincos True maxDegree ix u (-ln)
maxDegree = erfnMaxDegree c
cos ix f@(ERFnIntervalAny c) =
ERFnInterval 1 1 c ((-1) RA.\/ 1)
- cos ix (ERFnInterval u ln c g) =
+ cos ix f@(ERFnInterval u ln c g) =
-- unsafePrint
-- (
-- "ERFnInterval: RAEL.cos: "
@@ -363,13 +418,14 @@ instance
-- ++ "\n uCos = " ++ show uCos
-- ++ "\n lCosNeg = " ++ show lCosNeg
-- ) $
+ normalise $
ERFnInterval uCos (- lCos) c (RAEL.cos ix g)
where
(lCos, uCos) = sincos False maxDegree ix u (-ln)
maxDegree = erfnMaxDegree c
atan ix f@(ERFnIntervalAny c) =
ERFnInterval 1 1 c ((-1) RA.\/ 1)
- atan ix (ERFnInterval u ln c g) =
+ atan ix f@(ERFnInterval u ln c g) =
-- unsafePrint
-- (
-- "ERFnInterval: RAEL.atan: "
@@ -378,6 +434,7 @@ instance
-- ++ "\n uAtan = " ++ show uAtan
-- ++ "\n lAtanNeg = " ++ show lAtanNeg
-- ) $
+ normalise $
ERFnInterval uAtan lAtanNeg c (RAEL.atan ix g)
where
maxDegree = erfnMaxDegree c
@@ -514,18 +571,25 @@ sincos isSine maxDegree ix u l
(pLO `plusDown` (- errPoly), pHI + errPoly)
where
errPoly = UFB.const errB
+
+normalise f@(ERFnIntervalAny c) = f
+normalise f@(ERFnInterval u ln c g)
+ | UFB.isValid u && UFB.isValid ln = f
+ | otherwise = ERFnIntervalAny c
+check callerLocation f@(ERFnIntervalAny c) = f
+check callerLocation f@(ERFnInterval u ln c g) =
+ ERFnInterval
+ (UFB.check (callerLocation ++ "upper: ") u)
+ (UFB.check (callerLocation ++ "neg lower: ") ln)
+ c g
+
instance
(UFB.ERUnitFnBase boxb boxra varid b ra fb) =>
FA.ERFnApprox boxra varid ra ra (ERFnInterval fb ra)
where
- check callerLocation f@(ERFnIntervalAny c) = f
- check callerLocation (ERFnInterval u ln c g) =
- ERFnInterval
- (UFB.check (callerLocation ++ "upper: ") u)
- (UFB.check (callerLocation ++ "neg lower: ") ln)
- c g
+ check = check
domra2ranra _ = id
ranra2domra _ = id
setMaxDegree maxDegr (ERFnIntervalAny c) =
@@ -536,16 +600,26 @@ instance
(UFB.reduceDegreeUp maxDegr ln)
(c { erfnMaxDegree = maxDegr } )
g
- getRangeApprox (ERFnIntervalAny _) = RA.bottomApprox
+ getMaxDegree (ERFnIntervalAny c) =
+ erfnMaxDegree c
+ getMaxDegree (ERFnInterval _ _ c _) =
+ erfnMaxDegree c
+ getRangeApprox (ERFnIntervalAny _) =
+ RA.bottomApprox
getRangeApprox (ERFnInterval u ln c g) =
UFB.raFromEndpoints u
(
- (UFB.upperBound 10 u)
- ,
(- (UFB.upperBound 10 ln))
+ ,
+ (UFB.upperBound 10 u)
)
- scale ratio f@(ERFnIntervalAny c) = f
+ scale ratio f@(ERFnIntervalAny c) =
+ f
scale ratio f@(ERFnInterval u ln c g) =
+---- #ifdef RUNTIME_CHECKS
+---- FA.check ("ERFnInterval: scale:\n before:\n" ++ show f ++ "\n after:\n") $
+---- #endif
+ normalise $
case RA.compareReals ratio 0 of
Just GT ->
ERFnInterval (UFB.scaleApproxUp ratio u) (UFB.scaleApproxUp ratio ln) c g
@@ -560,7 +634,8 @@ instance
up = UFB.evalApprox ptBox u
lo = negate $ UFB.evalApprox ptBox ln
partialEval substitutions f@(ERFnIntervalAny c) = f
- partialEval substitutions (ERFnInterval u ln c g) =
+ partialEval substitutions f@(ERFnInterval u ln c g) =
+ normalise $
(ERFnInterval uP lnP c g)
where
uP = UFB.partialEvalApproxUp substitutions u
@@ -589,6 +664,10 @@ instance
ERFnIntervalAny erfnContextDefault
const [val]
| RA.isBounded val =
+---- #ifdef RUNTIME_CHECKS
+---- check ("ERFnInterval: const:\n") $
+---- #endif
+ normalise $
ERFnInterval
{
erfnUpper = fbH,
@@ -609,6 +688,10 @@ instance
}
affine [val] coeffsSingletons
| RA.isBounded val && (and $ map (RA.isBounded . head) $ Map.elems coeffsSingletons) =
+---- #ifdef RUNTIME_CHECKS
+---- check ("ERFnInterval: affine:\n") $
+---- #endif
+ normalise $
ERFnInterval
{
erfnUpper = fbH,
@@ -665,7 +748,12 @@ instance
Just LT -> (f1, 1) -- intersection made it worse, keep original
_ -> (intersection, improvementRA)
where
- intersection = f1 RA./\ f2
+ intersection =
+---- #ifdef RUNTIME_CHECKS
+---- check ("ERFnInterval: intersectMeasureImprovement:\n f1:\n" ++ show f1 ++ "\n f2:\n" ++ show f2 ++ "\n intersection:\n") $
+---- #endif
+ normalise $
+ f1 RA./\ f2
improvementRA
| 0 `RA.refines` intersectionVolume && 0 `RA.refines` f1Volume = 1
-- error $
@@ -682,9 +770,10 @@ instance
-- where
-- result =
UFB.raFromEndpoints u $ UFB.volumeAboveZero vars (u + ln)
+ integrate _ f@(ERFnIntervalAny c) _ _ _ = f
integrate
- ix (ERFnInterval u ln c g) x
- origin (ERFnInterval uInit lnInit cInit gInit) =
+ ix fD@(ERFnInterval u ln c g) x
+ origin fI@(ERFnInterval uInit lnInit cInit gInit) =
-- unsafePrint
-- (
-- "ERFnInterval: integrate: "
@@ -704,6 +793,10 @@ instance
-- ++ "\n uIov = " ++ show uIov
-- ++ "\n lnIov = " ++ show lnIov
-- )
+---- #ifdef RUNTIME_CHECKS
+---- check ("ERFnInterval: integrate:\n fD:\n" ++ show fD ++ "\n fI:\n" ++ show fI ++ "\n result:\n") $
+---- #endif
+ normalise $
(ERFnInterval uIov lnIov c gIov)
where
-- perform raw integration of both bounds:
diff --git a/src/Data/Number/ER/RnToRm/UnitDom/Base.hs b/src/Data/Number/ER/RnToRm/UnitDom/Base.hs
index 293c2d1..46d4752 100644
--- a/src/Data/Number/ER/RnToRm/UnitDom/Base.hs
+++ b/src/Data/Number/ER/RnToRm/UnitDom/Base.hs
@@ -39,7 +39,15 @@ class
where
initialiseBaseArithmetic :: ufb -> IO ()
initialiseBaseArithmetic _ =
- B.initialiseBaseArithmetic (0 :: b)
+ B.initialiseBaseArithmetic (0 :: b)
+ {-|
+ Check internal consistency, typically absence of NaN.
+ -}
+ isValid :: ufb -> Bool
+ {-|
+ A linear ordering, which can be syntactic and rather arbitrary.
+ -}
+ compareApprox :: ufb -> ufb -> Ordering
{-|
Check internal consistency of the function and report problem if any.
-}
diff --git a/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom.hs b/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom.hs
index 9304cf4..a100bab 100644
--- a/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom.hs
+++ b/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom.hs
@@ -59,7 +59,9 @@ instance
DomainIntBox boxra varid (ERInterval rb)) =>
(UFB.ERUnitFnBase boxb boxra varid rb (ERInterval rb) (ERChebPoly box rb))
where
+ isValid = chplHasNoNaNOrInfty
check = chplCheck
+ compareApprox = chplCompareApprox
getGranularity = chplGetGranularity
setMinGranularity = chplSetMinGranularity
setGranularity = chplSetGranularity
diff --git a/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Basic.hs b/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Basic.hs
index 4ede2ec..e720865 100644
--- a/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Basic.hs
+++ b/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Basic.hs
@@ -47,23 +47,36 @@ data ERChebPoly box b =
chplCoeffs :: (Map.Map (TermKey box) b)
}
deriving (Eq, Typeable, Data)
+
+type TermKey box = box
instance (Ord a, Binary a, Binary b) => Binary (ERChebPoly a b) where
put (ERChebPoly a) = put a
get = get >>= \a -> return (ERChebPoly a)
-
-chplCheck prgLocation p@(ERChebPoly coeffs)
- | ok = p
- | otherwise =
- unsafePrint (prgLocation ++ " problem with p = \n" ++ show p) p
+chplHasNoNaN p@(ERChebPoly coeffs) =
+ Map.fold (&&) True $ Map.map coeffOK coeffs
where
- ok =
- Map.fold (&&) True $ Map.map coeffOK coeffs
coeffOK c =
not $ B.isERNaN c
-
-type TermKey box = box
+
+chplHasNoNaNOrInfty p@(ERChebPoly coeffs) =
+ Map.fold (&&) True $ Map.map coeffOK coeffs
+ where
+ coeffOK = isBounded
+ isBounded c
+ | B.isERNaN c = False
+ | B.isPlusInfinity c = False
+ | B.isPlusInfinity (- c) = False
+ | otherwise = True
+
+chplCheck prgLocation p
+ | chplHasNoNaNOrInfty p = p
+ | otherwise =
+ unsafePrint (prgLocation ++ " problem with p = \n" ++ show p) p
+
+chplCompareApprox (ERChebPoly coeffs1) (ERChebPoly coeffs2) =
+ compare coeffs1 coeffs2
chplConstTermKey :: (DomainBox box varid d) => box
chplConstTermKey = DBox.noinfo
@@ -88,6 +101,20 @@ chplGetDegree ::
chplGetDegree (ERChebPoly coeffs) =
foldl max 0 $ map chplTermOrder $ Map.keys coeffs
+{-|
+ If the polynomial is constant, return the constant,
+ otherwise return Nothing.
+-}
+chplGetConst ::
+ (B.ERRealBase b, DomainBox box varid d, Num d, Ord d) =>
+ (ERChebPoly box b) ->
+ Maybe b
+chplGetConst (ERChebPoly coeffs) =
+ case Map.keys coeffs of
+ [key] | chplIsConstTermKey key ->
+ Just $ head $ Map.elems coeffs
+ _ -> Nothing
+
-- chplGetArity = length . chplGetVars
diff --git a/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Elementary.hs b/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Elementary.hs
index b3ad053..cbebc67 100644
--- a/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Elementary.hs
+++ b/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Elementary.hs
@@ -57,6 +57,17 @@ chplExp ::
ERChebPoly box b ->
(ERChebPoly box b, ERChebPoly box b)
chplExp maxDegree ix p =
+-- unsafePrint
+-- (
+-- "chplExp:" ++
+-- "\n expM = " ++ show expM ++
+-- "\n upperB = " ++ show upperB ++
+-- "\n lowerB = " ++ show lowerB ++
+-- "\n a_int = " ++ show a_int ++
+-- "\n expNear0Dn pNear0Dn = " ++ show (expNear0Dn pNear0Dn) ++
+-- "\n chplPow maxDegree (expNear0Up pNear0Up) 2000 = " ++ show (chplPow maxDegree (expNear0Up pNear0Up) 2000)
+-- )
+-- $
(expDownwards, expUpwards + valueAtRDnNeg + (chplConst expRUp))
where
expUpwards =
@@ -69,7 +80,8 @@ chplExp maxDegree ix p =
r = (upperB - lowerB) / 2
expMUp = erintv_right expM
expMDn = erintv_left expM
- expM = erExp_R ix (ERInterval m m)
+ expM =
+ erExp_R ix (ERInterval m m)
pNear0Up = (p - (chplConst m)) * (chplConst $ recip a_base)
pNear0Dn = - (((-p) + (chplConst m)) * (chplConst $ recip a_base))
a_base = fromInteger a_int
@@ -449,7 +461,9 @@ chplRecip ::
EffortIndex {-^ minimum approx degree -} ->
ERChebPoly box b ->
(ERChebPoly box b, ERChebPoly box b)
-chplRecip maxDegree ix p
+chplRecip maxDegree ix p@(ERChebPoly pCoeffs)
+ | pIsConst =
+ (chplConst $ - (recip (- pConst)), chplConst $ recip pConst)
| upperB < 0 = -- range negative
(\(lo, hi) -> (-hi, -lo)) $ chplRecip maxDegree ix (negate p)
| lowerB > 0 = -- range positive
@@ -470,7 +484,12 @@ chplRecip maxDegree ix p
ranp = ERInterval lowerB upperB
lowerB = - (chplUpperBoundAffine ix (- p))
upperB = chplUpperBoundAffine ix p
-
+
+ (pIsConst, pConst) =
+ case chplGetConst p of
+ Nothing -> (False, error "ChebyshevBase.Polynom.Elementary.chplRecip")
+ Just coeff -> (True, coeff)
+
tauDegree = effIx2int (max 2 $ ix `div` 3)
coeffGr = effIx2gran $ ix
diff --git a/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Field.hs b/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Field.hs
index 8fd15f1..0603f9b 100644
--- a/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Field.hs
+++ b/src/Data/Number/ER/RnToRm/UnitDom/ChebyshevBase/Polynom/Field.hs
@@ -86,7 +86,7 @@ instance (B.ERRealBase b, DomainBox box varid Int, Ord box) => Num (ERChebPoly b
(Map.unionWith (+) coeffs1 coeffs2)
-- point-wise sum of polynomials with coeffs rounded upwards
coeffsDown =
- (Map.unionWith (\c1 c2 -> - ((- c1) + (- c2))) coeffs1 coeffs2)
+ (Map.unionWith plusDown coeffs1 coeffs2)
-- point-wise sum of polynomials with coeffs rounded upwards
maxError =
Map.fold (+) 0 $
@@ -199,6 +199,8 @@ chplScaleApprox (ratioDown, ratioUp) (ERChebPoly coeffs) =
(ratioDown `timesDown` coeff, ratioUp `timesUp` coeff)
| coeff < 0 =
(ratioUp `timesDown` coeff, ratioDown `timesUp` coeff)
+ | otherwise =
+ error $ "chplScaleApprox: " ++ show coeff
instance (B.ERRealBase b, DomainBox box varid Int, Ord box) => Fractional (ERChebPoly box b)
@@ -212,7 +214,7 @@ instance (B.ERRealBase b, DomainBox box varid Int, Ord box) => Fractional (ERChe
instance (B.ERRealBase b, DomainBox box varid Int, Ord box) => Ord (ERChebPoly box b)
where
compare _ _ =
- errorModule "cannot compare polynomials, consider using RA.leqReals instead"
+ errorModule "cannot compare polynomials, consider using leqReals or compareApprox instead"
--instance (B.ERRealBase b, DomainBox box varid Int, Ord box) => Real (ERChebPoly box b)
-- where
diff --git a/tests/ISin3.hs b/tests/ISin3.hs
index 2413708..baaa75a 100644
--- a/tests/ISin3.hs
+++ b/tests/ISin3.hs
@@ -35,7 +35,7 @@ main =
head $
AERNFunc.eval (AERNFunc.unary 1) $
AERNFunc.integrateUnary 0
- (sin3 60 0 21 60)
+ (sin3 120 0 42 160)
0 (0 AERN.\/ 1) [0]
sin3 ix depth deg gran =
AERN.sin ix $