packages feed

limp-cbc (empty) → 0.2.8.6

raw patch · 9 files changed

+470/−0 lines, 9 filesdep +basedep +containersdep +limpsetup-changed

Dependencies added: base, containers, limp, vector

Files

+ LICENSE view
@@ -0,0 +1,19 @@+Copyright (c) 2014 <copyright holders>++Permission is hereby granted, free of charge, to any person obtaining a copy+of this software and associated documentation files (the "Software"), to deal+in the Software without restriction, including without limitation the rights+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell+copies of the Software, and to permit persons to whom the Software is+furnished to do so, subject to the following conditions:++The above copyright notice and this permission notice shall be included in+all copies or substantial portions of the Software.++THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN+THE SOFTWARE.
+ Setup.hs view
@@ -0,0 +1,2 @@+import Distribution.Simple+main = defaultMain
+ cbits/Cbc.cpp view
@@ -0,0 +1,81 @@+#include <coin/CbcModel.hpp>+#include <coin/OsiClpSolverInterface.hpp>++#define FROM_CPP+#include "Cbc.h"+extern "C" {++CbcModel* newModel()+{+    OsiClpSolverInterface solver;+    CbcModel* model = new CbcModel(solver);+    return model;+}+++void    freeModel(CbcModel* model)+{+    delete model;+}+++void    loadProblem(CbcModel* model, int ncols,    int nrows,+                    int* starts,     int* indices, double* values,+                    double* collb,   double* colub,+                    double* obj,+                    double* rowlb, double* rowub)+{+    model->solver()->loadProblem+            ( ncols, nrows+            , starts, indices, values+            , collb, colub+            , obj+            , rowlb, rowub);+}+++void    setInteger (CbcModel* model, int col)+{+    model->solver()->setInteger(col);+}+++void    branchAndBound (CbcModel* model)+{+    model->branchAndBound();+}+++int getNumCols (CbcModel* model)+{+    return model->solver()->getNumCols();+}+++const double* getBestSolution (CbcModel* model)+{+    return model->bestSolution();+}+++void setObjSense(CbcModel* model, double dir)+{+    model->solver()->setObjSense(dir);+}++void setLogLevel(CbcModel* model, int level)+{+    model->setLogLevel(level);+}++int isProvenInfeasible(CbcModel* model)+{+    return model->isProvenInfeasible();+}++double getCoinDblMax()+{+    return COIN_DBL_MAX;+}++}
+ limp-cbc.cabal view
@@ -0,0 +1,46 @@+name:                limp-cbc+version:             0.2.8.6+synopsis:            bindings for integer linear programming solver Coin/CBC+description:         very simple binding to external solver, CBC.+                     CBC is somewhat faster than GLPK, and also has a more permissive licence.++license:             MIT+license-file:        LICENSE+author:              Amos Robinson+maintainer:          amos.robinson@gmail.com+category:            numeric+build-type:          Simple+cabal-version:       >=1.10+homepage:            https://github.com/amosr/limp-cbc+++source-repository head+    type: git+    location: git://github.com/amosr/limp-cbc.git++library+  hs-source-dirs: src+  exposed-modules:+        Numeric.Limp.Solvers.Cbc+        Numeric.Limp.Solvers.Cbc.Solve+        Numeric.Limp.Solvers.Cbc.MatrixRepr++  other-modules:       +        Numeric.Limp.Solvers.Cbc.Internal.Foreign+        Numeric.Limp.Solvers.Cbc.Internal.Wrapper+  build-depends:+        base        < 5,+        containers  == 0.5.*,+        vector      == 0.10.*,+        limp++  ghc-options: -Wall -fno-warn-orphans+  default-language: Haskell2010+  default-extensions:       TemplateHaskell TypeFamilies FlexibleContexts GeneralizedNewtypeDeriving DataKinds GADTs RankNTypes++  build-tools:      c2hs++  extra-libraries:  Cbc Clp CbcSolver Cgl Osi OsiCbc OsiClp OsiCommonTests CoinUtils CoinMP stdc+++  Include-dirs:     cbits+  includes:         Cbc.h coin/Cbc_C_Interface.h+  C-sources:        cbits/Cbc.cpp
+ src/Numeric/Limp/Solvers/Cbc.hs view
@@ -0,0 +1,12 @@+module Numeric.Limp.Solvers.Cbc where+import qualified Numeric.Limp.Canon as C+import Numeric.Limp.Program+import Numeric.Limp.Rep++import qualified Numeric.Limp.Solvers.Cbc.Solve as S++solve :: (Ord z, Ord r) => Program z r IntDouble -> Either S.Error (Assignment z r IntDouble)+solve p+ = let pc = C.program p+   in  S.solve pc+
+ src/Numeric/Limp/Solvers/Cbc/Internal/Foreign.chs view
@@ -0,0 +1,102 @@+{-# OPTIONS_GHC -fno-warn-unused-binds -fno-warn-unused-do-bind -fno-warn-unused-matches #-}+module Numeric.Limp.Solvers.Cbc.Internal.Foreign where++import Foreign+import Foreign.C++import Control.Applicative++import qualified Data.Vector.Storable as V+import Data.Vector.Storable (Vector)++import Unsafe.Coerce++#include "Cbc.h"+++newtype CbcModel = CbcModel (ForeignPtr CbcModel)+{#pointer *CbcModel as CbcModelPtr -> CbcModel #}+++foreign import ccall "&freeModel"+   cbcDeleteModel_funptr :: FunPtr (Ptr CbcModel -> IO ())++mkCbcModel :: CbcModelPtr -> IO CbcModel+mkCbcModel p+    = CbcModel <$> newForeignPtr cbcDeleteModel_funptr p++withCbcModel :: CbcModel -> (CbcModelPtr -> IO b) -> IO b+withCbcModel (CbcModel ptr) f+    = do    withForeignPtr ptr f++{#fun newModel as ^+    { } -> `CbcModel' mkCbcModel* #}++withVecI    :: Vector Int -> (Ptr CInt -> IO b) -> IO b+withVecI v f+ = V.unsafeWith (V.map fromIntegral v) f++withVecD    :: Vector Double -> (Ptr CDouble -> IO b) -> IO b+withVecD v f+ = V.unsafeWith v (f . castPtr)+++{#fun loadProblem as ^+    { withCbcModel* `CbcModel'++    ,               `Int'               -- # rows+    ,               `Int'               -- # columns++    , withVecI*     `Vector Int'        -- starts+    , withVecI*     `Vector Int'        -- indices+    , withVecD*     `Vector Double'     -- values++    , withVecD*     `Vector Double'     -- col lower bounds+    , withVecD*     `Vector Double'     -- col upper bounds++    , withVecD*     `Vector Double'     -- objective coefficients+    , withVecD*     `Vector Double'     -- row lower bounds+    , withVecD*     `Vector Double'     -- row upper bounds+    }+    -> `()' #}+++{#fun setInteger as ^+    { withCbcModel* `CbcModel'+    ,               `Int' }+    -> `()' #}+++{#fun branchAndBound as ^+    { withCbcModel* `CbcModel' }+    -> `()' #}+++getSolution :: CbcModel -> IO (Vector Double)+getSolution m+ = do   ncols <- fromIntegral <$> withCbcModel m {#call getNumCols#}+        vd    <- unsafeCoerce <$> withCbcModel m {#call getBestSolution#}++        V.generateM ncols (peekElemOff vd)+++{#fun setObjSense as ^+    { withCbcModel* `CbcModel'+    ,               `Double' }+    -> `()' #}+++{#fun setLogLevel as ^+    { withCbcModel* `CbcModel'+    ,               `Int' }+    -> `()' #}+++{#fun isProvenInfeasible as ^+    { withCbcModel* `CbcModel' }+    -> `Bool' #}+++{#fun pure getCoinDblMax as ^+    { } -> `Double' #}+
+ src/Numeric/Limp/Solvers/Cbc/Internal/Wrapper.hs view
@@ -0,0 +1,49 @@+-- One small step above the Foreign interface+module Numeric.Limp.Solvers.Cbc.Internal.Wrapper+    ( F.CbcModel+    , F.newModel+    , F.branchAndBound+    , F.setInteger+    , F.setObjSense+    , F.setLogLevel+    , F.getSolution++    , F.isProvenInfeasible+    , F.getCoinDblMax++    , setQuiet+    , loadProblem+    ) where++import qualified Numeric.Limp.Solvers.Cbc.Internal.Foreign as F++import qualified Data.Vector.Storable as V+import Data.Vector.Storable (Vector)++setQuiet :: F.CbcModel -> IO ()+setQuiet = flip F.setLogLevel 0++loadProblem+ :: F.CbcModel+ -> Vector Int -> Vector Int -> Vector Double+ -> Vector Double -> Vector Double+ -> Vector Double -> Vector Double -> Vector Double+ -> IO ()+loadProblem model starts indices values collb colub obj rowlb rowub+ = let ncols = V.length collb+       nrows = V.length rowlb+   in  do   check_len "colub = ncols" colub ncols+            check_len "rowub = nrows" rowub nrows++            check_len "obj   = ncols" obj   ncols++            -- check_len "indices = values" indices (V.length values)+            F.loadProblem model ncols nrows starts indices values collb colub obj rowlb rowub+ where+  check_len e v l+   | V.length v == l+   = return ()+   | otherwise+   = putStrLn+   $ "Numeric.Limp.Solvers.Cbc.Internal.Foreign:loadProblem_check wrong length! " ++ e ++ ": " ++ show (V.length v, l)+
+ src/Numeric/Limp/Solvers/Cbc/MatrixRepr.hs view
@@ -0,0 +1,117 @@+module Numeric.Limp.Solvers.Cbc.MatrixRepr where+import Numeric.Limp.Canon+import Numeric.Limp.Rep++import Data.Either (isLeft)++import qualified Data.Map as M+import qualified Data.Set as S++import qualified Data.Vector.Storable as V+import Data.Vector.Storable (Vector)+import Numeric.Limp.Solvers.Cbc.Internal.Wrapper+++-- | A half-sparse, half-dense matrix (tableau?) representation of a linear program.+-- Columns are variables, and constraints are rows.+data MatrixRepr+ = MatrixRepr+ { _starts   :: Vector Int+ -- ^ starting indices segment descriptor for indsvals.+ -- one per row (constraint).+ -- but for some reason, there's an extra one at the end.+ -- starts ! 0 == 0 && starts ! (length starts - 1) == length indsvals+ -- (I think)+ --+ , _inds     :: Vector Int+ , _vals     :: Vector Double+ -- ^ segmented array. each pair is an index of the column (variable), and the coefficient.+ -- the segment number is the row (constraint).+ -- ^ Coefficients of the constraints++ , _colLs    :: Vector Double+ , _colUs    :: Vector Double+ -- ^ Lower and upper bound for each variable++ , _obj      :: Vector Double+ -- ^ Coefficients of each variable in the objective function++ , _rowLs    :: Vector Double+ , _rowUs    :: Vector Double+ -- ^ Lower and upper bound for each constraint++ , _ints     :: Vector Int+ -- ^ Indices of each integer variable.+ }+ deriving Show+++matrixReprOfProgram :: (Ord z, Ord r) => Program z r IntDouble -> MatrixRepr+matrixReprOfProgram p+ = MatrixRepr starts inds vals colLs colUs obj rowLs rowUs ints+ where++  starts+   = V.fromList+   $ scanl (+) 0+   $ map length nestedindvals++  inds = V.fromList $ map fst $ concat nestedindvals+  vals = V.fromList $ map snd $ concat nestedindvals++  nestedindvals = map getIndVals $ S.toList vs++  colLs = V.fromList $ map (fst . boundsOfVar) $ S.toList vs+  colUs = V.fromList $ map (snd . boundsOfVar) $ S.toList vs++  obj+   = case _objective p of+     Linear m+      -> V.fromList+       $ map (maybe 0 unwrapR . flip M.lookup m) $ S.toList vs++  rowLs = V.fromList $ map (fst . boundsOfCon) cons+  rowUs = V.fromList $ map (snd . boundsOfCon) cons++  ints = V.fromList $ map snd $ filter (isLeft . fst) $ S.toList vs `zip` [0..]+++  vs = varsOfProgram p+  -- nvars = S.size vs+  -- ncons = length cons+  cons = case _constraints p of+         Constraint cs -> cs++  getIndVals v+   = concatMap (\(C1 _ (Linear m) _, ix) ->+                    case M.lookup v m of+                    Nothing -> []+                    Just mul -> [(ix, unwrapR mul)]) (cons `zip` [0..])+++  boundsOfVar v+   = case M.lookup v $ _bounds p of+     Nothing    -> (lower Nothing, upper Nothing)+     Just (l,u) -> (lower l,       upper u)++  boundsOfCon (C1 l _ u) = (lower l, upper u)++  upper = maybe   getCoinDblMax  unwrapR+  lower = maybe (-getCoinDblMax) unwrapR+++makeAssignment :: (Ord z, Ord r) => Program z r IntDouble -> Vector Double -> Assignment z r IntDouble+makeAssignment p values+ = go M.empty M.empty+ $ S.toList (varsOfProgram p) `zip` V.toList   values+ where+  go mz mr []+   = Assignment mz mr++  go mz mr ((Left z, val) : rest)+   = go (M.insert z (Z $ truncate val) mz) mr rest++  go mz mr ((Right r, val) : rest)+   = go mz (M.insert r (R $ val) mr) rest++
+ src/Numeric/Limp/Solvers/Cbc/Solve.hs view
@@ -0,0 +1,42 @@+module Numeric.Limp.Solvers.Cbc.Solve where+import Numeric.Limp.Canon+import Numeric.Limp.Rep++import Numeric.Limp.Solvers.Cbc.MatrixRepr++import qualified Data.Vector.Storable as V+import Numeric.Limp.Solvers.Cbc.Internal.Wrapper++import System.IO.Unsafe++data Error+ = Infeasible+ deriving Show++solve :: (Ord z, Ord r) => Program z r IntDouble -> Either Error (Assignment z r IntDouble)+solve p+ = let mr = matrixReprOfProgram p+   in  unsafePerformIO $ do+            m <- newModel++            setQuiet m++            loadProblem m (_starts mr) (_inds mr)  (_vals mr)+                          (_colLs mr)  (_colUs mr) (_obj mr)+                          (_rowLs mr)  (_rowUs mr)++            V.forM_ (_ints mr) $ setInteger m++            setObjSense m 1++            branchAndBound m++            infeasible <- isProvenInfeasible m+            -- TODO get other statuses+            case infeasible of+             True+              ->    return $ Left  $ Infeasible++             False+              -> do vs <- getSolution m+                    return $ Right $ makeAssignment p vs