Cannot figure out a weight balancing puzzle using SMT-LIB and Z3

I have tried for a week to model my problem using the SMT-LIB but have found it really troubling to figure out exactly how to use my logic.

(declare-const w1 Int)
(declare-const w2 Int)
(declare-const w3 Int)
(declare-const w4 Int)

(assert (and
    (<= 1 w1 40)
    (<= 1 w2 40)
    (<= 1 w3 40)
    (<= 1 w4 40)
    
    
    (<= w1 w2)
    (<= w2 w3)
    (<= w3 w4)
    
    
    (forall ((i Int))
        (implies
            
            (<= 1 i 40)
            (exists((j Int) (k Int) (l Int) (m Int))
                (and
                    (<= -1 j 1)
                    (<= -1 k 1)
                    (<= -1 l 1)
                    (<= -1 m 1)
                    (= i (+ (* w1 j) (* w2 k) (* w3 l) (* w4 m)))
            
                )
            )
        )
    
    )
))

(check-sat)
(get-model)

I have tried to loop for i to 40 and use the -1 0 or 1 to represent side of weight or 0 if its not used

There's nothing wrong with the way you encoded the problem. Alas, SMT solvers aren't good at reasoning with quantifiers. Even though the search space isn't very big for this particular problem, the internal heuristics fail to recognize that; causing it to take too long, or loop forever. Instead, you should "expand" all the possibilities yourself and present a quantifier-free version of your problem.

Unfortunately, programming these sorts of constraints is hard to do with SMTLib, because it's really not meant for this kind of programming. Instead, I'd recommend using a higher-level API for z3 to address the problem instead. That is, use one of the z3 bindings from Python, Haskell, C, C++, Java, etc.; whichever you're more comfortable with, instead of coding this directly in SMTLib.

Below, I'll code it using SBV (the Haskell interface for z3), and also using the Python interface for z3; as those are the APIs I'm most familiar with.

Haskell

Using SBV (See http://leventerkok.github.io/sbv/ and https://hackage.haskell.org/package/sbv) and a bit of Haskell magic, you can code your problem as:

import Data.SBV

allCombs :: [a] -> [([a], [a], [a])]
allCombs []     = [([], [], [])]
allCombs (x:xs) = concat [[(x:l, m, r), (l, x:m, r), (l, m, x:r)] | (l, m, r) <- allCombs xs]

weigh :: ([SInteger], [SInteger], [SInteger]) -> SInteger
weigh (left, _, right) = sum left - sum right

pick :: IO SatResult
pick = sat $ do
        w1 <- sInteger "w1"
        w2 <- sInteger "w2"
        w3 <- sInteger "w3"
        w4 <- sInteger "w4"

        -- All weights must be positive
        constrain $ w1 .> 0
        constrain $ w2 .> 0
        constrain $ w3 .> 0
        constrain $ w4 .> 0

        -- Symmetry breaking: Order the weights. Strictly speaking this isn't
        -- necessary, but helps with presentation
        constrain $ w1 .<= w2
        constrain $ w2 .<= w3
        constrain $ w3 .<= w4

        let combs  = map weigh (allCombs [w1, w2, w3, w4])
            find g = sAny (.== literal g) combs

        constrain $ sAnd $ map find [1..40]

This essentially takes all possible orderings of the four weights into three groups: Those that are added, those that are ignored, and those that are subtracted. Then it asks the solver to find a combination that can achieve the summation for any value between 1 and 40. When I run this, I get the following answer fairly quickly:

*Main> pick
Satisfiable. Model:
  w1 =  1 :: Integer
  w2 =  3 :: Integer
  w3 =  9 :: Integer
  w4 = 27 :: Integer

And you can convince yourself that this is a good setting for the weights. (In fact, you can show that this is the only solution, by changing the call sat to allSat, and you'll see that SBV confirms this is the only solution.)

Python

The idea is similar, though Python code looks uglier (of course, this is a subjective claim!):

from z3 import *

def allCombs(lst):
    if not lst:
        yield ([], [], [])
    else:
        cur = [lst[0]]
        for (add, ign, sub) in allCombs (lst[1:]):
          yield (cur + add, ign, sub)
          yield (add, cur + ign, sub)
          yield (add, ign, cur + sub)

def weigh(comb):
    return sum(comb[0]) - sum(comb[2])

s = Solver()
w1, w2, w3, w4 = Ints('w1 w2 w3 w4')
s.add(w1 > 0)
s.add(w2 > 0)
s.add(w3 > 0)
s.add(w4 > 0)
s.add(w1 <= w2)
s.add(w2 <= w3)
s.add(w3 <= w4)

allSums = [weigh(comb) for comb in list(allCombs([w1, w2, w3, w4]))]

def constrain(goal):
    s.add(Or([goal == s for s in allSums]))

for i in range(40):
    constrain(i+1)

print(s.check())
print(s.model())

And this prints:

sat
[w1 = 1, w2 = 3, w3 = 9, w4 = 27]

This is the same solution z3 found via the SBV API.