msakai / toysolver / 496

Committed 10 Nov 2024 11:05AM UTC coverage: 69.994% (-1.1%) from 71.113%

Build # 496

Build Type

push

github

Committed by

web-flow

Commit Message

Merge pull request #117 from msakai/update-coveralls-and-haddock

GitHub Actions: Update coveralls and haddock configuration

Run Details

9872 of 14104 relevant lines covered (69.99%)

0.7 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

67.25

/src/ToySolver/SAT/Solver/SLS/ProbSAT.hs

{-# LANGUAGE CPP #-}
{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE DeriveDataTypeable #-}
{-# LANGUAGE ScopedTypeVariables #-}
{-# OPTIONS_GHC -Wall #-}
{-# OPTIONS_HADDOCK show-extensions #-}
----------------------------------------------------------------------
-- |
-- Module      :  ToySolver.SAT.Solver.SLS.ProbSAT
-- Copyright   :  (c) Masahiro Sakai 2017
-- License     :  BSD-style
--
-- Maintainer  :  masahiro.sakai@gmail.com
-- Stability   :  provisional
-- Portability :  non-portable
--
-- References:
--
----------------------------------------------------------------------
module ToySolver.SAT.Solver.SLS.ProbSAT
  ( Solver
  , newSolver
  , newSolverWeighted
  , getNumVars
  , getRandomGen
  , setRandomGen
  , getBestSolution
  , getStatistics

  , Options (..)
  , Callbacks (..)
  , Statistics (..)

  , generateUniformRandomSolution

  , probsat
  , walksat
  ) where

import Prelude hiding (break)

import Control.Exception
import Control.Loop
import Control.Monad
import Control.Monad.Primitive
import Control.Monad.Trans
import Control.Monad.Trans.Except
import Data.Array.Base (unsafeRead, unsafeWrite, unsafeAt)
import Data.Array.IArray
import Data.Array.IO
import Data.Array.Unboxed
import Data.Array.Unsafe
import Data.Bits
import Data.Default.Class
import qualified Data.Foldable as F
import Data.Int
import Data.IORef
import Data.Maybe
import Data.Sequence ((|>))
import qualified Data.Sequence as Seq
import Data.Typeable
import Data.Word
import System.Clock
import qualified System.Random.MWC as Rand
import qualified System.Random.MWC.Distributions as Rand
import qualified ToySolver.FileFormat.CNF as CNF
import ToySolver.Internal.Data.IOURef
import qualified ToySolver.Internal.Data.Vec as Vec
import qualified ToySolver.SAT.Types as SAT

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

data Solver
  = Solver
  { svClauses                :: !(Array ClauseId PackedClause)
  , svClauseWeights          :: !(Array ClauseId CNF.Weight)
  , svClauseWeightsF         :: !(UArray ClauseId Double)
  , svClauseNumTrueLits      :: !(IOUArray ClauseId Int32)
  , svClauseUnsatClauseIndex :: !(IOUArray ClauseId Int)
  , svUnsatClauses           :: !(Vec.UVec ClauseId)

  , svVarOccurs         :: !(Array SAT.Var (UArray Int ClauseId))
  , svVarOccursState    :: !(Array SAT.Var (IOUArray Int Bool))
  , svSolution          :: !(IOUArray SAT.Var Bool)

  , svObj               :: !(IORef CNF.Weight)

  , svRandomGen         :: !(IORef Rand.GenIO)
  , svBestSolution      :: !(IORef (CNF.Weight, SAT.Model))
  , svStatistics        :: !(IORef Statistics)
  }

type ClauseId = Int

type PackedClause = Array Int SAT.Lit

newSolver :: CNF.CNF -> IO Solver
newSolver cnf = do
  let wcnf =
        CNF.WCNF
        { CNF.wcnfNumVars    = CNF.cnfNumVars cnf
        , CNF.wcnfNumClauses = CNF.cnfNumClauses cnf
        , CNF.wcnfTopCost    = fromIntegral (CNF.cnfNumClauses cnf) + 1
        , CNF.wcnfClauses    = [(1,c) | c <- CNF.cnfClauses cnf]
        }
  newSolverWeighted wcnf

newSolverWeighted :: CNF.WCNF -> IO Solver
newSolverWeighted wcnf = do
  let m :: SAT.Var -> Bool
      m _ = False
      nv = CNF.wcnfNumVars wcnf

  objRef <- newIORef (0::Integer)

  cs <- liftM catMaybes $ forM (CNF.wcnfClauses wcnf) $ \(w,pc) -> do
    case SAT.normalizeClause (SAT.unpackClause pc) of
      Nothing -> return Nothing
      Just [] -> modifyIORef' objRef (w+) >> return Nothing
      Just c  -> do
        let c' = listArray (0, length c - 1) c
        seq c' $ return (Just (w,c'))
  let len = length cs
      clauses  = listArray (0, len - 1) (map snd cs)
      weights  :: Array ClauseId CNF.Weight
      weights  = listArray (0, len - 1) (map fst cs)
      weightsF :: UArray ClauseId Double
      weightsF = listArray (0, len - 1) (map (fromIntegral . fst) cs)

  (varOccurs' :: IOArray SAT.Var (Seq.Seq (Int, Bool))) <- newArray (1, nv) Seq.empty

  clauseNumTrueLits <- newArray (bounds clauses) 0
  clauseUnsatClauseIndex <- newArray (bounds clauses) (-1)
  unsatClauses <- Vec.new

  forAssocsM_ clauses $ \(c,clause) -> do
    let n = sum [1 | lit <- elems clause, SAT.evalLit m lit]
    writeArray clauseNumTrueLits c n
    when (n == 0) $ do
      i <- Vec.getSize unsatClauses
      writeArray clauseUnsatClauseIndex c i
      Vec.push unsatClauses c
      modifyIORef objRef ((weights ! c) +)
    forM_ (elems clause) $ \lit -> do
      let v = SAT.litVar lit
      let b = SAT.evalLit m lit
      seq b $ modifyArray varOccurs' v (|> (c,b))

  varOccurs <- do
    (arr::IOArray SAT.Var (UArray Int ClauseId)) <- newArray_ (1, nv)
    forM_ [1 .. nv] $ \v -> do
      s <- readArray varOccurs' v
      writeArray arr v $ listArray (0, Seq.length s - 1) (map fst (F.toList s))
    unsafeFreeze arr

  varOccursState <- do
    (arr::IOArray SAT.Var (IOUArray Int Bool)) <- newArray_ (1, nv)
    forM_ [1 .. nv] $ \v -> do
      s <- readArray varOccurs' v
      ss <- newArray_ (0, Seq.length s - 1)
      forM_ (zip [0..] (F.toList s)) $ \(j,a) -> writeArray ss j (snd a)
      writeArray arr v ss
    unsafeFreeze arr

  solution <- newListArray (1, nv) $ [SAT.evalVar m v | v <- [1..nv]]

  bestObj <- readIORef objRef
  bestSol <- freeze solution
  bestSolution <- newIORef (bestObj, bestSol)

  randGen <- newIORef =<< Rand.create

  stat <- newIORef def

  return $
    Solver
    { svClauses = clauses
    , svClauseWeights          = weights
    , svClauseWeightsF         = weightsF
    , svClauseNumTrueLits      = clauseNumTrueLits
    , svClauseUnsatClauseIndex = clauseUnsatClauseIndex
    , svUnsatClauses           = unsatClauses

    , svVarOccurs         = varOccurs
    , svVarOccursState    = varOccursState
    , svSolution          = solution

    , svObj = objRef

    , svRandomGen         = randGen
    , svBestSolution      = bestSolution
    , svStatistics        = stat
    }


flipVar :: Solver -> SAT.Var -> IO ()
flipVar solver v = mask_ $ do
  let occurs = svVarOccurs solver ! v
      occursState = svVarOccursState solver ! v
  seq occurs $ seq occursState $ return ()
  modifyArray (svSolution solver) v not
  forAssocsM_ occurs $ \(j,!c) -> do
    b <- unsafeRead occursState j
    n <- unsafeRead (svClauseNumTrueLits solver) c
    unsafeWrite occursState j (not b)
    if b then do
      unsafeWrite (svClauseNumTrueLits solver) c (n-1)
      when (n==1) $ do
        i <- Vec.getSize (svUnsatClauses solver)
        Vec.push (svUnsatClauses solver) c
        unsafeWrite (svClauseUnsatClauseIndex solver) c i
        modifyIORef' (svObj solver) (+ unsafeAt (svClauseWeights solver) c)
    else do
      unsafeWrite (svClauseNumTrueLits solver) c (n+1)
      when (n==0) $ do
        s <- Vec.getSize (svUnsatClauses solver)
        i <- unsafeRead (svClauseUnsatClauseIndex solver) c
        unless (i == s-1) $ do
          let i2 = s-1
          c2 <- Vec.unsafeRead (svUnsatClauses solver) i2
          Vec.unsafeWrite (svUnsatClauses solver) i2 c
          Vec.unsafeWrite (svUnsatClauses solver) i c2
          unsafeWrite (svClauseUnsatClauseIndex solver) c2 i
        _ <- Vec.unsafePop (svUnsatClauses solver)
        modifyIORef' (svObj solver) (subtract (unsafeAt (svClauseWeights solver) c))
        return ()

setSolution :: SAT.IModel m => Solver -> m -> IO ()
setSolution solver m = do
  b <- getBounds (svSolution solver)
  forM_ (range b) $ \v -> do
    val <- readArray (svSolution solver) v
    let val' = SAT.evalVar m v
    unless (val == val') $ do
      flipVar solver v

getNumVars :: Solver -> IO Int
getNumVars solver = return $ rangeSize $ bounds (svVarOccurs solver)

getRandomGen :: Solver -> IO Rand.GenIO
getRandomGen solver = readIORef (svRandomGen solver)

setRandomGen :: Solver -> Rand.GenIO -> IO ()
setRandomGen solver gen = writeIORef (svRandomGen solver) gen

getBestSolution :: Solver -> IO (CNF.Weight, SAT.Model)
getBestSolution solver = readIORef (svBestSolution solver)

getStatistics :: Solver -> IO Statistics
getStatistics solver = readIORef (svStatistics solver)

{-# INLINE getMakeValue #-}
getMakeValue :: Solver -> SAT.Var -> IO Double
getMakeValue solver v = do
  let occurs = svVarOccurs solver ! v
      (lb,ub) = bounds occurs
  seq occurs $ seq lb $ seq ub $
    numLoopState lb ub 0 $ \ !r !i -> do
      let c = unsafeAt occurs i
      n <- unsafeRead (svClauseNumTrueLits solver) c
      return $! if n == 0 then (r + unsafeAt (svClauseWeightsF solver) c) else r

{-# INLINE getBreakValue #-}
getBreakValue :: Solver -> SAT.Var -> IO Double
getBreakValue solver v = do
  let occurs = svVarOccurs solver ! v
      occursState = svVarOccursState solver ! v
      (lb,ub) = bounds occurs
  seq occurs $ seq occursState $ seq lb $ seq ub $
    numLoopState lb ub 0 $ \ !r !i -> do
      b <- unsafeRead occursState i
      if b then do
        let c = unsafeAt occurs i
        n <- unsafeRead (svClauseNumTrueLits solver) c
        return $! if n==1 then (r + unsafeAt (svClauseWeightsF solver) c) else r
      else
        return r

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

data Options
  = Options
  { optTarget   :: !CNF.Weight
  , optMaxTries :: !Int
  , optMaxFlips :: !Int
  , optPickClauseWeighted :: Bool
  }
  deriving (Eq, Show)

instance Default Options where
  def =
    Options
    { optTarget   = 0
    , optMaxTries = 1
    , optMaxFlips = 100000
    , optPickClauseWeighted = False
    }

data Callbacks
  = Callbacks
  { cbGenerateInitialSolution :: Solver -> IO SAT.Model
  , cbOnUpdateBestSolution :: Solver -> CNF.Weight -> SAT.Model -> IO ()
  }

instance Default Callbacks where
  def =
    Callbacks
    { cbGenerateInitialSolution = generateUniformRandomSolution
    , cbOnUpdateBestSolution = \_ _ _ -> return ()
    }

data Statistics
  = Statistics
  { statTotalCPUTime   :: !TimeSpec
  , statFlips          :: !Int
  , statFlipsPerSecond :: !Double
  }
  deriving (Eq, Show)

instance Default Statistics where
  def =
    Statistics
    { statTotalCPUTime = 0
    , statFlips = 0
    , statFlipsPerSecond = 0
    }

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

generateUniformRandomSolution :: Solver -> IO SAT.Model
generateUniformRandomSolution solver = do
  gen <- getRandomGen solver
  n <- getNumVars solver
  (a :: IOUArray Int Bool) <- newArray_ (1,n)
  forM_ [1..n] $ \v -> do
    b <- Rand.uniform gen
    writeArray a v b
  unsafeFreeze a

checkCurrentSolution :: Solver -> Callbacks -> IO ()
checkCurrentSolution solver cb = do
  best <- readIORef (svBestSolution solver)
  obj <- readIORef (svObj solver)
  when (obj < fst best) $ do
    sol <- freeze (svSolution solver)
    writeIORef (svBestSolution solver) (obj, sol)
    cbOnUpdateBestSolution cb solver obj sol

pickClause :: Solver -> Options -> IO PackedClause
pickClause solver opt = do
  gen <- getRandomGen solver
  if optPickClauseWeighted opt then do
    obj <- readIORef (svObj solver)
    let f !j !x = do
          c <- Vec.read (svUnsatClauses solver) j
          let w = svClauseWeights solver ! c
          if x < w then
            return c
          else
            f (j + 1) (x - w)
    x <- rand obj gen
    c <- f 0 x
    return $ (svClauses solver ! c)
  else do
    s <- Vec.getSize (svUnsatClauses solver)
    j <- Rand.uniformR (0, s - 1) gen -- For integral types inclusive range is used
    liftM (svClauses solver !) $ Vec.read (svUnsatClauses solver) j

rand :: PrimMonad m => Integer -> Rand.Gen (PrimState m) -> m Integer
rand n gen
  | n <= toInteger (maxBound :: Word32) = liftM toInteger $ Rand.uniformR (0, fromIntegral n - 1 :: Word32) gen
  | otherwise = do
      a <- rand (n `shiftR` 32) gen
      (b::Word32) <- Rand.uniform gen
      return $ (a `shiftL` 32) .|. toInteger b

data Finished = Finished
  deriving (Show, Typeable)

instance Exception Finished

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

probsat :: Solver -> Options -> Callbacks -> (Double -> Double -> Double) -> IO ()
probsat solver opt cb f = do
  gen <- getRandomGen solver
  let maxClauseLen =
        if rangeSize (bounds (svClauses solver)) == 0
        then 0
        else maximum $ map (rangeSize . bounds) $ elems (svClauses solver)
  (wbuf :: IOUArray Int Double) <- newArray_ (0, maxClauseLen-1)
  wsumRef <- newIOURef (0 :: Double)

  let pickVar :: PackedClause -> IO SAT.Var
      pickVar c = do
        writeIOURef wsumRef 0
        forAssocsM_ c $ \(k,lit) -> do
          let v = SAT.litVar lit
          m <- getMakeValue solver v
          b <- getBreakValue solver v
          let w = f m b
          writeArray wbuf k w
          modifyIOURef wsumRef (+w)
        wsum <- readIOURef wsumRef

        let go :: Int -> Double -> IO Int
            go !k !a = do
              if not (inRange (bounds c) k) then do
                return $! snd (bounds c)
              else do
                w <- readArray wbuf k
                if a <= w then
                  return k
                else
                  go (k + 1) (a - w)
        k <- go 0 =<< Rand.uniformR (0, wsum) gen
        return $! SAT.litVar (c ! k)

  startCPUTime <- getTime ProcessCPUTime
  flipsRef <- newIOURef (0::Int)

  -- It's faster to use Control.Exception than using Control.Monad.Except
  let body = do
        replicateM_ (optMaxTries opt) $ do
          sol <- cbGenerateInitialSolution cb solver
          setSolution solver sol
          checkCurrentSolution solver cb
          replicateM_ (optMaxFlips opt) $ do
            s <- Vec.getSize (svUnsatClauses solver)
            when (s == 0) $ throw Finished
            obj <- readIORef (svObj solver)
            when (obj <= optTarget opt) $ throw Finished
            c <- pickClause solver opt
            v <- pickVar c
            flipVar solver v
            modifyIOURef flipsRef inc
            checkCurrentSolution solver cb
  body `catch` (\(_::Finished) -> return ())

  endCPUTime <- getTime ProcessCPUTime
  flips <- readIOURef flipsRef
  let totalCPUTime = endCPUTime `diffTimeSpec` startCPUTime
      totalCPUTimeSec = fromIntegral (toNanoSecs totalCPUTime) / 10^(9::Int)
  writeIORef (svStatistics solver) $
    Statistics
    { statTotalCPUTime = totalCPUTime
    , statFlips = flips
    , statFlipsPerSecond = fromIntegral flips / totalCPUTimeSec
    }

  return ()



walksat :: Solver -> Options -> Callbacks -> Double -> IO ()
walksat solver opt cb p = do
  gen <- getRandomGen solver
  (buf :: Vec.UVec SAT.Var) <- Vec.new

  let pickVar :: PackedClause -> IO SAT.Var
      pickVar c = do
        Vec.clear buf
        let (lb,ub) = bounds c
        r <- runExceptT $ do
          _ <- numLoopState lb ub (1.0/0.0) $ \ !b0 !i -> do
            let v = SAT.litVar (c ! i)
            b <- lift $ getBreakValue solver v
            if b <= 0 then
              throwE v -- freebie move
            else if b < b0 then do
              lift $ Vec.clear buf >> Vec.push buf v
              return b
            else if b == b0 then do
              lift $ Vec.push buf v
              return b0
            else do
              return b0
          return ()
        case r of
          Left v -> return v
          Right _ -> do
            flag <- Rand.bernoulli p gen
            if flag then do
              -- random walk move
              i <- Rand.uniformR (lb,ub) gen
              return $! SAT.litVar (c ! i)
            else do
              -- greedy move
              s <- Vec.getSize buf
              if s == 1 then
                Vec.unsafeRead buf 0
              else do
                i <- Rand.uniformR (0, s - 1) gen
                Vec.unsafeRead buf i

  startCPUTime <- getTime ProcessCPUTime
  flipsRef <- newIOURef (0::Int)

  -- It's faster to use Control.Exception than using Control.Monad.Except
  let body = do
        replicateM_ (optMaxTries opt) $ do
          sol <- cbGenerateInitialSolution cb solver
          setSolution solver sol
          checkCurrentSolution solver cb
          replicateM_ (optMaxFlips opt) $ do
            s <- Vec.getSize (svUnsatClauses solver)
            when (s == 0) $ throw Finished
            obj <- readIORef (svObj solver)
            when (obj <= optTarget opt) $ throw Finished
            c <- pickClause solver opt
            v <- pickVar c
            flipVar solver v
            modifyIOURef flipsRef inc
            checkCurrentSolution solver cb
  body `catch` (\(_::Finished) -> return ())

  endCPUTime <- getTime ProcessCPUTime
  flips <- readIOURef flipsRef
  let totalCPUTime = endCPUTime `diffTimeSpec` startCPUTime
      totalCPUTimeSec = fromIntegral (toNanoSecs totalCPUTime) / 10^(9::Int)
  writeIORef (svStatistics solver) $
    Statistics
    { statTotalCPUTime = totalCPUTime
    , statFlips = flips
    , statFlipsPerSecond = fromIntegral flips / totalCPUTimeSec
    }

  return ()

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

#if !MIN_VERSION_array(0,5,6)

{-# INLINE modifyArray #-}
modifyArray :: (MArray a e m, Ix i) => a i e -> i -> (e -> e) -> m ()
modifyArray a i f = do
  e <- readArray a i
  writeArray a i (f e)

#endif

{-# INLINE forAssocsM_ #-}
forAssocsM_ :: (IArray a e, Monad m) => a Int e -> ((Int,e) -> m ()) -> m ()
forAssocsM_ a f = do
  let (lb,ub) = bounds a
  numLoop lb ub $ \i ->
    f (i, unsafeAt a i)

{-# INLINE inc #-}
inc :: Integral a => a -> a
inc a = a+1

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

1	{-# LANGUAGE CPP #-}
2	{-# LANGUAGE BangPatterns #-}
3	{-# LANGUAGE DeriveDataTypeable #-}
4	{-# LANGUAGE ScopedTypeVariables #-}
5	{-# OPTIONS_GHC -Wall #-}
6	{-# OPTIONS_HADDOCK show-extensions #-}
7	----------------------------------------------------------------------
8	-- \|
9	-- Module : ToySolver.SAT.Solver.SLS.ProbSAT
10	-- Copyright : (c) Masahiro Sakai 2017
11	-- License : BSD-style
12	--
13	-- Maintainer : masahiro.sakai@gmail.com
14	-- Stability : provisional
15	-- Portability : non-portable
16	--
17	-- References:
18	--
19	----------------------------------------------------------------------
20	module ToySolver.SAT.Solver.SLS.ProbSAT
21	( Solver
22	, newSolver
23	, newSolverWeighted
24	, getNumVars
25	, getRandomGen
26	, setRandomGen
27	, getBestSolution
28	, getStatistics
29
30	, Options (..)
31	, Callbacks (..)
32	, Statistics (..)
33
34	, generateUniformRandomSolution
35
36	, probsat
37	, walksat
38	) where
39
40	import Prelude hiding (break)
41
42	import Control.Exception
43	import Control.Loop
44	import Control.Monad
45	import Control.Monad.Primitive
46	import Control.Monad.Trans
47	import Control.Monad.Trans.Except
48	import Data.Array.Base (unsafeRead, unsafeWrite, unsafeAt)
49	import Data.Array.IArray
50	import Data.Array.IO
51	import Data.Array.Unboxed
52	import Data.Array.Unsafe
53	import Data.Bits
54	import Data.Default.Class
55	import qualified Data.Foldable as F
56	import Data.Int
57	import Data.IORef
58	import Data.Maybe
59	import Data.Sequence ((\|>))
60	import qualified Data.Sequence as Seq
61	import Data.Typeable
62	import Data.Word
63	import System.Clock
64	import qualified System.Random.MWC as Rand
65	import qualified System.Random.MWC.Distributions as Rand
66	import qualified ToySolver.FileFormat.CNF as CNF
67	import ToySolver.Internal.Data.IOURef
68	import qualified ToySolver.Internal.Data.Vec as Vec
69	import qualified ToySolver.SAT.Types as SAT
70
71	-- -------------------------------------------------------------------
72
73	data Solver
74	= Solver
75	{ svClauses :: !(Array ClauseId PackedClause)	1✔
76	, svClauseWeights :: !(Array ClauseId CNF.Weight)	1✔
77	, svClauseWeightsF :: !(UArray ClauseId Double)	1✔
78	, svClauseNumTrueLits :: !(IOUArray ClauseId Int32)	1✔
79	, svClauseUnsatClauseIndex :: !(IOUArray ClauseId Int)	1✔
80	, svUnsatClauses :: !(Vec.UVec ClauseId)	1✔
81
82	, svVarOccurs :: !(Array SAT.Var (UArray Int ClauseId))	1✔
83	, svVarOccursState :: !(Array SAT.Var (IOUArray Int Bool))	1✔
84	, svSolution :: !(IOUArray SAT.Var Bool)	1✔
85
86	, svObj :: !(IORef CNF.Weight)	1✔
87
88	, svRandomGen :: !(IORef Rand.GenIO)	1✔
89	, svBestSolution :: !(IORef (CNF.Weight, SAT.Model))	1✔
90	, svStatistics :: !(IORef Statistics)	1✔
91	}
92
93	type ClauseId = Int
94
95	type PackedClause = Array Int SAT.Lit
96
97	newSolver :: CNF.CNF -> IO Solver
98	newSolver cnf = do	1✔
99	let wcnf =	1✔
100	CNF.WCNF	1✔
101	{ CNF.wcnfNumVars = CNF.cnfNumVars cnf	1✔
102	, CNF.wcnfNumClauses = CNF.cnfNumClauses cnf	1✔
103	, CNF.wcnfTopCost = fromIntegral (CNF.cnfNumClauses cnf) + 1	1✔
104	, CNF.wcnfClauses = [(1,c) \| c <- CNF.cnfClauses cnf]	1✔
105	}
106	newSolverWeighted wcnf	1✔
107
108	newSolverWeighted :: CNF.WCNF -> IO Solver
109	newSolverWeighted wcnf = do	1✔
110	let m :: SAT.Var -> Bool
111	m _ = False	1✔
112	nv = CNF.wcnfNumVars wcnf	1✔
113
114	objRef <- newIORef (0::Integer)	1✔
115
116	cs <- liftM catMaybes $ forM (CNF.wcnfClauses wcnf) $ \(w,pc) -> do	1✔
117	case SAT.normalizeClause (SAT.unpackClause pc) of	1✔
118	Nothing -> return Nothing	1✔
119	Just [] -> modifyIORef' objRef (w+) >> return Nothing	1✔
120	Just c -> do	1✔
121	let c' = listArray (0, length c - 1) c	1✔
122	seq c' $ return (Just (w,c'))	1✔
123	let len = length cs	1✔
124	clauses = listArray (0, len - 1) (map snd cs)	1✔
125	weights :: Array ClauseId CNF.Weight
126	weights = listArray (0, len - 1) (map fst cs)	1✔
127	weightsF :: UArray ClauseId Double
128	weightsF = listArray (0, len - 1) (map (fromIntegral . fst) cs)	1✔
129
130	(varOccurs' :: IOArray SAT.Var (Seq.Seq (Int, Bool))) <- newArray (1, nv) Seq.empty	1✔
131
132	clauseNumTrueLits <- newArray (bounds clauses) 0	1✔
133	clauseUnsatClauseIndex <- newArray (bounds clauses) (-1)	1✔
134	unsatClauses <- Vec.new	1✔
135
136	forAssocsM_ clauses $ \(c,clause) -> do	1✔
137	let n = sum [1 \| lit <- elems clause, SAT.evalLit m lit]	1✔
138	writeArray clauseNumTrueLits c n	1✔
139	when (n == 0) $ do	1✔
140	i <- Vec.getSize unsatClauses	1✔
141	writeArray clauseUnsatClauseIndex c i	1✔
142	Vec.push unsatClauses c	1✔
143	modifyIORef objRef ((weights ! c) +)	1✔
144	forM_ (elems clause) $ \lit -> do	1✔
145	let v = SAT.litVar lit	1✔
146	let b = SAT.evalLit m lit	1✔
147	seq b $ modifyArray varOccurs' v (\|> (c,b))	1✔
148
149	varOccurs <- do	1✔
150	(arr::IOArray SAT.Var (UArray Int ClauseId)) <- newArray_ (1, nv)	1✔
151	forM_ [1 .. nv] $ \v -> do	1✔
152	s <- readArray varOccurs' v	1✔
153	writeArray arr v $ listArray (0, Seq.length s - 1) (map fst (F.toList s))	1✔
154	unsafeFreeze arr	1✔
155
156	varOccursState <- do	1✔
157	(arr::IOArray SAT.Var (IOUArray Int Bool)) <- newArray_ (1, nv)	1✔
158	forM_ [1 .. nv] $ \v -> do	1✔
159	s <- readArray varOccurs' v	1✔
160	ss <- newArray_ (0, Seq.length s - 1)	1✔
161	forM_ (zip [0..] (F.toList s)) $ \(j,a) -> writeArray ss j (snd a)	1✔
162	writeArray arr v ss	1✔
163	unsafeFreeze arr	1✔
164
165	solution <- newListArray (1, nv) $ [SAT.evalVar m v \| v <- [1..nv]]	×
166
167	bestObj <- readIORef objRef	1✔
168	bestSol <- freeze solution	1✔
169	bestSolution <- newIORef (bestObj, bestSol)	1✔
170
171	randGen <- newIORef =<< Rand.create	1✔
172
173	stat <- newIORef def	×
174
175	return $	1✔
176	Solver	1✔
177	{ svClauses = clauses	1✔
178	, svClauseWeights = weights	1✔
179	, svClauseWeightsF = weightsF	1✔
180	, svClauseNumTrueLits = clauseNumTrueLits	1✔
181	, svClauseUnsatClauseIndex = clauseUnsatClauseIndex	1✔
182	, svUnsatClauses = unsatClauses	1✔
183
184	, svVarOccurs = varOccurs	1✔
185	, svVarOccursState = varOccursState	1✔
186	, svSolution = solution	1✔
187
188	, svObj = objRef	1✔
189
190	, svRandomGen = randGen	1✔
191	, svBestSolution = bestSolution	1✔
192	, svStatistics = stat	1✔
193	}
194
195
196	flipVar :: Solver -> SAT.Var -> IO ()
197	flipVar solver v = mask_ $ do	1✔
198	let occurs = svVarOccurs solver ! v	1✔
199	occursState = svVarOccursState solver ! v	1✔
200	seq occurs $ seq occursState $ return ()	×
201	modifyArray (svSolution solver) v not	1✔
202	forAssocsM_ occurs $ \(j,!c) -> do	1✔
203	b <- unsafeRead occursState j	1✔
204	n <- unsafeRead (svClauseNumTrueLits solver) c	1✔
205	unsafeWrite occursState j (not b)	1✔
206	if b then do	1✔
207	unsafeWrite (svClauseNumTrueLits solver) c (n-1)	1✔
208	when (n==1) $ do	1✔
209	i <- Vec.getSize (svUnsatClauses solver)	1✔
210	Vec.push (svUnsatClauses solver) c	1✔
211	unsafeWrite (svClauseUnsatClauseIndex solver) c i	1✔
212	modifyIORef' (svObj solver) (+ unsafeAt (svClauseWeights solver) c)	1✔
213	else do	1✔
214	unsafeWrite (svClauseNumTrueLits solver) c (n+1)	1✔
215	when (n==0) $ do	1✔
216	s <- Vec.getSize (svUnsatClauses solver)	1✔
217	i <- unsafeRead (svClauseUnsatClauseIndex solver) c	1✔
218	unless (i == s-1) $ do	1✔
219	let i2 = s-1	1✔
220	c2 <- Vec.unsafeRead (svUnsatClauses solver) i2	1✔
221	Vec.unsafeWrite (svUnsatClauses solver) i2 c	1✔
222	Vec.unsafeWrite (svUnsatClauses solver) i c2	1✔
223	unsafeWrite (svClauseUnsatClauseIndex solver) c2 i	1✔
224	_ <- Vec.unsafePop (svUnsatClauses solver)	1✔
225	modifyIORef' (svObj solver) (subtract (unsafeAt (svClauseWeights solver) c))	1✔
226	return ()	×
227
228	setSolution :: SAT.IModel m => Solver -> m -> IO ()
229	setSolution solver m = do	1✔
230	b <- getBounds (svSolution solver)	1✔
231	forM_ (range b) $ \v -> do	1✔
232	val <- readArray (svSolution solver) v	1✔
233	let val' = SAT.evalVar m v	1✔
234	unless (val == val') $ do	1✔
235	flipVar solver v	1✔
236
237	getNumVars :: Solver -> IO Int
238	getNumVars solver = return $ rangeSize $ bounds (svVarOccurs solver)	1✔
239
240	getRandomGen :: Solver -> IO Rand.GenIO
241	getRandomGen solver = readIORef (svRandomGen solver)	1✔
242
243	setRandomGen :: Solver -> Rand.GenIO -> IO ()
244	setRandomGen solver gen = writeIORef (svRandomGen solver) gen	×
245
246	getBestSolution :: Solver -> IO (CNF.Weight, SAT.Model)
247	getBestSolution solver = readIORef (svBestSolution solver)	1✔
248
249	getStatistics :: Solver -> IO Statistics
250	getStatistics solver = readIORef (svStatistics solver)	×
251
252	{-# INLINE getMakeValue #-}
253	getMakeValue :: Solver -> SAT.Var -> IO Double
254	getMakeValue solver v = do	1✔
255	let occurs = svVarOccurs solver ! v	1✔
256	(lb,ub) = bounds occurs	1✔
257	seq occurs $ seq lb $ seq ub $	1✔
258	numLoopState lb ub 0 $ \ !r !i -> do	1✔
259	let c = unsafeAt occurs i	1✔
260	n <- unsafeRead (svClauseNumTrueLits solver) c	1✔
261	return $! if n == 0 then (r + unsafeAt (svClauseWeightsF solver) c) else r	1✔
262
263	{-# INLINE getBreakValue #-}
264	getBreakValue :: Solver -> SAT.Var -> IO Double
265	getBreakValue solver v = do	1✔
266	let occurs = svVarOccurs solver ! v	1✔
267	occursState = svVarOccursState solver ! v	1✔
268	(lb,ub) = bounds occurs	1✔
269	seq occurs $ seq occursState $ seq lb $ seq ub $	1✔
270	numLoopState lb ub 0 $ \ !r !i -> do	1✔
271	b <- unsafeRead occursState i	1✔
272	if b then do	1✔
273	let c = unsafeAt occurs i	1✔
274	n <- unsafeRead (svClauseNumTrueLits solver) c	1✔
275	return $! if n==1 then (r + unsafeAt (svClauseWeightsF solver) c) else r	1✔
276	else
277	return r	1✔
278
279	-- -------------------------------------------------------------------
280
281	data Options
282	= Options
283	{ optTarget :: !CNF.Weight	1✔
284	, optMaxTries :: !Int	1✔
285	, optMaxFlips :: !Int	1✔
286	, optPickClauseWeighted :: Bool	1✔
287	}
288	deriving (Eq, Show)	×
289
290	instance Default Options where
291	def =	1✔
292	Options	1✔
293	{ optTarget = 0	1✔
294	, optMaxTries = 1	1✔
295	, optMaxFlips = 100000	1✔
296	, optPickClauseWeighted = False	1✔
297	}
298
299	data Callbacks
300	= Callbacks
301	{ cbGenerateInitialSolution :: Solver -> IO SAT.Model	1✔
302	, cbOnUpdateBestSolution :: Solver -> CNF.Weight -> SAT.Model -> IO ()	1✔
303	}
304
305	instance Default Callbacks where
306	def =	1✔
307	Callbacks	1✔
308	{ cbGenerateInitialSolution = generateUniformRandomSolution	1✔
309	, cbOnUpdateBestSolution = \_ _ _ -> return ()	×
310	}
311
312	data Statistics
313	= Statistics
314	{ statTotalCPUTime :: !TimeSpec	×
315	, statFlips :: !Int	×
316	, statFlipsPerSecond :: !Double	×
317	}
318	deriving (Eq, Show)	×
319
320	instance Default Statistics where
321	def =	×
322	Statistics	×
323	{ statTotalCPUTime = 0	×
324	, statFlips = 0	×
325	, statFlipsPerSecond = 0	×
326	}
327
328	-- -------------------------------------------------------------------
329
330	generateUniformRandomSolution :: Solver -> IO SAT.Model
331	generateUniformRandomSolution solver = do	1✔
332	gen <- getRandomGen solver	1✔
333	n <- getNumVars solver	1✔
334	(a :: IOUArray Int Bool) <- newArray_ (1,n)	1✔
335	forM_ [1..n] $ \v -> do	1✔
336	b <- Rand.uniform gen	1✔
337	writeArray a v b	1✔
338	unsafeFreeze a	1✔
339
340	checkCurrentSolution :: Solver -> Callbacks -> IO ()
341	checkCurrentSolution solver cb = do	1✔
342	best <- readIORef (svBestSolution solver)	1✔
343	obj <- readIORef (svObj solver)	1✔
344	when (obj < fst best) $ do	1✔
345	sol <- freeze (svSolution solver)	1✔
346	writeIORef (svBestSolution solver) (obj, sol)	1✔
347	cbOnUpdateBestSolution cb solver obj sol	×
348
349	pickClause :: Solver -> Options -> IO PackedClause
350	pickClause solver opt = do	1✔
351	gen <- getRandomGen solver	1✔
352	if optPickClauseWeighted opt then do	×
353	obj <- readIORef (svObj solver)	×
354	let f !j !x = do	×
355	c <- Vec.read (svUnsatClauses solver) j	×
356	let w = svClauseWeights solver ! c	×
357	if x < w then	×
358	return c	×
359	else
360	f (j + 1) (x - w)	×
361	x <- rand obj gen	×
362	c <- f 0 x	×
363	return $ (svClauses solver ! c)	×
364	else do	1✔
365	s <- Vec.getSize (svUnsatClauses solver)	1✔
366	j <- Rand.uniformR (0, s - 1) gen -- For integral types inclusive range is used	1✔
367	liftM (svClauses solver !) $ Vec.read (svUnsatClauses solver) j	1✔
368
369	rand :: PrimMonad m => Integer -> Rand.Gen (PrimState m) -> m Integer
370	rand n gen	×
371	\| n <= toInteger (maxBound :: Word32) = liftM toInteger $ Rand.uniformR (0, fromIntegral n - 1 :: Word32) gen	×
372	\| otherwise = do	×
373	a <- rand (n `shiftR` 32) gen	×
374	(b::Word32) <- Rand.uniform gen	×
375	return $ (a `shiftL` 32) .\|. toInteger b	×
376
377	data Finished = Finished
378	deriving (Show, Typeable)	×
379
380	instance Exception Finished	×
381
382	-- -------------------------------------------------------------------
383
384	probsat :: Solver -> Options -> Callbacks -> (Double -> Double -> Double) -> IO ()
385	probsat solver opt cb f = do	1✔
386	gen <- getRandomGen solver	1✔
387	let maxClauseLen =	1✔
388	if rangeSize (bounds (svClauses solver)) == 0	1✔
389	then 0	1✔
390	else maximum $ map (rangeSize . bounds) $ elems (svClauses solver)	1✔
391	(wbuf :: IOUArray Int Double) <- newArray_ (0, maxClauseLen-1)	1✔
392	wsumRef <- newIOURef (0 :: Double)	1✔
393
394	let pickVar :: PackedClause -> IO SAT.Var
395	pickVar c = do	1✔
396	writeIOURef wsumRef 0	1✔
397	forAssocsM_ c $ \(k,lit) -> do	1✔
398	let v = SAT.litVar lit	1✔
399	m <- getMakeValue solver v	1✔
400	b <- getBreakValue solver v	1✔
401	let w = f m b	1✔
402	writeArray wbuf k w	1✔
403	modifyIOURef wsumRef (+w)	1✔
404	wsum <- readIOURef wsumRef	1✔
405
406	let go :: Int -> Double -> IO Int
407	go !k !a = do	1✔
408	if not (inRange (bounds c) k) then do	×
409	return $! snd (bounds c)	×
410	else do	1✔
411	w <- readArray wbuf k	1✔
412	if a <= w then	1✔
413	return k	1✔
414	else
415	go (k + 1) (a - w)	1✔
416	k <- go 0 =<< Rand.uniformR (0, wsum) gen	1✔
417	return $! SAT.litVar (c ! k)	1✔
418
419	startCPUTime <- getTime ProcessCPUTime	1✔
420	flipsRef <- newIOURef (0::Int)	1✔
421
422	-- It's faster to use Control.Exception than using Control.Monad.Except
423	let body = do	1✔
424	replicateM_ (optMaxTries opt) $ do	1✔
425	sol <- cbGenerateInitialSolution cb solver	1✔
426	setSolution solver sol	1✔
427	checkCurrentSolution solver cb	1✔
428	replicateM_ (optMaxFlips opt) $ do	1✔
429	s <- Vec.getSize (svUnsatClauses solver)	1✔
430	when (s == 0) $ throw Finished	×
431	obj <- readIORef (svObj solver)	1✔
432	when (obj <= optTarget opt) $ throw Finished	×
433	c <- pickClause solver opt	1✔
434	v <- pickVar c	1✔
435	flipVar solver v	1✔
436	modifyIOURef flipsRef inc	1✔
437	checkCurrentSolution solver cb	1✔
438	body `catch` (\(_::Finished) -> return ())	×
439
440	endCPUTime <- getTime ProcessCPUTime	1✔
441	flips <- readIOURef flipsRef	1✔
442	let totalCPUTime = endCPUTime `diffTimeSpec` startCPUTime	×
443	totalCPUTimeSec = fromIntegral (toNanoSecs totalCPUTime) / 10^(9::Int)	×
444	writeIORef (svStatistics solver) $	1✔
445	Statistics	×
446	{ statTotalCPUTime = totalCPUTime	×
447	, statFlips = flips	×
448	, statFlipsPerSecond = fromIntegral flips / totalCPUTimeSec	×
449	}
450
451	return ()	×
452
453
454
455	walksat :: Solver -> Options -> Callbacks -> Double -> IO ()
456	walksat solver opt cb p = do	×
457	gen <- getRandomGen solver	×
458	(buf :: Vec.UVec SAT.Var) <- Vec.new	×
459
460	let pickVar :: PackedClause -> IO SAT.Var
461	pickVar c = do	×
462	Vec.clear buf	×
463	let (lb,ub) = bounds c	×
464	r <- runExceptT $ do	×
465	_ <- numLoopState lb ub (1.0/0.0) $ \ !b0 !i -> do	×
466	let v = SAT.litVar (c ! i)	×
467	b <- lift $ getBreakValue solver v	×
468	if b <= 0 then	×
469	throwE v -- freebie move	×
470	else if b < b0 then do	×
471	lift $ Vec.clear buf >> Vec.push buf v	×
472	return b	×
473	else if b == b0 then do	×
474	lift $ Vec.push buf v	×
475	return b0	×
476	else do	×
477	return b0	×
478	return ()	×
479	case r of	×
480	Left v -> return v	×
481	Right _ -> do	×
482	flag <- Rand.bernoulli p gen	×
483	if flag then do	×
484	-- random walk move
485	i <- Rand.uniformR (lb,ub) gen	×
486	return $! SAT.litVar (c ! i)	×
487	else do	×
488	-- greedy move
489	s <- Vec.getSize buf	×
490	if s == 1 then	×
491	Vec.unsafeRead buf 0	×
492	else do	×
493	i <- Rand.uniformR (0, s - 1) gen	×
494	Vec.unsafeRead buf i	×
495
496	startCPUTime <- getTime ProcessCPUTime	×
497	flipsRef <- newIOURef (0::Int)	×
498
499	-- It's faster to use Control.Exception than using Control.Monad.Except
500	let body = do	×
501	replicateM_ (optMaxTries opt) $ do	×
502	sol <- cbGenerateInitialSolution cb solver	×
503	setSolution solver sol	×
504	checkCurrentSolution solver cb	×
505	replicateM_ (optMaxFlips opt) $ do	×
506	s <- Vec.getSize (svUnsatClauses solver)	×
507	when (s == 0) $ throw Finished	×
508	obj <- readIORef (svObj solver)	×
509	when (obj <= optTarget opt) $ throw Finished	×
510	c <- pickClause solver opt	×
511	v <- pickVar c	×
512	flipVar solver v	×
513	modifyIOURef flipsRef inc	×
514	checkCurrentSolution solver cb	×
515	body `catch` (\(_::Finished) -> return ())	×
516
517	endCPUTime <- getTime ProcessCPUTime	×
518	flips <- readIOURef flipsRef	×
519	let totalCPUTime = endCPUTime `diffTimeSpec` startCPUTime	×
520	totalCPUTimeSec = fromIntegral (toNanoSecs totalCPUTime) / 10^(9::Int)	×
521	writeIORef (svStatistics solver) $	×
522	Statistics	×
523	{ statTotalCPUTime = totalCPUTime	×
524	, statFlips = flips	×
525	, statFlipsPerSecond = fromIntegral flips / totalCPUTimeSec	×
526	}
527
528	return ()	×
529
530	-- -------------------------------------------------------------------
531
532	#if !MIN_VERSION_array(0,5,6)
533
534	{-# INLINE modifyArray #-}
535	modifyArray :: (MArray a e m, Ix i) => a i e -> i -> (e -> e) -> m ()
536	modifyArray a i f = do	1✔
537	e <- readArray a i	1✔
538	writeArray a i (f e)	1✔
539
540	#endif
541
542	{-# INLINE forAssocsM_ #-}
543	forAssocsM_ :: (IArray a e, Monad m) => a Int e -> ((Int,e) -> m ()) -> m ()
544	forAssocsM_ a f = do	1✔
545	let (lb,ub) = bounds a	1✔
546	numLoop lb ub $ \i ->	1✔
547	f (i, unsafeAt a i)	1✔
548
549	{-# INLINE inc #-}
550	inc :: Integral a => a -> a
551	inc a = a+1	1✔
552
553	-- -------------------------------------------------------------------

msakai / toysolver / 496

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous