{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE MultiParamTypeClasses #-}
{-# LANGUAGE TypeSynonymInstances #-}
{-# OPTIONS_GHC -Wall #-}
{-# OPTIONS_HADDOCK show-extensions #-}
-----------------------------------------------------------------------------
-- |
-- Module      :  ToySolver.Arith.OmegaTest.Base
-- Copyright   :  (c) Masahiro Sakai 2011
-- License     :  BSD-style
--
-- Maintainer  :  masahiro.sakai@gmail.com
-- Stability   :  provisional
-- Portability :  non-portable
--
-- (incomplete) implementation of Omega Test
--
-- References:
--
-- * William Pugh. The Omega test: a fast and practical integer
--   programming algorithm for dependence analysis. In Proceedings of
--   the 1991 ACM/IEEE conference on Supercomputing (1991), pp. 4-13.
--
-- * <http://users.cecs.anu.edu.au/~michaeln/pubs/arithmetic-dps.pdf>
--
-- See also:
--
-- * <http://hackage.haskell.org/package/Omega>
--
-----------------------------------------------------------------------------
module ToySolver.Arith.OmegaTest.Base
    ( Model
    , solve
    , solveQFLIRAConj
    , Options (..)
    , checkRealNoCheck
    , checkRealByFM

    -- * Exported just for testing
    , zmod
    ) where

import Control.Exception (assert)
import Control.Monad
import Data.Default.Class
import Data.List
import Data.Maybe
import Data.Ord
import Data.Ratio
import qualified Data.IntMap as IM
import qualified Data.IntSet as IS
import Data.VectorSpace

import ToySolver.Data.OrdRel
import ToySolver.Data.Boolean
import ToySolver.Data.DNF
import qualified ToySolver.Data.LA as LA
import ToySolver.Data.IntVar
import ToySolver.Internal.Util (combineMaybe)
import qualified ToySolver.Arith.FourierMotzkin as FM

-- ---------------------------------------------------------------------------

-- | Options for solving.
--
-- The default option can be obtained by 'def'.
data Options
  = Options
  { Options -> VarSet -> [Atom (Ratio Integer)] -> Bool
optCheckReal :: VarSet -> [LA.Atom Rational] -> Bool
  -- ^ optCheckReal is used for checking whether real shadow is satisfiable.
  --
  -- * If it returns @True@, the real shadow may or may not be satisfiable.
  --
  -- * If it returns @False@, the real shadow must be unsatisfiable
  }

instance Default Options where
  def :: Options
def =
    Options
    { optCheckReal :: VarSet -> [Atom (Ratio Integer)] -> Bool
optCheckReal =
        -- checkRealNoCheck
        VarSet -> [Atom (Ratio Integer)] -> Bool
checkRealByFM
    }

checkRealNoCheck :: VarSet -> [LA.Atom Rational] -> Bool
checkRealNoCheck :: VarSet -> [Atom (Ratio Integer)] -> Bool
checkRealNoCheck VarSet
_ [Atom (Ratio Integer)]
_ = Bool
True

checkRealByFM :: VarSet -> [LA.Atom Rational] -> Bool
checkRealByFM :: VarSet -> [Atom (Ratio Integer)] -> Bool
checkRealByFM VarSet
vs [Atom (Ratio Integer)]
as = forall a. Maybe a -> Bool
isJust forall a b. (a -> b) -> a -> b
$ VarSet -> [Atom (Ratio Integer)] -> Maybe (Model (Ratio Integer))
FM.solve VarSet
vs [Atom (Ratio Integer)]
as

-- ---------------------------------------------------------------------------

type ExprZ = LA.Expr Integer

-- | Atomic constraint
data Constr
  = IsNonneg ExprZ
  -- ^ e ≥ 0
  | IsZero ExprZ
  -- ^ e = 0
  deriving (Int -> Constr -> ShowS
[Constr] -> ShowS
Constr -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [Constr] -> ShowS
$cshowList :: [Constr] -> ShowS
show :: Constr -> String
$cshow :: Constr -> String
showsPrec :: Int -> Constr -> ShowS
$cshowsPrec :: Int -> Constr -> ShowS
Show, Constr -> Constr -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: Constr -> Constr -> Bool
$c/= :: Constr -> Constr -> Bool
== :: Constr -> Constr -> Bool
$c== :: Constr -> Constr -> Bool
Eq, Eq Constr
Constr -> Constr -> Bool
Constr -> Constr -> Ordering
Constr -> Constr -> Constr
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
min :: Constr -> Constr -> Constr
$cmin :: Constr -> Constr -> Constr
max :: Constr -> Constr -> Constr
$cmax :: Constr -> Constr -> Constr
>= :: Constr -> Constr -> Bool
$c>= :: Constr -> Constr -> Bool
> :: Constr -> Constr -> Bool
$c> :: Constr -> Constr -> Bool
<= :: Constr -> Constr -> Bool
$c<= :: Constr -> Constr -> Bool
< :: Constr -> Constr -> Bool
$c< :: Constr -> Constr -> Bool
compare :: Constr -> Constr -> Ordering
$ccompare :: Constr -> Constr -> Ordering
Ord)

instance Variables Constr where
  vars :: Constr -> VarSet
vars (IsNonneg ExprZ
t) = forall a. Variables a => a -> VarSet
vars ExprZ
t
  vars (IsZero ExprZ
t) = forall a. Variables a => a -> VarSet
vars ExprZ
t

-- 制約集合の単純化
-- It returns Nothing when a inconsistency is detected.
simplify :: [Constr] -> Maybe [Constr]
simplify :: [Constr] -> Maybe [Constr]
simplify = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap forall (t :: * -> *) a. Foldable t => t [a] -> [a]
concat forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall (t :: * -> *) (m :: * -> *) a b.
(Traversable t, Monad m) =>
(a -> m b) -> t a -> m (t b)
mapM Constr -> Maybe [Constr]
f
  where
    f :: Constr -> Maybe [Constr]
    f :: Constr -> Maybe [Constr]
f constr :: Constr
constr@(IsNonneg ExprZ
e) =
      case forall r. Num r => Expr r -> Maybe r
LA.asConst ExprZ
e of
        Just Integer
x -> forall (f :: * -> *). Alternative f => Bool -> f ()
guard (Integer
x forall a. Ord a => a -> a -> Bool
>= Integer
0) forall (m :: * -> *) a b. Monad m => m a -> m b -> m b
>> forall (m :: * -> *) a. Monad m => a -> m a
return []
        Maybe Integer
Nothing -> forall (m :: * -> *) a. Monad m => a -> m a
return [Constr
constr]
    f constr :: Constr
constr@(IsZero ExprZ
e) =
      case forall r. Num r => Expr r -> Maybe r
LA.asConst ExprZ
e of
        Just Integer
x -> forall (f :: * -> *). Alternative f => Bool -> f ()
guard (Integer
x forall a. Eq a => a -> a -> Bool
== Integer
0) forall (m :: * -> *) a b. Monad m => m a -> m b -> m b
>> forall (m :: * -> *) a. Monad m => a -> m a
return []
        Maybe Integer
Nothing -> forall (m :: * -> *) a. Monad m => a -> m a
return [Constr
constr]

leZ, ltZ, geZ, gtZ, eqZ :: ExprZ -> ExprZ -> Constr
-- Note that constants may be floored by division
leZ :: ExprZ -> ExprZ -> Constr
leZ ExprZ
e1 ExprZ
e2 = ExprZ -> Constr
IsNonneg (forall b a. (Num b, Eq b) => (a -> b) -> Expr a -> Expr b
LA.mapCoeff (forall a. Integral a => a -> a -> a
`div` Integer
d) ExprZ
e)
  where
    e :: ExprZ
e = ExprZ
e2 forall v. AdditiveGroup v => v -> v -> v
^-^ ExprZ
e1
    d :: Integer
d = forall a. Num a => a -> a
abs forall a b. (a -> b) -> a -> b
$ [Integer] -> Integer
gcd' [Integer
c | (Integer
c,Int
v) <- forall r. Expr r -> [(r, Int)]
LA.terms ExprZ
e, Int
v forall a. Eq a => a -> a -> Bool
/= Int
LA.unitVar]
ltZ :: ExprZ -> ExprZ -> Constr
ltZ ExprZ
e1 ExprZ
e2 = (ExprZ
e1 forall v. AdditiveGroup v => v -> v -> v
^+^ forall r. (Num r, Eq r) => r -> Expr r
LA.constant Integer
1) ExprZ -> ExprZ -> Constr
`leZ` ExprZ
e2
geZ :: ExprZ -> ExprZ -> Constr
geZ = forall a b c. (a -> b -> c) -> b -> a -> c
flip ExprZ -> ExprZ -> Constr
leZ
gtZ :: ExprZ -> ExprZ -> Constr
gtZ = forall a b c. (a -> b -> c) -> b -> a -> c
flip ExprZ -> ExprZ -> Constr
ltZ
eqZ :: ExprZ -> ExprZ -> Constr
eqZ ExprZ
e1 ExprZ
e2 =
  case ExprZ -> Maybe Constr
isZero (ExprZ
e1 forall v. AdditiveGroup v => v -> v -> v
^-^ ExprZ
e2) of
    Just Constr
c -> Constr
c
    Maybe Constr
Nothing -> ExprZ -> Constr
IsZero (forall r. (Num r, Eq r) => r -> Expr r
LA.constant Integer
1) -- unsatisfiable

isZero :: ExprZ -> Maybe Constr
isZero :: ExprZ -> Maybe Constr
isZero ExprZ
e
  = if forall r. Num r => Int -> Expr r -> r
LA.coeff Int
LA.unitVar ExprZ
e forall a. Integral a => a -> a -> a
`mod` Integer
d forall a. Eq a => a -> a -> Bool
== Integer
0
    then forall a. a -> Maybe a
Just forall a b. (a -> b) -> a -> b
$ ExprZ -> Constr
IsZero forall a b. (a -> b) -> a -> b
$ forall b a. (Num b, Eq b) => (a -> b) -> Expr a -> Expr b
LA.mapCoeff (forall a. Integral a => a -> a -> a
`div` Integer
d) ExprZ
e
    else forall a. Maybe a
Nothing
  where
    d :: Integer
d = forall a. Num a => a -> a
abs forall a b. (a -> b) -> a -> b
$ [Integer] -> Integer
gcd' [Integer
c | (Integer
c,Int
v) <- forall r. Expr r -> [(r, Int)]
LA.terms ExprZ
e, Int
v forall a. Eq a => a -> a -> Bool
/= Int
LA.unitVar]

applySubst1Constr :: Var -> ExprZ -> Constr -> Constr
applySubst1Constr :: Int -> ExprZ -> Constr -> Constr
applySubst1Constr Int
v ExprZ
e (IsNonneg ExprZ
e2) = forall r. (Num r, Eq r) => Int -> Expr r -> Expr r -> Expr r
LA.applySubst1 Int
v ExprZ
e ExprZ
e2 ExprZ -> ExprZ -> Constr
`geZ` forall r. (Num r, Eq r) => r -> Expr r
LA.constant Integer
0
applySubst1Constr Int
v ExprZ
e (IsZero ExprZ
e2)   = forall r. (Num r, Eq r) => Int -> Expr r -> Expr r -> Expr r
LA.applySubst1 Int
v ExprZ
e ExprZ
e2 ExprZ -> ExprZ -> Constr
`eqZ` forall r. (Num r, Eq r) => r -> Expr r
LA.constant Integer
0

instance Eval (Model Integer) Constr Bool where
  eval :: Model Integer -> Constr -> Bool
eval Model Integer
m (IsNonneg ExprZ
t) = forall m e v. Eval m e v => m -> e -> v
LA.eval Model Integer
m ExprZ
t forall a. Ord a => a -> a -> Bool
>= Integer
0
  eval Model Integer
m (IsZero ExprZ
t)   = forall m e v. Eval m e v => m -> e -> v
LA.eval Model Integer
m ExprZ
t forall a. Eq a => a -> a -> Bool
== Integer
0

-- ---------------------------------------------------------------------------

-- | (t,c) represents t/c, and c must be >0.
type Rat = (ExprZ, Integer)

{-
(ls,us) represents
{ x | ∀(M,c)∈ls. M/c≤x, ∀(M,c)∈us. x≤M/c }
-}
type BoundsZ = ([Rat],[Rat])

evalBoundsZ :: Model Integer -> BoundsZ -> IntervalZ
evalBoundsZ :: Model Integer -> BoundsZ -> IntervalZ
evalBoundsZ Model Integer
model ([Rat]
ls,[Rat]
us) =
  forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl' IntervalZ -> IntervalZ -> IntervalZ
intersectZ IntervalZ
univZ forall a b. (a -> b) -> a -> b
$
    [ (forall a. a -> Maybe a
Just (forall a b. (RealFrac a, Integral b) => a -> b
ceiling (forall m e v. Eval m e v => m -> e -> v
LA.eval Model Integer
model ExprZ
x forall a. Integral a => a -> a -> Ratio a
% Integer
c)), forall a. Maybe a
Nothing) | (ExprZ
x,Integer
c) <- [Rat]
ls ] forall a. [a] -> [a] -> [a]
++
    [ (forall a. Maybe a
Nothing, forall a. a -> Maybe a
Just (forall a b. (RealFrac a, Integral b) => a -> b
floor (forall m e v. Eval m e v => m -> e -> v
LA.eval Model Integer
model ExprZ
x forall a. Integral a => a -> a -> Ratio a
% Integer
c))) | (ExprZ
x,Integer
c) <- [Rat]
us ]

collectBoundsZ :: Var -> [Constr] -> (BoundsZ, [Constr])
collectBoundsZ :: Int -> [Constr] -> (BoundsZ, [Constr])
collectBoundsZ Int
v = forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
foldr Constr -> (BoundsZ, [Constr]) -> (BoundsZ, [Constr])
phi (([],[]),[])
  where
    phi :: Constr -> (BoundsZ,[Constr]) -> (BoundsZ,[Constr])
    phi :: Constr -> (BoundsZ, [Constr]) -> (BoundsZ, [Constr])
phi constr :: Constr
constr@(IsNonneg ExprZ
t) (([Rat]
ls,[Rat]
us),[Constr]
xs) =
      case forall r. Num r => Int -> Expr r -> (r, Expr r)
LA.extract Int
v ExprZ
t of
        (Integer
c,ExprZ
t') ->
          case Integer
c forall a. Ord a => a -> a -> Ordering
`compare` Integer
0 of
            Ordering
EQ -> (([Rat]
ls, [Rat]
us), Constr
constr forall a. a -> [a] -> [a]
: [Constr]
xs)
            Ordering
GT -> (((forall v. AdditiveGroup v => v -> v
negateV ExprZ
t', Integer
c) forall a. a -> [a] -> [a]
: [Rat]
ls, [Rat]
us), [Constr]
xs) -- 0 ≤ cx + M ⇔ -M/c ≤ x
            Ordering
LT -> (([Rat]
ls, (ExprZ
t', forall a. Num a => a -> a
negate Integer
c) forall a. a -> [a] -> [a]
: [Rat]
us), [Constr]
xs)   -- 0 ≤ cx + M ⇔ x ≤ M/-c
    phi constr :: Constr
constr@(IsZero ExprZ
_) (BoundsZ
bnd,[Constr]
xs) = (BoundsZ
bnd, Constr
constr forall a. a -> [a] -> [a]
: [Constr]
xs) -- we assume v does not appear in constr

-- ---------------------------------------------------------------------------

solve' :: Options -> VarSet -> [Constr] -> Maybe (Model Integer)
solve' :: Options -> VarSet -> [Constr] -> Maybe (Model Integer)
solve' Options
opt VarSet
vs2 [Constr]
xs = [Constr] -> Maybe [Constr]
simplify [Constr]
xs forall (m :: * -> *) a b. Monad m => m a -> (a -> m b) -> m b
>>= VarSet -> [Constr] -> Maybe (Model Integer)
go VarSet
vs2
  where
    go :: VarSet -> [Constr] -> Maybe (Model Integer)
    go :: VarSet -> [Constr] -> Maybe (Model Integer)
go VarSet
vs [Constr]
cs | VarSet -> Bool
IS.null VarSet
vs = do
      let m :: Model Integer
          m :: Model Integer
m = forall a. IntMap a
IM.empty
      forall (f :: * -> *). Alternative f => Bool -> f ()
guard forall a b. (a -> b) -> a -> b
$ forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
all (forall m e v. Eval m e v => m -> e -> v
eval Model Integer
m) [Constr]
cs
      forall (m :: * -> *) a. Monad m => a -> m a
return Model Integer
m
    go VarSet
vs [Constr]
cs | Just (ExprZ
e,[Constr]
cs2) <- [Constr] -> Maybe (ExprZ, [Constr])
extractEq [Constr]
cs = do
      (VarSet
vs',[Constr]
cs3, Model Integer -> Model Integer
mt) <- ExprZ
-> VarSet
-> [Constr]
-> Maybe (VarSet, [Constr], Model Integer -> Model Integer)
eliminateEq ExprZ
e VarSet
vs [Constr]
cs2
      Model Integer
m <- VarSet -> [Constr] -> Maybe (Model Integer)
go VarSet
vs' forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< [Constr] -> Maybe [Constr]
simplify [Constr]
cs3
      forall (m :: * -> *) a. Monad m => a -> m a
return (Model Integer -> Model Integer
mt Model Integer
m)
    go VarSet
vs [Constr]
cs = do
      forall (f :: * -> *). Alternative f => Bool -> f ()
guard forall a b. (a -> b) -> a -> b
$ Options -> VarSet -> [Atom (Ratio Integer)] -> Bool
optCheckReal Options
opt VarSet
vs (forall a b. (a -> b) -> [a] -> [b]
map Constr -> Atom (Ratio Integer)
toLAAtom [Constr]
cs)
      if BoundsZ -> Bool
isExact BoundsZ
bnd then
        Maybe (Model Integer)
case1
      else
        Maybe (Model Integer)
case1 forall (m :: * -> *) a. MonadPlus m => m a -> m a -> m a
`mplus` Maybe (Model Integer)
case2

      where
        (Int
v,VarSet
vs',bnd :: BoundsZ
bnd@([Rat]
ls,[Rat]
us),[Constr]
rest) = VarSet -> [Constr] -> (Int, VarSet, BoundsZ, [Constr])
chooseVariable VarSet
vs [Constr]
cs

        case1 :: Maybe (Model Integer)
case1 = do
          let zs :: [Constr]
zs = [ forall r. (Num r, Eq r) => r -> Expr r
LA.constant ((Integer
aforall a. Num a => a -> a -> a
-Integer
1)forall a. Num a => a -> a -> a
*(Integer
bforall a. Num a => a -> a -> a
-Integer
1)) ExprZ -> ExprZ -> Constr
`leZ` (Integer
a forall v. VectorSpace v => Scalar v -> v -> v
*^ ExprZ
d forall v. AdditiveGroup v => v -> v -> v
^-^ Integer
b forall v. VectorSpace v => Scalar v -> v -> v
*^ ExprZ
c)
                   | (ExprZ
c,Integer
a)<-[Rat]
ls , (ExprZ
d,Integer
b)<-[Rat]
us ]
          Model Integer
model <- VarSet -> [Constr] -> Maybe (Model Integer)
go VarSet
vs' forall (m :: * -> *) a b. Monad m => (a -> m b) -> m a -> m b
=<< [Constr] -> Maybe [Constr]
simplify ([Constr]
zs forall a. [a] -> [a] -> [a]
++ [Constr]
rest)
          case IntervalZ -> Maybe Integer
pickupZ (Model Integer -> BoundsZ -> IntervalZ
evalBoundsZ Model Integer
model BoundsZ
bnd) of
            Maybe Integer
Nothing  -> forall a. (?callStack::CallStack) => String -> a
error String
"ToySolver.OmegaTest.solve': should not happen"
            Just Integer
val -> forall (m :: * -> *) a. Monad m => a -> m a
return forall a b. (a -> b) -> a -> b
$ forall a. Int -> a -> IntMap a -> IntMap a
IM.insert Int
v Integer
val Model Integer
model

        case2 :: Maybe (Model Integer)
case2 = forall (t :: * -> *) (m :: * -> *) a.
(Foldable t, MonadPlus m) =>
t (m a) -> m a
msum
          [ do Constr
c <- ExprZ -> Maybe Constr
isZero forall a b. (a -> b) -> a -> b
$ Integer
a' forall v. VectorSpace v => Scalar v -> v -> v
*^ forall r. Num r => Int -> Expr r
LA.var Int
v forall v. AdditiveGroup v => v -> v -> v
^-^ (ExprZ
c' forall v. AdditiveGroup v => v -> v -> v
^+^ forall r. (Num r, Eq r) => r -> Expr r
LA.constant Integer
k)
               Model Integer
model <- VarSet -> [Constr] -> Maybe (Model Integer)
go VarSet
vs (Constr
c forall a. a -> [a] -> [a]
: [Constr]
cs)
               forall (m :: * -> *) a. Monad m => a -> m a
return Model Integer
model
          | let m :: Integer
m = forall (t :: * -> *) a. (Foldable t, Ord a) => t a -> a
maximum [Integer
b | (ExprZ
_,Integer
b)<-[Rat]
us]
          , (ExprZ
c',Integer
a') <- [Rat]
ls
          , Integer
k <- [Integer
0 .. (Integer
mforall a. Num a => a -> a -> a
*Integer
a'forall a. Num a => a -> a -> a
-Integer
a'forall a. Num a => a -> a -> a
-Integer
m) forall a. Integral a => a -> a -> a
`div` Integer
m]
          ]

isExact :: BoundsZ -> Bool
isExact :: BoundsZ -> Bool
isExact ([Rat]
ls,[Rat]
us) = forall (t :: * -> *). Foldable t => t Bool -> Bool
and [Integer
aforall a. Eq a => a -> a -> Bool
==Integer
1 Bool -> Bool -> Bool
|| Integer
bforall a. Eq a => a -> a -> Bool
==Integer
1 | (ExprZ
_,Integer
a)<-[Rat]
ls , (ExprZ
_,Integer
b)<-[Rat]
us]

toLAAtom :: Constr -> LA.Atom Rational
toLAAtom :: Constr -> Atom (Ratio Integer)
toLAAtom (IsNonneg ExprZ
e) = forall b a. (Num b, Eq b) => (a -> b) -> Expr a -> Expr b
LA.mapCoeff forall a. Num a => Integer -> a
fromInteger ExprZ
e forall e r. IsOrdRel e r => e -> e -> r
.>=. forall r. (Num r, Eq r) => r -> Expr r
LA.constant Ratio Integer
0
toLAAtom (IsZero ExprZ
e)   = forall b a. (Num b, Eq b) => (a -> b) -> Expr a -> Expr b
LA.mapCoeff forall a. Num a => Integer -> a
fromInteger ExprZ
e forall e r. IsEqRel e r => e -> e -> r
.==. forall r. (Num r, Eq r) => r -> Expr r
LA.constant Ratio Integer
0

chooseVariable :: VarSet -> [Constr] -> (Var, VarSet, BoundsZ, [Constr])
chooseVariable :: VarSet -> [Constr] -> (Int, VarSet, BoundsZ, [Constr])
chooseVariable VarSet
vs [Constr]
xs = forall a. [a] -> a
head forall a b. (a -> b) -> a -> b
$ [(Int, VarSet, BoundsZ, [Constr])
e | e :: (Int, VarSet, BoundsZ, [Constr])
e@(Int
_,VarSet
_,BoundsZ
bnd,[Constr]
_) <- [(Int, VarSet, BoundsZ, [Constr])]
table, BoundsZ -> Bool
isExact BoundsZ
bnd] forall a. [a] -> [a] -> [a]
++ [(Int, VarSet, BoundsZ, [Constr])]
table
  where
    table :: [(Int, VarSet, BoundsZ, [Constr])]
table = [ (Int
v, [Int] -> VarSet
IS.fromAscList [Int]
vs', BoundsZ
bnd, [Constr]
rest)
            | (Int
v,[Int]
vs') <- forall a. [a] -> [(a, [a])]
pickup (VarSet -> [Int]
IS.toAscList VarSet
vs), let (BoundsZ
bnd, [Constr]
rest) = Int -> [Constr] -> (BoundsZ, [Constr])
collectBoundsZ Int
v [Constr]
xs
            ]

-- Find an equality constraint (e==0) and extract e with rest of constraints.
extractEq :: [Constr] -> Maybe (ExprZ, [Constr])
extractEq :: [Constr] -> Maybe (ExprZ, [Constr])
extractEq = forall (t :: * -> *) (m :: * -> *) a.
(Foldable t, MonadPlus m) =>
t (m a) -> m a
msum forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a b. (a -> b) -> [a] -> [b]
map (Constr, [Constr]) -> Maybe (ExprZ, [Constr])
f forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a. [a] -> [(a, [a])]
pickup
  where
    f :: (Constr, [Constr]) -> Maybe (ExprZ, [Constr])
    f :: (Constr, [Constr]) -> Maybe (ExprZ, [Constr])
f (IsZero ExprZ
e, [Constr]
ls) = forall a. a -> Maybe a
Just (ExprZ
e, [Constr]
ls)
    f (Constr, [Constr])
_ = forall a. Maybe a
Nothing

-- Eliminate an equality equality constraint (e==0).
eliminateEq :: ExprZ -> VarSet -> [Constr] -> Maybe (VarSet, [Constr], Model Integer -> Model Integer)
eliminateEq :: ExprZ
-> VarSet
-> [Constr]
-> Maybe (VarSet, [Constr], Model Integer -> Model Integer)
eliminateEq ExprZ
e VarSet
vs [Constr]
cs | Just Integer
c <- forall r. Num r => Expr r -> Maybe r
LA.asConst ExprZ
e = forall (f :: * -> *). Alternative f => Bool -> f ()
guard (Integer
cforall a. Eq a => a -> a -> Bool
==Integer
0) forall (m :: * -> *) a b. Monad m => m a -> m b -> m b
>> forall (m :: * -> *) a. Monad m => a -> m a
return (VarSet
vs, [Constr]
cs, forall a. a -> a
id)
eliminateEq ExprZ
e VarSet
vs [Constr]
cs = do
  let (Integer
ak,Int
xk) = forall (t :: * -> *) a.
Foldable t =>
(a -> a -> Ordering) -> t a -> a
minimumBy (forall a b. Ord a => (b -> a) -> b -> b -> Ordering
comparing (forall a. Num a => a -> a
abs forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a b. (a, b) -> a
fst)) [(Integer
a,Int
x) | (Integer
a,Int
x) <- forall r. Expr r -> [(r, Int)]
LA.terms ExprZ
e, Int
x forall a. Eq a => a -> a -> Bool
/= Int
LA.unitVar]
  if forall a. Num a => a -> a
abs Integer
ak forall a. Eq a => a -> a -> Bool
== Integer
1 then
    case forall r. Num r => Int -> Expr r -> (r, Expr r)
LA.extract Int
xk ExprZ
e of
      (Integer
_, ExprZ
e') -> do
        let xk_def :: ExprZ
xk_def = forall a. Num a => a -> a
signum Integer
ak forall v. VectorSpace v => Scalar v -> v -> v
*^ forall v. AdditiveGroup v => v -> v
negateV ExprZ
e'
        forall (m :: * -> *) a. Monad m => a -> m a
return forall a b. (a -> b) -> a -> b
$
           ( Int -> VarSet -> VarSet
IS.delete Int
xk VarSet
vs
           , [Int -> ExprZ -> Constr -> Constr
applySubst1Constr Int
xk ExprZ
xk_def Constr
c | Constr
c <- [Constr]
cs]
           , \Model Integer
model -> forall a. Int -> a -> IntMap a -> IntMap a
IM.insert Int
xk (forall m e v. Eval m e v => m -> e -> v
LA.eval Model Integer
model ExprZ
xk_def) Model Integer
model
           )
  else do
    let m :: Integer
m = forall a. Num a => a -> a
abs Integer
ak forall a. Num a => a -> a -> a
+ Integer
1
    forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Integer
ak Integer -> Integer -> Integer
`zmod` Integer
m forall a. Eq a => a -> a -> Bool
== - forall a. Num a => a -> a
signum Integer
ak) forall a b. (a -> b) -> a -> b
$ forall (m :: * -> *) a. Monad m => a -> m a
return ()
    let -- sigma is a fresh variable
        sigma :: Int
sigma = if VarSet -> Bool
IS.null VarSet
vs then Int
0 else VarSet -> Int
IS.findMax VarSet
vs forall a. Num a => a -> a -> a
+ Int
1
        -- m *^ LA.var sigma == LA.fromTerms [(a `zmod` m, x) | (a,x) <- LA.terms e]
        -- m *^ LA.var sigma == LA.fromTerms [(a `zmod` m, x) | (a,x) <- LA.terms e, x /= xk] ^+^ (ak `zmod` m) *^ LA.var xk
        -- m *^ LA.var sigma == LA.fromTerms [(a `zmod` m, x) | (a,x) <- LA.terms e, x /= xk] ^+^ (- signum ak) *^ LA.var xk
        -- LA.var xk == (- signum ak * m) *^ LA.var sigma ^+^ LA.fromTerms [(signum ak * (a `zmod` m), x) | (a,x) <- LA.terms e, x /= xk]
        xk_def :: ExprZ
xk_def = (- forall a. Num a => a -> a
signum Integer
ak forall a. Num a => a -> a -> a
* Integer
m) forall v. VectorSpace v => Scalar v -> v -> v
*^ forall r. Num r => Int -> Expr r
LA.var Int
sigma forall v. AdditiveGroup v => v -> v -> v
^+^
                   forall r. (Num r, Eq r) => [(r, Int)] -> Expr r
LA.fromTerms [(forall a. Num a => a -> a
signum Integer
ak forall a. Num a => a -> a -> a
* (Integer
a Integer -> Integer -> Integer
`zmod` Integer
m), Int
x) | (Integer
a,Int
x) <- forall r. Expr r -> [(r, Int)]
LA.terms ExprZ
e, Int
x forall a. Eq a => a -> a -> Bool
/= Int
xk]
        -- e2 is normalized version of (LA.applySubst1 xk xk_def e).
        e2 :: ExprZ
e2 = (- forall a. Num a => a -> a
abs Integer
ak) forall v. VectorSpace v => Scalar v -> v -> v
*^ forall r. Num r => Int -> Expr r
LA.var Int
sigma forall v. AdditiveGroup v => v -> v -> v
^+^
                forall r. (Num r, Eq r) => [(r, Int)] -> Expr r
LA.fromTerms [((forall a b. (RealFrac a, Integral b) => a -> b
floor (Integer
aforall a. Integral a => a -> a -> Ratio a
%Integer
m forall a. Num a => a -> a -> a
+ Ratio Integer
1forall a. Fractional a => a -> a -> a
/Ratio Integer
2) forall a. Num a => a -> a -> a
+ (Integer
a Integer -> Integer -> Integer
`zmod` Integer
m)), Int
x) | (Integer
a,Int
x) <- forall r. Expr r -> [(r, Int)]
LA.terms ExprZ
e, Int
x forall a. Eq a => a -> a -> Bool
/= Int
xk]
    forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Integer
m forall v. VectorSpace v => Scalar v -> v -> v
*^ ExprZ
e2 forall a. Eq a => a -> a -> Bool
== forall r. (Num r, Eq r) => Int -> Expr r -> Expr r -> Expr r
LA.applySubst1 Int
xk ExprZ
xk_def ExprZ
e) forall a b. (a -> b) -> a -> b
$ forall (m :: * -> *) a. Monad m => a -> m a
return ()
    let mt :: Model Integer -> Model Integer
        mt :: Model Integer -> Model Integer
mt Model Integer
model = forall a. Int -> IntMap a -> IntMap a
IM.delete Int
sigma forall a b. (a -> b) -> a -> b
$ forall a. Int -> a -> IntMap a -> IntMap a
IM.insert Int
xk (forall m e v. Eval m e v => m -> e -> v
LA.eval Model Integer
model ExprZ
xk_def) Model Integer
model
    Constr
c <- ExprZ -> Maybe Constr
isZero ExprZ
e2
    forall (m :: * -> *) a. Monad m => a -> m a
return (Int -> VarSet -> VarSet
IS.insert Int
sigma (Int -> VarSet -> VarSet
IS.delete Int
xk VarSet
vs), Constr
c forall a. a -> [a] -> [a]
: [Int -> ExprZ -> Constr -> Constr
applySubst1Constr Int
xk ExprZ
xk_def Constr
c | Constr
c <- [Constr]
cs], Model Integer -> Model Integer
mt)

-- ---------------------------------------------------------------------------

type IntervalZ = (Maybe Integer, Maybe Integer)

univZ :: IntervalZ
univZ :: IntervalZ
univZ = (forall a. Maybe a
Nothing, forall a. Maybe a
Nothing)

intersectZ :: IntervalZ -> IntervalZ -> IntervalZ
intersectZ :: IntervalZ -> IntervalZ -> IntervalZ
intersectZ (Maybe Integer
l1,Maybe Integer
u1) (Maybe Integer
l2,Maybe Integer
u2) = (forall a. (a -> a -> a) -> Maybe a -> Maybe a -> Maybe a
combineMaybe forall a. Ord a => a -> a -> a
max Maybe Integer
l1 Maybe Integer
l2, forall a. (a -> a -> a) -> Maybe a -> Maybe a -> Maybe a
combineMaybe forall a. Ord a => a -> a -> a
min Maybe Integer
u1 Maybe Integer
u2)

pickupZ :: IntervalZ -> Maybe Integer
pickupZ :: IntervalZ -> Maybe Integer
pickupZ (Maybe Integer
Nothing,Maybe Integer
Nothing) = forall (m :: * -> *) a. Monad m => a -> m a
return Integer
0
pickupZ (Just Integer
x, Maybe Integer
Nothing) = forall (m :: * -> *) a. Monad m => a -> m a
return Integer
x
pickupZ (Maybe Integer
Nothing, Just Integer
x) = forall (m :: * -> *) a. Monad m => a -> m a
return Integer
x
pickupZ (Just Integer
x, Just Integer
y) = if Integer
x forall a. Ord a => a -> a -> Bool
<= Integer
y then forall (m :: * -> *) a. Monad m => a -> m a
return Integer
x else forall (m :: * -> *) a. MonadPlus m => m a
mzero

-- ---------------------------------------------------------------------------

-- | @'solve' opt {x1,…,xn} φ@ returns @Just M@ that @M ⊧_LIA φ@ when
-- such @M@ exists, returns @Nothing@ otherwise.
--
-- @FV(φ)@ must be a subset of @{x1,…,xn}@.
--
solve :: Options -> VarSet -> [LA.Atom Rational] -> Maybe (Model Integer)
solve :: Options
-> VarSet -> [Atom (Ratio Integer)] -> Maybe (Model Integer)
solve Options
opt VarSet
vs [Atom (Ratio Integer)]
cs = forall (t :: * -> *) (m :: * -> *) a.
(Foldable t, MonadPlus m) =>
t (m a) -> m a
msum [Options -> VarSet -> [Constr] -> Maybe (Model Integer)
solve' Options
opt VarSet
vs [Constr]
cs | [Constr]
cs <- forall lit. DNF lit -> [[lit]]
unDNF DNF Constr
dnf]
  where
    dnf :: DNF Constr
dnf = forall a. MonotoneBoolean a => [a] -> a
andB (forall a b. (a -> b) -> [a] -> [b]
map Atom (Ratio Integer) -> DNF Constr
f [Atom (Ratio Integer)]
cs)
    f :: Atom (Ratio Integer) -> DNF Constr
f (OrdRel Expr (Ratio Integer)
lhs RelOp
op Expr (Ratio Integer)
rhs) =
      case RelOp
op of
        RelOp
Lt  -> forall lit. [[lit]] -> DNF lit
DNF [[ExprZ
a ExprZ -> ExprZ -> Constr
`ltZ` ExprZ
b]]
        RelOp
Le  -> forall lit. [[lit]] -> DNF lit
DNF [[ExprZ
a ExprZ -> ExprZ -> Constr
`leZ` ExprZ
b]]
        RelOp
Gt  -> forall lit. [[lit]] -> DNF lit
DNF [[ExprZ
a ExprZ -> ExprZ -> Constr
`gtZ` ExprZ
b]]
        RelOp
Ge  -> forall lit. [[lit]] -> DNF lit
DNF [[ExprZ
a ExprZ -> ExprZ -> Constr
`geZ` ExprZ
b]]
        RelOp
Eql -> forall lit. [[lit]] -> DNF lit
DNF [[ExprZ
a ExprZ -> ExprZ -> Constr
`eqZ` ExprZ
b]]
        RelOp
NEq -> forall lit. [[lit]] -> DNF lit
DNF [[ExprZ
a ExprZ -> ExprZ -> Constr
`ltZ` ExprZ
b], [ExprZ
a ExprZ -> ExprZ -> Constr
`gtZ` ExprZ
b]]
      where
        (ExprZ
e1,Integer
c1) = Expr (Ratio Integer) -> Rat
g Expr (Ratio Integer)
lhs
        (ExprZ
e2,Integer
c2) = Expr (Ratio Integer) -> Rat
g Expr (Ratio Integer)
rhs
        a :: ExprZ
a = Integer
c2 forall v. VectorSpace v => Scalar v -> v -> v
*^ ExprZ
e1
        b :: ExprZ
b = Integer
c1 forall v. VectorSpace v => Scalar v -> v -> v
*^ ExprZ
e2
    g :: LA.Expr Rational -> (ExprZ, Integer)
    g :: Expr (Ratio Integer) -> Rat
g Expr (Ratio Integer)
a = (forall b a. (Num b, Eq b) => (a -> b) -> Expr a -> Expr b
LA.mapCoeff (\Ratio Integer
c -> forall a b. (RealFrac a, Integral b) => a -> b
floor (Ratio Integer
c forall a. Num a => a -> a -> a
* forall a. Num a => Integer -> a
fromInteger Integer
d)) Expr (Ratio Integer)
a, Integer
d)
      where
        d :: Integer
d = forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl' forall a. Integral a => a -> a -> a
lcm Integer
1 [forall a. Ratio a -> a
denominator Ratio Integer
c | (Ratio Integer
c,Int
_) <- forall r. Expr r -> [(r, Int)]
LA.terms Expr (Ratio Integer)
a]

-- | @'solveQFLIRAConj' {x1,…,xn} φ I@ returns @Just M@ that @M ⊧_LIRA φ@ when
-- such @M@ exists, returns @Nothing@ otherwise.
--
-- * @FV(φ)@ must be a subset of @{x1,…,xn}@.
--
-- * @I@ is a set of integer variables and must be a subset of @{x1,…,xn}@.
--
solveQFLIRAConj :: Options -> VarSet -> [LA.Atom Rational] -> VarSet -> Maybe (Model Rational)
solveQFLIRAConj :: Options
-> VarSet
-> [Atom (Ratio Integer)]
-> VarSet
-> Maybe (Model (Ratio Integer))
solveQFLIRAConj Options
opt VarSet
vs [Atom (Ratio Integer)]
cs VarSet
ivs = forall a. [a] -> Maybe a
listToMaybe forall a b. (a -> b) -> a -> b
$ do
  ([Atom (Ratio Integer)]
cs2, Model (Ratio Integer) -> Model (Ratio Integer)
mt) <- VarSet
-> [Atom (Ratio Integer)]
-> [([Atom (Ratio Integer)],
     Model (Ratio Integer) -> Model (Ratio Integer))]
FM.projectN VarSet
rvs [Atom (Ratio Integer)]
cs
  Model Integer
m <- forall a. Maybe a -> [a]
maybeToList forall a b. (a -> b) -> a -> b
$ Options
-> VarSet -> [Atom (Ratio Integer)] -> Maybe (Model Integer)
solve Options
opt VarSet
ivs [Atom (Ratio Integer)]
cs2
  forall (m :: * -> *) a. Monad m => a -> m a
return forall a b. (a -> b) -> a -> b
$ Model (Ratio Integer) -> Model (Ratio Integer)
mt forall a b. (a -> b) -> a -> b
$ forall a b. (a -> b) -> IntMap a -> IntMap b
IM.map forall a. Num a => Integer -> a
fromInteger Model Integer
m
  where
    rvs :: VarSet
rvs = VarSet
vs VarSet -> VarSet -> VarSet
`IS.difference` VarSet
ivs

-- ---------------------------------------------------------------------------

zmod :: Integer -> Integer -> Integer
Integer
a zmod :: Integer -> Integer -> Integer
`zmod` Integer
b = Integer
a forall a. Num a => a -> a -> a
- Integer
b forall a. Num a => a -> a -> a
* forall a b. (RealFrac a, Integral b) => a -> b
floor (Integer
a forall a. Integral a => a -> a -> Ratio a
% Integer
b forall a. Num a => a -> a -> a
+ Ratio Integer
1 forall a. Fractional a => a -> a -> a
/ Ratio Integer
2)

gcd' :: [Integer] -> Integer
gcd' :: [Integer] -> Integer
gcd' [] = Integer
1
gcd' [Integer]
xs = forall a. (a -> a -> a) -> [a] -> a
foldl1' forall a. Integral a => a -> a -> a
gcd [Integer]
xs

pickup :: [a] -> [(a,[a])]
pickup :: forall a. [a] -> [(a, [a])]
pickup [] = []
pickup (a
x:[a]
xs) = (a
x,[a]
xs) forall a. a -> [a] -> [a]
: [(a
y,a
xforall a. a -> [a] -> [a]
:[a]
ys) | (a
y,[a]
ys) <- forall a. [a] -> [(a, [a])]
pickup [a]
xs]

-- ---------------------------------------------------------------------------