{-|
    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)