Implementing SAT in OCaml Part 5

2026-09-18 · 8 min

Everything through Part 4 is a real, if unoptimized, SMT solver - boolean CDCL plus a difference-logic theory, sharing one conflict machinery. This binds Z3's own C API by hand, the same raw-stub technique as the c-interop-in-ocaml post, runs the identical difference-logic instances through both solvers, and times the gap.

§ 01

The binding surface

Z3 ships a plain C API (z3.h) behind its C++/Python/etc. wrappers - Z3_context, Z3_ast, Z3_solver are all opaque pointers, exactly the shape raw OCaml C stubs are good at wrapping. Only the handful of calls needed for one linear-arithmetic conjunction.

(* z3_stubs.c *)
#include <caml/mlvalues.h>
#include <caml/alloc.h>
#include <caml/memory.h>
#include <caml/fail.h>
#include <z3.h>
#define Context_val(v) (*((Z3_context *) Data_custom_val(v)))

Wrapping the opaque Z3_context pointer as an OCaml custom block, so the GC can call a finalizer to del the context when the OCaml value is collected - the same pattern any FFI wrapping a C handle needs.

static void finalize_context(value v) {
Z3_context ctx = Context_val(v);
if (ctx != NULL) Z3_del_context(ctx);
}
static struct custom_operations context_ops = {
"z3.context",
finalize_context,
custom_compare_default,
custom_hash_default,
custom_serialize_default,
custom_deserialize_default
};

Creating a context. Z3_mk_context (not the _rc variant) manages AST reference counts internally, which keeps the rest of this binding simpler at the cost of never manually decrementing a ref count.

CAMLprim value caml_z3_mk_context(value unit) {
CAMLparam1(unit);
CAMLlocal1(v_ctx);
Z3_config cfg = Z3_mk_config();
Z3_context ctx = Z3_mk_context(cfg);
Z3_del_config(cfg);
v_ctx = caml_alloc_custom(&context_ops, sizeof(Z3_context), 0, 1);
Context_val(v_ctx) = ctx;
CAMLreturn(v_ctx);
}

Integer constants and the sub/le/ge building blocks needed for x - y <= k atoms - each Z3_mk_* call builds an AST node, which Z3_mk_context's automatic ref counting keeps alive as long as the context lives.

CAMLprim value caml_z3_mk_int_var(value v_ctx, value v_name) {
CAMLparam2(v_ctx, v_name);
Z3_context ctx = Context_val(v_ctx);
Z3_symbol sym = Z3_mk_string_symbol(ctx, String_val(v_name));
Z3_sort int_sort = Z3_mk_int_sort(ctx);
Z3_ast var = Z3_mk_const(ctx, sym, int_sort);
CAMLreturn((value) var);
}
CAMLprim value caml_z3_mk_int(value v_ctx, value v_n) {
CAMLparam2(v_ctx, v_n);
Z3_context ctx = Context_val(v_ctx);
Z3_ast n = Z3_mk_int(ctx, Int_val(v_n), Z3_mk_int_sort(ctx));
CAMLreturn((value) n);
}

A caveat worth stating plainly, rather than glossing over: returning a raw Z3_ast cast to value sidesteps the GC-safety machinery every earlier part of this series (and the C-interop post) insisted on for OCaml-managed memory - it is safe here specifically because Z3_ast is a pointer Z3 itself owns and keeps alive via its own internal ref counting, never an OCaml heap block, so the OCaml GC never tries to move or collect it. Treating an opaque foreign pointer as an immediate value this way is standard for stable-address C handles; it would be wrong for anything the OCaml GC is responsible for.

(* Z3_ast is a stable pointer into memory Z3 manages, not the OCaml heap -
safe to smuggle through `value` as long as it never crosses back into
code that expects an OCaml-shaped value, e.g. never pattern-matched on *)

x - y <= k as one AST node - Z3_mk_sub takes an array of operands (n-ary subtraction in the API, used here with exactly two), Z3_mk_le compares it against the constant.

CAMLprim value caml_z3_mk_le_diff(
value v_ctx, value v_x, value v_y, value v_k) {
CAMLparam4(v_ctx, v_x, v_y, v_k);
Z3_context ctx = Context_val(v_ctx);
Z3_ast args[2] = { (Z3_ast) v_x, (Z3_ast) v_y };
Z3_ast diff = Z3_mk_sub(ctx, 2, args);
Z3_ast k = Z3_mk_int(ctx, Int_val(v_k), Z3_mk_int_sort(ctx));
Z3_ast le = Z3_mk_le(ctx, diff, k);
CAMLreturn((value) le);
}

Solver creation, asserting a formula, and checking - Z3_L_TRUE/Z3_L_FALSE/Z3_L_UNDEF are Z3's own three-valued result, mapped onto an OCaml variant at the boundary rather than leaking Z3's int encoding into calling code.

CAMLprim value caml_z3_mk_solver(value v_ctx) {
CAMLparam1(v_ctx);
Z3_context ctx = Context_val(v_ctx);
Z3_solver s = Z3_mk_solver(ctx);
Z3_solver_inc_ref(ctx, s);
CAMLreturn((value) s);
}
CAMLprim value caml_z3_solver_assert(value v_ctx, value v_solver, value v_ast) {
CAMLparam3(v_ctx, v_solver, v_ast);
Z3_solver_assert(Context_val(v_ctx), (Z3_solver) v_solver, (Z3_ast) v_ast);
CAMLreturn(Val_unit);
}
CAMLprim value caml_z3_solver_check(value v_ctx, value v_solver) {
CAMLparam2(v_ctx, v_solver);
Z3_lbool r = Z3_solver_check(Context_val(v_ctx), (Z3_solver) v_solver);
CAMLreturn(Val_int(r == Z3_L_TRUE ? 0 : r == Z3_L_FALSE ? 1 : 2));
}
§ 02

The OCaml side

dune - linking against libz3 directly. On most package managers z3-dev/z3 installs both the header and the shared library at a findable path; adjust -I/-L if yours differs.

(executable
(name main)
(libraries sat_lib cdcl dpll_t diff_logic dl_theory)
(foreign_stubs
(language c)
(names z3_stubs)
(flags (-I/usr/include)))
(c_library_flags (-lz3)))

lib/z3_binding.ml - the external declarations, and an OCaml-shaped result type replacing Z3's raw int encoding.

type context
type ast
type solver
type result = Sat | Unsat | Unknown
external mk_context : unit -> context = "caml_z3_mk_context"
external mk_int_var : context -> string -> ast = "caml_z3_mk_int_var"
external mk_int : context -> int -> ast = "caml_z3_mk_int"
external mk_le_diff : context -> ast -> ast -> int -> ast = "caml_z3_mk_le_diff"
external mk_solver : context -> solver = "caml_z3_mk_solver"
external solver_assert : context -> solver -> ast -> unit = "caml_z3_solver_assert"
external solver_check_raw : context -> solver -> int = "caml_z3_solver_check"
let solver_check ctx s =
match solver_check_raw ctx s with
| 0 -> Sat
| 1 -> Unsat
| _ -> Unknown

A thin wrapper turning the same Diff_logic.dl_atom list from Part 4 directly into Z3 assertions - the two solvers now consume literally the same input type, which is what makes a fair benchmark possible.

(* lib/z3_dl.ml *)
let check (edges : Diff_logic.dl_atom list) : Z3_binding.result =
let open Z3_binding in
let ctx = mk_context () in
let s = mk_solver ctx in
let vars = Hashtbl.create 16 in
let var_of i =
match Hashtbl.find_opt vars i with
| Some a -> a
| None ->
let a = mk_int_var ctx (Printf.sprintf "x%d" i) in
Hashtbl.replace vars i a;
a
in
List.iter
(fun (e : Diff_logic.dl_atom) ->
let ast = mk_le_diff ctx (var_of e.x) (var_of e.y) (int_of_float e.k) in
solver_assert ctx s ast)
edges;
solver_check ctx s

Confirming it against the two scheduling examples from Part 4 - the four-unit deadline (infeasible) and the six-unit one (feasible).

# Z3_dl.check
[ { x=1; y=0; k=(-3.0) }; { x=2; y=1; k=(-2.0) }; { x=0; y=2; k=4.0 } ];;
- : Z3_binding.result = Unsat
# Z3_dl.check
[ { x=1; y=0; k=(-3.0) }; { x=2; y=1; k=(-2.0) }; { x=0; y=2; k=6.0 } ];;
- : Z3_binding.result = Sat
§ 03

A shared benchmark generator

A random negative-cycle-free chain of difference constraints, plus one closing edge whose weight decides feasibility - large enough to actually take measurable time, small enough to fit in memory for the from-scratch solver too.

(* lib/dl_gen.ml *)
let random_instance ~(size : int) ~(feasible : bool) : Diff_logic.dl_atom list =
Random.self_init ();
let chain =
List.init (size - 1) (fun i ->
({ Diff_logic.x = i + 1; y = i; k = -.(1.0 +. Random.float 3.0) }
: Diff_logic.dl_atom))
in
let total_min = List.fold_left (fun acc (e : Diff_logic.dl_atom) -> acc -. e.k) 0.0 chain in
let closing_k = if feasible then total_min +. 1.0 else total_min -. 1.0 in
chain @ [ { Diff_logic.x = 0; y = size - 1; k = closing_k } ]

Boolean-satisfiable-but-theory-relevant instances for the from-scratch solver: one unit clause per atom, exactly as Part 4's worked examples, generated at whatever size is being benchmarked.

let to_cnf (atoms : Dl_theory.dl_predicate list) : Cdcl.clause list =
List.mapi
(fun i _ -> ({ Cdcl.lits = [| i + 1 |]; learnt = false } : Cdcl.clause))
atoms
§ 04

Timing both

A tiny timing harness - Unix.gettimeofday around each call, big enough a difference that clock resolution does not matter.

(* bin/bench.ml *)
let time label f =
let t0 = Unix.gettimeofday () in
let r = f () in
Printf.printf "%-28s %.4fs\n" label (Unix.gettimeofday () -. t0);
r
let () =
let size = int_of_string Sys.argv.(1) in
let edges = Dl_gen.random_instance ~size ~feasible:false in
let predicates =
List.map (fun (e : Diff_logic.dl_atom) -> Dl_theory.Le (e.x, e.y, e.k)) edges
in
let az = Atomize.make size predicates in
let f = Dl_gen.to_cnf predicates in
ignore
(time "from-scratch DPLL(T)" (fun () ->
Dpll_t.solve_loop (Cdcl.make_solver (List.length predicates) f) az));
ignore (time "z3" (fun () -> Z3_dl.check edges))

Both agree on small instances - the from-scratch solver is slower in absolute terms even here, since it rescans every clause on every propagation step (Part 2's named limitation), but the gap is not yet the story.

$ dune exec ./bin/bench.exe -- 50
from-scratch DPLL(T) 0.0140s
z3 0.0021s

At a few hundred variables the gap stops being a constant factor. Z3's incremental Simplex and watched-literal propagation scale roughly linearly on chain-shaped difference constraints; this solver's O(clauses) rescanning on every propagation step, times the number of decisions, does not.

$ dune exec ./bin/bench.exe -- 500
from-scratch DPLL(T) 1.8210s
z3 0.0087s
$ dune exec ./bin/bench.exe -- 2000
from-scratch DPLL(T) (still running after 60s - not shown)
z3 0.0340s
§ 05

What specifically causes the blowup

Named against Part 2 and Part 4's own honest limitations, now with a number attached: propagate rescans every clause on every single unit assignment - on a chain instance of size n, that is clauses scanned times per decision level, and up to n decision levels, for worst-case work before a single theory check even runs. Z3's two-watched-literal scheme makes each propagation step amortized instead.

this solver: O(n) clauses x O(n) propagation steps x O(n) decision levels
watched literals: O(1) amortized per propagation step, independent of clause count

Second: Diff_logic.check reruns Bellman-Ford from scratch - - on every single call, rather than incrementally updating shortest-path distances as edges are added and removed the way Z3's theory solvers do. On the chain benchmark, that is another per theory call, called at every boolean-complete point.

this solver: O(n^2) per theory check, from scratch every time
incremental: update only the paths touched by the newest edge, amortized far below that

Neither gap is a bug - both are exactly the two optimizations named as out of scope back in Part 2 and Part 4, now visible as actual wall-clock numbers instead of a bullet list. Watched literals and incremental theory propagation are what separate a solver that is correct from one that is also fast.

(* the from-scratch solver in this series was never meant to compete with Z3 -
it was meant to make every one of Z3's internal design decisions legible,
by building a smaller version of each one and feeling where it hurts *)
§ 06

The opam z3 package, for anything real

For actual production use, the official Microsoft-maintained bindings on opam wrap the full API (including incremental push/pop, unsat cores, and model extraction) far beyond the handful of calls bound by hand above - reach for these, not the hand-rolled stub file, once the point is solving problems rather than understanding the solver.

$ opam install z3

The same difference-logic check through the official bindings, for comparison against the by-hand version in section 2 - considerably less code, at the cost of depending on a much larger library surface.

let check_official (edges : Diff_logic.dl_atom list) =
let open Z3 in
let ctx = mk_context [] in
let vars = Hashtbl.create 16 in
let var_of i =
match Hashtbl.find_opt vars i with
| Some v -> v
| None ->
let v = Arithmetic.Integer.mk_const_s ctx (Printf.sprintf "x%d" i) in
Hashtbl.replace vars i v; v
in
let s = Solver.mk_solver ctx None in
List.iter
(fun (e : Diff_logic.dl_atom) ->
let lhs = Arithmetic.mk_sub ctx [ var_of e.x; var_of e.y ] in
let rhs = Arithmetic.Integer.mk_numeral_i ctx (int_of_float e.k) in
Solver.add s [ Arithmetic.mk_le ctx lhs rhs ])
edges;
Solver.check s []