Response to https://www.reddit.com/r/haskellquestions/comments/a8f5kv/combinations_that_satisfy_some_property_a_faster/

bundle.cabal

name:                bundle
version:             0.1.0.0
synopsis:            Solve bundling problem
-- description:
homepage:            https://github.com/minkless/sudoku#readme
license:             BSD3
license-file:        LICENSE
author:              Koji Miyazato
maintainer:          [email protected]
copyright:           Koji Miyazato
category:            Game
build-type:          Simple
cabal-version:       >=1.10
extra-source-files:  README.md

executable bundle
  hs-source-dirs:      src
  main-is:             bundle.hs
  -- other-extensions:    
  build-depends:       base         >=4.11 && <5,
                       containers,
                       glpk-hs
  ghc-options:         -Wall -threaded -O2
  default-language:    Haskell2010

bundle.hs

{-# LANGUAGE TupleSections #-}
module Main(main) where

import qualified Data.LinearProgram as LP

import           Data.Foldable      (foldl')

import           System.IO

import           Data.IntMap        (IntMap)
import qualified Data.IntMap        as IntMap
import qualified Data.Map.Strict    as Map
import           Data.Maybe         (fromMaybe)

--------------------------------------------------------------
-- Solving as MIP

data Var = X Int Int | P Int
    deriving (Show, Read, Eq, Ord)

type Bag = IntMap Int
type Input = Bag
type Answer = (Int, [Bag])
type Problem = LP.LP Var Int

bagFromList :: [Int] -> Bag
bagFromList = IntMap.fromListWith (+) . map (,1)

readInput :: String -> Input
readInput = bagFromList . validateInput . map read . words

qualityMax :: Int
qualityMax = 20

valueThreshold :: Int
valueThreshold = 40

validateInput :: [Int] -> [Int]
validateInput xs =
  if all (\x -> 1 <= x && x <= qualityMax) xs
    then xs
    else error "Out-of-bounds input"

prettyBags :: [Bag] -> String
prettyBags = unlines . map prettyBag

prettyBag :: Bag -> String
prettyBag bag = value ++ contents
  where
    value = "[value=" ++ show (bagTotal bag) ++ "]"
    contents = show bag

prettyAnswer :: Answer -> String
prettyAnswer (numValuableBags, bags) =
  "**Number of valuable bags: " ++ show numValuableBags ++ "**\n" ++
  prettyBags bags

printAnswer :: Maybe Answer -> IO ()
printAnswer = putStrLn . maybe "!NO_SOLUTION" prettyAnswer

bagTotal :: Bag -> Int
bagTotal = sum . map (uncurry (*)) . IntMap.toList

makeProblem :: Input -> Problem
makeProblem input =
  let kmax = bagTotal input `div` valueThreshold

      qs = [q | q <- [1..qualityMax], IntMap.member q input ]
      ks = [1..kmax]

      xs = X <$> qs <*> ks
      ps = P <$> ks

      countOf q = fromMaybe 0 $ IntMap.lookup q input

      goal = Map.fromList (zip ps (repeat 1))

      mkCountConstr q =
        let name = "Count(" ++ show q ++ ")"
            fun = Map.fromList [ (X q k, 1) | k <- ks ]
            bound = LP.UBound (countOf q)
        in LP.Constr (Just name) fun bound

      mkQualityConstr k =
        let name = "Quality(" ++ show k ++ ")"
            fun = Map.fromList $
              [ (P k, negate valueThreshold) ] ++
              [ (X q k, q) | q <- qs ]
            bound = LP.LBound 0
        in LP.Constr (Just name) fun bound

      constrs = (mkCountConstr <$> qs) ++ (mkQualityConstr <$> ks)

      bounds = Map.fromList $
        [ (X q k, LP.Bound 0 (countOf q)) | q <- qs, k <- ks ] ++
        [ (P k, LP.Free) | k <- ks ]

      varTypes = Map.fromList $
        zip xs (repeat LP.IntVar) ++
        zip ps (repeat LP.BinVar)
  in LP.LP LP.Max goal constrs bounds varTypes

solveProblem :: Problem -> IO (Maybe Answer)
solveProblem prob =
  do (retcode, maySolution) <- LP.glpSolveVars mipOption prob
     case retcode of
       LP.Success -> case maySolution of
         Nothing       -> fail "Succeeded but no solution???"
         Just (a, sol) -> return . Just $ (floor a, fromSolution sol)
       _ -> do hPrint stderr retcode
               return Nothing
  where
    mipOption = LP.mipDefaults{ LP.msgLev = LP.MsgErr }

fromSolution :: Map.Map Var Double -> [Bag]
fromSolution = IntMap.elems . foldl' step IntMap.empty . Map.toList
  where
    step accum (x, n) = case x of
      X q k | n > 0 ->
        let newElem = IntMap.singleton q (floor n)
        in IntMap.insertWith IntMap.union k newElem accum
      _ -> accum

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

test :: String -> IO ()
test probStr =
  do let input = readInput probStr
     putStrLn "-- Problem -----------"
     putStrLn (prettyBag input)
     putStrLn "-- Answer  -----------"
     maybeAns <- solveProblem (makeProblem input)
     printAnswer maybeAns

-- | Cheap PRNG
randomSeq :: Int -> [Int]
randomSeq n = map (\x -> x `mod` 20 + 1) $ drop 50 $ go n
  where
    go x = x `seq` x : go ((x * 10007 + 1) `mod` 999983)

-- | Trivial
input0 :: String
input0 = "1 2 3 4 16 17 18 19"

-- | Small (random)
input1 :: String
input1 = unwords . map show . take 20 $ randomSeq 1

-- | Medium (random)
input2 :: String
input2 = unwords . map show . take 40 $ randomSeq 2

-- | Big (random)
input3 :: String
input3 = unwords . map show . take 60 $ randomSeq 3

-- | Hard
--
-- This input have smaller total than @input3@, but solver runs
-- much slower on it.
--
-- (IDK the cause, but increasing the number of 19s beyond 46
-- made solving it abysmally slow.)
input4 :: String
input4 = unwords . map show . take 46 $ repeat (19 :: Int)

main :: IO ()
main = do
  test input0
  test input1
  test input2
  test input3
  test input4

stack.yaml

resolver: lts-12.21
packages:
- '.'

# For both glpk-hs and gasp, the versions on Hackage
# do not compiles with ghc-8.4
# Latest version on GitHub does work.
extra-deps:
- git: https://github.com/jyp/glpk-hs
  commit: e32f85cc7911c66fd3cd7bf70294ffb45d163824
- git: https://github.com/jyp/gasp
  commit: e63a9c1d91f59a2379cddf084655d0d6973f78d9

flags: {}

extra-package-dbs: []