Skip to content

Instantly share code, notes, and snippets.

@msakai
Created Mar 2, 2021
Embed
What would you like to do?
{-# LANGUAGE BangPatterns, ScopedTypeVariables #-}
import Test.QuickCheck
import qualified Test.QuickCheck.Monadic as QM
import qualified Z3.Monad as Z3
{-
a x + b y = n d where d = gcd(a,b) の解は必ず
(x, y) = (n x0 + b/d k, n y0 - a/d k)
の形で書けることを示したい。
quickCheck prop_bezout2 する限りはそういうのはなさそう。
-}
prop_bezout2 = QM.monadicIO $ do
NonZero (a :: Integer) <- QM.pick arbitrary
NonZero (b :: Integer) <- QM.pick arbitrary
let (d, x0, y0) = exgcd a b
QM.monitor (counterexample (show (d,x0,y0)))
(n :: Integer) <- QM.pick arbitrary
ret <- QM.run $ Z3.evalZ3 $ do
a_ <- Z3.mkIntNum a
b_ <- Z3.mkIntNum b
n_ <- Z3.mkIntNum n
x <- Z3.mkIntVar =<< Z3.mkStringSymbol "x"
y <- Z3.mkIntVar =<< Z3.mkStringSymbol "y"
e <- do
ax <- Z3.mkMul [a_, x]
by <- Z3.mkMul [b_, y]
Z3.mkAdd [ax, by]
Z3.solverAssertCnstr =<< Z3.mkEq e =<< Z3.mkIntNum (n * d)
x0_ <- Z3.mkIntNum x0
y0_ <- Z3.mkIntNum y0
tmp1 <- Z3.mkMul [n_, x0_] >>= \nx0 -> Z3.mkSub [x, nx0]
tmp2 <- Z3.mkMul [n_, y0_] >>= \ny0 -> Z3.mkSub [y, ny0]
zero <- Z3.mkIntNum 0
cond1 <- Z3.mkEq zero =<< Z3.mkMod tmp1 =<< Z3.mkIntNum (b `div` d)
cond2 <- Z3.mkEq zero =<< Z3.mkMod tmp2 =<< Z3.mkIntNum (a `div` d)
k1 <- Z3.mkDiv tmp1 =<< Z3.mkIntNum (b `div` d)
k2 <- Z3.mkUnaryMinus =<< Z3.mkDiv tmp2 =<< Z3.mkIntNum (a `div` d)
cond3 <- Z3.mkEq k1 k2
Z3.solverAssertCnstr =<< Z3.mkNot =<< Z3.mkAnd [cond1, cond2, cond3]
ret <- Z3.solverCheck
case ret of
Z3.Sat -> do
model <- Z3.solverGetModel
Just xVal <- Z3.evalInt model x
Just yVal <- Z3.evalInt model y
return $ Just (xVal, yVal)
_ -> return Nothing
case ret of
Nothing -> return ()
Just (x,y) -> do
QM.monitor (counterexample (show (x,y)))
QM.assert False
exgcd :: Integral a => a -> a -> (a, a, a)
exgcd f1 f2 = f $ go f1 f2 1 0 0 1
where
go !r0 !r1 !s0 !s1 !t0 !t1
| r1 == 0 = (r0, s0, t0)
| otherwise = go r1 r2 s1 s2 t1 t2
where
(q, r2) = r0 `divMod` r1
s2 = s0 - q*s1
t2 = t0 - q*t1
f (g,u,v)
| g < 0 = (-g, -u, -v)
| otherwise = (g,u,v)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment