comfort-glpk 0.0.1 → 0.1

raw patch · 11 files changed

+370/−609 lines, 11 filesdep +linear-programmingdep ~deepseqdep ~randomdep ~transformersPVP ok

version bump matches the API change (PVP)

Dependencies added: linear-programming

Dependency ranges changed: deepseq, random, transformers

API changes (from Hackage documentation)

Files

• Makefile +11/−0
• comfort-glpk.cabal +20/−4
• src/Numeric/GLPK.hs +27/−300
• src/Numeric/GLPK/Monad.hs +85/−6
• src/Numeric/GLPK/Private.hs +13/−33
• src/debug-off/Numeric/GLPK/Debug.hs +8/−0
• src/debug-on/Numeric/GLPK/Debug.hs +5/−0
• test/Test/Numeric/GLPK.hs +57/−52
• test/Test/Numeric/GLPK/Generator.hs +0/−208
• test/Test/Numeric/GLPK/Monad.hs +89/−6
• test/Test/Numeric/GLPK/Utility.hs +55/−0

@@ -0,0 +1,11 @@+run-test:	update-test+	runhaskell Setup configure --user --enable-tests+	runhaskell Setup build+	runhaskell Setup haddock+	runhaskell Setup test comfort-glpk-test --show-details=streaming++	runhaskell Setup configure --user -fdebug+	runhaskell Setup build++update-test:+	doctest-extract-0.1 -i src/ -o test/ --executable-main=Main.hs Numeric.GLPK Numeric.GLPK.Monad

@@ -1,6 +1,6 @@ Cabal-Version:    2.2 Name:             comfort-glpk-Version:          0.0.1+Version:          0.1 License:          BSD-3-Clause License-File:     LICENSE Author:           Henning Thielemann <haskell@henning-thielemann.de>@@ -32,8 +32,16 @@   .   Alternatives: @coinor-clp@, @hmatrix-glpk@, @glpk-hs@ +Extra-Source-Files:+  Makefile++Flag debug+  Description: Enable debug output+  Default:     False+  Manual:      True+ Source-Repository this-  Tag:         0.0.1+  Tag:         0.1   Type:        darcs   Location:    https://hub.darcs.net/thielema/comfort-glpk/ @@ -44,6 +52,7 @@ Library   Build-Depends:     glpk-headers >=0.4.1 && <0.6,+    linear-programming >=0.0 && <0.1,     comfort-array >=0.4 && <0.6,     deepseq >=1.3 && <1.5,     transformers >=0.3 && <0.7,@@ -52,20 +61,27 @@     base >=4.5 && <5    GHC-Options:      -Wall-  Hs-Source-Dirs:   src   Default-Language: Haskell98+  Hs-Source-Dirs:   src+  If flag(debug)+    Hs-Source-Dirs: src/debug-on+  Else+    Hs-Source-Dirs: src/debug-off   Extra-Libraries: glpk   Exposed-modules:     Numeric.GLPK     Numeric.GLPK.Monad   Other-Modules:     Numeric.GLPK.Private+    Numeric.GLPK.Debug  Test-Suite comfort-glpk-test   Type:             exitcode-stdio-1.0   Build-Depends:     comfort-glpk,+    linear-programming,     comfort-array >=0.5.2,+    transformers,     non-empty,     utility-ht >=0.0.17,     doctest-exitcode-stdio >=0.0 && <0.1,@@ -79,6 +95,6 @@   Default-Language: Haskell98   Main-Is: Main.hs   Other-Modules:-    Test.Numeric.GLPK.Generator+    Test.Numeric.GLPK.Utility     Test.Numeric.GLPK.Monad     Test.Numeric.GLPK

@@ -23,53 +23,38 @@ instead of @[i >=. a, i <=. b]@ use @i >=<. (a,b)@. -} module Numeric.GLPK (-   Term(..),+   Term,    Bound(..),    Inequality(..),-   free, (<=.), (>=.), (==.), (>=<.),-   NoSolutionType(..),+   LP.free, (LP.<=.), (LP.>=.), (LP.==.), (LP.>=<.),+   FailureType(..),    SolutionType(..),-   Solution,+   Result,    Constraints,    Direction(..),    Objective,    Bounds,    (.*),-   objectiveFromTerms,+   LP.objectiveFromTerms,    simplex,-   simplexMulti,-   simplexSuccessive,    exact,-   exactMulti,-   exactSuccessive,    interior,-   interiorMulti,-   interiorSuccessive,--   solveSuccessive,--   FormatIdentifier,-   formatMathProg,    ) where  import qualified Math.Programming.Glpk.Header as FFI+import qualified Numeric.GLPK.Debug as Debug+import qualified Numeric.LinearProgramming.Common as LP import Numeric.GLPK.Private+import Numeric.LinearProgramming.Common+         (Bound(..), Inequality(Inequality),+          Bounds, Direction(..), Objective, (.*))  import qualified Data.Array.Comfort.Storable.Mutable as Mutable import qualified Data.Array.Comfort.Storable as Array import qualified Data.Array.Comfort.Shape as Shape-import qualified Data.NonEmpty as NonEmpty-import qualified Data.List as List-import Data.Array.Comfort.Storable (Array)-import Data.Tuple.HT (mapFst, mapSnd)-import Data.Traversable (for) import Data.Foldable (for_) -import Text.Printf (printf)--import qualified Control.Monad.Trans.Except as ME-import qualified Control.Monad.Trans.State as MS-import Control.Monad (void, when)+import Control.Monad (void) import Control.Applicative (liftA2) import Control.Exception (bracket) @@ -80,7 +65,7 @@   {- $setup->>> import qualified Test.Numeric.GLPK.Generator as TestLP+>>> import qualified Numeric.LinearProgramming.Test as TestLP >>> import qualified Numeric.GLPK as LP >>> import Numeric.GLPK ((.*), (<=.), (==.)) >>>@@ -107,176 +92,37 @@ -}  -infix 7 .*--(.*) :: Double -> ix -> Term ix-(.*) = Term---infix 4 <=., >=., >=<., ==.--(<=.), (>=.), (==.) :: x -> Double -> Inequality x-x <=. bnd = Inequality x $ LessEqual bnd-x >=. bnd = Inequality x $ GreaterEqual bnd-x ==. bnd = Inequality x $ Equal bnd--(>=<.) :: x -> (Double,Double) -> Inequality x-x >=<. bnd = Inequality x $ uncurry Between bnd--free :: x -> Inequality x-free x = Inequality x Free----objectiveFromTerms ::-   (Shape.Indexed sh, Shape.Index sh ~ ix) => sh -> [Term ix] -> Objective sh-objectiveFromTerms sh =-   Array.fromAssociations 0 sh . map (\(Term x ix) -> (ix,x))-- {- | >>> case Shape.indexTupleFromShape tripletShape of (x1,x2,x3) -> mapSnd (mapSnd Array.toTuple) <$> LP.simplex [] [[2.*x1, 1.*x2] <=. 10, [1.*x2, 5.*x3] <=. 20] (LP.Maximize, Array.fromTuple (4,-3,2) :: Array.Array TripletShape Double) Right (Optimal,(28.0,(5.0,0.0,4.0)))  prop> \target -> case Shape.indexTupleFromShape pairShape of (pos,neg) -> case mapSnd (mapSnd Array.toTuple) <$> LP.simplex [] [[1.*pos, (-1).*neg] ==. target] (LP.Minimize, Array.fromTuple (1,1) :: Array.Array PairShape Double) of (Right (LP.Optimal,(absol,(posResult,negResult)))) -> QC.property (TestLP.approxReal 0.001 absol (abs target)) .&&. (posResult === 0 .||. negResult === 0); _ -> QC.property False prop> \(QC.Positive posWeight) (QC.Positive negWeight) target -> case Shape.indexTupleFromShape pairShape of (pos,neg) -> case mapSnd (mapSnd Array.toTuple) <$> LP.simplex [] [[1.*pos, (-1).*neg] ==. target] (LP.Minimize, Array.fromTuple (posWeight,negWeight) :: Array.Array PairShape Double) of (Right (LP.Optimal,(absol,(posResult,negResult)))) -> QC.property (absol>=0) .&&. (posResult === 0 .||. negResult === 0); _ -> QC.property False-prop> QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, _) -> True; _ -> False+prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, _) -> True; _ -> False+prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, (_,sol)) -> TestLP.checkFeasibility 0.1 bounds constrs sol; _ -> QC.property False+prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, (_,sol)) -> QC.forAll (QC.choose (0,1)) $ \lambda -> TestLP.checkFeasibility 0.1 bounds constrs $ TestLP.affineCombination lambda sol (Array.map fromIntegral origin); _ -> QC.property False+prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, (opt,sol)) -> QC.forAll (QC.choose (0,1)) $ \lambda -> let val = TestLP.scalarProduct obj $ TestLP.affineCombination lambda sol (Array.map fromIntegral origin) in (case dir of LP.Minimize -> opt-0.01 <= val; LP.Maximize -> opt+0.01 >= val); _ -> QC.property False -} simplex ::    (Shape.Indexed sh, Shape.Index sh ~ ix) =>    Bounds ix -> Constraints ix ->-   (Direction, Objective sh) -> Solution sh+   (Direction, Objective sh) -> Result sh simplex = solve (flip FFI.glp_simplex nullPtr) -{-# DEPRECATED simplexMulti "use GLPK.Monad instead" #-}-{-# DEPRECATED exactMulti "use GLPK.Monad instead" #-}-{-# DEPRECATED interiorMulti "run 'interior' in Either monad instead" #-}-{- |-Optimize for one objective after another.-That is, if the first optimization succeeds-then the optimum is fixed as constraint-and the optimization is continued with respect to the second objective and so on.-The iteration fails if one optimization fails.-The obtained objective values are returned as well.-Their number equals the number of attempted optimizations. -The last objective value is included in the Solution value.-This is a bit inconsistent,-but this way you have a warranty that there is an objective value-if the optimization is successful.--The objectives are expected as 'Term's-because after successful optimization step-they are used as (sparse) constraints.-It's also easy to assert that the same array shape-is used for all objectives.--The function does not work reliably,-because an added objective can make the system infeasible-due to rounding errors.-E.g. a non-negative objective can become very small but negative.---prop> QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAllShrink (TestLP.genProblem origin) TestLP.shrinkProblem $ \(bounds, constrs) -> QC.forAllShrink (TestLP.genObjectives origin) TestLP.shrinkObjectives $ \objs -> case LP.simplexMulti bounds constrs (Array.shape origin) objs of (_, Right (LP.Optimal, _)) -> QC.property True; result -> QC.counterexample (show result) False--The same property fails for 'exactMulti' and 'interiorMulti'.-I guess, due to rounding errors.--}-simplexMulti, exactMulti, interiorMulti ::-   (Shape.Indexed sh, Shape.Index sh ~ ix) =>-   Bounds ix -> Constraints ix ->-   sh -> NonEmpty.T [] (Direction, [Term ix]) -> ([Double], Solution sh)-simplexMulti = solveMulti . simplex-exactMulti = solveMulti . exact-interiorMulti = solveMulti . interior--solveMulti ::-   (Shape.Indexed sh, Shape.Index sh ~ ix) =>-   (Constraints ix -> (Direction, Objective sh) -> Solution sh) ->-   Constraints ix ->-   sh -> NonEmpty.T [] (Direction, [Term ix]) -> ([Double], Solution sh)-solveMulti solver constrs0 sh (NonEmpty.Cons obj0 objs0) =-   let go constrs curObj ((dir,obj):objs) (Right (Optimal, (opt,_))) =-         mapFst (opt:) $-         let extConstrs = (curObj==.opt) : constrs in-         go extConstrs obj objs $-         solver extConstrs (dir, objectiveFromTerms sh obj)-       go _ _ _ sol = ([], sol)-   in go constrs0 (snd obj0) objs0 $-      solver constrs0 $ mapSnd (objectiveFromTerms sh) obj0---{-# DEPRECATED simplexSuccessive "use GLPK.Monad instead" #-}-{-# DEPRECATED exactSuccessive "use GLPK.Monad instead" #-}-{-# DEPRECATED interiorSuccessive "run 'interior' in Either monad instead" #-} {- |-Like the @Multi@ functions,-but allows not only to fix the previously-found optimal solution as constraint,-but allows constraints with a tolerance.-This is necessary in the presence of rounding errors.--prop> QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAllShrink (TestLP.genProblem origin) TestLP.shrinkProblem $ \(bounds, constrs) -> QC.forAllShrink (TestLP.genObjectives origin) TestLP.shrinkObjectives $ \objs -> case uncurry (LP.simplexSuccessive bounds constrs) $ TestLP.successiveObjectives origin 0.01 objs of result -> QC.counterexample (show result) $ case result of Right results -> all (\r -> case r of (LP.Optimal, _) -> True; _ -> False) results; _ -> False-prop> QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAllShrink (TestLP.genProblem origin) TestLP.shrinkProblem $ \(bounds, constrs) -> QC.forAllShrink (TestLP.genObjectives origin) TestLP.shrinkObjectives $ \objs -> case uncurry (LP.exactSuccessive bounds constrs) $ TestLP.successiveObjectives origin 0.01 objs of result -> QC.counterexample (show result) $ case result of Right results -> all (\r -> case r of (LP.Optimal, _) -> True; _ -> False) results; _ -> False--}-simplexSuccessive, exactSuccessive, interiorSuccessive ::-   (Traversable f, Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix) =>-   Bounds ix -> Constraints ix ->-   (Direction, Objective sh) ->-   f ((SolutionType, (Double, Array sh Double)) -> Constraints ix,-      (Direction, Objective sh)) ->-   Either NoSolutionType-      (NonEmpty.T f (SolutionType, (Double, Array sh Double)))-simplexSuccessive = solveSuccessiveInPlace (flip FFI.glp_simplex nullPtr)-exactSuccessive = solveSuccessiveInPlace (flip FFI.glp_exact nullPtr)-interiorSuccessive = solveSuccessive . interior--{-# DEPRECATED solveSuccessive-      "run simple solvers in GLPK.Monad or Either monad instead" #-}-{- |-Allows for generic implementation of 'simplexSuccessive' et.al.-without reuse of interim results.--prop> QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjectives origin) $ (. TestLP.successiveObjectives origin 0.01) $ \(obj,objs) -> case (LP.simplexSuccessive bounds constrs obj objs, LP.solveSuccessive (LP.simplex bounds) constrs obj objs) of (resultA,resultB) -> TestLP.approxSuccession 0.01 resultA resultB-prop> QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjectives origin) $ (. TestLP.successiveObjectives origin 0.01) $ \(obj,objs) -> case (LP.exactSuccessive bounds constrs obj objs, LP.solveSuccessive (LP.exact bounds) constrs obj objs) of (resultA,resultB) -> TestLP.approxSuccession 0.01 resultA resultB--}-solveSuccessive ::-   (Traversable f, Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix) =>-   (Constraints ix -> (Direction, Objective sh) -> Solution sh) ->-   Constraints ix ->-   (Direction, Objective sh) ->-   f ((SolutionType, (Double, Array sh Double)) -> Constraints ix,-      (Direction, Objective sh)) ->-   Either NoSolutionType-      (NonEmpty.T f (SolutionType, (Double, Array sh Double)))-solveSuccessive solver constrs0 obj0 objs = do-   let checkShape obj =-         if Array.shape (snd obj0) == Array.shape obj-            then obj-            else error "GLPK.solveSuccessive: objective shapes mismatch"-   let solveWithConstraints constrs problem =-         (\sol -> (sol, (constrs,sol))) <$> solver constrs problem-   (sol0,state0) <- solveWithConstraints constrs0 obj0-   NonEmpty.cons sol0 <$>-      MS.evalStateT-         (for objs $-            \(newConstrs,(dir,obj)) -> MS.StateT $ \(constrs,sol) ->-               solveWithConstraints-                  (newConstrs sol ++ constrs)-                  (dir, checkShape obj))-         state0---{- | >>> case Shape.indexTupleFromShape tripletShape of (x1,x2,x3) -> mapSnd (mapSnd Array.toTuple) <$> LP.exact [] [[2.*x1, 1.*x2] <=. 10, [1.*x2, 5.*x3] <=. 20] (LP.Maximize, Array.fromTuple (4,-3,2) :: Array.Array TripletShape Double) Right (Optimal,(28.0,(5.0,0.0,4.0))) -prop> QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (LP.simplex bounds constrs (dir,obj), LP.exact bounds constrs (dir,obj)) of (Right (LP.Optimal, (optSimplex,_)), Right (LP.Optimal, (optExact,_))) -> TestLP.approx "optimum" 0.001 optSimplex optExact; _ -> QC.property False+prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (LP.simplex bounds constrs (dir,obj), LP.exact bounds constrs (dir,obj)) of (Right (LP.Optimal, (optSimplex,_)), Right (LP.Optimal, (optExact,_))) -> TestLP.approx "optimum" 0.001 optSimplex optExact; _ -> QC.property False+prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.exact bounds constrs (dir,obj) of Right (LP.Optimal, (_,sol)) -> TestLP.checkFeasibility 0.1 bounds constrs sol; _ -> QC.property False+prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.exact bounds constrs (dir,obj) of Right (LP.Optimal, (_,sol)) -> QC.forAll (QC.choose (0,1)) $ \lambda -> TestLP.checkFeasibility 0.01 bounds constrs $ TestLP.affineCombination lambda sol (Array.map fromIntegral origin); _ -> QC.property False+prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.exact bounds constrs (dir,obj) of Right (LP.Optimal, (opt,sol)) -> QC.forAll (QC.choose (0,1)) $ \lambda -> let val = TestLP.scalarProduct obj $ TestLP.affineCombination lambda sol (Array.map fromIntegral origin) in (case dir of LP.Minimize -> opt-0.01 <= val; LP.Maximize -> opt+0.01 >= val); _ -> QC.property False -} exact ::    (Shape.Indexed sh, Shape.Index sh ~ ix) =>    Bounds ix -> Constraints ix ->-   (Direction, Objective sh) -> Solution sh+   (Direction, Objective sh) -> Result sh exact = solve (flip FFI.glp_exact nullPtr)  @@ -285,70 +131,25 @@    (Shape.Indexed sh, Shape.Index sh ~ ix) =>    (Foreign.Ptr FFI.Problem -> IO FFI.GlpkSimplexStatus) ->    Bounds ix -> Constraints ix ->-   (Direction, Objective sh) -> Solution sh+   (Direction, Objective sh) -> Result sh solve solver bounds constrs (dir,obj) = unsafePerformIO $    bracket FFI.glp_create_prob FFI.glp_delete_prob $ \lp -> do    storeProblem bounds constrs (dir,obj) lp    void $ solver lp    peekSimplexSolution (Array.shape obj) lp -{-# INLINE solveSuccessiveInPlace #-}-solveSuccessiveInPlace ::-   (Traversable f, Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix) =>-   (Foreign.Ptr FFI.Problem -> IO FFI.GlpkSimplexStatus) ->-   Bounds ix -> Constraints ix ->-   (Direction, Objective sh) ->-   f ((SolutionType, (Double, Array sh Double)) -> Constraints ix,-      (Direction, Objective sh)) ->-   Either NoSolutionType-      (NonEmpty.T f (SolutionType, (Double, Array sh Double)))-solveSuccessiveInPlace solver bounds constrs0 (dir0,obj0) objs =-      unsafePerformIO $-   bracket FFI.glp_create_prob FFI.glp_delete_prob $ \lp -> ME.runExceptT $ do-   let shape = Array.shape obj0-   sol0 <- ME.ExceptT $ do-      storeProblem bounds constrs0 (dir0,obj0) lp-      void $ solver lp-      peekSimplexSolution shape lp-   NonEmpty.cons sol0 <$>-      MS.evalStateT-         (for objs $-            \(makeNewConstrs,(dir,obj)) -> MS.StateT $ \sol ->-                  fmap (\sol1 -> (sol1, sol1)) $ ME.ExceptT $ do-               setDirection lp dir-               when (shape /= Array.shape obj) $-                  error "GLPK.solveSuccessiveInplace: objective shapes mismatch"-               setObjective lp obj-               let newConstrs = makeNewConstrs sol-               newRow <- FFI.glp_add_rows lp $ fromIntegral $ length newConstrs-               for_ (zip [newRow..] (map prepareBounds newConstrs)) $-                     \(row, (terms,(bnd,lo,up))) -> do-                  FFI.glp_set_row_bnds lp row bnd lo up-                  let numTerms = length terms-                  allocaArray numTerms $ \indicesPtr ->-                     allocaArray numTerms $ \coeffsPtr -> do-                     for_ (zip [1..] terms) $-                        \(k, Term c x) -> do-                           pokeElem indicesPtr k (columnIndex shape x)-                           pokeElem coeffsPtr k (realToFrac c)-                     FFI.glp_set_mat_row lp row-                        (fromIntegral numTerms) indicesPtr coeffsPtr-               void $ solver lp-               peekSimplexSolution shape lp)-         sol0  - {- | >>> case Shape.indexTupleFromShape tripletShape of (x1,x2,x3) -> mapSnd (mapPair (round3, Array.toTuple . Array.map round3)) <$> LP.interior [] [[2.*x1, 1.*x2] <=. 10, [1.*x2, 5.*x3] <=. 20] (LP.Maximize, Array.fromTuple (4,-3,2) :: Array.Array TripletShape Double) Right (Optimal,(28.0,(5.0,0.0,4.0))) -prop> QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (LP.simplex bounds constrs (dir,obj), LP.interior bounds constrs (dir,obj)) of (Right (LP.Optimal, (optSimplex,_)), Right (LP.Optimal, (optExact,_))) -> TestLP.approx "optimum" 0.001 optSimplex optExact; _ -> QC.property False+prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (LP.simplex bounds constrs (dir,obj), LP.interior bounds constrs (dir,obj)) of (Right (LP.Optimal, (optSimplex,_)), Right (LP.Optimal, (optExact,_))) -> TestLP.approx "optimum" 0.001 optSimplex optExact; _ -> QC.property False -} interior ::    (Shape.Indexed sh, Shape.Index sh ~ ix) =>    Bounds ix -> Constraints ix ->-   (Direction, Objective sh) -> Solution sh+   (Direction, Objective sh) -> Result sh interior bounds constrs (dir,obj) = unsafePerformIO $    bracket FFI.glp_create_prob FFI.glp_delete_prob $ \lp -> do    storeProblem bounds constrs (dir,obj) lp@@ -369,7 +170,7 @@    Bounds ix -> Constraints ix ->    (Direction, Objective sh) -> Foreign.Ptr FFI.Problem -> IO () storeProblem bounds constrs (dir,obj) lp = do-   void $ FFI.glp_term_out FFI.glpkOff+   Debug.initLog    let shape = Array.shape obj    setDirection lp dir    firstRow <- FFI.glp_add_rows lp $ fromIntegral $ length constrs@@ -385,82 +186,8 @@       for_ (zip [1..] $ concat $             zipWith (map . (,)) [firstRow..] $             map (fst . prepareBounds) constrs) $-         \(k, (row, Term c x)) -> do+         \(k, (row, LP.Term c x)) -> do             pokeElem ia k row             pokeElem ja k (columnIndex shape x)             pokeElem ar k (realToFrac c)       FFI.glp_load_matrix lp (fromIntegral numTerms) ia ja ar-----class FormatIdentifier ix where-   formatIdentifier :: ix -> String--instance FormatIdentifier Char where-   formatIdentifier x = [x]--instance FormatIdentifier c => FormatIdentifier [c] where-   formatIdentifier = concatMap formatIdentifier--instance FormatIdentifier Int where-   formatIdentifier = printf "x%d"--instance FormatIdentifier Integer where-   formatIdentifier = printf "x%d"---formatBound :: (FormatIdentifier ix) => Inequality ix -> String-formatBound (Inequality ix bnd) =-   printf "var %s%s;" (formatIdentifier ix) $-   case bnd of-      LessEqual up -> printf ", <=%f" up-      GreaterEqual lo -> printf ", >=%f" lo-      Between lo up -> printf ", >=%f, <=%f" lo up-      Equal x -> printf ", =%f" x-      Free -> ""---formatSum :: (FormatIdentifier ix) => [Term ix] -> String-formatSum [] = "0"-formatSum xs =-   let formatTerm (Term c ix) = printf "%f*%s" c (formatIdentifier ix) in-   List.intercalate "+" $ map formatTerm xs--formatConstraint :: (FormatIdentifier ix) => Inequality [Term ix] -> String-formatConstraint (Inequality terms bnd) =-   let sumStr = formatSum terms in-   case bnd of-      LessEqual up -> printf "%s <= %f" sumStr up-      GreaterEqual lo -> printf "%f <= %s" lo sumStr-      Between lo up -> printf "%f <= %s <= %f" lo sumStr up-      Equal x -> printf "%s = %f" sumStr x-      Free -> sumStr--formatDirection :: Direction -> String-formatDirection Minimize = "minimize"-formatDirection Maximize = "maximize"--formatObjective ::-   (Shape.Indexed sh, Shape.Index sh ~ ix, FormatIdentifier ix) =>-   Objective sh -> String-formatObjective =-   formatSum . map (\(ix,c) -> Term c ix) . Array.toAssociations--formatMathProg ::-   (Shape.Indexed sh, Shape.Index sh ~ ix, FormatIdentifier ix) =>-   Bounds ix -> Constraints ix ->-   (Direction, Objective sh) -> [String]-formatMathProg bounds constrs (dir,obj) =-   map formatBound bounds ++-   "" :-   formatDirection dir :-   printf "value: %s;" (formatObjective obj) :-   "" :-   "subject to" :-   zipWith-      (\k constr -> printf "constr%d: %s;" k $ formatConstraint constr)-      [(0::Int)..] constrs ++-   "" :-   "end;" :-   []

@@ -4,6 +4,9 @@ {- | The monadic interface to GLPK allows to optimize with respect to multiple objectives, successively.++It is not intended to provide a general imperative interface to GLPK+that manipulates GLPK's state by setters and getters. -} module Numeric.GLPK.Monad (    T,@@ -14,11 +17,13 @@    ) where  import qualified Math.Programming.Glpk.Header as FFI+import qualified Numeric.GLPK.Debug as Debug import Numeric.GLPK.Private-         (Term(Term), Constraints, Solution,+         (Constraints, Result,           allocaArray, pokeElem, columnIndex, prepareBounds, storeBounds,           setDirection, setObjective, peekSimplexSolution) import Numeric.GLPK (Bounds, Direction(..), Objective)+import Numeric.LinearProgramming.Common (Term(Term))  import qualified Data.Array.Comfort.Storable as Array import qualified Data.Array.Comfort.Shape as Shape@@ -35,14 +40,82 @@   {- $setup+>>> :set -XTypeFamilies+>>> :set -XTypeOperators+>>> import qualified Numeric.LinearProgramming.Monad as LPMonad+>>> import qualified Numeric.LinearProgramming.Test as TestLP >>> import qualified Numeric.GLPK.Monad as LP+>>> import qualified Numeric.GLPK as GLPK+>>> import Test.Numeric.GLPK.Utility+>>>    (traverseLag, traverse_Lag, approxSuccession) >>> import Test.Numeric.GLPK (TripletShape, tripletShape)->>> import Numeric.GLPK ((.*), (<=.))+>>> import Numeric.GLPK (Bounds, Constraints, Objective, (.*), (<=.)) >>> >>> import qualified Data.Array.Comfort.Storable as Array >>> import qualified Data.Array.Comfort.Shape as Shape+>>> import qualified Data.NonEmpty as NonEmpty+>>> import Data.Array.Comfort.Storable (Array)+>>> import Data.Traversable (Traversable)+>>> import Data.Foldable (Foldable) >>>+>>> import qualified Control.Monad.Trans.Except as ME+>>> >>> import Data.Tuple.HT (mapSnd)+>>>+>>> import qualified Test.QuickCheck as QC+>>>+>>>+>>> runSuccessive ::+>>>    (Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix, Foldable t) =>+>>>    sh ->+>>>    Bounds ix ->+>>>    (Constraints ix, (GLPK.Direction, Objective sh)) ->+>>>    t (Double -> Constraints ix, (GLPK.Direction, Objective sh)) ->+>>>    Either GLPK.FailureType ()+>>> runSuccessive shape bounds (constrs,dirObj) objs =+>>>    LP.run shape bounds $ ME.runExceptT $ do+>>>       (_solType, (opt, _xs)) <- ME.ExceptT $ LP.simplex constrs dirObj+>>>       traverse_Lag opt+>>>          (\prevResult (newConstr, dirObjI) -> do+>>>              (_solType, (optI, _xs)) <-+>>>                 ME.ExceptT $ LP.simplex (newConstr prevResult) dirObjI+>>>              return optI)+>>>          objs+>>>+>>> solveSuccessiveWarm ::+>>>    (Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix, Traversable t) =>+>>>    sh ->+>>>    Bounds ix ->+>>>    (Constraints ix, (GLPK.Direction, Objective sh)) ->+>>>    t (Double -> Constraints ix, (GLPK.Direction, Objective sh)) ->+>>>    Either GLPK.FailureType+>>>       (NonEmpty.T t (GLPK.SolutionType, (Double, Array sh Double)))+>>> solveSuccessiveWarm shape bounds (constrs,dirObj) objs =+>>>    LP.run shape bounds $ ME.runExceptT $ do+>>>       result <- ME.ExceptT $ LP.simplex constrs dirObj+>>>       NonEmpty.Cons result <$>+>>>          traverseLag result+>>>             (\(_solType, (prevOpt, _xs)) (newConstr, dirObjI) ->+>>>                 ME.ExceptT $ LP.simplex (newConstr prevOpt) dirObjI)+>>>             objs+>>>+>>> solveSuccessiveGen ::+>>>    (Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix, Traversable t) =>+>>>    sh ->+>>>    Bounds ix ->+>>>    (Constraints ix, (GLPK.Direction, Objective sh)) ->+>>>    t (Double -> Constraints ix, (GLPK.Direction, Objective sh)) ->+>>>    Either GLPK.FailureType+>>>       (NonEmpty.T t (GLPK.SolutionType, (Double, Array sh Double)))+>>> solveSuccessiveGen shape bounds (constrs,dirObj) objs =+>>>    LPMonad.run shape bounds $ ME.runExceptT $ do+>>>       result <- ME.ExceptT $ LPMonad.lift GLPK.simplex constrs dirObj+>>>       NonEmpty.Cons result <$>+>>>          traverseLag result+>>>             (\(_solType, (prevOpt, _xs)) (newConstr, dirObjI) ->+>>>                 ME.ExceptT $+>>>                    LPMonad.lift GLPK.simplex (newConstr prevOpt) dirObjI)+>>>             objs -}  @@ -55,7 +128,7 @@ run shape bounds (Cons act) =    unsafePerformIO $    bracket FFI.glp_create_prob FFI.glp_delete_prob $ \lp -> do-      void $ FFI.glp_term_out FFI.glpkOff+      Debug.initLog       storeBounds lp shape bounds       liftIO $ MR.runReaderT act (shape, lp) @@ -66,17 +139,23 @@  >>> case Shape.indexTupleFromShape tripletShape of (x,y,z) -> mapSnd (mapSnd Array.toTuple) <$> LP.run tripletShape [] (LP.simplex [[2.*x, 1.*y] <=. 10, [1.*y, (5::Double).*z] <=. 20] (LP.Maximize, Array.fromTuple (4,-3,2) :: Array.Array TripletShape Double)) Right (Optimal,(28.0,(5.0,0.0,4.0)))++prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (GLPK.simplex bounds constrs (dir,obj), LP.run (Array.shape origin) bounds $ LP.simplex constrs (dir,obj)) of (Right (GLPK.Optimal, (optA,_)), Right (GLPK.Optimal, (optB,_))) -> TestLP.approxReal 0.1 optA optB; _ -> False++prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> TestLP.forAllObjectives origin $ \objs_ -> case TestLP.successiveObjectives origin 0.01 objs_ of (dirObj, objs) -> either (\msg -> QC.counterexample (show msg) False) (const $ QC.property True) $ runSuccessive (Array.shape origin) bounds (constrs,dirObj) objs++prop> TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> TestLP.forAllObjectives origin $ \objs_ -> case TestLP.successiveObjectives origin 0.01 objs_ of (dirObj, objs) -> approxSuccession 0.01 (solveSuccessiveWarm (Array.shape origin) bounds (constrs,dirObj) objs) (solveSuccessiveGen (Array.shape origin) bounds (constrs,dirObj) objs) -} simplex ::    (Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix) =>    Constraints ix ->-   (Direction, Objective sh) -> T sh (Solution sh)+   (Direction, Objective sh) -> T sh (Result sh) simplex = solve (flip FFI.glp_simplex nullPtr)  exact ::    (Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix) =>    Constraints ix ->-   (Direction, Objective sh) -> T sh (Solution sh)+   (Direction, Objective sh) -> T sh (Result sh) exact = solve (flip FFI.glp_exact nullPtr)  solve ::@@ -84,7 +163,7 @@    (Ptr FFI.Problem -> IO FFI.GlpkSimplexStatus) ->    Constraints ix ->    (Direction, Objective sh) ->-   T sh (Solution sh)+   T sh (Result sh) solve method constrs (dir,obj) = Cons $ do    (shape, lp) <- MR.ask    when (shape /= Array.shape obj) $

@@ -2,6 +2,11 @@ {-# LANGUAGE TypeOperators #-} module Numeric.GLPK.Private where +import qualified Numeric.LinearProgramming.Common as LP+import Numeric.LinearProgramming.Common+         (Bound(..), Inequality(Inequality),+          Bounds, Direction(..), Objective)+ import qualified Math.Programming.Glpk.Header as FFI  import qualified Data.Array.Comfort.Storable.Mutable as Mutable@@ -18,36 +23,11 @@ import Foreign.C.Types (CDouble)  --data Term ix = Term Double ix-   deriving (Show)---data Inequality x = Inequality x Bound-   deriving Show--data Bound =-     LessEqual Double-   | GreaterEqual Double-   | Between Double Double-   | Equal Double-   | Free-   deriving Show--instance Functor Inequality where-   fmap f (Inequality x bnd)  =  Inequality (f x) bnd--type Bounds ix = [Inequality ix]--type Constraints ix = [Inequality [Term ix]]--data Direction = Minimize | Maximize-   deriving (Eq, Enum, Bounded, Show)--type Objective sh = Array sh Double+type Term = LP.Term Double +type Constraints ix = LP.Constraints Double ix -data NoSolutionType =+data FailureType =      Undefined    | NoFeasible    | Unbounded@@ -59,7 +39,7 @@    | Optimal    deriving (Eq, Show) -instance NFData NoSolutionType where+instance NFData FailureType where     rnf NoFeasible = ()     rnf _ = () @@ -67,8 +47,8 @@     rnf Optimal = ()     rnf _ = () -type Solution sh =-      Either NoSolutionType (SolutionType, (Double, Array sh Double))+type Result sh =+      Either FailureType (SolutionType, (Double, Array sh Double))   {- |@@ -151,7 +131,7 @@    for_ (Shape.indices defShape) (act arr)    Array.reshape shape <$> Mutable.freeze arr -analyzeStatus :: FFI.GlpkSolutionStatus -> Either NoSolutionType SolutionType+analyzeStatus :: FFI.GlpkSolutionStatus -> Either FailureType SolutionType analyzeStatus status =    fromMaybe (error "glpk-solver: solution type unknown") $ lookup status $       (FFI.glpkUndefined,  Left Undefined) :@@ -164,7 +144,7 @@   peekSimplexSolution ::-   (Shape.C sh) => sh -> Foreign.Ptr FFI.Problem -> IO (Solution sh)+   (Shape.C sh) => sh -> Foreign.Ptr FFI.Problem -> IO (Result sh) peekSimplexSolution shape lp = do    let examine =          liftA2 (,)

@@ -0,0 +1,8 @@+module Numeric.GLPK.Debug where++import qualified Math.Programming.Glpk.Header as FFI+import Control.Monad (void)+++initLog :: IO ()+initLog = void $ FFI.glp_term_out FFI.glpkOff

@@ -0,0 +1,5 @@+module Numeric.GLPK.Debug where+++initLog :: IO ()+initLog = return ()

@@ -1,13 +1,13 @@ -- Do not edit! Automatically created with doctest-extract from src/Numeric/GLPK.hs-{-# LINE 82 "src/Numeric/GLPK.hs" #-}+{-# LINE 67 "src/Numeric/GLPK.hs" #-}  module Test.Numeric.GLPK where  import Test.DocTest.Base import qualified Test.DocTest.Driver as DocTest -{-# LINE 83 "src/Numeric/GLPK.hs" #-}-import     qualified Test.Numeric.GLPK.Generator as TestLP+{-# LINE 68 "src/Numeric/GLPK.hs" #-}+import     qualified Numeric.LinearProgramming.Test as TestLP import     qualified Numeric.GLPK as LP import     Numeric.GLPK ((.*), (<=.), (==.)) @@ -34,71 +34,76 @@  test :: DocTest.T () test = do- DocTest.printPrefix "Numeric.GLPK:141: "-{-# LINE 141 "src/Numeric/GLPK.hs" #-}+ DocTest.printPrefix "Numeric.GLPK:99: "+{-# LINE 99 "src/Numeric/GLPK.hs" #-}  DocTest.property-{-# LINE 141 "src/Numeric/GLPK.hs" #-}+{-# LINE 99 "src/Numeric/GLPK.hs" #-}      (\target -> case Shape.indexTupleFromShape pairShape of (pos,neg) -> case mapSnd (mapSnd Array.toTuple) <$> LP.simplex [] [[1.*pos, (-1).*neg] ==. target] (LP.Minimize, Array.fromTuple (1,1) :: Array.Array PairShape Double) of (Right (LP.Optimal,(absol,(posResult,negResult)))) -> QC.property (TestLP.approxReal 0.001 absol (abs target)) .&&. (posResult === 0 .||. negResult === 0); _ -> QC.property False)- DocTest.printPrefix "Numeric.GLPK:142: "-{-# LINE 142 "src/Numeric/GLPK.hs" #-}+ DocTest.printPrefix "Numeric.GLPK:100: "+{-# LINE 100 "src/Numeric/GLPK.hs" #-}  DocTest.property-{-# LINE 142 "src/Numeric/GLPK.hs" #-}+{-# LINE 100 "src/Numeric/GLPK.hs" #-}      (\(QC.Positive posWeight) (QC.Positive negWeight) target -> case Shape.indexTupleFromShape pairShape of (pos,neg) -> case mapSnd (mapSnd Array.toTuple) <$> LP.simplex [] [[1.*pos, (-1).*neg] ==. target] (LP.Minimize, Array.fromTuple (posWeight,negWeight) :: Array.Array PairShape Double) of (Right (LP.Optimal,(absol,(posResult,negResult)))) -> QC.property (absol>=0) .&&. (posResult === 0 .||. negResult === 0); _ -> QC.property False)- DocTest.printPrefix "Numeric.GLPK:143: "-{-# LINE 143 "src/Numeric/GLPK.hs" #-}+ DocTest.printPrefix "Numeric.GLPK:101: "+{-# LINE 101 "src/Numeric/GLPK.hs" #-}  DocTest.property-{-# LINE 143 "src/Numeric/GLPK.hs" #-}-     (QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, _) -> True; _ -> False)- DocTest.printPrefix "Numeric.GLPK:138: "-{-# LINE 138 "src/Numeric/GLPK.hs" #-}+{-# LINE 101 "src/Numeric/GLPK.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, _) -> True; _ -> False)+ DocTest.printPrefix "Numeric.GLPK:102: "+{-# LINE 102 "src/Numeric/GLPK.hs" #-}+ DocTest.property+{-# LINE 102 "src/Numeric/GLPK.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, (_,sol)) -> TestLP.checkFeasibility 0.1 bounds constrs sol; _ -> QC.property False)+ DocTest.printPrefix "Numeric.GLPK:103: "+{-# LINE 103 "src/Numeric/GLPK.hs" #-}+ DocTest.property+{-# LINE 103 "src/Numeric/GLPK.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, (_,sol)) -> QC.forAll (QC.choose (0,1)) $ \lambda -> TestLP.checkFeasibility 0.1 bounds constrs $ TestLP.affineCombination lambda sol (Array.map fromIntegral origin); _ -> QC.property False)+ DocTest.printPrefix "Numeric.GLPK:104: "+{-# LINE 104 "src/Numeric/GLPK.hs" #-}+ DocTest.property+{-# LINE 104 "src/Numeric/GLPK.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.simplex bounds constrs (dir,obj) of Right (LP.Optimal, (opt,sol)) -> QC.forAll (QC.choose (0,1)) $ \lambda -> let val = TestLP.scalarProduct obj $ TestLP.affineCombination lambda sol (Array.map fromIntegral origin) in (case dir of LP.Minimize -> opt-0.01 <= val; LP.Maximize -> opt+0.01 >= val); _ -> QC.property False)+ DocTest.printPrefix "Numeric.GLPK:96: "+{-# LINE 96 "src/Numeric/GLPK.hs" #-}  DocTest.example-{-# LINE 138 "src/Numeric/GLPK.hs" #-}+{-# LINE 96 "src/Numeric/GLPK.hs" #-}    (case Shape.indexTupleFromShape tripletShape of (x1,x2,x3) -> mapSnd (mapSnd Array.toTuple) <$> LP.simplex [] [[2.*x1, 1.*x2] <=. 10, [1.*x2, 5.*x3] <=. 20] (LP.Maximize, Array.fromTuple (4,-3,2) :: Array.Array TripletShape Double))   [ExpectedLine [LineChunk "Right (Optimal,(28.0,(5.0,0.0,4.0)))"]]- DocTest.printPrefix "Numeric.GLPK:180: "-{-# LINE 180 "src/Numeric/GLPK.hs" #-}- DocTest.property-{-# LINE 180 "src/Numeric/GLPK.hs" #-}-     (QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAllShrink (TestLP.genProblem origin) TestLP.shrinkProblem $ \(bounds, constrs) -> QC.forAllShrink (TestLP.genObjectives origin) TestLP.shrinkObjectives $ \objs -> case LP.simplexMulti bounds constrs (Array.shape origin) objs of (_, Right (LP.Optimal, _)) -> QC.property True; result -> QC.counterexample (show result) False)- DocTest.printPrefix "Numeric.GLPK:219: "-{-# LINE 219 "src/Numeric/GLPK.hs" #-}- DocTest.property-{-# LINE 219 "src/Numeric/GLPK.hs" #-}-     (QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAllShrink (TestLP.genProblem origin) TestLP.shrinkProblem $ \(bounds, constrs) -> QC.forAllShrink (TestLP.genObjectives origin) TestLP.shrinkObjectives $ \objs -> case uncurry (LP.simplexSuccessive bounds constrs) $ TestLP.successiveObjectives origin 0.01 objs of result -> QC.counterexample (show result) $ case result of Right results -> all (\r -> case r of (LP.Optimal, _) -> True; _ -> False) results; _ -> False)- DocTest.printPrefix "Numeric.GLPK:220: "-{-# LINE 220 "src/Numeric/GLPK.hs" #-}+ DocTest.printPrefix "Numeric.GLPK:117: "+{-# LINE 117 "src/Numeric/GLPK.hs" #-}  DocTest.property-{-# LINE 220 "src/Numeric/GLPK.hs" #-}-     (QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAllShrink (TestLP.genProblem origin) TestLP.shrinkProblem $ \(bounds, constrs) -> QC.forAllShrink (TestLP.genObjectives origin) TestLP.shrinkObjectives $ \objs -> case uncurry (LP.exactSuccessive bounds constrs) $ TestLP.successiveObjectives origin 0.01 objs of result -> QC.counterexample (show result) $ case result of Right results -> all (\r -> case r of (LP.Optimal, _) -> True; _ -> False) results; _ -> False)- DocTest.printPrefix "Numeric.GLPK:240: "-{-# LINE 240 "src/Numeric/GLPK.hs" #-}+{-# LINE 117 "src/Numeric/GLPK.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (LP.simplex bounds constrs (dir,obj), LP.exact bounds constrs (dir,obj)) of (Right (LP.Optimal, (optSimplex,_)), Right (LP.Optimal, (optExact,_))) -> TestLP.approx "optimum" 0.001 optSimplex optExact; _ -> QC.property False)+ DocTest.printPrefix "Numeric.GLPK:118: "+{-# LINE 118 "src/Numeric/GLPK.hs" #-}  DocTest.property-{-# LINE 240 "src/Numeric/GLPK.hs" #-}-     (QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjectives origin) $ (. TestLP.successiveObjectives origin 0.01) $ \(obj,objs) -> case (LP.simplexSuccessive bounds constrs obj objs, LP.solveSuccessive (LP.simplex bounds) constrs obj objs) of (resultA,resultB) -> TestLP.approxSuccession 0.01 resultA resultB)- DocTest.printPrefix "Numeric.GLPK:241: "-{-# LINE 241 "src/Numeric/GLPK.hs" #-}+{-# LINE 118 "src/Numeric/GLPK.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.exact bounds constrs (dir,obj) of Right (LP.Optimal, (_,sol)) -> TestLP.checkFeasibility 0.1 bounds constrs sol; _ -> QC.property False)+ DocTest.printPrefix "Numeric.GLPK:119: "+{-# LINE 119 "src/Numeric/GLPK.hs" #-}  DocTest.property-{-# LINE 241 "src/Numeric/GLPK.hs" #-}-     (QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjectives origin) $ (. TestLP.successiveObjectives origin 0.01) $ \(obj,objs) -> case (LP.exactSuccessive bounds constrs obj objs, LP.solveSuccessive (LP.exact bounds) constrs obj objs) of (resultA,resultB) -> TestLP.approxSuccession 0.01 resultA resultB)- DocTest.printPrefix "Numeric.GLPK:274: "-{-# LINE 274 "src/Numeric/GLPK.hs" #-}+{-# LINE 119 "src/Numeric/GLPK.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.exact bounds constrs (dir,obj) of Right (LP.Optimal, (_,sol)) -> QC.forAll (QC.choose (0,1)) $ \lambda -> TestLP.checkFeasibility 0.01 bounds constrs $ TestLP.affineCombination lambda sol (Array.map fromIntegral origin); _ -> QC.property False)+ DocTest.printPrefix "Numeric.GLPK:120: "+{-# LINE 120 "src/Numeric/GLPK.hs" #-}  DocTest.property-{-# LINE 274 "src/Numeric/GLPK.hs" #-}-     (QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (LP.simplex bounds constrs (dir,obj), LP.exact bounds constrs (dir,obj)) of (Right (LP.Optimal, (optSimplex,_)), Right (LP.Optimal, (optExact,_))) -> TestLP.approx "optimum" 0.001 optSimplex optExact; _ -> QC.property False)- DocTest.printPrefix "Numeric.GLPK:271: "-{-# LINE 271 "src/Numeric/GLPK.hs" #-}+{-# LINE 120 "src/Numeric/GLPK.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case LP.exact bounds constrs (dir,obj) of Right (LP.Optimal, (opt,sol)) -> QC.forAll (QC.choose (0,1)) $ \lambda -> let val = TestLP.scalarProduct obj $ TestLP.affineCombination lambda sol (Array.map fromIntegral origin) in (case dir of LP.Minimize -> opt-0.01 <= val; LP.Maximize -> opt+0.01 >= val); _ -> QC.property False)+ DocTest.printPrefix "Numeric.GLPK:114: "+{-# LINE 114 "src/Numeric/GLPK.hs" #-}  DocTest.example-{-# LINE 271 "src/Numeric/GLPK.hs" #-}+{-# LINE 114 "src/Numeric/GLPK.hs" #-}    (case Shape.indexTupleFromShape tripletShape of (x1,x2,x3) -> mapSnd (mapSnd Array.toTuple) <$> LP.exact [] [[2.*x1, 1.*x2] <=. 10, [1.*x2, 5.*x3] <=. 20] (LP.Maximize, Array.fromTuple (4,-3,2) :: Array.Array TripletShape Double))   [ExpectedLine [LineChunk "Right (Optimal,(28.0,(5.0,0.0,4.0)))"]]- DocTest.printPrefix "Numeric.GLPK:346: "-{-# LINE 346 "src/Numeric/GLPK.hs" #-}+ DocTest.printPrefix "Numeric.GLPK:147: "+{-# LINE 147 "src/Numeric/GLPK.hs" #-}  DocTest.property-{-# LINE 346 "src/Numeric/GLPK.hs" #-}-     (QC.forAllShrink TestLP.genOrigin TestLP.shrinkOrigin $ \origin -> QC.forAll (TestLP.genProblem origin) $ \(bounds, constrs) -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (LP.simplex bounds constrs (dir,obj), LP.interior bounds constrs (dir,obj)) of (Right (LP.Optimal, (optSimplex,_)), Right (LP.Optimal, (optExact,_))) -> TestLP.approx "optimum" 0.001 optSimplex optExact; _ -> QC.property False)- DocTest.printPrefix "Numeric.GLPK:343: "-{-# LINE 343 "src/Numeric/GLPK.hs" #-}+{-# LINE 147 "src/Numeric/GLPK.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (LP.simplex bounds constrs (dir,obj), LP.interior bounds constrs (dir,obj)) of (Right (LP.Optimal, (optSimplex,_)), Right (LP.Optimal, (optExact,_))) -> TestLP.approx "optimum" 0.001 optSimplex optExact; _ -> QC.property False)+ DocTest.printPrefix "Numeric.GLPK:144: "+{-# LINE 144 "src/Numeric/GLPK.hs" #-}  DocTest.example-{-# LINE 343 "src/Numeric/GLPK.hs" #-}+{-# LINE 144 "src/Numeric/GLPK.hs" #-}    (case Shape.indexTupleFromShape tripletShape of (x1,x2,x3) -> mapSnd (mapPair (round3, Array.toTuple . Array.map round3)) <$> LP.interior [] [[2.*x1, 1.*x2] <=. 10, [1.*x2, 5.*x3] <=. 20] (LP.Maximize, Array.fromTuple (4,-3,2) :: Array.Array TripletShape Double))   [ExpectedLine [LineChunk "Right (Optimal,(28.0,(5.0,0.0,4.0)))"]]

@@ -1,208 +0,0 @@-{-# LANGUAGE TypeFamilies #-}-module Test.Numeric.GLPK.Generator where--import qualified Numeric.GLPK as LP-import Numeric.GLPK ((<=.), (>=.))--import qualified Test.QuickCheck as QC-import System.Random (Random)--import qualified Data.Array.Comfort.Boxed as BoxedArray-import qualified Data.Array.Comfort.Storable as Array-import qualified Data.Array.Comfort.Shape as Shape-import qualified Data.NonEmpty as NonEmpty-import qualified Data.List.HT as ListHT-import qualified Data.Ix as Ix-import Data.Array.Comfort.Storable (Array, (!))-import Data.Traversable (sequenceA, for)-import Data.Tuple.HT (mapSnd)-import Data.Maybe (fromMaybe)-import Data.Int (Int64)--import Control.Applicative (liftA2)--import Text.Printf (PrintfArg, printf)--import Foreign.Storable (Storable)---{- |-Generate constraints in the form of a polyhedron-which contains warrantedly the zero vector.-That is, there is an admissible solution.-In order to assert that the polyhedron is closed,-we bound all variables by a hypercube.--}-genProblem ::-   (Shape.Indexed sh, Shape.Index sh ~ ix, Element a) =>-   Array sh a -> QC.Gen (LP.Bounds ix, LP.Constraints ix)-genProblem origin =-   liftA2 (,)-      (for (Array.toAssociations origin) $ \(ix,x) ->-         LP.Inequality ix <$>-         liftA2 LP.Between-            (doubleFromElement . (x+) <$> QC.choose (-100,-50))-            (doubleFromElement . (x+) <$> QC.choose (50,100)))-      (do-         numConstraints <- QC.choose (1,20)-         QC.vectorOf numConstraints $ do-            ixs <- QC.sublistOf $ Shape.indices $ Array.shape origin-            terms <- for ixs $ \ix -> do-               coeff <- QC.choose (-10,10)-               return (coeff, ix)-            let offset = scalarProduct terms origin-            let deviation = 25-            LP.Inequality-               (map (uncurry (LP.Term . doubleFromElement)) terms)-               <$>-               QC.oneof (-                  (do bound <- QC.choose (offset-deviation, offset+deviation)-                      return $-                         if bound > offset-                            then LP.LessEqual    $ doubleFromElement bound-                            else LP.GreaterEqual $ doubleFromElement bound) :-                  (liftA2 LP.Between-                     (doubleFromElement <$>-                        QC.choose (offset-deviation, offset))-                     (doubleFromElement <$>-                        QC.choose (offset, offset+deviation))) :-                  []))--scalarProduct ::-   (Shape.Indexed sh, Shape.Index sh ~ ix, Storable a, Num a) =>-   [(a,ix)] -> Array sh a -> a-scalarProduct terms origin =-   sum $ map (\(coeff, ix) -> coeff * origin!ix) terms--genVarShape :: QC.Gen (Shape.Range Char)-genVarShape = Shape.Range 'a' <$> QC.choose ('a','j')--genOrigin :: QC.Gen (Array (Shape.Range Char) Int64)-genOrigin = genVector =<< genVarShape--_genOrigin :: QC.Gen (Array (Shape.Range Char) Double)-_genOrigin = genVector =<< genVarShape---shrinkVarShape :: Shape.Range Char -> [Shape.Range Char]-shrinkVarShape (Shape.Range from to) =-   if from<to then [Shape.Range from (pred to)] else []--shrinkOrigin ::-   (Storable a) => Array (Shape.Range Char) a -> [Array (Shape.Range Char) a]-shrinkOrigin vec =-   case Array.shape vec of-      Shape.Range from to ->-         if from<to-            then [Array.sample (Shape.Range from (pred to)) (vec!)]-            else []---class (Storable a, Random a, Num a, Ord a) => Element a where-   doubleFromElement :: a -> Double--instance Element Double where-   doubleFromElement = id--instance Element Int64 where-   doubleFromElement = fromIntegral--genObjective ::-   (Shape.Indexed sh, Shape.Index sh ~ ix, Element a) =>-   Array sh a -> QC.Gen (LP.Direction, LP.Objective sh)-genObjective origin =-   liftA2 (,) QC.arbitraryBoundedEnum-      (fmap (Array.map doubleFromElement . flip asTypeOf origin) $-       genVector $ Array.shape origin)--genVector :: (Shape.Indexed sh, Element a) => sh -> QC.Gen (Array sh a)-genVector shape =-   fmap Array.fromBoxed $ sequenceA $-   BoxedArray.fromAssociations (QC.choose (-10,10)) shape []--shrinkProblem ::-   (LP.Bounds ix, LP.Constraints ix) ->-   [(LP.Bounds ix, LP.Constraints ix)]-shrinkProblem (bounds, constraints) =-   map (\shrinked -> (bounds, shrinked)) $-   filter (not . null) $ QC.shrinkList (const []) constraints--genObjectives ::-   (Shape.Indexed sh, Shape.Index sh ~ ix, Element a) =>-   Array sh a -> QC.Gen (NonEmpty.T [] (LP.Direction, [LP.Term ix]))-genObjectives origin = do-   let shape = Array.shape origin-   let stageRange :: (Int,Int)-       stageRange = (0,3)-   stages <- for (Shape.indices shape) $ \ix -> (,) ix <$> QC.choose stageRange-   let varSets =-         fromMaybe (error "there should be at least one stage") $-         NonEmpty.fetch $-         filter (not . null) $-         map (\k -> map fst $ filter ((k==) . snd) stages) $-         Ix.range stageRange-   let asTypeOfElement :: a -> f a -> a-       asTypeOfElement = const-   for varSets $ \varSet ->-      liftA2 (,)-         QC.arbitraryBoundedEnum-         (for varSet $ \ix ->-            flip LP.Term ix . doubleFromElement-               <$> QC.choose (-10, 10 `asTypeOfElement` origin))--shrinkObjectives ::-   NonEmpty.T [] (LP.Direction, [LP.Term ix]) ->-   [NonEmpty.T [] (LP.Direction, [LP.Term ix])]-shrinkObjectives (NonEmpty.Cons obj objs) =-   map (NonEmpty.Cons obj) $-   QC.shrinkList-      (\(dir,terms) ->-         map ((,) dir) $ filter (not . null) $-         QC.shrinkList (const []) terms)-      objs--successiveObjectives ::-   (Shape.Indexed sh, Shape.Index sh ~ ix) =>-   Array sh a -> Double ->-   NonEmpty.T [] (LP.Direction, [LP.Term ix]) ->-   ((LP.Direction, LP.Objective sh),-    [((LP.SolutionType, (Double, Array sh Double)) -> LP.Constraints ix,-      (LP.Direction, LP.Objective sh))])-successiveObjectives origin tol xs =-   let shape = Array.shape origin in-   (mapSnd (LP.objectiveFromTerms shape) $ NonEmpty.head xs,-    NonEmpty.mapAdjacent-      (\(dir0,obj0) y1 ->-         (\(_sol,(opt,_vec)) ->-            case dir0 of-               LP.Minimize -> [obj0 <=. opt + tol]-               LP.Maximize -> [obj0 >=. opt - tol],-          mapSnd (LP.objectiveFromTerms shape) y1))-      xs)---approxReal :: (Ord a, Num a) => a -> a -> a -> Bool-approxReal tol x y = abs (x-y) <= tol--approx :: (PrintfArg a, Ord a, Num a) => String -> a -> a -> a -> QC.Property-approx name tol x y =-   QC.counterexample (printf "%s: %f - %f" name x y) (approxReal tol x y)--approxSuccession ::-   (Shape.C sh, Show sh, Show a, Ord a, Num a, Storable a) =>-   a ->-   Either LP.NoSolutionType-      (NonEmpty.T [] (LP.SolutionType, (a, Array sh a))) ->-   Either LP.NoSolutionType-      (NonEmpty.T [] (LP.SolutionType, (a, Array sh a))) ->-   QC.Property-approxSuccession tol x y =-   QC.counterexample (show x) $-   QC.counterexample (show y) $-   case (x,y) of-      (Left sx, Left sy) -> sx==sy-      (Right (NonEmpty.Cons xh xs), Right (NonEmpty.Cons yh ys)) ->-         let equalSol (solX, (optX, _)) (solY, (optY, _)) =-               solX == solY && approxReal tol optX optY-         in equalSol xh yh  &&  ListHT.equalWith equalSol xs ys-      _ -> False

@@ -1,26 +1,109 @@ -- Do not edit! Automatically created with doctest-extract from src/Numeric/GLPK/Monad.hs-{-# LINE 37 "src/Numeric/GLPK/Monad.hs" #-}+{-# LINE 42 "src/Numeric/GLPK/Monad.hs" #-} +{-# OPTIONS_GHC -XTypeFamilies #-}+{-# OPTIONS_GHC -XTypeOperators #-} module Test.Numeric.GLPK.Monad where  import Test.DocTest.Base import qualified Test.DocTest.Driver as DocTest -{-# LINE 38 "src/Numeric/GLPK/Monad.hs" #-}+{-# LINE 45 "src/Numeric/GLPK/Monad.hs" #-}+import     qualified Numeric.LinearProgramming.Monad as LPMonad+import     qualified Numeric.LinearProgramming.Test as TestLP import     qualified Numeric.GLPK.Monad as LP+import     qualified Numeric.GLPK as GLPK+import     Test.Numeric.GLPK.Utility+       (traverseLag, traverse_Lag, approxSuccession) import     Test.Numeric.GLPK (TripletShape, tripletShape)-import     Numeric.GLPK ((.*), (<=.))+import     Numeric.GLPK (Bounds, Constraints, Objective, (.*), (<=.))  import     qualified Data.Array.Comfort.Storable as Array import     qualified Data.Array.Comfort.Shape as Shape+import     qualified Data.NonEmpty as NonEmpty+import     Data.Array.Comfort.Storable (Array)+import     Data.Traversable (Traversable)+import     Data.Foldable (Foldable) +import     qualified Control.Monad.Trans.Except as ME+ import     Data.Tuple.HT (mapSnd) +import     qualified Test.QuickCheck as QC+++runSuccessive     ::+       (Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix, Foldable t) =>+       sh ->+       Bounds ix ->+       (Constraints ix, (GLPK.Direction, Objective sh)) ->+       t (Double -> Constraints ix, (GLPK.Direction, Objective sh)) ->+       Either GLPK.FailureType ()+runSuccessive     shape bounds (constrs,dirObj) objs =+       LP.run shape bounds $ ME.runExceptT $ do+          (_solType, (opt, _xs)) <- ME.ExceptT $ LP.simplex constrs dirObj+          traverse_Lag opt+             (\prevResult (newConstr, dirObjI) -> do+                 (_solType, (optI, _xs)) <-+                    ME.ExceptT $ LP.simplex (newConstr prevResult) dirObjI+                 return optI)+             objs++solveSuccessiveWarm     ::+       (Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix, Traversable t) =>+       sh ->+       Bounds ix ->+       (Constraints ix, (GLPK.Direction, Objective sh)) ->+       t (Double -> Constraints ix, (GLPK.Direction, Objective sh)) ->+       Either GLPK.FailureType+          (NonEmpty.T t (GLPK.SolutionType, (Double, Array sh Double)))+solveSuccessiveWarm     shape bounds (constrs,dirObj) objs =+       LP.run shape bounds $ ME.runExceptT $ do+          result <- ME.ExceptT $ LP.simplex constrs dirObj+          NonEmpty.Cons result <$>+             traverseLag result+                (\(_solType, (prevOpt, _xs)) (newConstr, dirObjI) ->+                    ME.ExceptT $ LP.simplex (newConstr prevOpt) dirObjI)+                objs++solveSuccessiveGen     ::+       (Eq sh, Shape.Indexed sh, Shape.Index sh ~ ix, Traversable t) =>+       sh ->+       Bounds ix ->+       (Constraints ix, (GLPK.Direction, Objective sh)) ->+       t (Double -> Constraints ix, (GLPK.Direction, Objective sh)) ->+       Either GLPK.FailureType+          (NonEmpty.T t (GLPK.SolutionType, (Double, Array sh Double)))+solveSuccessiveGen     shape bounds (constrs,dirObj) objs =+       LPMonad.run shape bounds $ ME.runExceptT $ do+          result <- ME.ExceptT $ LPMonad.lift GLPK.simplex constrs dirObj+          NonEmpty.Cons result <$>+             traverseLag result+                (\(_solType, (prevOpt, _xs)) (newConstr, dirObjI) ->+                    ME.ExceptT $+                       LPMonad.lift GLPK.simplex (newConstr prevOpt) dirObjI)+                objs+ test :: DocTest.T () test = do- DocTest.printPrefix "Numeric.GLPK.Monad:67: "-{-# LINE 67 "src/Numeric/GLPK/Monad.hs" #-}+ DocTest.printPrefix "Numeric.GLPK.Monad:143: "+{-# LINE 143 "src/Numeric/GLPK/Monad.hs" #-}+ DocTest.property+{-# LINE 143 "src/Numeric/GLPK/Monad.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> QC.forAll (TestLP.genObjective origin) $ \(dir,obj) -> case (GLPK.simplex bounds constrs (dir,obj), LP.run (Array.shape origin) bounds $ LP.simplex constrs (dir,obj)) of (Right (GLPK.Optimal, (optA,_)), Right (GLPK.Optimal, (optB,_))) -> TestLP.approxReal 0.1 optA optB; _ -> False)+ DocTest.printPrefix "Numeric.GLPK.Monad:145: "+{-# LINE 145 "src/Numeric/GLPK/Monad.hs" #-}+ DocTest.property+{-# LINE 145 "src/Numeric/GLPK/Monad.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> TestLP.forAllObjectives origin $ \objs_ -> case TestLP.successiveObjectives origin 0.01 objs_ of (dirObj, objs) -> either (\msg -> QC.counterexample (show msg) False) (const $ QC.property True) $ runSuccessive (Array.shape origin) bounds (constrs,dirObj) objs)+ DocTest.printPrefix "Numeric.GLPK.Monad:147: "+{-# LINE 147 "src/Numeric/GLPK/Monad.hs" #-}+ DocTest.property+{-# LINE 147 "src/Numeric/GLPK/Monad.hs" #-}+     (TestLP.forAllOrigin $ \origin -> TestLP.forAllProblem origin $ \bounds constrs -> TestLP.forAllObjectives origin $ \objs_ -> case TestLP.successiveObjectives origin 0.01 objs_ of (dirObj, objs) -> approxSuccession 0.01 (solveSuccessiveWarm (Array.shape origin) bounds (constrs,dirObj) objs) (solveSuccessiveGen (Array.shape origin) bounds (constrs,dirObj) objs))+ DocTest.printPrefix "Numeric.GLPK.Monad:140: "+{-# LINE 140 "src/Numeric/GLPK/Monad.hs" #-}  DocTest.example-{-# LINE 67 "src/Numeric/GLPK/Monad.hs" #-}+{-# LINE 140 "src/Numeric/GLPK/Monad.hs" #-}    (case Shape.indexTupleFromShape tripletShape of (x,y,z) -> mapSnd (mapSnd Array.toTuple) <$> LP.run tripletShape [] (LP.simplex [[2.*x, 1.*y] <=. 10, [1.*y, (5::Double).*z] <=. 20] (LP.Maximize, Array.fromTuple (4,-3,2) :: Array.Array TripletShape Double)))   [ExpectedLine [LineChunk "Right (Optimal,(28.0,(5.0,0.0,4.0)))"]]

@@ -0,0 +1,55 @@+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeOperators #-}+module Test.Numeric.GLPK.Utility where++import qualified Numeric.GLPK as LP+import Numeric.LinearProgramming.Test (approxReal)++import qualified Test.QuickCheck as QC++import qualified Data.Array.Comfort.Shape as Shape+import qualified Data.NonEmpty as NonEmpty+import qualified Data.List.HT as ListHT+import Data.Array.Comfort.Storable (Array)+import Data.Tuple.HT (double)+import Data.Traversable (Traversable, traverse)+import Data.Foldable (Foldable, traverse_)++import qualified Control.Monad.Trans.State as MS++import Foreign.Storable (Storable)+++approxSuccession ::+   (Shape.C sh, Show sh, Show a, Ord a, Num a, Storable a) =>+   a ->+   Either LP.FailureType+      (NonEmpty.T [] (LP.SolutionType, (a, Array sh a))) ->+   Either LP.FailureType+      (NonEmpty.T [] (LP.SolutionType, (a, Array sh a))) ->+   QC.Property+approxSuccession tol x y =+   QC.counterexample (show x) $+   QC.counterexample (show y) $+   case (x,y) of+      (Left sx, Left sy) -> sx==sy+      (Right (NonEmpty.Cons xh xs), Right (NonEmpty.Cons yh ys)) ->+         let equalSol (solX, (optX, _)) (solY, (optY, _)) =+               solX == solY && approxReal tol optX optY+         in equalSol xh yh  &&  ListHT.equalWith equalSol xs ys+      _ -> False+++traverse_Lag ::+   (Foldable t, Monad m) =>+   b -> (b -> a -> m b) -> t a -> m ()+traverse_Lag b0 f =+   flip MS.evalStateT b0 .+   traverse_ (\a -> MS.StateT $ \b -> fmap double $ f b a)++traverseLag ::+   (Traversable t, Monad m) =>+   b -> (b -> a -> m b) -> t a -> m (t b)+traverseLag b0 f =+   flip MS.evalStateT b0 .+   traverse (\a -> MS.StateT $ \b -> fmap double $ f b a)