{-# LANGUAGE ScopedTypeVariables #-}
{-# OPTIONS_GHC -Wall #-}
{-# OPTIONS_HADDOCK show-extensions #-}
-----------------------------------------------------------------------------
-- |
-- Module      :  ToySolver.Arith.Simplex.Textbook
-- Copyright   :  (c) Masahiro Sakai 2011
-- License     :  BSD-style
--
-- Maintainer  :  masahiro.sakai@gmail.com
-- Stability   :  provisional
-- Portability :  non-portable
--
-- Naïve implementation of Simplex method
--
-- Reference:
--
-- * <http://www.math.cuhk.edu.hk/~wei/lpch3.pdf>
--
-----------------------------------------------------------------------------
module ToySolver.Arith.Simplex.Textbook
  (
  -- * The @Tableau@ type
    Tableau
  , RowIndex
  , ColIndex
  , Row
  , emptyTableau
  , objRowIndex
  , pivot
  , lookupRow
  , addRow
  , setObjFun

  -- * Optimization direction
  , module Data.OptDir

  -- * Reading status
  , currentValue
  , currentObjValue
  , isFeasible
  , isOptimal

  -- * High-level solving functions
  , simplex
  , dualSimplex
  , phaseI
  , primalDualSimplex

  -- * For debugging
  , isValidTableau
  , toCSV
  ) where

import Data.Ord
import Data.List
import qualified Data.IntMap as IM
import qualified Data.IntSet as IS
import Data.OptDir
import Data.VectorSpace
import Control.Exception

import qualified ToySolver.Data.LA as LA
import ToySolver.Data.IntVar

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

type Tableau r = VarMap (Row r)
{-
tbl ! v == (m, val)
==>
var v .+. m .==. constant val
-}

-- | Basic variables
type RowIndex = Int

-- | Non-basic variables
type ColIndex = Int

type Row r = (VarMap r, r)

data PivotResult r = PivotUnbounded | PivotFinished | PivotSuccess (Tableau r)
  deriving (Int -> PivotResult r -> ShowS
[PivotResult r] -> ShowS
PivotResult r -> String
(Int -> PivotResult r -> ShowS)
-> (PivotResult r -> String)
-> ([PivotResult r] -> ShowS)
-> Show (PivotResult r)
forall r. Show r => Int -> PivotResult r -> ShowS
forall r. Show r => [PivotResult r] -> ShowS
forall r. Show r => PivotResult r -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [PivotResult r] -> ShowS
$cshowList :: forall r. Show r => [PivotResult r] -> ShowS
show :: PivotResult r -> String
$cshow :: forall r. Show r => PivotResult r -> String
showsPrec :: Int -> PivotResult r -> ShowS
$cshowsPrec :: forall r. Show r => Int -> PivotResult r -> ShowS
Show, PivotResult r -> PivotResult r -> Bool
(PivotResult r -> PivotResult r -> Bool)
-> (PivotResult r -> PivotResult r -> Bool) -> Eq (PivotResult r)
forall r. Eq r => PivotResult r -> PivotResult r -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: PivotResult r -> PivotResult r -> Bool
$c/= :: forall r. Eq r => PivotResult r -> PivotResult r -> Bool
== :: PivotResult r -> PivotResult r -> Bool
$c== :: forall r. Eq r => PivotResult r -> PivotResult r -> Bool
Eq, Eq (PivotResult r)
Eq (PivotResult r)
-> (PivotResult r -> PivotResult r -> Ordering)
-> (PivotResult r -> PivotResult r -> Bool)
-> (PivotResult r -> PivotResult r -> Bool)
-> (PivotResult r -> PivotResult r -> Bool)
-> (PivotResult r -> PivotResult r -> Bool)
-> (PivotResult r -> PivotResult r -> PivotResult r)
-> (PivotResult r -> PivotResult r -> PivotResult r)
-> Ord (PivotResult r)
PivotResult r -> PivotResult r -> Bool
PivotResult r -> PivotResult r -> Ordering
PivotResult r -> PivotResult r -> PivotResult r
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
forall r. Ord r => Eq (PivotResult r)
forall r. Ord r => PivotResult r -> PivotResult r -> Bool
forall r. Ord r => PivotResult r -> PivotResult r -> Ordering
forall r. Ord r => PivotResult r -> PivotResult r -> PivotResult r
min :: PivotResult r -> PivotResult r -> PivotResult r
$cmin :: forall r. Ord r => PivotResult r -> PivotResult r -> PivotResult r
max :: PivotResult r -> PivotResult r -> PivotResult r
$cmax :: forall r. Ord r => PivotResult r -> PivotResult r -> PivotResult r
>= :: PivotResult r -> PivotResult r -> Bool
$c>= :: forall r. Ord r => PivotResult r -> PivotResult r -> Bool
> :: PivotResult r -> PivotResult r -> Bool
$c> :: forall r. Ord r => PivotResult r -> PivotResult r -> Bool
<= :: PivotResult r -> PivotResult r -> Bool
$c<= :: forall r. Ord r => PivotResult r -> PivotResult r -> Bool
< :: PivotResult r -> PivotResult r -> Bool
$c< :: forall r. Ord r => PivotResult r -> PivotResult r -> Bool
compare :: PivotResult r -> PivotResult r -> Ordering
$ccompare :: forall r. Ord r => PivotResult r -> PivotResult r -> Ordering
$cp1Ord :: forall r. Ord r => Eq (PivotResult r)
Ord)

emptyTableau :: Tableau r
emptyTableau :: Tableau r
emptyTableau = Tableau r
forall a. IntMap a
IM.empty

objRowIndex :: RowIndex
objRowIndex :: Int
objRowIndex = -Int
1

pivot :: (Fractional r, Eq r) => RowIndex -> ColIndex -> Tableau r -> Tableau r
{-# INLINE pivot #-}
{-# SPECIALIZE pivot :: RowIndex -> ColIndex -> Tableau Rational -> Tableau Rational #-}
{-# SPECIALIZE pivot :: RowIndex -> ColIndex -> Tableau Double -> Tableau Double #-}
pivot :: Int -> Int -> Tableau r -> Tableau r
pivot Int
r Int
s Tableau r
tbl =
    Bool -> Tableau r -> Tableau r
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Tableau r -> Bool
forall r. Tableau r -> Bool
isValidTableau Tableau r
tbl) (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$  -- precondition
    Bool -> Tableau r -> Tableau r
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Tableau r -> Bool
forall r. Tableau r -> Bool
isValidTableau Tableau r
tbl') (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ -- postcondition
      Tableau r
tbl'
  where
    tbl' :: Tableau r
tbl' = Int -> (IntMap r, r) -> Tableau r -> Tableau r
forall a. Int -> a -> IntMap a -> IntMap a
IM.insert Int
s (IntMap r, r)
row_s (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ ((IntMap r, r) -> (IntMap r, r)) -> Tableau r -> Tableau r
forall a b. (a -> b) -> IntMap a -> IntMap b
IM.map (IntMap r, r) -> (IntMap r, r)
f (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ Int -> Tableau r -> Tableau r
forall a. Int -> IntMap a -> IntMap a
IM.delete Int
r (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ Tableau r
tbl
    f :: (IntMap r, r) -> (IntMap r, r)
f orig :: (IntMap r, r)
orig@(IntMap r
row_i, r
row_i_val) =
      case Int -> IntMap r -> Maybe r
forall a. Int -> IntMap a -> Maybe a
IM.lookup Int
s IntMap r
row_i of
        Maybe r
Nothing -> (IntMap r, r)
orig
        Just r
c ->
          ( (r -> Bool) -> IntMap r -> IntMap r
forall a. (a -> Bool) -> IntMap a -> IntMap a
IM.filter (r
0r -> r -> Bool
forall a. Eq a => a -> a -> Bool
/=) (IntMap r -> IntMap r) -> IntMap r -> IntMap r
forall a b. (a -> b) -> a -> b
$ (r -> r -> r) -> IntMap r -> IntMap r -> IntMap r
forall a. (a -> a -> a) -> IntMap a -> IntMap a -> IntMap a
IM.unionWith r -> r -> r
forall a. Num a => a -> a -> a
(+) (Int -> IntMap r -> IntMap r
forall a. Int -> IntMap a -> IntMap a
IM.delete Int
s IntMap r
row_i) ((r -> r) -> IntMap r -> IntMap r
forall a b. (a -> b) -> IntMap a -> IntMap b
IM.map (r -> r
forall a. Num a => a -> a
negate r
c r -> r -> r
forall a. Num a => a -> a -> a
*) IntMap r
row_r')
          , r
row_i_val r -> r -> r
forall a. Num a => a -> a -> a
- r
cr -> r -> r
forall a. Num a => a -> a -> a
*r
row_r_val'
          )
    (IntMap r
row_r, r
row_r_val) = Int -> Tableau r -> (IntMap r, r)
forall r. Int -> Tableau r -> Row r
lookupRow Int
r Tableau r
tbl
    a_rs :: r
a_rs = IntMap r
row_r IntMap r -> Int -> r
forall a. IntMap a -> Int -> a
IM.! Int
s
    row_r' :: IntMap r
row_r' = (r -> r) -> IntMap r -> IntMap r
forall a b. (a -> b) -> IntMap a -> IntMap b
IM.map (r -> r -> r
forall a. Fractional a => a -> a -> a
/ r
a_rs) (IntMap r -> IntMap r) -> IntMap r -> IntMap r
forall a b. (a -> b) -> a -> b
$ Int -> r -> IntMap r -> IntMap r
forall a. Int -> a -> IntMap a -> IntMap a
IM.insert Int
r r
1 (IntMap r -> IntMap r) -> IntMap r -> IntMap r
forall a b. (a -> b) -> a -> b
$ Int -> IntMap r -> IntMap r
forall a. Int -> IntMap a -> IntMap a
IM.delete Int
s IntMap r
row_r
    row_r_val' :: r
row_r_val' = r
row_r_val r -> r -> r
forall a. Fractional a => a -> a -> a
/ r
a_rs
    row_s :: (IntMap r, r)
row_s = (IntMap r
row_r', r
row_r_val')

-- | Lookup a row by basic variable
lookupRow :: RowIndex -> Tableau r -> Row r
lookupRow :: Int -> Tableau r -> Row r
lookupRow Int
r Tableau r
m = Tableau r
m Tableau r -> Int -> Row r
forall a. IntMap a -> Int -> a
IM.! Int
r

-- 行の基底変数の列が0になるように変形
normalizeRow :: (Num r, Eq r) => Tableau r -> Row r -> Row r
normalizeRow :: Tableau r -> Row r -> Row r
normalizeRow Tableau r
a (VarMap r
row0,r
val0) = Row r
obj'
  where
    obj' :: Row r
obj' = Row r -> Row r
forall a b. (Eq a, Num a) => (IntMap a, b) -> (IntMap a, b)
g (Row r -> Row r) -> Row r -> Row r
forall a b. (a -> b) -> a -> b
$ (Row r -> Row r -> Row r) -> Row r -> [Row r] -> Row r
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl' Row r -> Row r -> Row r
forall a b.
(Num a, Num b) =>
(IntMap a, b) -> (IntMap a, b) -> (IntMap a, b)
f (VarMap r
forall a. IntMap a
IM.empty, r
val0) ([Row r] -> Row r) -> [Row r] -> Row r
forall a b. (a -> b) -> a -> b
$
           [ case Int -> Tableau r -> Maybe (Row r)
forall a. Int -> IntMap a -> Maybe a
IM.lookup Int
j Tableau r
a of
               Maybe (Row r)
Nothing -> (Int -> r -> VarMap r
forall a. Int -> a -> IntMap a
IM.singleton Int
j r
x, r
0)
               Just (VarMap r
row,r
val) -> ((r -> r) -> VarMap r -> VarMap r
forall a b. (a -> b) -> IntMap a -> IntMap b
IM.map ((-r
x)r -> r -> r
forall a. Num a => a -> a -> a
*) (Int -> VarMap r -> VarMap r
forall a. Int -> IntMap a -> IntMap a
IM.delete Int
j VarMap r
row), -r
xr -> r -> r
forall a. Num a => a -> a -> a
*r
val)
           | (Int
j,r
x) <- VarMap r -> [(Int, r)]
forall a. IntMap a -> [(Int, a)]
IM.toList VarMap r
row0 ]
    f :: (IntMap a, b) -> (IntMap a, b) -> (IntMap a, b)
f (IntMap a
m1,b
v1) (IntMap a
m2,b
v2) = ((a -> a -> a) -> IntMap a -> IntMap a -> IntMap a
forall a. (a -> a -> a) -> IntMap a -> IntMap a -> IntMap a
IM.unionWith a -> a -> a
forall a. Num a => a -> a -> a
(+) IntMap a
m1 IntMap a
m2, b
v1b -> b -> b
forall a. Num a => a -> a -> a
+b
v2)
    g :: (IntMap a, b) -> (IntMap a, b)
g (IntMap a
m,b
v) = ((a -> Bool) -> IntMap a -> IntMap a
forall a. (a -> Bool) -> IntMap a -> IntMap a
IM.filter (a
0a -> a -> Bool
forall a. Eq a => a -> a -> Bool
/=) IntMap a
m, b
v)

setRow :: (Num r, Eq r) => Tableau r -> RowIndex -> Row r -> Tableau r
setRow :: Tableau r -> Int -> Row r -> Tableau r
setRow Tableau r
tbl Int
i Row r
row = Bool -> Tableau r -> Tableau r
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Tableau r -> Bool
forall r. Tableau r -> Bool
isValidTableau Tableau r
tbl) (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ Bool -> Tableau r -> Tableau r
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Tableau r -> Bool
forall r. Tableau r -> Bool
isValidTableau Tableau r
tbl') (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ Tableau r
tbl'
  where
    tbl' :: Tableau r
tbl' = Int -> Row r -> Tableau r -> Tableau r
forall a. Int -> a -> IntMap a -> IntMap a
IM.insert Int
i (Tableau r -> Row r -> Row r
forall r. (Num r, Eq r) => Tableau r -> Row r -> Row r
normalizeRow Tableau r
tbl Row r
row) Tableau r
tbl

addRow :: (Num r, Eq r) => Tableau r -> RowIndex -> Row r -> Tableau r
addRow :: Tableau r -> Int -> Row r -> Tableau r
addRow Tableau r
tbl Int
i Row r
row = Bool -> Tableau r -> Tableau r
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Int
i Int -> Tableau r -> Bool
forall a. Int -> IntMap a -> Bool
`IM.notMember` Tableau r
tbl) (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ Tableau r -> Int -> Row r -> Tableau r
forall r. (Num r, Eq r) => Tableau r -> Int -> Row r -> Tableau r
setRow Tableau r
tbl Int
i Row r
row

setObjFun :: (Num r, Eq r) => Tableau r -> LA.Expr r -> Tableau r
setObjFun :: Tableau r -> Expr r -> Tableau r
setObjFun Tableau r
tbl Expr r
e = Tableau r -> Int -> Row r -> Tableau r
forall r. (Num r, Eq r) => Tableau r -> Int -> Row r -> Tableau r
addRow Tableau r
tbl Int
objRowIndex Row r
row
  where
    row :: Row r
row =
      case Int -> Expr r -> (r, Expr r)
forall r. Num r => Int -> Expr r -> (r, Expr r)
LA.extract Int
LA.unitVar Expr r
e of
        (r
c, Expr r
e') -> (Expr r -> IntMap r
forall r. Expr r -> IntMap r
LA.coeffMap (Expr r -> Expr r
forall v. AdditiveGroup v => v -> v
negateV Expr r
e'), r
c)

copyObjRow :: (Num r, Eq r) => Tableau r -> Tableau r -> Tableau r
copyObjRow :: Tableau r -> Tableau r -> Tableau r
copyObjRow Tableau r
from Tableau r
to =
  case Int -> Tableau r -> Maybe (Row r)
forall a. Int -> IntMap a -> Maybe a
IM.lookup Int
objRowIndex Tableau r
from of
    Maybe (Row r)
Nothing -> Int -> Tableau r -> Tableau r
forall a. Int -> IntMap a -> IntMap a
IM.delete Int
objRowIndex Tableau r
to
    Just Row r
row -> Tableau r -> Int -> Row r -> Tableau r
forall r. (Num r, Eq r) => Tableau r -> Int -> Row r -> Tableau r
addRow Tableau r
to Int
objRowIndex Row r
row

currentValue :: Num r => Tableau r -> Var -> r
currentValue :: Tableau r -> Int -> r
currentValue Tableau r
tbl Int
v =
  case Int -> Tableau r -> Maybe (Row r)
forall a. Int -> IntMap a -> Maybe a
IM.lookup Int
v Tableau r
tbl of
    Maybe (Row r)
Nothing -> r
0
    Just (VarMap r
_, r
val) -> r
val

currentObjValue :: Tableau r -> r
currentObjValue :: Tableau r -> r
currentObjValue = (VarMap r, r) -> r
forall a b. (a, b) -> b
snd ((VarMap r, r) -> r)
-> (Tableau r -> (VarMap r, r)) -> Tableau r -> r
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Int -> Tableau r -> (VarMap r, r)
forall r. Int -> Tableau r -> Row r
lookupRow Int
objRowIndex

isValidTableau :: Tableau r -> Bool
isValidTableau :: Tableau r -> Bool
isValidTableau Tableau r
tbl =
  [Bool] -> Bool
forall (t :: * -> *). Foldable t => t Bool -> Bool
and [Int
v Int -> IntMap r -> Bool
forall a. Int -> IntMap a -> Bool
`IM.notMember` IntMap r
m | (Int
v, (IntMap r
m,r
_)) <- Tableau r -> [(Int, Row r)]
forall a. IntMap a -> [(Int, a)]
IM.toList Tableau r
tbl, Int
v Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= Int
objRowIndex] Bool -> Bool -> Bool
&&
  [Bool] -> Bool
forall (t :: * -> *). Foldable t => t Bool -> Bool
and [IntSet -> Bool
IS.null (IntMap r -> IntSet
forall a. IntMap a -> IntSet
IM.keysSet IntMap r
m IntSet -> IntSet -> IntSet
`IS.intersection` IntSet
vs) | (IntMap r
m,r
_) <- Tableau r -> [Row r]
forall a. IntMap a -> [a]
IM.elems Tableau r
tbl']
  where
    tbl' :: Tableau r
tbl' = Int -> Tableau r -> Tableau r
forall a. Int -> IntMap a -> IntMap a
IM.delete Int
objRowIndex Tableau r
tbl
    vs :: IntSet
vs = Tableau r -> IntSet
forall a. IntMap a -> IntSet
IM.keysSet Tableau r
tbl'

isFeasible :: Real r => Tableau r -> Bool
isFeasible :: Tableau r -> Bool
isFeasible Tableau r
tbl =
  [Bool] -> Bool
forall (t :: * -> *). Foldable t => t Bool -> Bool
and [r
val r -> r -> Bool
forall a. Ord a => a -> a -> Bool
>= r
0 | (Int
v, (VarMap r
_,r
val)) <- Tableau r -> [(Int, Row r)]
forall a. IntMap a -> [(Int, a)]
IM.toList Tableau r
tbl, Int
v Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= Int
objRowIndex]

isOptimal :: Real r => OptDir -> Tableau r -> Bool
isOptimal :: OptDir -> Tableau r -> Bool
isOptimal OptDir
optdir Tableau r
tbl =
  [Bool] -> Bool
forall (t :: * -> *). Foldable t => t Bool -> Bool
and [Bool -> Bool
not (r -> Bool
cmp r
cj) | r
cj <- IntMap r -> [r]
forall a. IntMap a -> [a]
IM.elems ((IntMap r, r) -> IntMap r
forall a b. (a, b) -> a
fst (Int -> Tableau r -> (IntMap r, r)
forall r. Int -> Tableau r -> Row r
lookupRow Int
objRowIndex Tableau r
tbl))]
  where
    cmp :: r -> Bool
cmp = case OptDir
optdir of
            OptDir
OptMin -> (r
0r -> r -> Bool
forall a. Ord a => a -> a -> Bool
<)
            OptDir
OptMax -> (r
0r -> r -> Bool
forall a. Ord a => a -> a -> Bool
>)

isImproving :: Real r => OptDir -> Tableau r -> Tableau r -> Bool
isImproving :: OptDir -> Tableau r -> Tableau r -> Bool
isImproving OptDir
OptMin Tableau r
from Tableau r
to = Tableau r -> r
forall r. Tableau r -> r
currentObjValue Tableau r
to r -> r -> Bool
forall a. Ord a => a -> a -> Bool
<= Tableau r -> r
forall r. Tableau r -> r
currentObjValue Tableau r
from
isImproving OptDir
OptMax Tableau r
from Tableau r
to = Tableau r -> r
forall r. Tableau r -> r
currentObjValue Tableau r
to r -> r -> Bool
forall a. Ord a => a -> a -> Bool
>= Tableau r -> r
forall r. Tableau r -> r
currentObjValue Tableau r
from

-- ---------------------------------------------------------------------------
-- primal simplex

simplex :: (Real r, Fractional r) => OptDir -> Tableau r -> (Bool, Tableau r)
{-# SPECIALIZE simplex :: OptDir -> Tableau Rational -> (Bool, Tableau Rational) #-}
{-# SPECIALIZE simplex :: OptDir -> Tableau Double -> (Bool, Tableau Double) #-}
simplex :: OptDir -> Tableau r -> (Bool, Tableau r)
simplex OptDir
optdir = Tableau r -> (Bool, Tableau r)
forall r. (Real r, Fractional r) => Tableau r -> (Bool, Tableau r)
go
  where
    go :: Tableau r -> (Bool, Tableau r)
go Tableau r
tbl = Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Tableau r -> Bool
forall r. Real r => Tableau r -> Bool
isFeasible Tableau r
tbl) ((Bool, Tableau r) -> (Bool, Tableau r))
-> (Bool, Tableau r) -> (Bool, Tableau r)
forall a b. (a -> b) -> a -> b
$
      case OptDir -> Tableau r -> PivotResult r
forall r.
(Real r, Fractional r) =>
OptDir -> Tableau r -> PivotResult r
primalPivot OptDir
optdir Tableau r
tbl of
        PivotResult r
PivotFinished  -> Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (OptDir -> Tableau r -> Bool
forall r. Real r => OptDir -> Tableau r -> Bool
isOptimal OptDir
optdir Tableau r
tbl) (Bool
True, Tableau r
tbl)
        PivotResult r
PivotUnbounded -> (Bool
False, Tableau r
tbl)
        PivotSuccess Tableau r
tbl' -> Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (OptDir -> Tableau r -> Tableau r -> Bool
forall r. Real r => OptDir -> Tableau r -> Tableau r -> Bool
isImproving OptDir
optdir Tableau r
tbl Tableau r
tbl') ((Bool, Tableau r) -> (Bool, Tableau r))
-> (Bool, Tableau r) -> (Bool, Tableau r)
forall a b. (a -> b) -> a -> b
$ Tableau r -> (Bool, Tableau r)
go Tableau r
tbl'

primalPivot :: (Real r, Fractional r) => OptDir -> Tableau r -> PivotResult r
{-# INLINE primalPivot #-}
primalPivot :: OptDir -> Tableau r -> PivotResult r
primalPivot OptDir
optdir Tableau r
tbl
  | [(Int, r)] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [(Int, r)]
cs   = PivotResult r
forall r. PivotResult r
PivotFinished
  | [(Int, r)] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [(Int, r)]
rs   = PivotResult r
forall r. PivotResult r
PivotUnbounded
  | Bool
otherwise = Tableau r -> PivotResult r
forall r. Tableau r -> PivotResult r
PivotSuccess (Int -> Int -> Tableau r -> Tableau r
forall r.
(Fractional r, Eq r) =>
Int -> Int -> Tableau r -> Tableau r
pivot Int
r Int
s Tableau r
tbl)
  where
    cmp :: r -> Bool
cmp = case OptDir
optdir of
            OptDir
OptMin -> (r
0r -> r -> Bool
forall a. Ord a => a -> a -> Bool
<)
            OptDir
OptMax -> (r
0r -> r -> Bool
forall a. Ord a => a -> a -> Bool
>)
    cs :: [(Int, r)]
cs = [(Int
j,r
cj) | (Int
j,r
cj) <- IntMap r -> [(Int, r)]
forall a. IntMap a -> [(Int, a)]
IM.toList ((IntMap r, r) -> IntMap r
forall a b. (a, b) -> a
fst (Int -> Tableau r -> (IntMap r, r)
forall r. Int -> Tableau r -> Row r
lookupRow Int
objRowIndex Tableau r
tbl)), r -> Bool
cmp r
cj]
    -- smallest subscript rule
    s :: Int
s = (Int, r) -> Int
forall a b. (a, b) -> a
fst ((Int, r) -> Int) -> (Int, r) -> Int
forall a b. (a -> b) -> a -> b
$ [(Int, r)] -> (Int, r)
forall a. [a] -> a
head [(Int, r)]
cs
    -- classical rule
    --s = fst $ (if optdir==OptMin then maximumBy else minimumBy) (compare `on` snd) cs
    rs :: [(Int, r)]
rs = [ (Int
i, r
y_i0 r -> r -> r
forall a. Fractional a => a -> a -> a
/ r
y_is)
         | (Int
i, (IntMap r
row_i, r
y_i0)) <- Tableau r -> [(Int, (IntMap r, r))]
forall a. IntMap a -> [(Int, a)]
IM.toList Tableau r
tbl, Int
i Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= Int
objRowIndex
         , let y_is :: r
y_is = r -> Int -> IntMap r -> r
forall a. a -> Int -> IntMap a -> a
IM.findWithDefault r
0 Int
s IntMap r
row_i, r
y_is r -> r -> Bool
forall a. Ord a => a -> a -> Bool
> r
0
         ]
    r :: Int
r = (Int, r) -> Int
forall a b. (a, b) -> a
fst ((Int, r) -> Int) -> (Int, r) -> Int
forall a b. (a -> b) -> a -> b
$ ((Int, r) -> (Int, r) -> Ordering) -> [(Int, r)] -> (Int, r)
forall (t :: * -> *) a.
Foldable t =>
(a -> a -> Ordering) -> t a -> a
minimumBy (((Int, r) -> r) -> (Int, r) -> (Int, r) -> Ordering
forall a b. Ord a => (b -> a) -> b -> b -> Ordering
comparing (Int, r) -> r
forall a b. (a, b) -> b
snd) [(Int, r)]
rs

-- ---------------------------------------------------------------------------
-- dual simplex

dualSimplex :: (Real r, Fractional r) => OptDir -> Tableau r -> (Bool, Tableau r)
{-# SPECIALIZE dualSimplex :: OptDir -> Tableau Rational -> (Bool, Tableau Rational) #-}
{-# SPECIALIZE dualSimplex :: OptDir -> Tableau Double -> (Bool, Tableau Double) #-}
dualSimplex :: OptDir -> Tableau r -> (Bool, Tableau r)
dualSimplex OptDir
optdir = Tableau r -> (Bool, Tableau r)
forall r. (Real r, Fractional r) => Tableau r -> (Bool, Tableau r)
go
  where
    go :: Tableau r -> (Bool, Tableau r)
go Tableau r
tbl = Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (OptDir -> Tableau r -> Bool
forall r. Real r => OptDir -> Tableau r -> Bool
isOptimal OptDir
optdir Tableau r
tbl) ((Bool, Tableau r) -> (Bool, Tableau r))
-> (Bool, Tableau r) -> (Bool, Tableau r)
forall a b. (a -> b) -> a -> b
$
      case OptDir -> Tableau r -> PivotResult r
forall r.
(Real r, Fractional r) =>
OptDir -> Tableau r -> PivotResult r
dualPivot OptDir
optdir Tableau r
tbl of
        PivotResult r
PivotFinished  -> Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Tableau r -> Bool
forall r. Real r => Tableau r -> Bool
isFeasible Tableau r
tbl) ((Bool, Tableau r) -> (Bool, Tableau r))
-> (Bool, Tableau r) -> (Bool, Tableau r)
forall a b. (a -> b) -> a -> b
$ (Bool
True, Tableau r
tbl)
        PivotResult r
PivotUnbounded -> (Bool
False, Tableau r
tbl)
        PivotSuccess Tableau r
tbl' -> Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (OptDir -> Tableau r -> Tableau r -> Bool
forall r. Real r => OptDir -> Tableau r -> Tableau r -> Bool
isImproving OptDir
optdir Tableau r
tbl' Tableau r
tbl) ((Bool, Tableau r) -> (Bool, Tableau r))
-> (Bool, Tableau r) -> (Bool, Tableau r)
forall a b. (a -> b) -> a -> b
$ Tableau r -> (Bool, Tableau r)
go Tableau r
tbl'

dualPivot :: (Real r, Fractional r) => OptDir -> Tableau r -> PivotResult r
{-# INLINE dualPivot #-}
dualPivot :: OptDir -> Tableau r -> PivotResult r
dualPivot OptDir
optdir Tableau r
tbl
  | [(Int, VarMap r)] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [(Int, VarMap r)]
rs   = PivotResult r
forall r. PivotResult r
PivotFinished
  | [(Int, r)] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [(Int, r)]
cs   = PivotResult r
forall r. PivotResult r
PivotUnbounded
  | Bool
otherwise = Tableau r -> PivotResult r
forall r. Tableau r -> PivotResult r
PivotSuccess (Int -> Int -> Tableau r -> Tableau r
forall r.
(Fractional r, Eq r) =>
Int -> Int -> Tableau r -> Tableau r
pivot Int
r Int
s Tableau r
tbl)
  where
    rs :: [(Int, VarMap r)]
rs = [(Int
i, VarMap r
row_i) | (Int
i, (VarMap r
row_i, r
y_i0)) <- Tableau r -> [(Int, Row r)]
forall a. IntMap a -> [(Int, a)]
IM.toList Tableau r
tbl, Int
i Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= Int
objRowIndex, r
0 r -> r -> Bool
forall a. Ord a => a -> a -> Bool
> r
y_i0]
    (Int
r, VarMap r
row_r) = [(Int, VarMap r)] -> (Int, VarMap r)
forall a. [a] -> a
head [(Int, VarMap r)]
rs
    cs :: [(Int, r)]
cs = [ (Int
j, if OptDir
optdirOptDir -> OptDir -> Bool
forall a. Eq a => a -> a -> Bool
==OptDir
OptMin then r
y_0j r -> r -> r
forall a. Fractional a => a -> a -> a
/ r
y_rj else - r
y_0j r -> r -> r
forall a. Fractional a => a -> a -> a
/ r
y_rj)
         | (Int
j, r
y_rj) <- VarMap r -> [(Int, r)]
forall a. IntMap a -> [(Int, a)]
IM.toList VarMap r
row_r
         , r
y_rj r -> r -> Bool
forall a. Ord a => a -> a -> Bool
< r
0
         , let y_0j :: r
y_0j = r -> Int -> VarMap r -> r
forall a. a -> Int -> IntMap a -> a
IM.findWithDefault r
0 Int
j VarMap r
obj
         ]
    s :: Int
s = (Int, r) -> Int
forall a b. (a, b) -> a
fst ((Int, r) -> Int) -> (Int, r) -> Int
forall a b. (a -> b) -> a -> b
$ ((Int, r) -> (Int, r) -> Ordering) -> [(Int, r)] -> (Int, r)
forall (t :: * -> *) a.
Foldable t =>
(a -> a -> Ordering) -> t a -> a
minimumBy (((Int, r) -> r) -> (Int, r) -> (Int, r) -> Ordering
forall a b. Ord a => (b -> a) -> b -> b -> Ordering
comparing (Int, r) -> r
forall a b. (a, b) -> b
snd) [(Int, r)]
cs
    (VarMap r
obj,r
_) = Int -> Tableau r -> Row r
forall r. Int -> Tableau r -> Row r
lookupRow Int
objRowIndex Tableau r
tbl

-- ---------------------------------------------------------------------------
-- phase I of the two-phased method

phaseI :: (Real r, Fractional r) => Tableau r -> VarSet -> (Bool, Tableau r)
{-# SPECIALIZE phaseI :: Tableau Rational -> VarSet -> (Bool, Tableau Rational) #-}
{-# SPECIALIZE phaseI :: Tableau Double -> VarSet -> (Bool, Tableau Double) #-}
phaseI :: Tableau r -> IntSet -> (Bool, Tableau r)
phaseI Tableau r
tbl IntSet
avs
  | Tableau r -> r
forall r. Tableau r -> r
currentObjValue Tableau r
tbl1' r -> r -> Bool
forall a. Eq a => a -> a -> Bool
/= r
0 = (Bool
False, Tableau r
tbl1')
  | Bool
otherwise = (Bool
True, Tableau r -> Tableau r -> Tableau r
forall r. (Num r, Eq r) => Tableau r -> Tableau r -> Tableau r
copyObjRow Tableau r
tbl (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ IntSet -> Tableau r -> Tableau r
forall r.
(Real r, Fractional r) =>
IntSet -> Tableau r -> Tableau r
removeArtificialVariables IntSet
avs (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ Tableau r
tbl1')
  where
    optdir :: OptDir
optdir = OptDir
OptMax
    tbl1 :: Tableau r
tbl1 = Tableau r -> Expr r -> Tableau r
forall r. (Num r, Eq r) => Tableau r -> Expr r -> Tableau r
setObjFun Tableau r
tbl (Expr r -> Tableau r) -> Expr r -> Tableau r
forall a b. (a -> b) -> a -> b
$ Expr r -> Expr r
forall v. AdditiveGroup v => v -> v
negateV (Expr r -> Expr r) -> Expr r -> Expr r
forall a b. (a -> b) -> a -> b
$ [Expr r] -> Expr r
forall (f :: * -> *) v. (Foldable f, AdditiveGroup v) => f v -> v
sumV [Int -> Expr r
forall r. Num r => Int -> Expr r
LA.var Int
v | Int
v <- IntSet -> [Int]
IS.toList IntSet
avs]
    tbl1' :: Tableau r
tbl1' = Tableau r -> Tableau r
forall r. (Real r, Fractional r) => Tableau r -> Tableau r
go Tableau r
tbl1
    go :: Tableau r -> Tableau r
go Tableau r
tbl2
      | Tableau r -> r
forall r. Tableau r -> r
currentObjValue Tableau r
tbl2 r -> r -> Bool
forall a. Eq a => a -> a -> Bool
== r
0 = Tableau r
tbl2
      | Bool
otherwise =
        case OptDir -> Tableau r -> PivotResult r
forall r.
(Real r, Fractional r) =>
OptDir -> Tableau r -> PivotResult r
primalPivot OptDir
optdir Tableau r
tbl2 of
          PivotSuccess Tableau r
tbl2' -> Bool -> Tableau r -> Tableau r
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (OptDir -> Tableau r -> Tableau r -> Bool
forall r. Real r => OptDir -> Tableau r -> Tableau r -> Bool
isImproving OptDir
optdir Tableau r
tbl2 Tableau r
tbl2') (Tableau r -> Tableau r) -> Tableau r -> Tableau r
forall a b. (a -> b) -> a -> b
$ Tableau r -> Tableau r
go Tableau r
tbl2'
          PivotResult r
PivotFinished -> Bool -> Tableau r -> Tableau r
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (OptDir -> Tableau r -> Bool
forall r. Real r => OptDir -> Tableau r -> Bool
isOptimal OptDir
optdir Tableau r
tbl2) Tableau r
tbl2
          PivotResult r
PivotUnbounded -> String -> Tableau r
forall a. (?callStack::CallStack) => String -> a
error String
"phaseI: should not happen"

-- post-processing of phaseI
-- pivotを使ってartificial variablesを基底から除いて、削除
removeArtificialVariables :: (Real r, Fractional r) => VarSet -> Tableau r -> Tableau r
removeArtificialVariables :: IntSet -> Tableau r -> Tableau r
removeArtificialVariables IntSet
avs Tableau r
tbl0 = Tableau r
tbl2
  where
    tbl1 :: Tableau r
tbl1 = (Tableau r -> Int -> Tableau r) -> Tableau r -> [Int] -> Tableau r
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl' Tableau r -> Int -> Tableau r
forall r.
(Eq r, Fractional r) =>
IntMap (IntMap r, r) -> Int -> IntMap (IntMap r, r)
f (Int -> Tableau r -> Tableau r
forall a. Int -> IntMap a -> IntMap a
IM.delete Int
objRowIndex Tableau r
tbl0) (IntSet -> [Int]
IS.toList IntSet
avs)
    tbl2 :: Tableau r
tbl2 = ((IntMap r, r) -> (IntMap r, r)) -> Tableau r -> Tableau r
forall a b. (a -> b) -> IntMap a -> IntMap b
IM.map (\(IntMap r
row,r
val) ->  ((Int -> r -> Bool) -> IntMap r -> IntMap r
forall a. (Int -> a -> Bool) -> IntMap a -> IntMap a
IM.filterWithKey (\Int
j r
_ -> Int
j Int -> IntSet -> Bool
`IS.notMember` IntSet
avs) IntMap r
row, r
val)) Tableau r
tbl1
    f :: IntMap (IntMap r, r) -> Int -> IntMap (IntMap r, r)
f IntMap (IntMap r, r)
tbl Int
i =
      case Int -> IntMap (IntMap r, r) -> Maybe (IntMap r, r)
forall a. Int -> IntMap a -> Maybe a
IM.lookup Int
i IntMap (IntMap r, r)
tbl of
        Maybe (IntMap r, r)
Nothing -> IntMap (IntMap r, r)
tbl
        Just (IntMap r, r)
row ->
          case [Int
j | (Int
j,r
c) <- IntMap r -> [(Int, r)]
forall a. IntMap a -> [(Int, a)]
IM.toList ((IntMap r, r) -> IntMap r
forall a b. (a, b) -> a
fst (IntMap r, r)
row), r
c r -> r -> Bool
forall a. Eq a => a -> a -> Bool
/= r
0, Int
j Int -> IntSet -> Bool
`IS.notMember` IntSet
avs] of
            [] -> Int -> IntMap (IntMap r, r) -> IntMap (IntMap r, r)
forall a. Int -> IntMap a -> IntMap a
IM.delete Int
i IntMap (IntMap r, r)
tbl
            Int
j:[Int]
_ -> Int -> Int -> IntMap (IntMap r, r) -> IntMap (IntMap r, r)
forall r.
(Fractional r, Eq r) =>
Int -> Int -> Tableau r -> Tableau r
pivot Int
i Int
j IntMap (IntMap r, r)
tbl

-- ---------------------------------------------------------------------------
-- primal-dual simplex

data PDResult = PDUnsat | PDOptimal | PDUnbounded

primalDualSimplex :: (Real r, Fractional r) => OptDir -> Tableau r -> (Bool, Tableau r)
{-# SPECIALIZE primalDualSimplex :: OptDir -> Tableau Rational -> (Bool, Tableau Rational) #-}
{-# SPECIALIZE primalDualSimplex :: OptDir -> Tableau Double -> (Bool, Tableau Double) #-}
primalDualSimplex :: OptDir -> Tableau r -> (Bool, Tableau r)
primalDualSimplex OptDir
optdir = Tableau r -> (Bool, Tableau r)
forall r. (Real r, Fractional r) => Tableau r -> (Bool, Tableau r)
go
  where
    go :: Tableau r -> (Bool, Tableau r)
go Tableau r
tbl =
      case OptDir -> Tableau r -> Either PDResult (Tableau r)
forall r.
(Real r, Fractional r) =>
OptDir -> Tableau r -> Either PDResult (Tableau r)
pdPivot OptDir
optdir Tableau r
tbl of
        Left PDResult
PDOptimal   -> Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Tableau r -> Bool
forall r. Real r => Tableau r -> Bool
isFeasible Tableau r
tbl) ((Bool, Tableau r) -> (Bool, Tableau r))
-> (Bool, Tableau r) -> (Bool, Tableau r)
forall a b. (a -> b) -> a -> b
$ Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (OptDir -> Tableau r -> Bool
forall r. Real r => OptDir -> Tableau r -> Bool
isOptimal OptDir
optdir Tableau r
tbl) ((Bool, Tableau r) -> (Bool, Tableau r))
-> (Bool, Tableau r) -> (Bool, Tableau r)
forall a b. (a -> b) -> a -> b
$ (Bool
True, Tableau r
tbl)
        Left PDResult
PDUnsat     -> Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Bool -> Bool
not (Tableau r -> Bool
forall r. Real r => Tableau r -> Bool
isFeasible Tableau r
tbl)) ((Bool, Tableau r) -> (Bool, Tableau r))
-> (Bool, Tableau r) -> (Bool, Tableau r)
forall a b. (a -> b) -> a -> b
$ (Bool
False, Tableau r
tbl)
        Left PDResult
PDUnbounded -> Bool -> (Bool, Tableau r) -> (Bool, Tableau r)
forall a. (?callStack::CallStack) => Bool -> a -> a
assert (Bool -> Bool
not (OptDir -> Tableau r -> Bool
forall r. Real r => OptDir -> Tableau r -> Bool
isOptimal OptDir
optdir Tableau r
tbl)) ((Bool, Tableau r) -> (Bool, Tableau r))
-> (Bool, Tableau r) -> (Bool, Tableau r)
forall a b. (a -> b) -> a -> b
$ (Bool
False, Tableau r
tbl)
        Right Tableau r
tbl' -> Tableau r -> (Bool, Tableau r)
go Tableau r
tbl'

pdPivot :: (Real r, Fractional r) => OptDir -> Tableau r -> Either PDResult (Tableau r)
pdPivot :: OptDir -> Tableau r -> Either PDResult (Tableau r)
pdPivot OptDir
optdir Tableau r
tbl
  | [(Either Int Any, r)] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [(Either Int Any, r)]
forall b. [(Either Int b, r)]
ps Bool -> Bool -> Bool
&& [(Either Any Int, r)] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [(Either Any Int, r)]
forall a. [(Either a Int, r)]
qs = PDResult -> Either PDResult (Tableau r)
forall a b. a -> Either a b
Left PDResult
PDOptimal -- optimal
  | Bool
otherwise =
      case Either Int Int
ret of
        Left Int
p -> -- primal update
          let rs :: [(Int, r)]
rs = [ (Int
i, (r
bi r -> r -> r
forall a. Num a => a -> a -> a
- r
t) r -> r -> r
forall a. Fractional a => a -> a -> a
/ r
y_ip)
                   | (Int
i, (VarMap r
row_i, r
bi)) <- Tableau r -> [(Int, Row r)]
forall a. IntMap a -> [(Int, a)]
IM.toList Tableau r
tbl, Int
i Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= Int
objRowIndex
                   , let y_ip :: r
y_ip = r -> Int -> VarMap r -> r
forall a. a -> Int -> IntMap a -> a
IM.findWithDefault r
0 Int
p VarMap r
row_i, r
y_ip r -> r -> Bool
forall a. Ord a => a -> a -> Bool
> r
0
                   ]
              q :: Int
q = (Int, r) -> Int
forall a b. (a, b) -> a
fst ((Int, r) -> Int) -> (Int, r) -> Int
forall a b. (a -> b) -> a -> b
$ ((Int, r) -> (Int, r) -> Ordering) -> [(Int, r)] -> (Int, r)
forall (t :: * -> *) a.
Foldable t =>
(a -> a -> Ordering) -> t a -> a
minimumBy (((Int, r) -> r) -> (Int, r) -> (Int, r) -> Ordering
forall a b. Ord a => (b -> a) -> b -> b -> Ordering
comparing (Int, r) -> r
forall a b. (a, b) -> b
snd) [(Int, r)]
rs
          in if [(Int, r)] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [(Int, r)]
rs
             then PDResult -> Either PDResult (Tableau r)
forall a b. a -> Either a b
Left PDResult
PDUnsat
             else Tableau r -> Either PDResult (Tableau r)
forall a b. b -> Either a b
Right (Int -> Int -> Tableau r -> Tableau r
forall r.
(Fractional r, Eq r) =>
Int -> Int -> Tableau r -> Tableau r
pivot Int
q Int
p Tableau r
tbl)
        Right Int
q -> -- dual update
          let (VarMap r
row_q, r
_bq) = (Tableau r
tbl Tableau r -> Int -> Row r
forall a. IntMap a -> Int -> a
IM.! Int
q)
              cs :: [(Int, r)]
cs = [ (Int
j, (r
cj'r -> r -> r
forall a. Num a => a -> a -> a
-r
t) r -> r -> r
forall a. Fractional a => a -> a -> a
/ (-r
y_qj))
                   | (Int
j, r
y_qj) <- VarMap r -> [(Int, r)]
forall a. IntMap a -> [(Int, a)]
IM.toList VarMap r
row_q
                   , r
y_qj r -> r -> Bool
forall a. Ord a => a -> a -> Bool
< r
0
                   , let cj :: r
cj  = r -> Int -> VarMap r -> r
forall a. a -> Int -> IntMap a -> a
IM.findWithDefault r
0 Int
j VarMap r
obj
                   , let cj' :: r
cj' = if OptDir
optdirOptDir -> OptDir -> Bool
forall a. Eq a => a -> a -> Bool
==OptDir
OptMax then r
cj else -r
cj
                   ]
              p :: Int
p = (Int, r) -> Int
forall a b. (a, b) -> a
fst ((Int, r) -> Int) -> (Int, r) -> Int
forall a b. (a -> b) -> a -> b
$ ((Int, r) -> (Int, r) -> Ordering) -> [(Int, r)] -> (Int, r)
forall (t :: * -> *) a.
Foldable t =>
(a -> a -> Ordering) -> t a -> a
maximumBy (((Int, r) -> r) -> (Int, r) -> (Int, r) -> Ordering
forall a b. Ord a => (b -> a) -> b -> b -> Ordering
comparing (Int, r) -> r
forall a b. (a, b) -> b
snd) [(Int, r)]
cs
              (VarMap r
obj,r
_) = Int -> Tableau r -> Row r
forall r. Int -> Tableau r -> Row r
lookupRow Int
objRowIndex Tableau r
tbl
          in if [(Int, r)] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [(Int, r)]
cs
             then PDResult -> Either PDResult (Tableau r)
forall a b. a -> Either a b
Left PDResult
PDUnbounded -- dual infeasible
             else Tableau r -> Either PDResult (Tableau r)
forall a b. b -> Either a b
Right (Int -> Int -> Tableau r -> Tableau r
forall r.
(Fractional r, Eq r) =>
Int -> Int -> Tableau r -> Tableau r
pivot Int
q Int
p Tableau r
tbl)
  where
    qs :: [(Either a Int, r)]
qs = [ (Int -> Either a Int
forall a b. b -> Either a b
Right Int
i, r
bi) | (Int
i, (VarMap r
_row_i, r
bi)) <- Tableau r -> [(Int, Row r)]
forall a. IntMap a -> [(Int, a)]
IM.toList Tableau r
tbl, Int
i Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= Int
objRowIndex, r
0 r -> r -> Bool
forall a. Ord a => a -> a -> Bool
> r
bi ]
    ps :: [(Either Int b, r)]
ps = [ (Int -> Either Int b
forall a b. a -> Either a b
Left Int
j, r
cj')
         | (Int
j,r
cj) <- VarMap r -> [(Int, r)]
forall a. IntMap a -> [(Int, a)]
IM.toList (Row r -> VarMap r
forall a b. (a, b) -> a
fst (Int -> Tableau r -> Row r
forall r. Int -> Tableau r -> Row r
lookupRow Int
objRowIndex Tableau r
tbl))
         , let cj' :: r
cj' = if OptDir
optdirOptDir -> OptDir -> Bool
forall a. Eq a => a -> a -> Bool
==OptDir
OptMax then r
cj else -r
cj
         , r
0 r -> r -> Bool
forall a. Ord a => a -> a -> Bool
> r
cj' ]
    (Either Int Int
ret, r
t) = ((Either Int Int, r) -> (Either Int Int, r) -> Ordering)
-> [(Either Int Int, r)] -> (Either Int Int, r)
forall (t :: * -> *) a.
Foldable t =>
(a -> a -> Ordering) -> t a -> a
minimumBy (((Either Int Int, r) -> r)
-> (Either Int Int, r) -> (Either Int Int, r) -> Ordering
forall a b. Ord a => (b -> a) -> b -> b -> Ordering
comparing (Either Int Int, r) -> r
forall a b. (a, b) -> b
snd) ([(Either Int Int, r)]
forall a. [(Either a Int, r)]
qs [(Either Int Int, r)]
-> [(Either Int Int, r)] -> [(Either Int Int, r)]
forall a. [a] -> [a] -> [a]
++ [(Either Int Int, r)]
forall b. [(Either Int b, r)]
ps)

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

toCSV :: (Num r) => (r -> String) -> Tableau r -> String
toCSV :: (r -> String) -> Tableau r -> String
toCSV r -> String
showCell Tableau r
tbl = [String] -> String
unlines ([String] -> String)
-> ([[String]] -> [String]) -> [[String]] -> String
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ([String] -> String) -> [[String]] -> [String]
forall a b. (a -> b) -> [a] -> [b]
map ([String] -> String
forall (t :: * -> *) a. Foldable t => t [a] -> [a]
concat ([String] -> String)
-> ([String] -> [String]) -> [String] -> String
forall b c a. (b -> c) -> (a -> b) -> a -> c
. String -> [String] -> [String]
forall a. a -> [a] -> [a]
intersperse String
",") ([[String]] -> String) -> [[String]] -> String
forall a b. (a -> b) -> a -> b
$ [String]
header [String] -> [[String]] -> [[String]]
forall a. a -> [a] -> [a]
: [[String]]
body
  where
    header :: [String]
    header :: [String]
header = String
"" String -> [String] -> [String]
forall a. a -> [a] -> [a]
: (Int -> String) -> [Int] -> [String]
forall a b. (a -> b) -> [a] -> [b]
map Int -> String
colName [Int]
cols [String] -> [String] -> [String]
forall a. [a] -> [a] -> [a]
++ [String
""]

    body :: [[String]]
    body :: [[String]]
body = [Int -> (IntMap r, r) -> [String]
showRow Int
i (Int -> Tableau r -> (IntMap r, r)
forall r. Int -> Tableau r -> Row r
lookupRow Int
i Tableau r
tbl) | Int
i <- [Int]
rows]

    rows :: [RowIndex]
    rows :: [Int]
rows = Tableau r -> [Int]
forall a. IntMap a -> [Int]
IM.keys (Int -> Tableau r -> Tableau r
forall a. Int -> IntMap a -> IntMap a
IM.delete Int
objRowIndex Tableau r
tbl) [Int] -> [Int] -> [Int]
forall a. [a] -> [a] -> [a]
++ [Int
objRowIndex]

    cols :: [ColIndex]
    cols :: [Int]
cols = [Int
0..Int
colMax]
      where
        colMax :: Int
colMax = [Int] -> Int
forall (t :: * -> *) a. (Foldable t, Ord a) => t a -> a
maximum (-Int
1 Int -> [Int] -> [Int]
forall a. a -> [a] -> [a]
: [Int
c | (IntMap r
row, r
_) <- Tableau r -> [(IntMap r, r)]
forall a. IntMap a -> [a]
IM.elems Tableau r
tbl, Int
c <- IntMap r -> [Int]
forall a. IntMap a -> [Int]
IM.keys IntMap r
row])

    rowName :: RowIndex -> String
    rowName :: Int -> String
rowName Int
i = if Int
iInt -> Int -> Bool
forall a. Eq a => a -> a -> Bool
==Int
objRowIndex then String
"obj" else String
"x" String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show Int
i

    colName :: ColIndex -> String
    colName :: Int -> String
colName Int
j = String
"x" String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show Int
j

    showRow :: Int -> (IntMap r, r) -> [String]
showRow Int
i (IntMap r
row, r
row_val) = Int -> String
rowName Int
i String -> [String] -> [String]
forall a. a -> [a] -> [a]
: [r -> String
showCell (r -> Int -> IntMap r -> r
forall a. a -> Int -> IntMap a -> a
IM.findWithDefault r
0 Int
j IntMap r
row') | Int
j <- [Int]
cols] [String] -> [String] -> [String]
forall a. [a] -> [a] -> [a]
++ [r -> String
showCell r
row_val]
      where row' :: IntMap r
row' = Int -> r -> IntMap r -> IntMap r
forall a. Int -> a -> IntMap a -> IntMap a
IM.insert Int
i r
1 IntMap r
row

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