AERN-RnToRm-0.3.0: src/Data/Number/ER/RnToRm/BisectionTree/Integration.hs
{-|
Module : Data.Number.ER.RnToRm.BisectionTree.Integration
Description : abstract zipping of domain partitions used for integration
Copyright : (c) 2007-2008 Michal Konecny
License : BSD3
Maintainer : mik@konecny.aow.cz
Stability : experimental
Portability : portable
To be imported qualified, usually with prefix BTINTEG.
-}
module Data.Number.ER.RnToRm.BisectionTree.Integration
(
zipFromOrigin, zipOnSubdomain
)
where
import qualified Data.Number.ER.RnToRm.BisectionTree as BISTR
import qualified Data.Number.ER.Real.Approx as RA
import qualified Data.Number.ER.Real.DomainBox as DBox
import Data.Number.ER.Real.DomainBox (VariableID(..), DomainBox, DomainIntBox)
import Data.Number.ER.BasicTypes
import Data.Number.ER.Misc
--import qualified Data.Sequence as Seq
--import qualified Data.Map as Map
import Data.Maybe
{-|
Transform a bunch of bisection trees over the same domain
by "integrating" them in a very abstract sense.
The trees are unified in their splitting patterns in the process.
By supplying certain parameters, this function can in fact
perform numerical integration of piece-wise polynomial functions.
It can be also viewed as a "zipping+folding" operator over bisection trees that
generates another bunch of bisection trees, synchronously traversing the original trees
from a certain point on a selected axis outwards in both directions,
carrying some data along.
-}
zipFromOrigin ::
(RA.ERIntApprox d, DomainIntBox box varid d, Show v1, Show v2, Show valPass) =>
BISTR.ValueSplitter box varid d v1 ->
BISTR.ValueCombiner box varid d v1 ->
EffortIndex ->
varid
{-^ variable @x@ (ie axis or direction) to integrate in -} ->
d
{-^ origin in terms of variable @x@ -} ->
(Maybe (box))
{-^ support, ie the domain on which to zip
(automatically extended to include origin when projected to @x@) -} ->
(Maybe valPass -> Maybe valPass -> [BISTR.BisectionTree box varid d v1] -> [BISTR.BisectionTree box varid d v2])
{-^ what to do outside the support,
possibly being passed values from left/right
when leaving the support -} ->
(EffortIndex -> BISTR.Depth -> (box) -> [v1] -> [v2] -> Bool)
{-^ should a leaf be split? -} ->
(EffortIndex -> BISTR.Depth -> (box) -> [v1] -> (valPass,[v2],valPass))
{-^ integrator for a leaf containing the origin -} ->
(EffortIndex -> BISTR.Depth -> (box) -> valPass -> [v1] -> ([v2], valPass))
{-^ integrator over a leaf that sees the origin towards -infinity -} ->
(EffortIndex -> BISTR.Depth -> (box) -> [v1] -> valPass -> (valPass, [v2]))
{-^ integrator over a leaf that sees the origin towards +infinity -} ->
[BISTR.BisectionTree box varid d v1]
{-^ input functions -} ->
[BISTR.BisectionTree box varid d v2]
{-^ output functions
The number of output functions does not have to be
the same as the number of input functions.
-}
zipFromOrigin
valSplitter valCombiner ix
ivar origin maybeResultSupport outerValTransformer
decideShouldSplit integrLeafOH integrLeafOL integrLeafOR
bistrs =
resultBistrs
where
(_, resultBistrs, _) =
integrateBistrOriginHere $ BISTR.syncMany valSplitter ix bistrs
maybeSupport = -- extend resultSupport to cover the origin
fmap extendToOrigin maybeResultSupport
where
extendToOrigin domB =
case DBox.member ivar domB of
True -> DBox.insertWith (RA.\/) ivar origin domB
False -> domB
-- the following function is used when we know the origin is within the current sub-domain:
integrateBistrOriginHere bistrs@((BISTR.Leaf depth dom _) : _)
| decideShouldSplit ix depth dom vals integrVals = -- must descend
integrateBistrOriginHere $
map (BISTR.split valSplitter ix var pt dom) bistrs
| otherwise =
(Just lVal, map (\v -> BISTR.Leaf depth dom v) integrVals, Just rVal)
where
(var, pt) = DBox.bestSplit dom
vals = map BISTR.bistrVal bistrs
(lVal, integrVals, rVal) =
integrLeafOH ix depth dom vals
integrateBistrOriginHere bistrs@((BISTR.Node depth dom var pt lBounds rBounds):_)
| origin `RA.refines` rDom =
-- unsafePrint
-- ("BTINTEG: integrateBistrOriginHere: rDom = " ++ show rDom ++
-- " origin = " ++ show origin ++
-- " lValHI = " ++ show lValHI ++
-- " rValHI = " ++ show rValHI)
(lValHI, bistrsIntgHI, rValHI)
| origin `RA.refines` lDom =
-- unsafePrint
-- ("BTINTEG: integrateBistrOriginHere: lDom = " ++ show lDom ++
-- " origin = " ++ show origin ++
-- " lValLO = " ++ show lValLO ++
-- " rValLO = " ++ show rValLO)
(lValLO, bistrsIntgLO, rValLO)
| otherwise = -- origin overlaps both sides
-- have to amalgamate these trees:
integrateBistrOriginHere $
map (\b -> BISTR.Leaf depth dom (valCombiner ix depth b)) bistrs
where
lDom = DBox.lookup "BTINTEG: zipFromOrigin: Here: L: " var (BISTR.bistrDom lBounds)
rDom = DBox.lookup "BTINTEG: zipFromOrigin: Here: R: " var (BISTR.bistrDom rBounds)
-- recursion when origin is entirely to the right of the centre:
bistrsIntgHI =
zipWith
(\lo hi -> BISTR.Node depth dom var pt lo hi)
lBoundsIntgHI rBoundsIntgHI
(lValHIHI, rBoundsIntgHI, rValHI) =
integrateBistrOriginHere $
BISTR.syncMany valSplitter ix $ map BISTR.bistrHI bistrs
(lValHI, lBoundsIntgHI) =
integrateBistrOriginRight
(BISTR.syncMany valSplitter ix $ map BISTR.bistrLO bistrs)
lValHIHI
-- recursion when origin is entirely to the left of the centre:
bistrsIntgLO =
zipWith
(\lo hi -> BISTR.Node depth dom var pt lo hi)
lBoundsIntgLO rBoundsIntgLO
(lValLO, lBoundsIntgLO, rValLOLO) =
integrateBistrOriginHere $
BISTR.syncMany valSplitter ix $ map BISTR.bistrLO bistrs
(rBoundsIntgLO, rValLO) =
integrateBistrOriginLeft
rValLOLO
(BISTR.syncMany valSplitter ix $ map BISTR.bistrHI bistrs)
-- the following function is used when we know
-- the origin is to the LEFT of the current sub-domain:
integrateBistrOriginLeft Nothing bistrs =
-- previously detected as being outside the support
(outerValTransformer Nothing Nothing bistrs, Nothing)
integrateBistrOriginLeft (Just lVal) bistrs@(bistr:_)
| (isJust maybeSupport) &&
(and $ Prelude.map snd $
DBox.zipWithDefaultSecond RA.bottomApprox RA.isInteriorDisjoint
(BISTR.bistrDom bistr)
(fromJust maybeSupport)) =
-- outside the integration domain
(outerValTransformer (Just lVal) Nothing bistrs, Nothing)
integrateBistrOriginLeft (Just lVal) bistrs@((BISTR.Leaf depth dom _) : _)
| decideShouldSplit ix depth dom vals integrVals = -- improve granularity by splitting
integrateBistrOriginLeft (Just lVal) $
map (BISTR.split valSplitter ix var pt dom) bistrs
| otherwise =
(map (\v -> BISTR.Leaf depth dom v) integrVals,
Just rVal)
where
(var, pt) = DBox.bestSplit dom
vals = map BISTR.bistrVal bistrs
(integrVals, rVal) =
integrLeafOL ix depth dom lVal vals
integrateBistrOriginLeft mlVal bistrs@((BISTR.Node depth dom var pt _ _):_) =
(bistrsIntg, mrVal2)
where
bistrsIntg =
zipWith (\lo hi -> BISTR.Node depth dom var pt lo hi) lBoundsINT rBoundsINT
(lBoundsINT, mrVal1) =
integrateBistrOriginLeft mlVal $
BISTR.syncMany valSplitter ix $ map BISTR.bistrLO bistrs
(rBoundsINT, mrVal2) =
integrateBistrOriginLeft mrVal1 $
BISTR.syncMany valSplitter ix $ map BISTR.bistrHI bistrs
-- -- the following function is used when we know
-- -- the origin is to the RIGHT of the current sub-domain:
integrateBistrOriginRight bistrs Nothing =
-- previously detected as being outside the support
(Nothing, outerValTransformer Nothing Nothing bistrs)
integrateBistrOriginRight bistrs@(bistr:_) (Just rVal)
| (isJust maybeSupport) &&
(and $ Prelude.map snd $
DBox.zipWithDefaultSecond RA.bottomApprox RA.isInteriorDisjoint
(BISTR.bistrDom bistr)
(fromJust maybeSupport)) =
-- outside the integration domain
(Nothing, outerValTransformer Nothing (Just rVal) bistrs)
integrateBistrOriginRight bistrs@((BISTR.Leaf depth dom _) : _) (Just rVal)
| decideShouldSplit ix depth dom vals integrVals = -- improve granularity by splitting
integrateBistrOriginRight
(map (BISTR.split valSplitter ix var pt dom) bistrs)
(Just rVal)
| otherwise =
(Just lVal,
map (\v -> BISTR.Leaf depth dom v) integrVals)
where
(var, pt) = DBox.bestSplit dom
vals = map BISTR.bistrVal bistrs
(lVal, integrVals) =
integrLeafOR ix depth dom vals rVal
integrateBistrOriginRight bistrs@((BISTR.Node depth dom var pt _ _):_) mrVal =
(mlVal2, bistrsIntg)
where
bistrsIntg =
zipWith (\lo hi -> BISTR.Node depth dom var pt lo hi) lBoundsINT rBoundsINT
(mlVal2, lBoundsINT) =
integrateBistrOriginRight
(BISTR.syncMany valSplitter ix $ map BISTR.bistrLO bistrs) mlVal1
(mlVal1, rBoundsINT) =
integrateBistrOriginRight
(BISTR.syncMany valSplitter ix $ map BISTR.bistrHI bistrs) mrVal
{-|
Zip a list of bisection trees in synchrony but do something
else inside and not inside a given subdomain.
Further splitting at default points will be done up to the given depth
in an attempt to separate the subdomain as well as possible.
If the subdomain is not properly isolated by the splitting at the
maximum depth, splits are made at irregular points to ensure full isolation
of the subdomain.
-}
zipOnSubdomain ::
(RA.ERIntApprox d, DomainIntBox box varid d) =>
BISTR.ValueSplitter box varid d v1 ->
EffortIndex ->
BISTR.Depth
{-^ depth limit -} ->
box
{-^ subdomain @sd@ -} ->
(box -> [v1] -> [v2])
{-^ what to do with values /inside/ @sd@ -} ->
(box -> [v1] -> [v2])
{-^ what to do with values /outside/ @sd@ but /touching/ it -} ->
(box -> [v1] -> [v2])
{-^ what to do with values /outside/ @sd@ -} ->
[BISTR.BisectionTree box varid d v1] ->
[BISTR.BisectionTree box varid d v2]
zipOnSubdomain valSplitter ix maxDepth sdom updateInside updateTouch updateAway bistrs =
resultBistrs
where
resultBistrs =
zz $ BISTR.syncMany valSplitter ix bistrs
zz bistrs@(BISTR.Leaf depth dom _ : _)
| intersect =
case depth < maxDepth of
True ->
zz $ map (BISTR.split valSplitter ix var pt dom) bistrs
False ->
error "BTINTEG: zipOnSubdomain: maxDepth reached but irregular splitting not implemented yet"
| away = lift updateAway
| touch = lift updateTouch
| inside = lift updateInside
where
(var, pt) = DBox.bestSplit dom
lift updateFn =
map (BISTR.Leaf depth dom) $
updateFn dom $
map BISTR.bistrVal bistrs
(away, touch, intersect, inside) =
DBox.classifyPosition dom sdom
zz bistrs@(BISTR.Node depth dom var pt _ _ : _) =
zipWith
(\bLO bHI -> BISTR.Node depth dom var pt bLO bHI)
(zz $ map BISTR.bistrLO bistrs)
(zz $ map BISTR.bistrHI bistrs)