Why does this SBV code stop before hitting the limit I set?

Question

I have this theorem (not sure if that's the right word), and I want to get all the solutions.

pairCube limit = do
    m <- natural exists "m"
    n <- natural exists "n"
    a <- natural exists "a"
    constrain $ m^3 .== n^2
    constrain $ m .< limit
    return $ m + n .== a^2

res <- allSat (pairCube 1000)

-- Run from ghci
extractModels res :: [[Integer]]

This is trying to solve the problem:

There are infinite pairs of integers (m, n) such that m^3 = n^2 and m + n is a perfect square. What is the pair with the greatest m less than 1000?

I know the actual answer, just through brute forcing, but I want to do with SBV.

However, when I run the code it gives only the following values (in the form [m, n, a]): [[9,27,6],[64,512,24],[]]

However, there are several other solutions with an m value less than 1000 that aren't included.

Answer 1

It's always good to give a full program:

{-# LANGUAGE ScopedTypeVariables #-}

import Data.SBV

pairCube :: SInteger -> Symbolic SBool
pairCube limit = do
        (m :: SInteger) <- exists "m"
        (n :: SInteger) <- exists "n"
        (a :: SInteger) <- exists "a"
        constrain $ m^(3::Integer) .== n^(2::Integer)
        constrain $ m .< limit
        return $ m + n .== a^(2::Integer)

main :: IO ()
main = print =<< allSat (pairCube 1000)

When I run it, I get:

Main> main
Solution #1:
  m = 0 :: Integer
  n = 0 :: Integer
  a = 0 :: Integer
Solution #2:
  m =  9 :: Integer
  n = 27 :: Integer
  a = -6 :: Integer
Solution #3:
  m =  1 :: Integer
  n = -1 :: Integer
  a =  0 :: Integer
Solution #4:
  m =  9 :: Integer
  n = 27 :: Integer
  a =  6 :: Integer
Solution #5:
  m =  64 :: Integer
  n = 512 :: Integer
  a = -24 :: Integer
Solution #6:
  m =  64 :: Integer
  n = 512 :: Integer
  a =  24 :: Integer
Unknown

Note the final Unknown.

Essentially, SBV queried Z3, and got 6 solutions; when SBV asked for the 7th, Z3 said "I don't know if there's any other solution." With non-linear arithmetic, this behavior is expected.

To answer the original question (i.e., find the max m), I changed the constraint to read:

constrain $ m .== limit

and coupled it with the following "driver:"

main :: IO ()
main = loop 1000
  where loop (-1) = putStrLn "Can't find the largest m!"
        loop m    = do putStrLn $ "Trying: " ++ show m
                       mbModel <- extractModel `fmap` sat (pairCube m)
                       case mbModel of
                         Nothing -> loop (m-1)
                         Just r  -> print (r :: (Integer, Integer, Integer))

After running about 50 minutes on my machine, Z3 produced:

(576,13824,-120)

So, clearly the allSat based approach is causing Z3 to give-up way earlier than what it can actually achieve if we fix m and iterate ourself. With a non-linear problem, expecting anything faster/better would be too much to ask of a general purpose SMT solver..