sbv-14.3: Documentation/SBV/Examples/TP/Numeric.hs
-----------------------------------------------------------------------------
-- |
-- Module : Documentation.SBV.Examples.TP.Numeric
-- Copyright : (c) Levent Erkok
-- License : BSD3
-- Maintainer: erkokl@gmail.com
-- Stability : experimental
--
-- Example use of inductive TP proofs, over integers.
-----------------------------------------------------------------------------
{-# LANGUAGE CPP #-}
{-# LANGUAGE DataKinds #-}
{-# LANGUAGE QuasiQuotes #-}
{-# LANGUAGE ScopedTypeVariables #-}
{-# LANGUAGE TypeAbstractions #-}
{-# LANGUAGE TypeApplications #-}
{-# OPTIONS_GHC -Wall -Werror #-}
module Documentation.SBV.Examples.TP.Numeric where
import Prelude hiding (sum, map, product, length, (^), replicate, elem)
import Data.SBV
import Data.SBV.TP
import Data.SBV.List
#ifdef DOCTEST
-- $setup
-- >>> :set -XScopedTypeVariables
-- >>> import Data.SBV
-- >>> import Data.SBV.TP
-- >>> import Control.Exception
#endif
-- * Sum of constants
-- | \(\sum_{i=1}^{n} c = c \cdot n\)
--
-- >>> runTP $ sumConstProof (uninterpret "c")
-- Inductive lemma: sumConst_correct
-- Step: Base Q.E.D.
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Step: 3 Q.E.D.
-- Step: 4 Q.E.D.
-- Result: Q.E.D.
-- Functions proven terminating: sbv.foldr, sbv.replicate
-- [Proven] sumConst_correct :: Ɐn ∷ Integer → Bool
sumConstProof :: SInteger -> TP (Proof (Forall "n" Integer -> SBool))
sumConstProof c = induct "sumConst_correct"
(\(Forall n) -> n .>= 0 .=> sum (replicate n c) .== c * n) $
\ih n -> [n .>= 0] |- sum (replicate (n+1) c)
=: sum (c .: replicate n c)
=: c + sum (replicate n c)
?? ih
=: c + c*n
=: c*(n+1)
=: qed
-- * Sum of numbers
-- | \(\sum_{i=0}^{n} i = \frac{n(n+1)}{2}\)
--
-- NB. We define the sum of numbers from @0@ to @n@ as @sum [sEnum|n, n-1 .. 0|]@, i.e., we
-- construct the list starting from @n@ going down to @0@. Contrast this to the perhaps more natural
-- definition of @sum [sEnum|0 .. n]@, i.e., going up. While the latter is equivalent functionality, the former
-- works much better with the proof-structure: Since we induct on @n@, in each step we strip of one
-- layer, and the recursion in the down-to construction matches the inductive schema.
--
-- >>> runTP sumProof
-- Inductive lemma: sum_correct
-- Step: Base Q.E.D.
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Step: 3 Q.E.D.
-- Result: Q.E.D.
-- Functions proven terminating: EnumSymbolic.Integer.enumFromThenTo.down, sbv.foldr
-- [Proven] sum_correct :: Ɐn ∷ Integer → Bool
sumProof :: TP (Proof (Forall "n" Integer -> SBool))
sumProof = induct "sum_correct"
(\(Forall n) -> n .>= 0 .=> sum [sEnum|n, n-1 .. 0|] .== (n * (n+1)) `sEDiv` 2) $
\ih n -> [n .>= 0] |- sum [sEnum|n+1, n .. 0|]
=: n+1 + sum [sEnum|n, n-1 .. 0|]
?? ih
=: n+1 + (n * (n+1)) `sEDiv` 2
=: ((n+1) * (n+2)) `sEDiv` 2
=: qed
-- * Sum of squares of numbers
--
-- | \(\sum_{i=0}^{n} i^2 = \frac{n(n+1)(2n+1)}{6}\)
--
-- >>> runTP sumSquareProof
-- Inductive lemma: sumSquare_correct
-- Step: Base Q.E.D.
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Step: 3 Q.E.D.
-- Step: 4 Q.E.D.
-- Step: 5 Q.E.D.
-- Step: 6 Q.E.D.
-- Result: Q.E.D.
-- Functions proven terminating: EnumSymbolic.Integer.enumFromThenTo.down, sbv.foldr, sbv.map
-- [Proven] sumSquare_correct :: Ɐn ∷ Integer → Bool
sumSquareProof :: TP (Proof (Forall "n" Integer -> SBool))
sumSquareProof = do
let sq :: SInteger -> SInteger
sq k = k * k
sumSquare n = sum $ map sq [sEnum|n, n-1 .. 0|]
induct "sumSquare_correct"
(\(Forall n) -> n .>= 0 .=> sumSquare n .== (n*(n+1)*(2*n+1)) `sEDiv` 6) $
\ih n -> [n .>= 0] |- sumSquare (n+1)
=: sum (map sq [sEnum|n+1, n .. 0|])
=: sum (map sq (n+1 .: [sEnum|n, n-1 .. 0|]))
=: sum ((n+1)*(n+1) .: map sq [sEnum|n, n-1 .. 0|])
=: (n+1)*(n+1) + sum (map sq [sEnum|n, n-1 .. 0|])
?? ih
=: (n+1)*(n+1) + (n*(n+1)*(2*n+1)) `sEDiv` 6
=: ((n+1)*(n+2)*(2*n+3)) `sEDiv` 6
=: qed
-- * Sum of cubes of numbers
-- | \(\sum_{i=0}^{n} i^3 = \left( \sum_{i=0}^{n} i \right)^2 = \left( \frac{n(n+1)}{2} \right)^2\)
--
-- This is attributed to Nicomachus, hence the name.
--
-- >>> runTP nicomachus
-- Inductive lemma: sum_correct
-- Step: Base Q.E.D.
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Step: 3 Q.E.D.
-- Result: Q.E.D.
-- Lemma: evenHalfSquared Q.E.D.
-- Inductive lemma: nn1IsEven
-- Step: Base Q.E.D.
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Step: 3 Q.E.D.
-- Result: Q.E.D.
-- Lemma: sum_squared
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Result: Q.E.D.
-- Inductive lemma: nicomachus
-- Step: Base Q.E.D.
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Result: Q.E.D.
-- Functions proven terminating: EnumSymbolic.Integer.enumFromThenTo.down, sbv.foldr, sumCubed
-- [Proven] nicomachus :: Ɐn ∷ Integer → Bool
nicomachus :: TP (Proof (Forall "n" Integer -> SBool))
nicomachus = do
let (^) :: SInteger -> Integer -> SInteger
_ ^ 0 = 1
b ^ n = b * b ^ (n-1)
infixr 8 ^
sumCubed :: SInteger -> SInteger
sumCubed = smtFunction "sumCubed" $ \n -> [sCase| n of
_ | n .<= 0 -> 0
_ -> n^3 + sumCubed (n - 1)
|]
-- Grab the proof of regular summation formula
sp <- sumProof
-- Square of the summation result. This is a trivial lemma for humans, but there are lots
-- of multiplications involved making the problem non-linear and we need to spell it out.
ssp <- do
-- Squaring half of an even number? You can square the number and divide by 4 instead:
-- z3 can prove this out of the box, but without it being explicitly expressed, the
-- following proof doesn't go through.
evenHalfSquared <- lemma "evenHalfSquared"
(\(Forall n) -> 2 `sDivides` n .=> (n `sEDiv` 2) ^ 2 .== (n ^ 2) `sEDiv` 4)
[]
-- The multiplication @n * (n+1)@ is always even. It's surprising that I had to use induction here
-- but neither z3 nor cvc5 can converge on this out-of-the-box.
nn1IsEven <- induct "nn1IsEven"
(\(Forall n) -> n .>= 0 .=> 2 `sDivides` (n * (n+1))) $
\ih n -> [n .>= 0] |- 2 `sDivides` ((n+1) * (n+2))
=: 2 `sDivides` (n*(n+1) + 2*(n+1))
=: 2 `sDivides` (n*(n+1))
?? ih
=: sTrue
=: qed
calc "sum_squared"
(\(Forall @"n" n) -> n .>= 0 .=> sum [sEnum|n, n-1 .. 0|] ^ 2 .== (n^2 * (n+1)^2) `sEDiv` 4) $
\n -> [n .>= 0] |- sum [sEnum|n, n-1 .. 0|] ^ 2
?? sp `at` Inst @"n" n
=: ((n * (n+1)) `sEDiv` 2)^2
?? nn1IsEven `at` Inst @"n" n
?? evenHalfSquared `at` Inst @"n" (n * (n+1))
=: ((n * (n+1))^2) `sEDiv` 4
=: qed
-- We can finally put it together:
induct "nicomachus"
(\(Forall n) -> n .>= 0 .=> sumCubed n .== sum [sEnum|n, n-1 .. 0|] ^ 2) $
\ih n -> [n .>= 0]
|- sumCubed (n+1)
=: (n+1)^3 + sumCubed n
?? ih
?? ssp
=: sum [sEnum|n+1, n .. 0|] ^ 2
=: qed
-- * Exponents and divisibility by 7
-- | \(7 \mid \left(11^n - 4^n\right)\)
--
-- NB. As of Feb 2025, z3 struggles with the inductive step in this proof, but cvc5 performs just fine.
--
-- >>> runTP elevenMinusFour
-- Lemma: powN Q.E.D.
-- Inductive lemma: elevenMinusFour
-- Step: Base Q.E.D.
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Step: 3 Q.E.D.
-- Step: 4 Q.E.D.
-- Step: 5 Q.E.D.
-- Step: 6 Q.E.D.
-- Step: 7 Q.E.D.
-- Step: 8 Q.E.D.
-- Result: Q.E.D.
-- Functions proven terminating: pow
-- [Proven] elevenMinusFour :: Ɐn ∷ Integer → Bool
elevenMinusFour :: TP (Proof (Forall "n" Integer -> SBool))
elevenMinusFour = do
let pow :: SInteger -> SInteger -> SInteger
pow = smtFunction "pow" $ \x y -> [sCase| y of
_ | y .<= 0 -> 1
_ -> x * pow x (y - 1)
|]
emf :: SInteger -> SBool
emf n = 7 `sDivides` (11 `pow` n - 4 `pow` n)
-- helper
powN <- lemma "powN" (\(Forall x) (Forall n) -> n .>= 0 .=> x `pow` (n+1) .== x * x `pow` n) []
inductWith cvc5 "elevenMinusFour"
(\(Forall n) -> n .>= 0 .=> emf n) $
\ih n -> [n .>= 0]
|- emf (n+1)
=: 7 `sDivides` (11 `pow` (n+1) - 4 `pow` (n+1))
?? powN `at` (Inst @"x" 11, Inst @"n" n)
=: 7 `sDivides` (11 * 11 `pow` n - 4 `pow` (n+1))
?? powN `at` (Inst @"x" 4, Inst @"n" n)
=: 7 `sDivides` (11 * 11 `pow` n - 4 * 4 `pow` n)
=: 7 `sDivides` (7 * 11 `pow` n + 4 * 11 `pow` n - 4 * 4 `pow` n)
=: 7 `sDivides` (7 * 11 `pow` n + 4 * (11 `pow` n - 4 `pow` n))
?? ih
=: let x = some "x" (\v -> 7*v .== 11 `pow` n - 4 `pow` n) -- Apply the IH and grab the witness for it
in 7 `sDivides` (7 * 11 `pow` n + 4 * 7 * x)
=: 7 `sDivides` (7 * (11 `pow` n + 4 * x))
=: sTrue
=: qed
-- * A proof about factorials
-- | \(\sum_{k=0}^{n} k \cdot k! = (n+1)! - 1\)
--
-- >>> runTP sumMulFactorial
-- Lemma: fact (n+1)
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Step: 3 Q.E.D.
-- Step: 4 Q.E.D.
-- Result: Q.E.D.
-- Inductive lemma: sumMulFactorial
-- Step: Base Q.E.D.
-- Step: 1 Q.E.D.
-- Step: 2 Q.E.D.
-- Step: 3 Q.E.D.
-- Step: 4 Q.E.D.
-- Step: 5 Q.E.D.
-- Step: 6 Q.E.D.
-- Step: 7 Q.E.D.
-- Result: Q.E.D.
-- Functions proven terminating: EnumSymbolic.Integer.enumFromThenTo.down, sbv.foldr, sbv.map
-- [Proven] sumMulFactorial :: Ɐn ∷ Integer → Bool
sumMulFactorial :: TP (Proof (Forall "n" Integer -> SBool))
sumMulFactorial = do
let fact :: SInteger -> SInteger
fact n = product [sEnum|n, n-1 .. 1|]
-- This is pure expansion, but without it z3 struggles in the next lemma.
helper <- calc "fact (n+1)"
(\(Forall n) -> n .>= 0 .=> fact (n+1) .== (n+1) * fact n) $
\n -> [n .>= 0] |- fact (n+1)
=: product [sEnum|n+1, n .. 1|]
=: product (n+1 .: [sEnum|n, n-1 .. 1|])
=: (n+1) * product [sEnum|n, n-1 .. 1|]
=: (n+1) * fact n
=: qed
induct "sumMulFactorial"
(\(Forall n) -> n .>= 0 .=> sum (map (\k -> k * fact k) [sEnum|n, n-1 .. 0|]) .== fact (n+1) - 1) $
\ih n -> [n .>= 0] |- sum (map (\k -> k * fact k) [sEnum|n+1, n .. 0|])
=: sum (map (\k -> k * fact k) (n+1 .: [sEnum|n, n-1 .. 0|]))
=: sum ((n+1) * fact (n+1) .: map (\k -> k * fact k) [sEnum|n, n-1 .. 0|])
=: (n+1) * fact (n+1) + sum (map (\k -> k * fact k) [sEnum|n, n-1 .. 0|])
?? ih
=: (n+1) * fact (n+1) + fact (n+1) - 1
=: ((n+1) + 1) * fact (n+1) - 1
=: (n+2) * fact (n+1) - 1
?? helper `at` Inst @"n" (n+1)
=: fact (n+2) - 1
=: qed
-- * Product with 0
-- | \(\prod_{x \in xs} x = 0 \iff 0 \in xs\)
--
-- >>> runTP product0
-- Inductive lemma: product0
-- Step: Base Q.E.D.
-- Step: 1 Q.E.D.
-- Step: 2 (2 way case split)
-- Step: 2.1 Q.E.D.
-- Step: 2.2.1 Q.E.D.
-- Step: 2.2.2 Q.E.D.
-- Step: 2.Completeness Q.E.D.
-- Result: Q.E.D.
-- Functions proven terminating: sbv.foldr
-- [Proven] product0 :: Ɐxs ∷ [Integer] → Bool
product0 :: TP (Proof (Forall "xs" [Integer] -> SBool))
product0 =
induct "product0"
(\(Forall @"xs" (xs :: SList Integer)) -> product xs .== 0 .<=> 0 `elem` xs) $
\ih (x, xs) -> [] |- (product (x .: xs) .== 0 .<=> 0 `elem` (x .: xs))
=: (x * product xs .== 0 .<=> x .== 0 .|| 0 `elem` xs)
=: cases [ x .== 0 ==> trivial
, x ./= 0 ==> (x * product xs .== 0 .<=> 0 `elem` xs)
?? ih
=: sTrue
=: qed
]
-- * A negative example
-- | The regular inductive proof on integers (i.e., proving at @0@, assuming at @n@ and proving at
-- @n+1@ will not allow you to conclude things when @n < 0@. The following example demonstrates this with the most
-- obvious example:
--
-- >>> badNonNegative `catch` (\(_ :: SomeException) -> pure ())
-- Inductive lemma: badNonNegative
-- Step: Base Q.E.D.
-- Step: 1
-- *** Failed to prove badNonNegative.1.
-- Falsifiable. Counter-example:
-- n = -2 :: Integer
badNonNegative :: IO ()
badNonNegative = runTP $ do
_ <- induct "badNonNegative"
(\(Forall @"n" (n :: SInteger)) -> n .>= 0) $
\ih n -> [] |- n + 1 .>= 0
?? ih
=: sTrue
=: qed
pure ()