From a908f2b3f28d306d8836cf62baf9f8763cd23a57 Mon Sep 17 00:00:00 2001 From: Simon Cruanes Date: Mon, 15 Feb 2021 16:19:13 -0500 Subject: [PATCH] feat(arith): integrate simplex2 into sidekick; remove old simplex --- src/arith/base-term/Base_types.ml | 28 +- src/arith/lra/linear_expr.ml | 33 -- src/arith/lra/linear_expr.mli | 8 - src/arith/lra/linear_expr_intf.ml | 45 -- src/arith/lra/sidekick_arith_lra.ml | 352 +++++++------ src/arith/lra/simplex.ml | 765 ---------------------------- src/arith/lra/simplex.mli | 33 -- src/arith/tests/run_tests.ml | 1 - src/arith/tests/test_simplex.ml | 180 ------- src/smtlib/Process.ml | 3 + 10 files changed, 220 insertions(+), 1228 deletions(-) delete mode 100644 src/arith/lra/simplex.ml delete mode 100644 src/arith/lra/simplex.mli delete mode 100644 src/arith/tests/test_simplex.ml diff --git a/src/arith/base-term/Base_types.ml b/src/arith/base-term/Base_types.ml index affedece..8d32a736 100644 --- a/src/arith/base-term/Base_types.ml +++ b/src/arith/base-term/Base_types.ml @@ -12,6 +12,8 @@ type 'a lra_view = 'a Sidekick_arith_lra.lra_view = | LRA_op of lra_op * 'a * 'a | LRA_mult of Q.t * 'a | LRA_const of Q.t + | LRA_simplex_var of 'a + | LRA_simplex_pred of 'a * Sidekick_arith_lra.S_op.t * Q.t | LRA_other of 'a (* main term cell. *) @@ -238,6 +240,10 @@ let pp_term_view_gen ~pp_id ~pp_t out = function | LRA_mult (n, x) -> Fmt.fprintf out "(@[*@ %a@ %a@])" Q.pp_print n pp_t x | LRA_const q -> Q.pp_print out q + | LRA_simplex_var v -> pp_t out v + | LRA_simplex_pred (v, op, q) -> + Fmt.fprintf out "(@[%a@ %s %a@])" + pp_t v (Sidekick_arith_lra.S_op.to_string op) Q.pp_print q | LRA_other x -> pp_t out x end @@ -593,6 +599,8 @@ end = struct let sub_hash = A.hash let sub_eq = A.equal + let hash_q q = Hash.string (Q.to_string q) + let hash (t:A.t view) : int = match t with | Bool b -> Hash.bool b | App_fun (f,l) -> @@ -607,8 +615,11 @@ end = struct | LRA_op (p, a, b) -> Hash.combine4 82 (Hash.poly p) (sub_hash a) (sub_hash b) | LRA_mult (n, x) -> - Hash.combine3 83 (Hash.string @@ Q.to_string n) (sub_hash x) - | LRA_const q -> Hash.combine2 84 (Hash.string @@ Q.to_string q) + Hash.combine3 83 (hash_q n) (sub_hash x) + | LRA_const q -> Hash.combine2 84 (hash_q q) + | LRA_simplex_var v -> Hash.combine2 99 (sub_hash v) + | LRA_simplex_pred (v,op,q) -> + Hash.combine4 120 (sub_hash v) (Hash.poly op) (hash_q q) | LRA_other x -> sub_hash x end @@ -630,8 +641,11 @@ end = struct | LRA_const a1, LRA_const a2 -> Q.equal a1 a2 | LRA_mult (n1,x1), LRA_mult (n2,x2) -> Q.equal n1 n2 && sub_eq x1 x2 | LRA_other x1, LRA_other x2 -> sub_eq x1 x2 - | (LRA_pred _ | LRA_op _ | LRA_const _ - | LRA_mult _ | LRA_other _), _ -> false + | LRA_simplex_var v1, LRA_simplex_var v2 -> sub_eq v1 v2 + | LRA_simplex_pred (v1,op1,q1), LRA_simplex_pred (v2,op2,q2) -> + sub_eq v1 v2 && (op1==op2) && Q.equal q1 q2 + | (LRA_pred _ | LRA_op _ | LRA_const _ | LRA_simplex_var _ + | LRA_mult _ | LRA_other _ | LRA_simplex_pred _), _ -> false end | (Bool _ | App_fun _ | Eq _ | Not _ | Ite _ | LRA _), _ -> false @@ -700,8 +714,8 @@ end = struct end | LRA l -> begin match l with - | LRA_pred _ -> Ty.bool - | LRA_op _ | LRA_const _ | LRA_mult _ -> Ty.real + | LRA_pred _ | LRA_simplex_pred _ -> Ty.bool + | LRA_op _ | LRA_const _ | LRA_mult _ | LRA_simplex_var _ -> Ty.real | LRA_other x -> x.term_ty end @@ -716,6 +730,8 @@ end = struct begin match l with | LRA_pred (_, a, b)|LRA_op (_, a, b) -> f a; f b | LRA_mult (_,x) | LRA_other x -> f x + | LRA_simplex_var x -> f x + | LRA_simplex_pred (x,_,_) -> f x | LRA_const _ -> () end diff --git a/src/arith/lra/linear_expr.ml b/src/arith/lra/linear_expr.ml index 306fcd87..8691f6cc 100644 --- a/src/arith/lra/linear_expr.ml +++ b/src/arith/lra/linear_expr.ml @@ -4,9 +4,6 @@ module type COEFF = Linear_expr_intf.COEFF module type VAR = Linear_expr_intf.VAR -module type FRESH = Linear_expr_intf.FRESH -module type VAR_GEN = Linear_expr_intf.VAR_GEN -module type VAR_EXTENDED = Linear_expr_intf.VAR_EXTENDED module type S = Linear_expr_intf.S @@ -195,33 +192,3 @@ module Make(C : COEFF)(Var : VAR) = struct end end[@@inline] -module Make_var_gen(Var : VAR) - : VAR_EXTENDED with type user_var = Var.t - and type lit = Var.lit -= struct - type user_var = Var.t - - type t = - | User of user_var - | Internal of int - - let compare (a:t) b : int = match a, b with - | User a, User b -> Var.compare a b - | User _, Internal _ -> -1 - | Internal _, User _ -> 1 - | Internal i, Internal j -> CCInt.compare i j - - let pp out = function - | User v -> Var.pp out v - | Internal i -> Format.fprintf out "internal_v_%d" i - - type lit = Var.lit - let pp_lit = Var.pp_lit - - module Fresh = struct - type t = int ref - let create() = ref 0 - let copy r = ref !r - let fresh r = Internal (CCRef.get_then_incr r) - end -end[@@inline] diff --git a/src/arith/lra/linear_expr.mli b/src/arith/lra/linear_expr.mli index 3720cb85..35d18228 100644 --- a/src/arith/lra/linear_expr.mli +++ b/src/arith/lra/linear_expr.mli @@ -7,9 +7,6 @@ module type COEFF = Linear_expr_intf.COEFF module type VAR = Linear_expr_intf.VAR -module type FRESH = Linear_expr_intf.FRESH -module type VAR_GEN = Linear_expr_intf.VAR_GEN -module type VAR_EXTENDED = Linear_expr_intf.VAR_EXTENDED module type S = Linear_expr_intf.S @@ -19,8 +16,3 @@ module Make(C : COEFF)(Var : VAR) : S with module C = C and module Var = Var and module Var_map = CCMap.Make(Var) - -module Make_var_gen(Var : VAR) - : VAR_EXTENDED - with type user_var = Var.t - and type lit = Var.lit diff --git a/src/arith/lra/linear_expr_intf.ml b/src/arith/lra/linear_expr_intf.ml index a301ea90..672818c9 100644 --- a/src/arith/lra/linear_expr_intf.ml +++ b/src/arith/lra/linear_expr_intf.ml @@ -57,51 +57,6 @@ module type VAR = sig val pp_lit : lit Fmt.printer end -(** {2 Fresh variables} - - Standard interface for variables with an infinite number - of 'fresh' variables. A 'fresh' variable should be distinct - from any other. -*) -module type FRESH = sig - type var - (** The type of variables. *) - - type t - (** A type of state for creating fresh variables. *) - - val copy : t -> t - (** Copy state *) - - val fresh : t -> var - (** Create a fresh variable using an existing variable as base. - TODO: need some explaining, about the difference with {!create}. *) -end - -(** {2 Generative Variable interface} - - Standard interface for variables that are meant to be used - in expressions. Furthermore, fresh variables can be generated - (which is useful to refactor and/or put problems in specific - formats used by algorithms). -*) -module type VAR_GEN = sig - include VAR - - (** Generate fresh variables on demand *) - module Fresh : FRESH with type var := t -end - -module type VAR_EXTENDED = sig - type user_var (** original variables *) - - type t = - | User of user_var - | Internal of int - - include VAR_GEN with type t := t -end - type bool_op = Predicate.t = Leq | Geq | Lt | Gt | Eq | Neq (** {2 Linear expressions & formulas} *) diff --git a/src/arith/lra/sidekick_arith_lra.ml b/src/arith/lra/sidekick_arith_lra.ml index 7d48baf8..09d23ceb 100644 --- a/src/arith/lra/sidekick_arith_lra.ml +++ b/src/arith/lra/sidekick_arith_lra.ml @@ -6,11 +6,12 @@ open Sidekick_core -module Simplex = Simplex module Simplex2 = Simplex2 module Predicate = Predicate module Linear_expr = Linear_expr +module S_op = Simplex2.Op + type pred = Linear_expr_intf.bool_op = Leq | Geq | Lt | Gt | Eq | Neq type op = Plus | Minus @@ -19,6 +20,8 @@ type 'a lra_view = | LRA_op of op * 'a * 'a | LRA_mult of Q.t * 'a | LRA_const of Q.t + | LRA_simplex_var of 'a (* an opaque variable *) + | LRA_simplex_pred of 'a * S_op.t * Q.t (* an atomic constraint *) | LRA_other of 'a let map_view f (l:_ lra_view) : _ lra_view = @@ -27,6 +30,8 @@ let map_view f (l:_ lra_view) : _ lra_view = | LRA_op (p, a, b) -> LRA_op (p, f a, f b) | LRA_mult (n,a) -> LRA_mult (n, f a) | LRA_const q -> LRA_const q + | LRA_simplex_var v -> LRA_simplex_var (f v) + | LRA_simplex_pred (v, op, q) -> LRA_simplex_pred (f v, op, q) | LRA_other x -> LRA_other (f x) end @@ -44,6 +49,9 @@ module type ARG = sig val ty_lra : S.T.Term.state -> ty + val mk_and : S.T.Term.state -> term -> term -> term + val mk_or : S.T.Term.state -> term -> term -> term + val has_ty_real : term -> bool (** Does this term have the type [Real] *) @@ -98,9 +106,8 @@ module Make(A : ARG) : S with module A = A = struct end module SimpVar - : Linear_expr.VAR_GEN + : Linear_expr.VAR with type t = A.term - and type Fresh.t = A.Gensym.t and type lit = Tag.t = struct type t = A.term @@ -108,21 +115,12 @@ module Make(A : ARG) : S with module A = A = struct let compare = A.S.T.Term.compare type lit = Tag.t let pp_lit = Tag.pp - module Fresh = struct - type t = A.Gensym.t - let copy = A.Gensym.copy - let fresh (st:t) = - let ty = A.ty_lra (A.Gensym.tst st) in - A.Gensym.fresh_term ~pre:"_lra" st ty - end end - module SimpSolver = Simplex.Make_full(SimpVar) - - (* linear expressions *) - module LComb = SimpSolver.L.Comb - module LE = SimpSolver.L.Expr - module LConstr = SimpSolver.L.Constr + module LE_ = Linear_expr.Make(struct include Q let pp=pp_print end)(SimpVar) + module LE = LE_.Expr + module SimpSolver = Simplex2.Make(SimpVar) + module LConstr = SimpSolver.Constraint type state = { tst: T.state; @@ -130,11 +128,13 @@ module Make(A : ARG) : S with module A = A = struct gensym: A.Gensym.t; neq_encoded: unit T.Tbl.t; (* if [a != b] asserted and not in this table, add clause [a = b \/ ab] *) - needs_th_combination: LE.t T.Tbl.t; (* terms that require theory combination *) + needs_th_combination: LE_.Comb.t T.Tbl.t; (* terms that require theory combination *) t_defs: LE.t T.Tbl.t; (* term definitions *) pred_defs: (pred * LE.t * LE.t * T.t * T.t) T.Tbl.t; (* predicate definitions *) local_eqs: (N.t * N.t) Backtrack_stack.t; (* inferred by the congruence closure *) + simplex: SimpSolver.t; } + (* TODO *) let create tst : state = { tst; @@ -145,13 +145,16 @@ module Make(A : ARG) : S with module A = A = struct t_defs=T.Tbl.create 8; pred_defs=T.Tbl.create 16; local_eqs = Backtrack_stack.create(); + simplex=SimpSolver.create (); } let push_level self = + SimpSolver.push_level self.simplex; Backtrack_stack.push_level self.local_eqs; () let pop_levels self n = + SimpSolver.pop_levels self.simplex n; Backtrack_stack.pop_levels self.local_eqs n ~f:(fun _ -> ()); () @@ -212,7 +215,7 @@ module Make(A : ARG) : S with module A = A = struct let open LE.Infix in match A.view_as_lra t with | LRA_other _ -> LE.monomial1 (f t) - | LRA_pred _ -> + | LRA_pred _ | LRA_simplex_pred _ -> Error.errorf "type error: in linexp, LRA predicate %a" T.pp t | LRA_op (op, t1, t2) -> let t1 = as_linexp ~f t1 in @@ -224,6 +227,7 @@ module Make(A : ARG) : S with module A = A = struct | LRA_mult (n, x) -> let t = as_linexp ~f x in LE.( n * t ) + | LRA_simplex_var v -> LE.monomial1 v | LRA_const q -> LE.of_const q let as_linexp_id = as_linexp ~f:CCFun.id @@ -233,98 +237,175 @@ module Make(A : ARG) : S with module A = A = struct *) (* preprocess linear expressions away *) - let preproc_lra (self:state) si ~recurse ~mk_lit:_ ~add_clause:_ (t:T.t) : T.t option = + let preproc_lra (self:state) si ~recurse ~mk_lit ~add_clause (t:T.t) : T.t option = Log.debugf 50 (fun k->k "lra.preprocess %a" T.pp t); let tst = SI.tst si in + + let mk_eq x q = + let t1 = A.mk_lra tst (LRA_simplex_pred (x, Leq, q)) in + let t2 = A.mk_lra tst (LRA_simplex_pred (x, Geq, q)) in + A.mk_and tst t1 t2 + and mk_neq x q = + let t1 = A.mk_lra tst (LRA_simplex_pred (x, Lt, q)) in + let t2 = A.mk_lra tst (LRA_simplex_pred (x, Gt, q)) in + A.mk_or tst t1 t2 + in + match A.view_as_lra t with - | LRA_pred ((Eq|Neq) as pred, t1, t2) -> - (* keep equality as is, needed for congruence closure *) - let t1 = recurse t1 in - let t2 = recurse t2 in - let u = A.mk_lra tst (LRA_pred (pred, t1, t2)) in - if T.equal t u then None else Some u | LRA_pred (pred, t1, t2) -> let l1 = as_linexp ~f:recurse t1 in let l2 = as_linexp ~f:recurse t2 in - let proxy = fresh_term self ~pre:"_pred_lra_" Ty.bool in - T.Tbl.add self.pred_defs proxy (pred, l1, l2, t1, t2); - Log.debugf 5 (fun k->k"@[lra.preprocess.step %a@ :into %a@ :def %a@]" - T.pp t T.pp proxy pp_pred_def (pred,l1,l2)); - Some proxy + let le = LE.(l1 - l2) in + let le_comb, le_const = LE.comb le, LE.const le in + let le_const = Q.neg le_const in + (* now we have [le_comb le_const] *) + + begin match LE_.Comb.as_singleton le_comb, pred with + | None, _ -> + (* non trivial linexp, give it a fresh name in the simplex *) + let proxy = fresh_term self ~pre:"_le" (T.ty t1) in + T.Tbl.replace self.needs_th_combination proxy le_comb; + + let le_comb = LE_.Comb.to_list le_comb in + List.iter (fun (_,v) -> SimpSolver.add_var self.simplex v) le_comb; + SimpSolver.define self.simplex proxy le_comb; + + let new_t = + match pred with + | Eq -> mk_eq proxy le_const + | Neq -> mk_neq proxy le_const + | Leq -> A.mk_lra tst (LRA_simplex_pred (proxy, S_op.Leq, le_const)) + | Lt -> A.mk_lra tst (LRA_simplex_pred (proxy, S_op.Lt, le_const)) + | Geq -> A.mk_lra tst (LRA_simplex_pred (proxy, S_op.Geq, le_const)) + | Gt -> A.mk_lra tst (LRA_simplex_pred (proxy, S_op.Gt, le_const)) + in + + Log.debugf 10 (fun k->k "lra.preprocess@ :%a@ :into %a" T.pp t T.pp new_t); + Some new_t + + | Some (coeff, v), Eq -> + let q = Q.(le_const / coeff) in + Some (mk_eq v q) (* turn into [c.v <= const /\ … >= ..] *) + + | Some (coeff, v), Neq -> + let q = Q.(le_const / coeff) in + Some (mk_neq v q) (* turn into [c.v < const \/ … > ..] *) + + | Some (coeff, v), pred -> + (* [c . v <= const] becomes a direct simplex constraint [v <= const/c] *) + let negate = Q.sign coeff < 0 in + let q = Q.div le_const coeff in + + let op = match pred with + | Leq -> S_op.Leq + | Lt -> S_op.Lt + | Geq -> S_op.Geq + | Gt -> S_op.Gt + | Eq | Neq -> assert false + in + let op = if negate then S_op.neg_sign op else op in + + let new_t = A.mk_lra tst (LRA_simplex_pred (v, op, q)) in + Log.debugf 10 (fun k->k "lra.preprocess@ :%a@ :into %a" T.pp t T.pp new_t); + Some new_t + end + | LRA_op _ | LRA_mult _ -> let le = as_linexp ~f:recurse t in - let proxy = fresh_term self ~pre:"_e_lra_" (T.ty t) in - T.Tbl.add self.t_defs proxy le; - T.Tbl.add self.needs_th_combination proxy le; - Log.debugf 5 (fun k->k"@[lra.preprocess.step %a@ :into %a@ :def %a@]" - T.pp t T.pp proxy LE.pp le); - Some proxy + let le_comb, le_const = LE.comb le, LE.const le in + let le_comb = LE_.Comb.to_list le_comb in + List.iter (fun (_,v) -> SimpSolver.add_var self.simplex v) le_comb; + + let proxy = fresh_term self ~pre:"_le" (T.ty t) in + + if Q.(le_const = zero) then ( + (* if there is no constant, define [proxy] as [proxy := le_comb] and + return [proxy] *) + SimpSolver.define self.simplex proxy le_comb; + Some proxy + ) else ( + (* a bit more complicated: we cannot just define [proxy := le_comb] + because of the coefficient. + Instead we assert [proxy - le_comb = le_const] using a secondary + variable [proxy2 := le_comb - proxy] + and asserting [proxy2 = -le_const] *) + let proxy2 = fresh_term self ~pre:"_le_diff" (T.ty t) in + SimpSolver.define self.simplex proxy2 + ((Q.minus_one, proxy) :: le_comb); + + add_clause [ + mk_lit (A.mk_lra tst (LRA_simplex_pred (proxy2, Leq, Q.neg le_const))) + ]; + add_clause [ + mk_lit (A.mk_lra tst (LRA_simplex_pred (proxy2, Geq, Q.neg le_const))) + ]; + + Some proxy + ) + | LRA_other t when A.has_ty_real t -> - let le = LE.monomial1 t in + let le = LE_.Comb.monomial1 t in T.Tbl.replace self.needs_th_combination t le; None - | LRA_const _ | LRA_other _ -> None - - (* ensure that [a != b] triggers the clause - [a=b \/ ab] *) - let encode_neq self si acts trail : unit = - let tst = self.tst in - begin - trail - |> Iter.filter_map - (fun lit -> - let t = Lit.term lit in - Log.debugf 50 (fun k->k "@[lra: check lit %a@ :t %a@ :sign %B@]" - Lit.pp lit T.pp t (Lit.sign lit)); - - let check_pred pred a b = - let pred = if Lit.sign lit then pred else Predicate.neg pred in - Log.debugf 50 (fun k->k "pred = `%s`" (Predicate.to_string pred)); - if pred = Neq && not (T.Tbl.mem self.neq_encoded t) then ( - Some (lit, a, b) - ) else None - in - - begin match T.Tbl.find self.pred_defs t with - | (pred, _, _, ta, tb) -> check_pred pred ta tb - | exception Not_found -> - begin match A.view_as_lra t with - | LRA_pred (pred, a, b) -> check_pred pred a b - | _ -> None - end - end) - |> Iter.iter - (fun (lit,a,b) -> - Log.debugf 50 (fun k->k "encode neq in %a" Lit.pp lit); - let c = [ - Lit.neg lit; - SI.mk_lit si acts (A.mk_lra tst (LRA_pred (Lt, a, b))); - SI.mk_lit si acts (A.mk_lra tst (LRA_pred (Lt, b, a))); - ] in - SI.add_clause_permanent si acts c; - T.Tbl.add self.neq_encoded (Lit.term (Lit.abs lit)) (); - ) - end + | LRA_const _ | LRA_simplex_pred _ | LRA_simplex_var _ | LRA_other _ -> None module Q_map = CCMap.Make(Q) - let final_check_ (self:state) si (acts:SI.actions) (trail:_ Iter.t) : unit = + (* raise conflict from certificate *) + let fail_with_cert si acts cert : 'a = + (* TODO: check certificate *) + let confl = + SimpSolver.Unsat_cert.lits cert + |> CCList.flat_map (Tag.to_lits si) + |> List.rev_map SI.Lit.neg + in + SI.raise_conflict si acts confl SI.P.default + + let check_simplex_ self si acts : SimpSolver.Subst.t = + Log.debug 5 "lra: call arith solver"; + let res = Profile.with1 "simplex.solve" SimpSolver.check self.simplex in + begin match res with + | SimpSolver.Sat m -> m + | SimpSolver.Unsat cert -> + Log.debugf 10 + (fun k->k "(@[lra.check.unsat@ :cert %a@])" + SimpSolver.Unsat_cert.pp cert); + fail_with_cert si acts cert + end + + let partial_check_ self si acts trail : unit = + Profile.with_ "lra.partial-check" @@ fun () -> + trail + (fun lit -> + let sign = SI.Lit.sign lit in + let lit_t = SI.Lit.term lit in + match A.view_as_lra lit_t with + | LRA_simplex_pred (v, op, q) -> + let op = if sign then op else S_op.neg_sign op in + let constr = SimpSolver.Constraint.mk v op q in + Log.debugf 10 + (fun k->k "(@[lra.partial-check.assert@ %a@])" + SimpSolver.Constraint.pp constr); + begin + try + SimpSolver.add_var self.simplex v; + SimpSolver.add_constraint self.simplex constr (Tag.Lit lit) + with SimpSolver.E_unsat cert -> + Log.debugf 10 + (fun k->k "(@[lra.partial-check.unsat@ :cert %a@])" + SimpSolver.Unsat_cert.pp cert); + fail_with_cert si acts cert + end + | _ -> ()); + + (* incremental check *) + ignore (check_simplex_ self si acts : SimpSolver.Subst.t); + () + + let final_check_ (self:state) si (acts:SI.actions) (_trail:_ Iter.t) : unit = Log.debug 5 "(th-lra.final-check)"; Profile.with_ "lra.final-check" @@ fun () -> - let simplex = SimpSolver.create self.gensym in - encode_neq self si acts trail; - (* first, add definitions *) - begin - T.Tbl.iter - (fun t le -> - let open LE.Infix in - let le = le - LE.monomial1 t in - let c = LConstr.eq0 le in - let lit = Tag.By_def in - SimpSolver.add_constr simplex c lit - ) - self.t_defs - end; + (* FIXME (* add congruence closure equalities *) Backtrack_stack.iter self.local_eqs ~f:(fun (n1,n2) -> @@ -333,50 +414,24 @@ module Make(A : ARG) : S with module A = A = struct let c = LConstr.eq0 LE.(t1 - t2) in let lit = Tag.CC_eq (n1,n2) in SimpSolver.add_constr simplex c lit); - (* add trail *) - begin - trail - |> Iter.iter - (fun lit -> - let sign = Lit.sign lit in - let t = Lit.term lit in - let assert_pred pred a b = - let pred = if sign then pred else Predicate.neg pred in - if pred = Neq then ( - Log.debugf 50 (fun k->k "(@[LRA.skip-neq@ :in %a@])" T.pp t); - ) else ( - let c = LConstr.of_expr LE.(a-b) pred in - SimpSolver.add_constr simplex c (Tag.Lit lit); - ) - in - begin match T.Tbl.find self.pred_defs t with - | (pred, a, b, _, _) -> assert_pred pred a b - | exception Not_found -> - begin match A.view_as_lra t with - | LRA_pred (pred, a, b) -> - let a = try T.Tbl.find self.t_defs a with _ -> as_linexp_id a in - let b = try T.Tbl.find self.t_defs b with _ -> as_linexp_id b in - assert_pred pred a b - | _ -> () - end - end) - end; + *) + Log.debug 5 "lra: call arith solver"; - let res = Profile.with1 "simplex.solve" SimpSolver.solve simplex in - begin match res with - | SimpSolver.Solution _m -> - Log.debug 5 "lra: solver returns SAT"; - let n_th_comb = - T.Tbl.keys self.needs_th_combination |> Iter.length - in - if n_th_comb > 0 then ( - Log.debugf 5 - (fun k->k "(@[LRA.needs-th-combination@ :n-lits %d@])" n_th_comb); - ); - Log.debugf 50 - (fun k->k "(@[LRA.needs-th-combination@ :lits %a@])" - (Util.pp_iter @@ Fmt.within "`" "`" T.pp) (T.Tbl.keys self.needs_th_combination)); - (* FIXME: theory combination + let model = check_simplex_ self si acts in + Log.debugf 20 (fun k->k "(@[lra.model@ %a@])" SimpSolver.Subst.pp model); + Log.debug 5 "lra: solver returns SAT"; + let n_th_comb = + T.Tbl.keys self.needs_th_combination |> Iter.length + in + if n_th_comb > 0 then ( + Log.debugf 5 + (fun k->k "(@[LRA.needs-th-combination@ :n-lits %d@])" n_th_comb); + ); + Log.debugf 50 + (fun k->k "(@[LRA.needs-th-combination@ :lits %a@])" + (Util.pp_iter @@ Fmt.within "`" "`" T.pp) (T.Tbl.keys self.needs_th_combination)); + + (* FIXME: theory combination let lazy model = model in Log.debugf 30 (fun k->k "(@[LRA.model@ %a@])" FM_A.pp_model model); @@ -418,26 +473,8 @@ module Make(A : ARG) : S with module A = A = struct by_val; () end; - *) - () - | SimpSolver.Unsatisfiable cert -> - let unsat_core = - match SimpSolver.check_cert simplex cert with - | `Ok unsat_core -> unsat_core (* TODO *) - | _ -> assert false (* some kind of fatal error ? *) - in - Log.debugf 5 (fun k->k"lra: solver returns UNSAT@ with cert %a" - (Fmt.Dump.list Tag.pp) unsat_core); - (* TODO: produce and store a proper LRA resolution proof *) - let confl = - unsat_core - |> Iter.of_list - |> Iter.flat_map_l (fun tag -> Tag.to_lits si tag) - |> Iter.map Lit.neg - |> Iter.to_rev_list - in - SI.raise_conflict si acts confl SI.P.default - end; + *) + () let create_and_setup si = @@ -446,6 +483,7 @@ module Make(A : ARG) : S with module A = A = struct (* TODO SI.add_simplifier si (simplify st); *) SI.add_preprocess si (preproc_lra st); SI.on_final_check si (final_check_ st); + SI.on_partial_check si (partial_check_ st); SI.on_cc_post_merge si (fun _ _ n1 n2 -> if A.has_ty_real (N.term n1) then ( diff --git a/src/arith/lra/simplex.ml b/src/arith/lra/simplex.ml deleted file mode 100644 index a9e66673..00000000 --- a/src/arith/lra/simplex.ml +++ /dev/null @@ -1,765 +0,0 @@ -(* - copyright (c) 2014-2018, Guillaume Bury, Simon Cruanes - *) - -(* OPTIMS: - * - distinguish separate systems (that do not interact), such as in { 1 <= 3x = 3y <= 2; z <= 3} ? - * - Implement gomorry cuts ? -*) - -module Fmt = CCFormat -module type VAR = Linear_expr_intf.VAR -module type FRESH = Linear_expr_intf.FRESH -module type VAR_GEN = Linear_expr_intf.VAR_GEN - -module type S = Simplex_intf.S -module type S_FULL = Simplex_intf.S_FULL - -module Vec = CCVector - -module Matrix : sig - type 'a t - - val create : unit -> 'a t - val get : 'a t -> int -> int -> 'a - val set : 'a t -> int -> int -> 'a -> unit - val get_row : 'a t -> int -> 'a Vec.vector - val copy : 'a t -> 'a t - val n_row : _ t -> int - val n_col : _ t -> int - val push_row : 'a t -> 'a -> unit (* new row, filled with element *) - val push_col : 'a t -> 'a -> unit (* new column, filled with element *) - - (**/**) - val check_invariants : _ t -> bool - (**/**) -end = struct - type 'a t = { - mutable n_col: int; (* num of columns *) - tab: 'a Vec.vector Vec.vector; - } - - let[@inline] create() : _ = {tab=Vec.create(); n_col=0} - - let[@inline] get m i j = Vec.get (Vec.get m.tab i) j - let[@inline] get_row m i = Vec.get m.tab i - let[@inline] set (m:_ t) i j x = Vec.set (Vec.get m.tab i) j x - let[@inline] copy m = {m with tab=Vec.map Vec.copy m.tab} - - let[@inline] n_row m = Vec.length m.tab - let[@inline] n_col m = m.n_col - - let push_row m x = Vec.push m.tab (Vec.make (n_col m) x) - let push_col m x = - m.n_col <- m.n_col + 1; - Vec.iter (fun row -> Vec.push row x) m.tab - - let check_invariants m = Vec.for_all (fun r -> Vec.length r = n_col m) m.tab -end - -(* use non-polymorphic comparison ops *) -open CCInt.Infix - -(* Simplex Implementation *) -module Make_inner - (Var: VAR) - (VMap : CCMap.S with type key=Var.t) - (Param: sig type t val copy : t -> t end) -= struct - module Var_map = VMap - module M = Var_map - - type param = Param.t - type var = Var.t - type lit = Var.lit - - type erat = { - base: Q.t; (* reference number *) - eps_factor: Q.t; (* coefficient for epsilon, the infinitesimal *) - } - - (** Epsilon-rationals, used for strict bounds *) - module Erat = struct - type t = erat - - let zero : t = {base=Q.zero; eps_factor=Q.zero} - - let[@inline] make base eps_factor : t = {base; eps_factor} - let[@inline] base t = t.base - let[@inline] eps_factor t = t.eps_factor - let[@inline] mul k e = make Q.(k * e.base) Q.(k * e.eps_factor) - let[@inline] sum e1 e2 = make Q.(e1.base + e2.base) Q.(e1.eps_factor + e2.eps_factor) - let[@inline] compare e1 e2 = match Q.compare e1.base e2.base with - | 0 -> Q.compare e1.eps_factor e2.eps_factor - | x -> x - - let lt a b = compare a b < 0 - let gt a b = compare a b > 0 - - let[@inline] min x y = if compare x y <= 0 then x else y - let[@inline] max x y = if compare x y >= 0 then x else y - - let[@inline] evaluate (epsilon:Q.t) (e:t) : Q.t = Q.(e.base + epsilon * e.eps_factor) - - let pp out e = - if Q.equal Q.zero (eps_factor e) - then Q.pp_print out (base e) - else - Fmt.fprintf out "(@[%a + @<1>ε * %a@])" - Q.pp_print (base e) Q.pp_print (eps_factor e) - end - - let str_of_var = Fmt.to_string Var.pp - let str_of_erat = Fmt.to_string Erat.pp - let str_of_q = Fmt.to_string Q.pp_print - - type bound = { - value : Erat.t; - reason : lit option; - } - - (* state associated with a variable *) - type var_state = { - var: var; - mutable assign: Erat.t option; (* current assignment *) - mutable l_bound: bound; (* lower bound *) - mutable u_bound: bound; (* upper bound *) - mutable idx_basic: int; (* index in [t.nbasic] *) - mutable idx_nbasic: int; (* index in [t.nbasic] *) - } - - (* Exceptions *) - exception Unsat of var_state - exception AbsurdBounds of var_state - exception NoneSuitable - - type basic_var = var_state - type nbasic_var = var_state - - type t = { - param: param; - mutable var_states: var_state M.t; (* var -> its state *) - tab : Q.t Matrix.t; (* the matrix of coefficients *) - basic : basic_var Vec.vector; (* basic variables *) - nbasic : nbasic_var Vec.vector; (* non basic variables *) - } - - type cert = { - cert_var: var; - cert_expr: (Q.t * var) list; - } - - type res = - | Solution of Q.t Var_map.t - | Unsatisfiable of cert - - let create param : t = { - param; - var_states = M.empty; - tab = Matrix.create (); - basic = Vec.create (); - nbasic = Vec.create (); - } - - let[@inline] index_basic (x:basic_var) : int = x.idx_basic - let[@inline] index_nbasic (x:nbasic_var) : int = x.idx_nbasic - - let[@inline] is_basic (x:var_state) : bool = x.idx_basic >= 0 - let[@inline] is_nbasic (x:var_state) : bool = x.idx_nbasic >= 0 - - (* check invariants, for test purposes *) - let check_invariants (t:t) : bool = - Matrix.check_invariants t.tab && - Vec.for_all (fun v -> is_basic v) t.basic && - Vec.for_all (fun v -> is_nbasic v) t.nbasic && - Vec.for_all (fun v -> not (is_nbasic v)) t.basic && - Vec.for_all (fun v -> not (is_basic v)) t.nbasic && - Vec.for_all (fun v -> CCOpt.is_some v.assign) t.nbasic && - Vec.for_all (fun v -> CCOpt.is_none v.assign) t.basic && - true - - (* find the definition of the basic variable [x], - as a linear combination of non basic variables *) - let find_expr_basic_opt t (x:var_state) : Q.t Vec.vector option = - begin match index_basic x with - | -1 -> None - | i -> Some (Matrix.get_row t.tab i) - end - - (* expression that defines a basic variable in terms of non-basic variables *) - let find_expr_basic t (x:basic_var) : Q.t Vec.vector = - let i = index_basic x in - assert (i >= 0); - Matrix.get_row t.tab i - - (* build the expression [y = \sum_i (if x_i=y then 1 else 0)·x_i] *) - let find_expr_nbasic t (x:nbasic_var) : Q.t Vec.vector = - Vec.map - (fun y -> if x == y then Q.one else Q.zero) - t.nbasic - - (* find expression of [x] *) - let find_expr_total (t:t) (x:var_state) : Q.t Vec.vector = - match find_expr_basic_opt t x with - | Some e -> e - | None -> - assert (is_nbasic x); - find_expr_nbasic t x - - (* compute value of basic variable. - It can be computed by using [x]'s definition - in terms of nbasic variables, which have values *) - let value_basic (t:t) (x:basic_var) : Erat.t = - assert (is_basic x); - let res = ref Erat.zero in - let expr = find_expr_basic t x in - for i = 0 to Vec.length expr - 1 do - let val_nbasic_i = - match (Vec.get t.nbasic i).assign with - | None -> assert false - | Some e -> e - in - res := Erat.sum !res (Erat.mul (Vec.get expr i) val_nbasic_i) - done; - !res - - (* extract a value for [x] *) - let[@inline] value (t:t) (x:var_state) : Erat.t = - match x.assign with - | Some e -> e - | None -> value_basic t x - - (* trivial bounds *) - let empty_bounds : bound * bound = - { value = Erat.make Q.minus_inf Q.zero; reason = None; }, - { value = Erat.make Q.inf Q.zero; reason = None; } - - (* find bounds of [x] *) - let[@inline] get_bounds (x:var_state) : bound * bound = - x.l_bound, x.u_bound - - let[@inline] get_bounds_values (x:var_state) : Erat.t * Erat.t = - let l, u = get_bounds x in - l.value, u.value - - (* is [value x] within the bounds for [x]? *) - let is_within_bounds (t:t) (x:var_state) : bool * Erat.t = - let v = value t x in - let low, upp = get_bounds_values x in - if Erat.compare v low < 0 then - false, low - else if Erat.compare v upp > 0 then - false, upp - else - true, v - - (* add [v] as a non-basic variable, or return its state if already mapped *) - let get_var_or_add_as_nbasic (t:t) (v:var) : var_state = - match M.get v t.var_states with - | Some v -> v - | None -> - let l_bound, u_bound = empty_bounds in - let idx_nbasic = Vec.length t.nbasic in - let vs = { - var=v; l_bound; u_bound; - assign=Some Erat.zero; - idx_nbasic; idx_basic=(-1); - } in - t.var_states <- M.add v vs t.var_states; - Vec.push t.nbasic vs; - Matrix.push_col t.tab Q.zero; (* new empty column *) - vs - - (* add new variables as nbasic variables, return them, ignore - the already existing variables *) - let add_vars_as_nbasic (t:t) (l:var list) : unit = - List.iter - (fun x -> - if not (M.mem x t.var_states) then ( - (* allocate new index for [x] *) - ignore (get_var_or_add_as_nbasic t x : var_state) - )) - l - - (* define basic variable [x] by [eq] in [t] *) - let add_eq (t:t) (x, eq : var * _ list) : unit = - let eq = List.map (fun (coeff,x) -> coeff, get_var_or_add_as_nbasic t x) eq in - (* add [x] as a basic var *) - begin match M.get x t.var_states with - | Some _ -> - invalid_arg (Fmt.sprintf "Variable `%a` already defined." Var.pp x); - | None -> - let l_bound, u_bound = empty_bounds in - let idx_basic = Vec.length t.basic in - let vs = { - var=x; l_bound; u_bound; assign=None; idx_basic; - idx_nbasic=(-1); - } in - Vec.push t.basic vs; - t.var_states <- M.add x vs t.var_states; - end; - (* add new row for defining [x] *) - assert (Matrix.n_col t.tab > 0); - Matrix.push_row t.tab Q.zero; - let row_i = Matrix.n_row t.tab - 1 in - assert (row_i >= 0); - - (* now put into the row the coefficients corresponding to [eq], - expanding basic variables to their definition *) - List.iter - (fun (c, x) -> - (* FIXME(perf): replace with a `idx -> Q.t` function, do not allocate vector *) - let expr = find_expr_total t x in - assert (Vec.length expr = Matrix.n_col t.tab); - Vec.iteri - (fun j c' -> - if not (Q.equal Q.zero c') then ( - Matrix.set t.tab row_i j Q.(Matrix.get t.tab row_i j + c * c') - )) - expr) - eq; - () - - (* add bounds to [x] in [t] *) - let add_bound_aux (x:var_state) - (low:Erat.t) (low_reason:lit option) - (upp:Erat.t) (upp_reason:lit option) : unit = - let l, u = get_bounds x in - let l' = if Erat.lt low l.value then l else { value = low; reason = low_reason; } in - let u' = if Erat.gt upp u.value then u else { value = upp; reason = upp_reason; } in - x.l_bound <- l'; - x.u_bound <- u'; - () - - let add_bounds (t:t) - ?strict_lower:(slow=false) ?strict_upper:(supp=false) - ?lower_reason ?upper_reason (x, l, u) : unit = - let x = get_var_or_add_as_nbasic t x in - let e1 = if slow then Q.one else Q.zero in - let e2 = if supp then Q.neg Q.one else Q.zero in - add_bound_aux x (Erat.make l e1) lower_reason (Erat.make u e2) upper_reason; - if is_nbasic x then ( - let b, v = is_within_bounds t x in - if not b then ( - x.assign <- Some v; - ) - ) - - let add_lower_bound t ?strict ~reason x l = - add_bounds t ?strict_lower:strict ~lower_reason:reason (x,l,Q.inf) - - let add_upper_bound t ?strict ~reason x u = - add_bounds t ?strict_upper:strict ~upper_reason:reason (x,Q.minus_inf,u) - - let iter_all_vars (t:t) : var_state Iter.t = - Iter.append (Vec.to_iter t.nbasic) (Vec.to_iter t.basic) - - (* full assignment *) - let full_assign (t:t) : (var * Erat.t) Iter.t = - Iter.append (Vec.to_iter t.nbasic) (Vec.to_iter t.basic) - |> Iter.map (fun x -> x.var, value t x) - - let[@inline] min x y = if Q.compare x y < 0 then x else y - - (* Find an epsilon that is small enough for finding a solution, yet - it must be positive. - - {!Erat.t} values are used to turn strict bounds ([X > 0]) into - non-strict bounds ([X >= 0 + ε]), because the simplex algorithm - only deals with non-strict bounds. - When a solution is found, we need to turn {!Erat.t} into {!Q.t} by - finding a rational value that is small enough that it will fit into - all the intervals of [t]. This rational will be the actual value of [ε]. - *) - let solve_epsilon (t:t) : Q.t = - let emax = - Iter.fold - (fun emax x -> - let { value = {base=low;eps_factor=e_low}; _} = x.l_bound in - let { value = {base=upp;eps_factor=e_upp}; _} = x.u_bound in - let {base=v; eps_factor=e_v} = value t x in - (* lower bound *) - let emax = - if Q.compare low Q.minus_inf > 0 && Q.compare e_v e_low < 0 - then min emax Q.((low - v) / (e_v - e_low)) - else emax - in - (* upper bound *) - if Q.compare upp Q.inf < 0 && Q.compare e_v e_upp > 0 - then min emax Q.((upp - v) / (e_v - e_upp)) - else emax) - Q.inf (iter_all_vars t) - in - if Q.compare emax Q.one >= 0 then Q.one else emax - - let get_full_assign_seq (t:t) : _ Iter.t = - let e = solve_epsilon t in - let f = Erat.evaluate e in - full_assign t - |> Iter.map (fun (x,v) -> x, f v) - - let get_full_assign t : Q.t Var_map.t = Var_map.of_iter (get_full_assign_seq t) - - (* Find nbasic variable suitable for pivoting with [x]. - A nbasic variable [y] is suitable if it "goes into the right direction" - (its coefficient in the definition of [x] is of the adequate sign) - and if it hasn't reached its bound in this direction. - - precondition: [x] is a basic variable whose value in current assignment - is outside its bounds - - We return the smallest (w.r.t Var.compare) suitable variable. - This is important for termination. - *) - let find_suitable_nbasic_for_pivot (t:t) (x:basic_var) : nbasic_var * Q.t = - Profile.with_ "simplex.find-pivot-var" @@ fun () -> - assert (is_basic x); - let _, v = is_within_bounds t x in - let b = Erat.compare (value t x) v < 0 in - (* is nbasic var [y], with coeff [a] in definition of [x], suitable? *) - let test (y:nbasic_var) (a:Q.t) : bool = - assert (is_nbasic y); - let v = value t y in - let low, upp = get_bounds_values y in - if b then ( - (Erat.lt v upp && Q.compare a Q.zero > 0) || - (Erat.gt v low && Q.compare a Q.zero < 0) - ) else ( - (Erat.gt v low && Q.compare a Q.zero > 0) || - (Erat.lt v upp && Q.compare a Q.zero < 0) - ) - in - let nbasic_vars = t.nbasic in - let expr = find_expr_basic t x in - (* find best suitable variable *) - let rec aux i = - if i = Vec.length nbasic_vars then ( - assert (i = Vec.length expr); - None - ) else ( - let y = Vec.get nbasic_vars i in - let a = Vec.get expr i in - if test y a then ( - (* see if other variables are better suited *) - begin match aux (i+1) with - | None -> Some (y,a) - | Some (z, _) as res_tail -> - if Var.compare y.var z.var <= 0 - then Some (y,a) - else res_tail - end - ) else ( - aux (i+1) - ) - ) - in - begin match aux 0 with - | Some res -> res - | None -> raise NoneSuitable - end - - (* pivot to exchange [x] and [y] *) - let pivot (t:t) (x:basic_var) (y:nbasic_var) (a:Q.t) : unit = - Profile.with_ "simplex.pivot" @@ fun () -> - (* swap values ([x] becomes assigned) *) - let val_x = value t x in - y.assign <- None; - x.assign <- Some val_x; - (* Matrix Pivot operation *) - let kx = index_basic x in - let ky = index_nbasic y in - for j = 0 to Vec.length t.nbasic - 1 do - if y == Vec.get t.nbasic j then ( - Matrix.set t.tab kx j Q.(inv a) - ) else ( - Matrix.set t.tab kx j Q.(neg (Matrix.get t.tab kx j) / a) - ) - done; - for i = 0 to Vec.length t.basic - 1 do - if i <> kx then ( - let c = Matrix.get t.tab i ky in - Matrix.set t.tab i ky Q.zero; - for j = 0 to Vec.length t.nbasic - 1 do - Matrix.set t.tab i j Q.(Matrix.get t.tab i j + c * Matrix.get t.tab kx j) - done - ) - done; - (* Switch x and y in basic and nbasic vars *) - Vec.set t.basic kx y; - Vec.set t.nbasic ky x; - x.idx_basic <- -1; - y.idx_basic <- kx; - x.idx_nbasic <- ky; - y.idx_nbasic <- -1; - () - - (* find minimum element of [arr] (wrt [cmp]) that satisfies predicate [f] *) - let find_min_filter ~cmp (f:'a -> bool) (arr:('a,_) Vec.t) : 'a option = - (* find the first element that satisfies [f] *) - let rec aux_find_first i = - if i = Vec.length arr then None - else ( - let x = Vec.get arr i in - if f x - then aux_compare_with x (i+1) - else aux_find_first (i+1) - ) - (* find if any element of [l] satisfies [f] and is smaller than [x] *) - and aux_compare_with x i = - if i = Vec.length arr then Some x - else ( - let y = Vec.get arr i in - let best = if f y && cmp y x < 0 then y else x in - aux_compare_with best (i+1) - ) - in - aux_find_first 0 - - (* check bounds *) - let check_bounds (t:t) : unit = - iter_all_vars t - (fun x -> - let l = x.l_bound in - let u = x.u_bound in - if Erat.gt l.value u.value then raise (AbsurdBounds x)) - - let[@inline] compare_by_var x y = Var.compare x.var y.var - - (* actual solving algorithm *) - let solve_aux (t:t) : unit = - Profile.instant - (Printf.sprintf "(simplex.solve :basic %d :non-basic %d)" - (Vec.length t.basic) (Vec.length t.nbasic)); - check_bounds t; - (* select the smallest basic variable that is not satisfied in the current - assignment. *) - let rec aux_select_basic_var () = - match - Profile.with_ "simplex.select-basic-var" @@ fun () -> - find_min_filter ~cmp:compare_by_var - (fun x -> not (fst (is_within_bounds t x))) - t.basic - with - | Some x -> aux_pivot_on_basic x - | None -> () - (* remove the basic variable *) - and aux_pivot_on_basic x = - let _b, v = is_within_bounds t x in - assert (not _b); - match find_suitable_nbasic_for_pivot t x with - | y, a -> - (* exchange [x] and [y] by pivoting *) - pivot t x y a; - (* assign [x], now a nbasic variable, to the faulty bound [v] *) - x.assign <- Some v; - (* next iteration *) - aux_select_basic_var () - | exception NoneSuitable -> - raise (Unsat x) - in - aux_select_basic_var (); - () - - (* main method for the user to call *) - let solve (t:t) : res = - try - solve_aux t; - Solution (get_full_assign t) - with - | Unsat x -> - let cert_expr = - List.combine - (Vec.to_list (find_expr_basic t x)) - (Vec.to_list t.nbasic |> CCList.map (fun x -> x.var)) - in - Unsatisfiable { cert_var=x.var; cert_expr; } (* FIXME *) - | AbsurdBounds x -> - Unsatisfiable { cert_var=x.var; cert_expr=[]; } - - (* add [c·x] to [m] *) - let add_expr_ (x:var) (c:Q.t) (m:Q.t M.t) = - let c' = M.get_or ~default:Q.zero x m in - let c' = Q.(c + c') in - if Q.equal Q.zero c' then M.remove x m else M.add x c' m - - (* dereference basic variables from [c·x], and add the result to [m] *) - let rec deref_var_ t x c m = - match find_expr_basic_opt t x with - | None -> add_expr_ x.var c m - | Some expr_x -> - let m = ref m in - Vec.iteri - (fun i c_i -> - let y_i = Vec.get t.nbasic i in - m := deref_var_ t y_i Q.(c * c_i) !m) - expr_x; - !m - - (* maybe invert bounds, if [c < 0] *) - let scale_bounds c (l,u) : bound * bound = - match Q.compare c Q.zero with - | 0 -> - let b = { value = Erat.zero; reason = None; } in - b, b - | n when n<0 -> - { u with value = Erat.mul c u.value; }, - { l with value = Erat.mul c l.value; } - | _ -> - { l with value = Erat.mul c l.value; }, - { u with value = Erat.mul c u.value; } - - let add_to_unsat_core acc = function - | None -> acc - | Some reason -> reason :: acc - - let check_cert (t:t) (c:cert) = - let x = M.get c.cert_var t.var_states |> CCOpt.get_lazy (fun()->assert false) in - let { value = low_x; reason = low_x_reason; } = x.l_bound in - let { value = up_x; reason = upp_x_reason; } = x.u_bound in - begin match c.cert_expr with - | [] -> - if Erat.compare low_x up_x > 0 - then `Ok (add_to_unsat_core (add_to_unsat_core [] low_x_reason) upp_x_reason) - else `Bad_bounds (str_of_erat low_x, str_of_erat up_x) - | expr -> - let e0 = deref_var_ t x (Q.neg Q.one) M.empty in - (* compute bounds for the expression [c.cert_expr], - and also compute [c.cert_expr - x] to check if it's 0] *) - let low, low_unsat_core, up, up_unsat_core, expr_minus_x = - List.fold_left - (fun (l, luc, u, uuc, expr_minus_x) (c, y) -> - let y = M.get y t.var_states |> CCOpt.get_lazy (fun ()->assert false) in - let ly, uy = scale_bounds c (get_bounds y) in - assert (Erat.compare ly.value uy.value <= 0); - let expr_minus_x = deref_var_ t y c expr_minus_x in - let luc = add_to_unsat_core luc ly.reason in - let uuc = add_to_unsat_core uuc uy.reason in - Erat.sum l ly.value, luc, Erat.sum u uy.value, uuc, expr_minus_x) - (Erat.zero, [], Erat.zero, [], e0) - expr - in - (* check that the expanded expression is [x], and that - one of the bounds on [x] is incompatible with bounds of [c.cert_expr] *) - if M.is_empty expr_minus_x then ( - if Erat.compare low_x up > 0 - then `Ok (add_to_unsat_core up_unsat_core low_x_reason) - else if Erat.compare up_x low < 0 - then `Ok (add_to_unsat_core low_unsat_core upp_x_reason) - else `Bad_bounds (str_of_erat low, str_of_erat up) - ) else `Diff_not_0 expr_minus_x - end - - (* printer *) - - let matrix_pp_width = ref 8 - - let fmt_head = format_of_string "|%*s|| " - let fmt_cell = format_of_string "%*s| " - - let pp_cert out (c:cert) = match c.cert_expr with - | [] -> Fmt.fprintf out "(@[inconsistent-bounds %a@])" Var.pp c.cert_var - | _ -> - let pp_pair = Fmt.(hvbox ~i:2 @@ pair ~sep:(return "@ * ") Q.pp_print Var.pp) in - Fmt.fprintf out "(@[cert@ :var %a@ :linexp %a@])" - Var.pp c.cert_var - Fmt.(within "[" "]" @@ hvbox @@ list ~sep:(return "@ + ") pp_pair) - c.cert_expr - - let pp_mat out t = - let open Fmt in - fprintf out "@["; - (* header *) - fprintf out fmt_head !matrix_pp_width ""; - Vec.iter (fun x -> fprintf out fmt_cell !matrix_pp_width (str_of_var x.var)) t.nbasic; - fprintf out "@,"; - (* rows *) - for i=0 to Matrix.n_row t.tab-1 do - if i>0 then fprintf out "@,"; - let v = Vec.get t.basic i in - fprintf out fmt_head !matrix_pp_width (str_of_var v.var); - let row = Matrix.get_row t.tab i in - Vec.iter (fun q -> fprintf out fmt_cell !matrix_pp_width (str_of_q q)) row; - done; - fprintf out "@]" - - let pp_vars = - let ppv out v = - Fmt.fprintf out "(@[var %a@ :assign %a@ :lbound %a@ :ubound %a@])" - Var.pp v.var (Fmt.Dump.option Erat.pp) v.assign - Erat.pp v.l_bound.value Erat.pp v.u_bound.value - in - Fmt.(within "(" ")" @@ hvbox @@ iter ppv) - - let pp_full_state out (t:t) : unit = - (* print main matrix *) - Fmt.fprintf out - "(@[simplex@ :n-row %d :n-col %d@ :mat %a@ :vars %a @])" - (Matrix.n_row t.tab) (Matrix.n_col t.tab) pp_mat t - pp_vars (iter_all_vars t) -end - -module Make(Var:VAR) = - Make_inner(Var)(CCMap.Make(Var))(struct - type t = unit - let copy ()=() - end) - -module Make_full_for_expr(V : VAR_GEN) - (L : Linear_expr.S - with type Var.t = V.t - and type C.t = Q.t - and type Var.lit = V.lit) - : S_FULL with type var = V.t - and type lit = V.lit - and module L = L - and module Var_map = L.Var_map - and type L.var = V.t - and type L.Comb.t = L.Comb.t - and type param = V.Fresh.t -= struct - include Make_inner(V)(L.Var_map)(V.Fresh) - module L = L - - type op = Predicate.t = Leq | Geq | Lt | Gt | Eq | Neq - - type constr = L.Constr.t - - (* add a constraint *) - let add_constr (t:t) (c:constr) (reason:lit) : unit = - let e, op, q = L.Constr.split c in - match L.Comb.as_singleton e with - | Some (c0, x0) -> - (* no need for a fresh variable, just add constraint on [x0] *) - let q = Q.div q c0 in - let op = if Q.sign c0 < 0 then Predicate.neg_sign op else op in - begin match op with - | Leq -> add_upper_bound t ~strict:false ~reason x0 q - | Geq -> add_lower_bound t ~strict:false ~reason x0 q - | Lt -> add_upper_bound t ~strict:true ~reason x0 q - | Gt -> add_lower_bound t ~strict:true ~reason x0 q - | Eq -> add_bounds t (x0,q,q) - ~strict_lower:false ~strict_upper:false - ~lower_reason:reason ~upper_reason:reason - | Neq -> assert false - end - | None -> - let (x:var) = V.Fresh.fresh t.param in - add_eq t (x, L.Comb.to_list e); - begin match op with - | Leq -> add_upper_bound t ~strict:false ~reason x q - | Geq -> add_lower_bound t ~strict:false ~reason x q - | Lt -> add_upper_bound t ~strict:true ~reason x q - | Gt -> add_lower_bound t ~strict:true ~reason x q - | Eq -> add_bounds t (x,q,q) - ~strict_lower:false ~strict_upper:false - ~lower_reason:reason ~upper_reason:reason - | Neq -> assert false - end -end - -module Make_full(V : VAR_GEN) - : S_FULL with type var = V.t - and type lit = V.lit - and type L.var = V.t - and type param = V.Fresh.t - = Make_full_for_expr(V)(Linear_expr.Make(struct include Q let pp = pp_print end)(V)) diff --git a/src/arith/lra/simplex.mli b/src/arith/lra/simplex.mli deleted file mode 100644 index 8bbcef5d..00000000 --- a/src/arith/lra/simplex.mli +++ /dev/null @@ -1,33 +0,0 @@ - -(** Solving Linear systems of rational equations. *) - -module type VAR = Linear_expr_intf.VAR -module type FRESH = Linear_expr_intf.FRESH -module type VAR_GEN = Linear_expr_intf.VAR_GEN - -module type S = Simplex_intf.S -module type S_FULL = Simplex_intf.S_FULL - -(** Low level simplex interface *) -module Make(V : VAR) : - S with type var = V.t - and type lit = V.lit - and type param = unit - and module Var_map = CCMap.Make(V) - -(** High-level simplex interface *) -module Make_full_for_expr(V : VAR_GEN) - (L : Linear_expr.S with type Var.t = V.t and type Var.lit = V.lit and type C.t = Q.t) - : S_FULL with type var = V.t - and type lit = V.lit - and module L = L - and module Var_map = L.Var_map - and type L.var = V.t - and type L.Comb.t = L.Comb.t - and type param = V.Fresh.t - -module Make_full(V : VAR_GEN) - : S_FULL with type var = V.t - and type lit = V.lit - and type L.var = V.t - and type param = V.Fresh.t diff --git a/src/arith/tests/run_tests.ml b/src/arith/tests/run_tests.ml index 2a766dd7..2fd6a4f3 100644 --- a/src/arith/tests/run_tests.ml +++ b/src/arith/tests/run_tests.ml @@ -6,7 +6,6 @@ let tests : unit Alcotest.test list = [ let props = List.flatten [ Test_simplex2.props; - Test_simplex.props; ] let () = diff --git a/src/arith/tests/test_simplex.ml b/src/arith/tests/test_simplex.ml deleted file mode 100644 index 44c34523..00000000 --- a/src/arith/tests/test_simplex.ml +++ /dev/null @@ -1,180 +0,0 @@ -module Fmt = CCFormat -module QC = QCheck - -module Var = struct - include CCInt - - let pp out x = Format.fprintf out "X_%d" x - - let rand n : t QC.arbitrary = QC.make ~print:(Fmt.to_string pp) @@ QC.Gen.(0--n) - type lit = int - let pp_lit = Fmt.int - - module Fresh = struct - type var = t - type t = int ref - let copy r = ref !r - let create() = ref ~-1 - let fresh r = decr r; !r - end -end - -module L = Sidekick_arith_lra.Linear_expr.Make(struct include Q let pp=pp_print end)(Var) -module Spl = Sidekick_arith_lra.Simplex.Make_full_for_expr(Var)(L) -module Var_map = Spl.Var_map - -let rand_n low n : Z.t QC.arbitrary = - QC.map ~rev:Z.to_int Z.of_int QC.(low -- n) - -let rand_q : Q.t QC.arbitrary = - let n1 = rand_n ~-100_000 100_000 in - let n2 = rand_n 1 40_000 in - let qc = - QC.map ~rev:(fun q -> Q.num q, Q.den q) - (fun (x,y) -> Q.make x y) - (QC.pair n1 n2) - in - (* avoid [undef] when shrinking *) - let shrink q yield = - CCOpt.get_exn qc.QC.shrink q (fun x -> if Q.is_real x then yield x) - in - QC.set_shrink shrink qc - -type subst = Spl.L.subst - -(* NOTE: should arrive in qcheck at some point *) -let filter_shrink (f:'a->bool) (a:'a QC.arbitrary) : 'a QC.arbitrary = - match a.QC.shrink with - | None -> a - | Some shr -> QC.set_shrink (QC.Shrink.filter f shr) a - -module Comb = struct - include Spl.L.Comb - - let rand n : t QC.arbitrary = - let a = - QC.map_same_type (fun e -> if is_empty e then monomial1 0 else e) @@ - QC.map ~rev:to_list of_list @@ - QC.list_of_size QC.Gen.(1--n) @@ QC.pair rand_q (Var.rand 10) - in - filter_shrink (fun e -> not (is_empty e)) a -end - -module Expr = struct - include Spl.L.Expr - - let rand n : t QC.arbitrary = - QC.map ~rev:(fun e->comb e, const e) (CCFun.uncurry make) @@ - QC.pair (Comb.rand n) rand_q -end - -module Constr = struct - include Spl.L.Constr - - let shrink c : t QC.Iter.t = - let open QC.Iter in - CCOpt.map_or ~default:empty - (fun s -> s c.expr >|= fun expr -> {c with expr}) - (Expr.rand 5).QC.shrink - - let rand n : t QC.arbitrary = - let gen = - QC.Gen.( - return of_expr <*> - (Expr.rand n).QC.gen <*> - oneofl Sidekick_arith_lra.Predicate.([Leq;Geq;Lt;Gt;Eq]) - ) - in - QC.make ~print:(Fmt.to_string pp) ~shrink gen -end - -module Problem = struct - type t = Constr.t list - - module Infix = struct - let (&&) = List.rev_append - end - include Infix - - let eval subst = List.for_all (L.Constr.eval subst) - - let pp out pb = Fmt.(hvbox @@ list ~sep:(return "@ @<1>∧ ") L.Constr.pp) out pb - - let rand ?min:(m=3) n : t QC.arbitrary = - let n = max m (max n 6) in - QC.list_of_size QC.Gen.(m -- n) @@ Constr.rand 10 -end - -let add_problem (t:Spl.t) (pb:Problem.t) : unit = - let lit = 0 in - List.iter (fun constr -> Spl.add_constr t constr lit) pb - -let pp_subst : subst Fmt.printer = - Fmt.(map Spl.L.Var_map.to_iter @@ - within "{" "}" @@ hvbox @@ iter ~sep:(return ",@ ") @@ - pair ~sep:(return "@ @<1>→ ") Var.pp Q.pp_print - ) - -let check_invariants = - let prop pb = - let simplex = Spl.create (Var.Fresh.create()) in - add_problem simplex pb; - Spl.check_invariants simplex - in - QC.Test.make ~long_factor:10 ~count:50 ~name:"simplex_invariants" (Problem.rand 20) prop - -let check_invariants_after_solve = - let prop pb = - let simplex = Spl.create (Var.Fresh.create()) in - add_problem simplex pb; - ignore (Spl.solve simplex); - if Spl.check_invariants simplex then true - else ( - QC.Test.fail_reportf "(@[bad-invariants@ %a@])" Spl.pp_full_state simplex - ) - in - QC.Test.make ~long_factor:10 ~count:50 ~name:"simplex_invariants_after_solve" (Problem.rand 20) prop - -let check_sound = - let prop pb = - let simplex = Spl.create (Var.Fresh.create()) in - add_problem simplex pb; - begin match Spl.solve simplex with - | Spl.Solution subst -> - if Problem.eval subst pb then true - else ( - QC.Test.fail_reportf - "(@[bad-solution@ :problem %a@ :sol %a@ :simplex %a@])" - Problem.pp pb pp_subst subst Spl.pp_full_state simplex - ) - | Spl.Unsatisfiable cert -> - begin match Spl.check_cert simplex cert with - | `Ok _ -> true - | `Bad_bounds (low, up) -> - QC.Test.fail_reportf - "(@[bad-certificat@ :problem %a@ :cert %a@ :low %s :up %s@ :simplex %a@])" - Problem.pp pb Spl.pp_cert cert low up Spl.pp_full_state simplex - | `Diff_not_0 e -> - QC.Test.fail_reportf - "(@[bad-certificat@ :problem %a@ :cert %a@ :diff %a@ :simplex %a@])" - Problem.pp pb Spl.pp_cert cert Comb.pp (Comb.of_map e) Spl.pp_full_state simplex - end - end - in - QC.Test.make ~long_factor:10 ~count:300 ~name:"simplex_sound" (Problem.rand 20) prop - -let check_scalable = - let prop pb = - let simplex = Spl.create (Var.Fresh.create()) in - add_problem simplex pb; - ignore (Spl.solve simplex); - true - in - QC.Test.make ~long_factor:2 ~count:10 ~name:"simplex_scalable" (Problem.rand ~min:150 150) prop - - -let props = [ - check_invariants; - check_sound; - check_scalable; -] diff --git a/src/smtlib/Process.ml b/src/smtlib/Process.ml index a6e72523..fbecaa78 100644 --- a/src/smtlib/Process.ml +++ b/src/smtlib/Process.ml @@ -312,6 +312,9 @@ module Th_lra = Sidekick_arith_lra.Make(struct type term = S.T.Term.t type ty = S.T.Ty.t + let mk_and = Form.and_ + let mk_or = Form.or_ + let mk_lra = T.lra let view_as_lra t = match T.view t with | T.LRA l -> l