{-# OPTIONS_GHC -Wall #-}
-----------------------------------------------------------------------------
-- |
-- Module      :  ToySolver.Arith.FourierMotzkin.Optimization
-- Copyright   :  (c) Masahiro Sakai 2014-2015
-- License     :  BSD-style
--
-- Maintainer  :  masahiro.sakai@gmail.com
-- Stability   :  provisional
-- Portability :  portable
--
-- Naïve implementation of Fourier-Motzkin Variable Elimination
--
-- Reference:
--
-- * <http://users.cecs.anu.edu.au/~michaeln/pubs/arithmetic-dps.pdf>
--
-----------------------------------------------------------------------------
module ToySolver.Arith.FourierMotzkin.Optimization
  ( optimize
  ) where

import Control.Exception (assert)
import Control.Monad
import Data.Maybe
import qualified Data.IntMap as IntMap
import qualified Data.IntSet as IntSet
import Data.ExtendedReal
import qualified Data.Interval as Interval
import Data.Interval (Interval, (<=..<), (<..<=))
import Data.OptDir

import ToySolver.Data.DNF
import qualified ToySolver.Data.LA as LA
import ToySolver.Data.IntVar
import ToySolver.Arith.FourierMotzkin.Base

-- | @optimize dir obj φ@ returns @(I, lift)@ where
--
-- * @I@ is convex hull of feasible region, and
--
-- * @lift@ is a function, that takes @x ∈ I@ and returns the feasible solution with objective value /better than or equal to/ @x@.
--
-- Note:
--
-- * @'Interval.lowerBound' i@ (resp. @'Interval.upperBound' i@) is the optimal value in case of minimization (resp. maximization).
--
-- * If @I@ is empty, then the problem is INFEASIBLE.
--
optimize :: VarSet -> OptDir -> LA.Expr Rational -> [LA.Atom Rational] -> (Interval Rational, Rational -> Model Rational)
optimize :: VarSet
-> OptDir
-> Expr Rational
-> [Atom Rational]
-> (Interval Rational, Rational -> Model Rational)
optimize VarSet
vs OptDir
dir Expr Rational
obj [Atom Rational]
cs = (Interval Rational
ival, Rational -> Model Rational
m)
  where
    rs :: [(Interval Rational, Rational -> Model Rational)]
rs = [VarSet
-> Expr Rational
-> [Constr]
-> (Interval Rational, Rational -> Model Rational)
projectToObj' VarSet
vs Expr Rational
obj [Constr]
cs' | [Constr]
cs' <- DNF Constr -> [[Constr]]
forall lit. DNF lit -> [[lit]]
unDNF ([Atom Rational] -> DNF Constr
constraintsToDNF [Atom Rational]
cs)]
    ival :: Interval Rational
ival = [Interval Rational] -> Interval Rational
forall r. Ord r => [Interval r] -> Interval r
Interval.hulls (((Interval Rational, Rational -> Model Rational)
 -> Interval Rational)
-> [(Interval Rational, Rational -> Model Rational)]
-> [Interval Rational]
forall a b. (a -> b) -> [a] -> [b]
map (Interval Rational, Rational -> Model Rational)
-> Interval Rational
forall a b. (a, b) -> a
fst [(Interval Rational, Rational -> Model Rational)]
rs)

    m :: Rational -> Model Rational
    m :: Rational -> Model Rational
m Rational
x = Maybe (Model Rational) -> Model Rational
forall a. HasCallStack => Maybe a -> a
fromJust (Maybe (Model Rational) -> Model Rational)
-> Maybe (Model Rational) -> Model Rational
forall a b. (a -> b) -> a -> b
$ [Maybe (Model Rational)] -> Maybe (Model Rational)
forall (t :: * -> *) (m :: * -> *) a.
(Foldable t, MonadPlus m) =>
t (m a) -> m a
msum ([Maybe (Model Rational)] -> Maybe (Model Rational))
-> [Maybe (Model Rational)] -> Maybe (Model Rational)
forall a b. (a -> b) -> a -> b
$ ((Interval Rational, Rational -> Model Rational)
 -> Maybe (Model Rational))
-> [(Interval Rational, Rational -> Model Rational)]
-> [Maybe (Model Rational)]
forall a b. (a -> b) -> [a] -> [b]
map (Interval Rational, Rational -> Model Rational)
-> Maybe (Model Rational)
f [(Interval Rational, Rational -> Model Rational)]
rs
      where
        f :: (Interval Rational, Rational -> Model Rational) -> Maybe (Model Rational)
        f :: (Interval Rational, Rational -> Model Rational)
-> Maybe (Model Rational)
f (Interval Rational
i1,Rational -> Model Rational
m1) = do
          Rational
x' <- Interval Rational -> Maybe Rational
forall r. RealFrac r => Interval r -> Maybe Rational
Interval.simplestRationalWithin (Interval Rational -> Interval Rational -> Interval Rational
forall r. Ord r => Interval r -> Interval r -> Interval r
Interval.intersection Interval Rational
i1 Interval Rational
ib)
          Model Rational -> Maybe (Model Rational)
forall (m :: * -> *) a. Monad m => a -> m a
return (Rational -> Model Rational
m1 Rational
x')

        ib :: Interval Rational
ib = case OptDir
dir of
               OptDir
OptMin -> Extended Rational
forall r. Extended r
NegInf Extended Rational -> Extended Rational -> Interval Rational
forall r. Ord r => Extended r -> Extended r -> Interval r
<..<= Rational -> Extended Rational
forall r. r -> Extended r
Finite Rational
x
               OptDir
OptMax -> Rational -> Extended Rational
forall r. r -> Extended r
Finite Rational
x Extended Rational -> Extended Rational -> Interval Rational
forall r. Ord r => Extended r -> Extended r -> Interval r
<=..< Extended Rational
forall r. Extended r
PosInf

projectToObj' :: VarSet -> LA.Expr Rational -> [Constr] -> (Interval Rational, Rational -> Model Rational)
projectToObj' :: VarSet
-> Expr Rational
-> [Constr]
-> (Interval Rational, Rational -> Model Rational)
projectToObj' VarSet
vs Expr Rational
obj [Constr]
cs = VarSet
-> [Constr]
-> Var
-> (Interval Rational, Rational -> Model Rational)
projectTo' VarSet
vs (Rat -> Rat -> Constr
eqR (Expr Rational -> Rat
toRat (Var -> Expr Rational
forall r. Num r => Var -> Expr r
LA.var Var
z)) (Expr Rational -> Rat
toRat Expr Rational
obj) Constr -> [Constr] -> [Constr]
forall a. a -> [a] -> [a]
: [Constr]
cs) Var
z
  where
    z :: Var
z = Var -> Maybe Var -> Var
forall a. a -> Maybe a -> a
fromMaybe Var
0 (((Var, VarSet) -> Var) -> Maybe (Var, VarSet) -> Maybe Var
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap ((Var -> Var -> Var
forall a. Num a => a -> a -> a
+Var
1) (Var -> Var) -> ((Var, VarSet) -> Var) -> (Var, VarSet) -> Var
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (Var, VarSet) -> Var
forall a b. (a, b) -> a
fst) (VarSet -> Maybe (Var, VarSet)
IntSet.maxView VarSet
vs))

projectTo' :: VarSet -> [Constr] -> Var -> (Interval Rational, Rational -> Model Rational)
projectTo' :: VarSet
-> [Constr]
-> Var
-> (Interval Rational, Rational -> Model Rational)
projectTo' VarSet
vs [Constr]
cs Var
z = (Interval Rational, Rational -> Model Rational)
-> Maybe (Interval Rational, Rational -> Model Rational)
-> (Interval Rational, Rational -> Model Rational)
forall a. a -> Maybe a -> a
fromMaybe (Interval Rational
forall r. Ord r => Interval r
Interval.empty, \Rational
_ -> [Char] -> Model Rational
forall a. HasCallStack => [Char] -> a
error [Char]
"ToySolver.FourierMotzkin.projectTo': should not be called") (Maybe (Interval Rational, Rational -> Model Rational)
 -> (Interval Rational, Rational -> Model Rational))
-> Maybe (Interval Rational, Rational -> Model Rational)
-> (Interval Rational, Rational -> Model Rational)
forall a b. (a -> b) -> a -> b
$ do
  ([Constr]
ys,Model Rational -> Model Rational
mt) <- VarSet
-> [Constr] -> Maybe ([Constr], Model Rational -> Model Rational)
projectN' VarSet
vs ([Constr] -> Maybe ([Constr], Model Rational -> Model Rational))
-> Maybe [Constr]
-> Maybe ([Constr], Model Rational -> Model Rational)
forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< [Constr] -> Maybe [Constr]
simplify [Constr]
cs
  let (Bounds
bs,[Constr]
ws) = Var -> [Constr] -> (Bounds, [Constr])
collectBounds Var
z [Constr]
ys
  Bool -> Maybe () -> Maybe ()
forall a. HasCallStack => Bool -> a -> a
assert ([Constr] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [Constr]
ws) (Maybe () -> Maybe ()) -> Maybe () -> Maybe ()
forall a b. (a -> b) -> a -> b
$ () -> Maybe ()
forall (m :: * -> *) a. Monad m => a -> m a
return ()
  let ival :: Interval Rational
ival = Model Rational -> Bounds -> Interval Rational
evalBounds Model Rational
forall a. IntMap a
IntMap.empty Bounds
bs
  (Interval Rational, Rational -> Model Rational)
-> Maybe (Interval Rational, Rational -> Model Rational)
forall (m :: * -> *) a. Monad m => a -> m a
return (Interval Rational
ival, \Rational
v -> Model Rational -> Model Rational
mt (Var -> Rational -> Model Rational
forall a. Var -> a -> IntMap a
IntMap.singleton Var
z Rational
v))