[micromega] More pre-procesing

- Remove obviously redundant constraints
- Perform (partial) Fourier elimination to detect (easy) cutting-planes

Closes #13227

6 files changed, 501 insertions, 177 deletions

diff --git a/plugins/micromega/certificate.ml b/plugins/micromega/certificate.ml
index 9008691bca..dc36372d3d 100644
--- a/plugins/micromega/certificate.ml
+++ b/plugins/micromega/certificate.ml
@@ -385,6 +385,16 @@ let subst sys =
       sys';
   sys'
 
+let tr_sys str f sys =
+  let sys' = f sys in
+  if debug then (
+    Printf.fprintf stdout "[%s\n" str;
+    List.iter (fun s -> Printf.fprintf stdout "%a\n" WithProof.output s) sys;
+    Printf.fprintf stdout "\n => \n";
+    List.iter (fun s -> Printf.fprintf stdout "%a\n" WithProof.output s) sys';
+    Printf.fprintf stdout "]\n" );
+  sys'
+
 (** [saturate_linear_equality sys] generate new constraints
     obtained by eliminating linear equalities by pivoting.
     For integers, the obtained constraints are  sound but not complete.
@@ -392,11 +402,7 @@ let subst sys =
 let saturate_by_linear_equalities sys0 = WithProof.saturate_subst false sys0
 
 let saturate_by_linear_equalities sys =
-  let sys' = saturate_by_linear_equalities sys in
-  if debug then
-    Printf.fprintf stdout "[saturate_by_linear_equalities:\n%a\n==>\n%a\n]"
-      output_sys sys output_sys sys';
-  sys'
+  tr_sys "saturate_by_linear_equalities" saturate_by_linear_equalities sys
 
 let bound_monomials (sys : WithProof.t list) =
   let l =
@@ -497,10 +503,10 @@ let nlinear_prover prfdepth sys =
   let sys = List.map (fun ((p, o), prf) -> (cstr_of_poly (p, o), prf)) sys in
   let id =
     List.fold_left
-      (fun acc (_, r) -> max acc (ProofFormat.pr_rule_max_id r))
+      (fun acc (_, r) -> max acc (ProofFormat.pr_rule_max_hyp r))
       0 sys
   in
-  let env = CList.interval 0 id in
+  let env = List.map (fun i -> ProofFormat.Hyp i) (CList.interval 0 id) in
   match linear_prover_cstr sys with
   | None -> Unknown
   | Some cert -> Prf (ProofFormat.cmpl_prf_rule Mc.normQ CamlToCoq.q env cert)
@@ -514,7 +520,7 @@ let linear_prover_with_cert prfdepth sys =
   | Some cert ->
     Prf
       (ProofFormat.cmpl_prf_rule Mc.normQ CamlToCoq.q
-         (List.mapi (fun i e -> i) sys)
+         (List.mapi (fun i e -> ProofFormat.Hyp i) sys)
          cert)
 
 (* The prover is (probably) incomplete --
@@ -885,6 +891,11 @@ let check_sys sys =
 
 open ProofFormat
 
+let output_cstr_sys sys =
+  (pp_list ";" (fun o (c, wp) ->
+       Printf.fprintf o "%a by %a" output_cstr c ProofFormat.output_prf_rule wp))
+    sys
+
 let xlia (can_enum : bool) reduction_equations sys =
   let rec enum_proof (id : int) (sys : prf_sys) =
     if debug then (
@@ -922,16 +933,10 @@ let xlia (can_enum : bool) reduction_equations sys =
         | _ -> Unknown )
   and aux_lia (id : int) (sys : prf_sys) =
     assert (check_sys sys);
-    if debug then
-      Printf.printf "xlia:  %a \n"
-        (pp_list ";" (fun o (c, _) -> output_cstr o c))
-        sys;
+    if debug then Printf.printf "xlia:  %a \n" output_cstr_sys sys;
     try
       let sys = reduction_equations sys in
-      if debug then
-        Printf.printf "after reduction:  %a \n"
-          (pp_list ";" (fun o (c, _) -> output_cstr o c))
-          sys;
+      if debug then Printf.printf "after reduction:  %a \n" output_cstr_sys sys;
       match linear_prover_cstr sys with
       | Some prf -> Prf (Step (id, prf, Done))
       | None -> if can_enum then enum_proof id sys else Unknown
@@ -943,7 +948,7 @@ let xlia (can_enum : bool) reduction_equations sys =
   let id =
     1
     + List.fold_left
-        (fun acc (_, r) -> max acc (ProofFormat.pr_rule_max_id r))
+        (fun acc (_, r) -> max acc (ProofFormat.pr_rule_max_hyp r))
         0 sys
   in
   let orpf =
@@ -973,7 +978,7 @@ let xlia_simplex env red sys =
     let id =
       1
       + List.fold_left
-          (fun acc (_, r) -> max acc (ProofFormat.pr_rule_max_id r))
+          (fun acc (_, r) -> max acc (ProofFormat.pr_rule_max_hyp r))
           0 sys
     in
     let env = CList.interval 0 (id - 1) in
@@ -1007,6 +1012,143 @@ let gen_bench (tac, prover) can_enum prfdepth sys =
     flush o; close_out o );
   res
 
+let normalise sys =
+  List.fold_left
+    (fun acc s ->
+      match WithProof.cutting_plane s with
+      | None -> s :: acc
+      | Some s' -> s' :: acc)
+    [] sys
+
+let normalise = tr_sys "normalise" normalise
+
+let elim_redundant sys =
+  let module VectMap = Map.Make (Vect) in
+  let elim_eq sys =
+    List.fold_left
+      (fun acc (((v, o), prf) as wp) ->
+        match o with
+        | Gt -> assert false
+        | Ge -> wp :: acc
+        | Eq -> wp :: WithProof.neg wp :: acc)
+      [] sys
+  in
+  let of_list l =
+    List.fold_left
+      (fun m (((v, o), prf) as wp) ->
+        let q, v' = Vect.decomp_cst v in
+        try
+          let q', wp' = VectMap.find v' m in
+          match Q.compare q q' with
+          | 0 -> if o = Eq then VectMap.add v' (q, wp) m else m
+          | 1 -> m
+          | _ -> VectMap.add v' (q, wp) m
+        with Not_found -> VectMap.add v' (q, wp) m)
+      VectMap.empty l
+  in
+  let to_list m = VectMap.fold (fun _ (_, wp) sys -> wp :: sys) m [] in
+  to_list (of_list (elim_eq sys))
+
+let elim_redundant sys = tr_sys "elim_redundant" elim_redundant sys
+
+(** [fourier_small] performs some variable elimination and keeps the cutting planes.
+    To decide which elimination to perform, the constraints are sorted according to
+    1 - the number of variables
+    2 - the value of the smallest coefficient
+    Given the smallest constraint, we eliminate the variable with the smallest coefficient.
+    The rational is that a constraint with a single variable provides some bound information.
+    When there are several variables, we hope to eliminate all the variables.
+    A necessary condition is to take the variable with the smallest coefficient *)
+
+let fourier_small (sys : WithProof.t list) =
+  let size ((p, o), prf) =
+    let _, p' = Vect.decomp_cst p in
+    let (x, q), p' = Vect.decomp_fst p' in
+    Vect.fold
+      (fun (l, (q, x)) x' q' ->
+        let q' = Q.abs q' in
+        (l + 1, if q </ q then (q, x) else (q', x')))
+      (1, (Q.abs q, x))
+      p
+  in
+  let cmp ((l1, (q1, _)), ((_, o), _)) ((l2, (q2, _)), ((_, o'), _)) =
+    if l1 < l2 then -1 else if l1 = l2 then Q.compare q1 q2 else 1
+  in
+  let syss = List.sort cmp (List.rev_map (fun wp -> (size wp, wp)) sys) in
+  let gen_pivot acc (q, x) wp l =
+    List.fold_left
+      (fun acc (s, wp') ->
+        match WithProof.simple_pivot (q, x) wp wp' with
+        | None -> acc
+        | Some wp2 -> (
+          match WithProof.cutting_plane wp2 with
+          | Some wp2 -> (s, wp2) :: acc
+          | _ -> acc ))
+      acc l
+  in
+  let rec all_pivots acc l =
+    match l with
+    | [] -> acc
+    | ((_, qx), wp) :: l' -> all_pivots (gen_pivot acc qx wp (acc @ l')) l'
+  in
+  List.rev_map snd (all_pivots [] syss)
+
+let fourier_small = tr_sys "fourier_small" fourier_small
+
+(** [propagate_bounds sys] generate new constraints by exploiting bounds.
+    A bound is a constraint of the form c + a.x >= 0
+ *)
+
+(*let propagate_bounds sys =
+  let bounds, sys' =
+    List.fold_left
+      (fun (b, r) (((c, o), prf) as wp) ->
+        match Vect.Bound.of_vect c with
+        | None -> (b, wp :: r)
+        | Some b' -> ((b', wp) :: b, r))
+      ([], []) sys
+  in
+  let exploit_bound acc (b, wp) =
+    let cf = b.Vect.Bound.coeff in
+    let vr = b.Vect.Bound.var in
+    List.fold_left
+      (fun acc (((c, o), prf) as wp') ->
+        let cf' = Vect.get vr c in
+        if Q.sign (cf */ cf') = -1 then
+          WithProof.(
+            let wf2 =
+              addition
+                (mult (LinPoly.constant (Q.abs cf')) wp)
+                (mult (LinPoly.constant (Q.abs cf)) wp')
+            in
+            match cutting_plane wf2 with None -> acc | Some cp -> cp :: acc)
+        else acc)
+      acc sys'
+  in
+  List.fold_left exploit_bound [] bounds
+ *)
+
+let rev_concat l =
+  let rec conc acc l =
+    match l with [] -> acc | l1 :: lr -> conc (List.rev_append l1 acc) lr
+  in
+  conc [] l
+
+
+let pre_process sys =
+  let sys = normalise sys in
+  let bnd1 = bound_monomials sys in
+  let sys1 = normalise (subst sys) in
+  let pbnd1 = fourier_small sys1 in
+  let sys2 = elim_redundant (List.rev_append pbnd1 sys1) in
+  let bnd2 = bound_monomials sys2 in
+  let pbnd2 = [] (*fourier_small sys2*) in
+  (* Should iterate ? *)
+  let sys =
+    rev_concat [pbnd2; bnd1; bnd2; saturate_by_linear_equalities sys2; sys2]
+  in
+  sys
+
 let lia (can_enum : bool) (prfdepth : int) sys =
   let sys = develop_constraints prfdepth z_spec sys in
   if debug then begin
@@ -1020,11 +1162,7 @@ let lia (can_enum : bool) (prfdepth : int) sys =
           p)
       sys
   end;
-  let bnd1 = bound_monomials sys in
-  let sys = subst sys in
-  let bnd2 = bound_monomials sys in
-  (* To deal with non-linear monomials *)
-  let sys = bnd1 @ bnd2 @ saturate_by_linear_equalities sys @ sys in
+  let sys = pre_process sys in
   let sys' = List.map (fun ((p, o), prf) -> (cstr_of_poly (p, o), prf)) sys in
   xlia (List.map fst sys) can_enum reduction_equations sys'
 
@@ -1039,7 +1177,8 @@ let nlia enum prfdepth sys =
     List.iter (fun s -> Printf.fprintf stdout "%a\n" WithProof.output s) sys
   end;
   if is_linear then
-    xlia (List.map fst sys) enum reduction_equations (make_cstr_system sys)
+    xlia (List.map fst sys) enum reduction_equations
+      (make_cstr_system (pre_process sys))
   else
     (*
       let sys1 = elim_every_substitution sys in
diff --git a/plugins/micromega/polynomial.ml b/plugins/micromega/polynomial.ml
index 5c0aa9ef0d..b0317ef643 100644
--- a/plugins/micromega/polynomial.ml
+++ b/plugins/micromega/polynomial.ml
@@ -254,6 +254,16 @@ let is_strict c = c.op = Gt
 let eval_op = function Eq -> ( =/ ) | Ge -> ( >=/ ) | Gt -> ( >/ )
 let string_of_op = function Eq -> "=" | Ge -> ">=" | Gt -> ">"
 
+let compare_op o1 o2 =
+  match (o1, o2) with
+  | Eq, Eq -> 0
+  | Eq, _ -> -1
+  | _, Eq -> 1
+  | Ge, Ge -> 0
+  | Ge, _ -> -1
+  | _, Ge -> 1
+  | Gt, Gt -> 0
+
 let output_cstr o {coeffs; op; cst} =
   Printf.fprintf o "%a %s %s" Vect.pp coeffs (string_of_op op) (Q.to_string cst)
 
@@ -284,7 +294,11 @@ module LinPoly = struct
       if !fresh > vr then failwith (Printf.sprintf "Cannot reserve %i" vr)
       else fresh := vr + 1
 
-    let get_fresh () = !fresh
+    let safe_reserve vr = if !fresh > vr then () else fresh := vr + 1
+
+    let get_fresh () =
+      let vr = !fresh in
+      incr fresh; vr
 
     let register m =
       try MonoMap.find m !index_of_monomial
@@ -489,23 +503,35 @@ module ProofFormat = struct
     | CutPrf p -> pr_size p
     | MulC (v, p) -> pr_size p
 
-  let rec pr_rule_max_id = function
-    | Annot (_, p) -> pr_rule_max_id p
-    | Hyp i | Def i -> i
+  let rec pr_rule_max_hyp = function
+    | Annot (_, p) -> pr_rule_max_hyp p
+    | Hyp i -> i
+    | Def i -> -1
+    | Cst _ | Zero | Square _ -> -1
+    | MulC (_, p) | CutPrf p | Gcd (_, p) -> pr_rule_max_hyp p
+    | MulPrf (p1, p2) | AddPrf (p1, p2) ->
+      max (pr_rule_max_hyp p1) (pr_rule_max_hyp p2)
+
+  let rec pr_rule_max_def = function
+    | Annot (_, p) -> pr_rule_max_hyp p
+    | Hyp i -> -1
+    | Def i -> i
     | Cst _ | Zero | Square _ -> -1
-    | MulC (_, p) | CutPrf p | Gcd (_, p) -> pr_rule_max_id p
+    | MulC (_, p) | CutPrf p | Gcd (_, p) -> pr_rule_max_def p
     | MulPrf (p1, p2) | AddPrf (p1, p2) ->
-      max (pr_rule_max_id p1) (pr_rule_max_id p2)
+      max (pr_rule_max_def p1) (pr_rule_max_def p2)
 
-  let rec proof_max_id = function
+  let rec proof_max_def = function
     | Done -> -1
-    | Step (i, pr, prf) -> max i (max (pr_rule_max_id pr) (proof_max_id prf))
+    | Step (i, pr, prf) -> max i (max (pr_rule_max_def pr) (proof_max_def prf))
     | Enum (i, p1, _, p2, l) ->
-      let m = max (pr_rule_max_id p1) (pr_rule_max_id p2) in
-      List.fold_left (fun i prf -> max i (proof_max_id prf)) (max i m) l
+      let m = max (pr_rule_max_def p1) (pr_rule_max_def p2) in
+      List.fold_left (fun i prf -> max i (proof_max_def prf)) (max i m) l
     | ExProof (i, j, k, _, _, _, prf) ->
-      max (max (max i j) k) (proof_max_id prf)
+      max (max (max i j) k) (proof_max_def prf)
 
+  (** [pr_rule_def_cut id pr] gives an explicit [id] to cut rules.
+      This is because the Coq proof format only accept they as a proof-step *)
   let rec pr_rule_def_cut id = function
     | Annot (_, p) -> pr_rule_def_cut id p
     | MulC (p, prf) ->
@@ -536,33 +562,35 @@ module ProofFormat = struct
 
   let rec implicit_cut p = match p with CutPrf p -> implicit_cut p | _ -> p
 
-  let rec pr_rule_collect_hyps pr =
+  let rec pr_rule_collect_defs pr =
     match pr with
-    | Annot (_, pr) -> pr_rule_collect_hyps pr
-    | Hyp i | Def i -> ISet.add i ISet.empty
+    | Annot (_, pr) -> pr_rule_collect_defs pr
+    | Def i -> ISet.add i ISet.empty
+    | Hyp i -> ISet.empty
     | Cst _ | Zero | Square _ -> ISet.empty
-    | MulC (_, pr) | Gcd (_, pr) | CutPrf pr -> pr_rule_collect_hyps pr
+    | MulC (_, pr) | Gcd (_, pr) | CutPrf pr -> pr_rule_collect_defs pr
     | MulPrf (p1, p2) | AddPrf (p1, p2) ->
-      ISet.union (pr_rule_collect_hyps p1) (pr_rule_collect_hyps p2)
+      ISet.union (pr_rule_collect_defs p1) (pr_rule_collect_defs p2)
 
+  (** [simplify_proof p] removes proof steps that are never re-used. *)
   let simplify_proof p =
     let rec simplify_proof p =
       match p with
       | Done -> (Done, ISet.empty)
-      | Step (i, pr, Done) -> (p, ISet.add i (pr_rule_collect_hyps pr))
+      | Step (i, pr, Done) -> (p, ISet.add i (pr_rule_collect_defs pr))
       | Step (i, pr, prf) ->
         let prf', hyps = simplify_proof prf in
         if not (ISet.mem i hyps) then (prf', hyps)
         else
           ( Step (i, pr, prf')
-          , ISet.add i (ISet.union (pr_rule_collect_hyps pr) hyps) )
+          , ISet.add i (ISet.union (pr_rule_collect_defs pr) hyps) )
       | Enum (i, p1, v, p2, pl) ->
         let pl, hl = List.split (List.map simplify_proof pl) in
         let hyps = List.fold_left ISet.union ISet.empty hl in
         ( Enum (i, p1, v, p2, pl)
         , ISet.add i
             (ISet.union
-               (ISet.union (pr_rule_collect_hyps p1) (pr_rule_collect_hyps p2))
+               (ISet.union (pr_rule_collect_defs p1) (pr_rule_collect_defs p2))
                hyps) )
       | ExProof (i, j, k, x, z, t, prf) ->
         let prf', hyps = simplify_proof prf in
@@ -652,9 +680,9 @@ module ProofFormat = struct
       | Gcd (b1, p1), Gcd (b2, p2) ->
         cmp_pair Z.compare compare (b1, p1) (b2, p2)
       | MulPrf (p1, q1), MulPrf (p2, q2) ->
-        cmp_pair compare compare (p1, q1) (p2, q2)
-      | AddPrf (p1, q1), MulPrf (p2, q2) ->
-        cmp_pair compare compare (p1, q1) (p2, q2)
+        cmp_pair compare compare (p1, p2) (q1, q2)
+      | AddPrf (p1, q1), AddPrf (p2, q2) ->
+        cmp_pair compare compare (p1, p2) (q1, q2)
       | CutPrf p, CutPrf p' -> compare p p'
       | _, _ -> Int.compare (id_of_constr p1) (id_of_constr p2)
   end
@@ -746,16 +774,23 @@ module ProofFormat = struct
       Zero vect
 
   module Env = struct
-    let rec string_of_int_list l =
+    let output_hyp_or_def o = function
+      | Hyp i -> Printf.fprintf o "Hyp %i" i
+      | Def i -> Printf.fprintf o "Def %i" i
+      | _ -> ()
+
+    let rec output_hyps o l =
       match l with
-      | [] -> ""
-      | i :: l -> Printf.sprintf "%i,%s" i (string_of_int_list l)
+      | [] -> ()
+      | i :: l -> Printf.fprintf o "%a,%a" output_hyp_or_def i output_hyps l
 
     let id_of_hyp hyp l =
       let rec xid_of_hyp i l' =
         match l' with
         | [] ->
-          failwith (Printf.sprintf "id_of_hyp %i %s" hyp (string_of_int_list l))
+          Printf.fprintf stdout "\nid_of_hyp: %a notin [%a]\n" output_hyp_or_def
+            hyp output_hyps l;
+          failwith "Cannot find hyp or def"
         | hyp' :: l' -> if hyp = hyp' then i else xid_of_hyp (i + 1) l'
       in
       xid_of_hyp 0 l
@@ -764,7 +799,7 @@ module ProofFormat = struct
   let cmpl_prf_rule norm (cst : Q.t -> 'a) env prf =
     let rec cmpl = function
       | Annot (s, p) -> cmpl p
-      | Hyp i | Def i -> Mc.PsatzIn (CamlToCoq.nat (Env.id_of_hyp i env))
+      | (Hyp _ | Def _) as h -> Mc.PsatzIn (CamlToCoq.nat (Env.id_of_hyp h env))
       | Cst i -> Mc.PsatzC (cst i)
       | Zero -> Mc.PsatzZ
       | MulPrf (p1, p2) -> Mc.PsatzMulE (cmpl p1, cmpl p2)
@@ -785,20 +820,22 @@ module ProofFormat = struct
     | Step (i, p, prf) -> (
       match p with
       | CutPrf p' ->
-        Mc.CutProof (cmpl_prf_rule_z env p', cmpl_proof (i :: env) prf)
-      | _ -> Mc.RatProof (cmpl_prf_rule_z env p, cmpl_proof (i :: env) prf) )
+        Mc.CutProof (cmpl_prf_rule_z env p', cmpl_proof (Def i :: env) prf)
+      | _ -> Mc.RatProof (cmpl_prf_rule_z env p, cmpl_proof (Def i :: env) prf)
+      )
     | Enum (i, p1, _, p2, l) ->
       Mc.EnumProof
         ( cmpl_prf_rule_z env p1
         , cmpl_prf_rule_z env p2
-        , List.map (cmpl_proof (i :: env)) l )
+        , List.map (cmpl_proof (Def i :: env)) l )
     | ExProof (i, j, k, x, _, _, prf) ->
-      Mc.ExProof (CamlToCoq.positive x, cmpl_proof (i :: j :: k :: env) prf)
+      Mc.ExProof
+        (CamlToCoq.positive x, cmpl_proof (Def i :: Def j :: Def k :: env) prf)
 
   let compile_proof env prf =
-    let id = 1 + proof_max_id prf in
+    let id = 1 + proof_max_def prf in
     let _, prf = normalise_proof id prf in
-    cmpl_proof env prf
+    cmpl_proof (List.map (fun i -> Hyp i) env) prf
 
   let rec eval_prf_rule env = function
     | Annot (s, p) -> eval_prf_rule env p
@@ -863,7 +900,7 @@ module WithProof = struct
   let compare : t -> t -> int =
    fun ((lp1, o1), _) ((lp2, o2), _) ->
     let c = Vect.compare lp1 lp2 in
-    if c = 0 then compare o1 o2 else c
+    if c = 0 then compare_op o1 o2 else c
 
   let annot s (p, prf) = (p, ProofFormat.Annot (s, prf))
 
@@ -887,6 +924,13 @@ module WithProof = struct
    fun ((p1, o1), prf1) ((p2, o2), prf2) ->
     ((Vect.add p1 p2, opAdd o1 o2), ProofFormat.add_proof prf1 prf2)
 
+  let neg : t -> t =
+   fun ((p1, o1), prf1) ->
+    match o1 with
+    | Eq ->
+      ((Vect.mul Q.minus_one p1, o1), ProofFormat.mul_cst_proof Q.minus_one prf1)
+    | _ -> failwith "neg: invalid proof"
+
   let mult p ((p1, o1), prf1) =
     match o1 with
     | Eq -> ((LinPoly.product p p1, o1), ProofFormat.sMulC p prf1)
@@ -912,13 +956,13 @@ module WithProof = struct
       else
         match o with
         | Eq ->
-          Some ((Vect.set 0 Q.minus_one Vect.null, Eq), ProofFormat.Gcd (g, prf))
+          Some ((Vect.set 0 Q.minus_one Vect.null, Eq), ProofFormat.CutPrf prf)
         | Gt -> failwith "cutting_plane ignore strict constraints"
         | Ge ->
           (* This is a non-trivial common divisor *)
           Some
             ( (Vect.set 0 c1' (Vect.div (Q.of_bigint g) p), o)
-            , ProofFormat.Gcd (g, prf) )
+            , ProofFormat.CutPrf prf )
 
   let construct_sign p =
     let c, p' = Vect.decomp_cst p in
@@ -1029,6 +1073,26 @@ module WithProof = struct
     in
     saturate select gen sys0
 
+  let simple_pivot (q1, x) ((v1, o1), prf1) ((v2, o2), prf2) =
+    let q2 = Vect.get x v2 in
+    if q2 =/ Q.zero then None
+    else
+      let cv1, cv2 =
+        if Q.sign q1 <> Q.sign q2 then (Q.abs q2, Q.abs q1)
+        else
+          match (o1, o2) with
+          | Eq, _ -> (q2, Q.abs q1)
+          | _, Eq -> (Q.abs q2, q2)
+          | _, _ -> (Q.zero, Q.zero)
+      in
+      if cv2 =/ Q.zero then None
+      else
+        Some
+          ( (Vect.mul_add cv1 v1 cv2 v2, opAdd o1 o2)
+          , ProofFormat.add_proof
+              (ProofFormat.mul_cst_proof cv1 prf1)
+              (ProofFormat.mul_cst_proof cv2 prf2) )
+
   open Vect.Bound
 
   let mul_bound w1 w2 =
diff --git a/plugins/micromega/polynomial.mli b/plugins/micromega/polynomial.mli
index 9c09f76691..2402a03e0d 100644
--- a/plugins/micromega/polynomial.mli
+++ b/plugins/micromega/polynomial.mli
@@ -120,6 +120,7 @@ type cstr = {coeffs : Vect.t; op : op; cst : Q.t}
 and op = Eq | Ge | Gt
 
 val eval_op : op -> Q.t -> Q.t -> bool
+val compare_op : op -> op -> int
 
 (*val opMult : op -> op -> op*)
 
@@ -153,6 +154,9 @@ module LinPoly : sig
     (** [reserve i] reserves the integer i *)
     val reserve : int -> unit
 
+    (** [safe_reserve i] reserves the integer i *)
+    val safe_reserve : int -> unit
+
     (** [get_fresh ()] return the first fresh variable *)
     val get_fresh : unit -> int
 
@@ -289,8 +293,9 @@ module ProofFormat : sig
   (* x = z - t, z >= 0, t >= 0 *)
 
   val pr_size : prf_rule -> Q.t
-  val pr_rule_max_id : prf_rule -> int
-  val proof_max_id : proof -> int
+  val pr_rule_max_def : prf_rule -> int
+  val pr_rule_max_hyp : prf_rule -> int
+  val proof_max_def : proof -> int
   val normalise_proof : int -> proof -> int * proof
   val output_prf_rule : out_channel -> prf_rule -> unit
   val output_proof : out_channel -> proof -> unit
@@ -302,13 +307,15 @@ module ProofFormat : sig
   val cmpl_prf_rule :
        ('a Micromega.pExpr -> 'a Micromega.pol)
     -> (Q.t -> 'a)
-    -> int list
+    -> prf_rule list
     -> prf_rule
     -> 'a Micromega.psatz
 
   val proof_of_farkas : prf_rule IMap.t -> Vect.t -> prf_rule
   val eval_prf_rule : (int -> LinPoly.t * op) -> prf_rule -> LinPoly.t * op
   val eval_proof : (LinPoly.t * op) IMap.t -> proof -> bool
+
+  module PrfRuleMap : Map.S with type key = prf_rule
 end
 
 val output_cstr : out_channel -> cstr -> unit
@@ -344,6 +351,12 @@ module WithProof : sig
       @return the polynomial p+q with its sign and proof *)
   val addition : t -> t -> t
 
+  (** [neg p]
+      @return the polynomial -p with its sign and proof
+      @raise an error if this not an equality
+ *)
+  val neg : t -> t
+
   (** [mult p q]
       @return the polynomial p*q with its sign and proof.
       @raise InvalidProof if p is not a constant and p  is not an equality *)
@@ -360,6 +373,10 @@ module WithProof : sig
    *)
   val linear_pivot : t list -> t -> Vect.var -> t -> t option
 
+  (** [simple_pivot (c,x) p q]  performs a pivoting over the variable [x] where
+      p = c+a1.x1+....+c.x+...an.xn and c <> 0 *)
+  val simple_pivot : Q.t * var -> t -> t -> t option
+
   (** [subst sys] performs the equivalent of the 'subst' tactic of Coq.
     For every p=0 \in sys such that p is linear in x with coefficient +/- 1
                                i.e. p = 0 <-> x = e and x \notin e.
diff --git a/plugins/micromega/simplex.ml b/plugins/micromega/simplex.ml
index bd47ccae8b..be9b990fe1 100644
--- a/plugins/micromega/simplex.ml
+++ b/plugins/micromega/simplex.ml
@@ -60,6 +60,77 @@ let get_profile_info () =
       ( try (p.success_pivots + p.failure_pivots) / p.average_pivots
         with Division_by_zero -> 0 ) }
 
+(* SMT output for debugging *)
+
+(*
+let pp_smt_row o (k, v) =
+  Printf.fprintf o "(assert (= x%i %a))\n" k Vect.pp_smt v
+
+let pp_smt_assert_tbl o tbl = IMap.iter (fun k v -> pp_smt_row o (k, v)) tbl
+
+let pp_smt_goal_tbl o tbl =
+  let pp_rows o tbl =
+    IMap.iter (fun k v -> Printf.fprintf o "(= x%i %a)" k Vect.pp_smt v) tbl
+  in
+  Printf.fprintf o "(assert (not (and %a)))\n" pp_rows tbl
+
+let pp_smt_vars s o var =
+  ISet.iter
+    (fun i ->
+      Printf.fprintf o "(declare-const x%i %s);%a\n" i s LinPoly.pp_var i)
+    (ISet.remove 0 var)
+
+let pp_smt_goal s o tbl1 tbl2 =
+  let set_of_row vr v = ISet.add vr (Vect.variables v) in
+  let var =
+    IMap.fold (fun k v acc -> ISet.union (set_of_row k v) acc) tbl1 ISet.empty
+  in
+  Printf.fprintf o "(echo \"%s\")\n(push) %a  %a %a (check-sat) (pop)\n" s
+    (pp_smt_vars "Real") var pp_smt_assert_tbl tbl1 pp_smt_goal_tbl tbl2;
+  flush stdout
+
+let pp_smt_cut o lp c =
+  let var =
+    ISet.remove 0
+      (List.fold_left
+         (fun acc ((c, o), _) -> ISet.union (Vect.variables c) acc)
+         ISet.empty lp)
+  in
+  let pp_list o l =
+    List.iter
+      (fun ((c, _), _) -> Printf.fprintf o "(assert (>= %a 0))\n" Vect.pp_smt c)
+      l
+  in
+  Printf.fprintf o
+    "(push) \n\
+     (echo \"new cut\")\n\
+     %a %a (assert (not (>= %a 0)))\n\
+     (check-sat) (pop)\n"
+    (pp_smt_vars "Int") var pp_list lp Vect.pp_smt c
+
+let pp_smt_sat o lp sol =
+  let var =
+    ISet.remove 0
+      (List.fold_left
+         (fun acc ((c, o), _) -> ISet.union (Vect.variables c) acc)
+         ISet.empty lp)
+  in
+  let pp_list o l =
+    List.iter
+      (fun ((c, _), _) -> Printf.fprintf o "(assert (>= %a 0))\n" Vect.pp_smt c)
+      l
+  in
+  let pp_model o v =
+    Vect.fold
+      (fun () v x ->
+        Printf.fprintf o "(assert (= x%i %a))\n" v Vect.pp_smt (Vect.cst x))
+      () v
+  in
+  Printf.fprintf o
+    "(push) \n(echo \"check base\")\n%a %a %a\n(check-sat) (pop)\n"
+    (pp_smt_vars "Real") var pp_list lp pp_model sol
+ *)
+
 type iset = unit IMap.t
 
 (** Mapping basic variables to their equation.
@@ -375,38 +446,6 @@ open Polynomial
 
 (*type varmap = (int * bool) IMap.t*)
 
-let make_certificate vm l =
-  Vect.normalise
-    (Vect.fold
-       (fun acc x n ->
-         let x', b = IMap.find x vm in
-         Vect.set x' (if b then n else Q.neg n) acc)
-       Vect.null l)
-
-(** [eliminate_equalities vr0 l]
-    represents an equality e = 0 of index idx in the list l
-    by 2 constraints (vr:e >= 0) and (vr+1:-e >= 0)
-    The mapping vm maps vr to idx
- *)
-
-let eliminate_equalities (vr0 : var) (l : Polynomial.cstr list) =
-  let rec elim idx vr vm l acc =
-    match l with
-    | [] -> (vr, vm, acc)
-    | c :: l -> (
-      match c.op with
-      | Ge ->
-        let v = Vect.set 0 (Q.neg c.cst) c.coeffs in
-        elim (idx + 1) (vr + 1) (IMap.add vr (idx, true) vm) l ((vr, v) :: acc)
-      | Eq ->
-        let v1 = Vect.set 0 (Q.neg c.cst) c.coeffs in
-        let v2 = Vect.mul Q.minus_one v1 in
-        let vm = IMap.add vr (idx, true) (IMap.add (vr + 1) (idx, false) vm) in
-        elim (idx + 1) (vr + 2) vm l ((vr, v1) :: (vr + 1, v2) :: acc)
-      | Gt -> raise Strict )
-  in
-  elim 0 vr0 IMap.empty l []
-
 let find_solution rst tbl =
   IMap.fold
     (fun vr v res ->
@@ -440,19 +479,9 @@ let rec solve opt l (rst : Restricted.t) (t : tableau) =
   | Some ((vr, v), l) -> (
     match push_real opt vr v (Restricted.set_exc vr rst) t with
     | Sat (t', x) -> (
-      (*        let t' = remove_redundant rst t' in*)
-      match l with
-      | [] -> Inl (rst, t', x)
-      | _ -> solve opt l rst t' )
+      match l with [] -> Inl (rst, t', x) | _ -> solve opt l rst t' )
     | Unsat c -> Inr c )
 
-let find_unsat_certificate (l : Polynomial.cstr list) =
-  let vr = LinPoly.MonT.get_fresh () in
-  let _, vm, l' = eliminate_equalities vr l in
-  match solve false l' (Restricted.make vr) IMap.empty with
-  | Inr c -> Some (make_certificate vm c)
-  | Inl _ -> None
-
 let fresh_var l =
   1
   +
@@ -463,64 +492,110 @@ let fresh_var l =
          ISet.empty l)
   with Not_found -> 0
 
+module PrfEnv = struct
+  type t = WithProof.t IMap.t
+
+  let empty = IMap.empty
+
+  let register prf env =
+    let fr = LinPoly.MonT.get_fresh () in
+    (fr, IMap.add fr prf env)
+
+  (* let register_def (v, op) {fresh; env} =
+     LinPoly.MonT.reserve fresh;
+     (fresh, {fresh = fresh + 1; env = IMap.add fresh ((v, op), Def fresh) env}) *)
+
+  let set_prf i prf env = IMap.add i prf env
+  let find idx env = IMap.find idx env
+
+  let rec of_list acc env l =
+    match l with
+    | [] -> (acc, env)
+    | (((lp, op), prf) as wp) :: l -> (
+      match op with
+      | Gt -> raise Strict (* Should be eliminated earlier *)
+      | Ge ->
+        (* Simply register *)
+        let f, env' = register wp env in
+        of_list ((f, lp) :: acc) env' l
+      | Eq ->
+        (* Generate two constraints *)
+        let f1, env = register wp env in
+        let wp' = WithProof.neg wp in
+        let f2, env = register wp' env in
+        of_list ((f1, lp) :: (f2, fst (fst wp')) :: acc) env l )
+
+  let map f env = IMap.map f env
+end
+
+let make_env (l : Polynomial.cstr list) =
+  PrfEnv.of_list [] PrfEnv.empty
+    (List.rev_map WithProof.of_cstr
+       (List.mapi (fun i x -> (x, ProofFormat.Hyp i)) l))
+
 let find_point (l : Polynomial.cstr list) =
   let vr = fresh_var l in
-  let _, vm, l' = eliminate_equalities vr l in
+  LinPoly.MonT.safe_reserve vr;
+  let l', _ = make_env l in
   match solve false l' (Restricted.make vr) IMap.empty with
   | Inl (rst, t, _) -> Some (find_solution rst t)
   | _ -> None
 
 let optimise obj l =
-  let vr0 = LinPoly.MonT.get_fresh () in
-  let _, vm, l' = eliminate_equalities (vr0 + 1) l in
+  let vr = fresh_var l in
+  LinPoly.MonT.safe_reserve vr;
+  let l', _ = make_env l in
   let bound pos res =
     match res with
     | Opt (_, Max n) -> Some (if pos then n else Q.neg n)
     | Opt (_, Ubnd _) -> None
     | Opt (_, Feas) -> None
   in
-  match solve false l' (Restricted.make vr0) IMap.empty with
+  match solve false l' (Restricted.make vr) IMap.empty with
   | Inl (rst, t, _) ->
     Some
-      ( bound false (simplex true vr0 rst (add_row vr0 t (Vect.uminus obj)))
-      , bound true (simplex true vr0 rst (add_row vr0 t obj)) )
+      ( bound false (simplex true vr rst (add_row vr t (Vect.uminus obj)))
+      , bound true (simplex true vr rst (add_row vr t obj)) )
   | _ -> None
 
-open Polynomial
+(** [make_certificate env l] makes very strong assumptions
+    about the form of the environment.
+    Each proof is assumed to be either:
+    - an hypothesis Hyp i
+    - or, the negation of an hypothesis (MulC(-1,Hyp i))
+ *)
+
+let make_certificate env l =
+  Vect.normalise
+    (Vect.fold
+       (fun acc x n ->
+         let _, prf = PrfEnv.find x env in
+         ProofFormat.(
+           match prf with
+           | Hyp i -> Vect.set i n acc
+           | MulC (_, Hyp i) -> Vect.set i (Q.neg n) acc
+           | _ -> failwith "make_certificate: invalid proof"))
+       Vect.null l)
 
-let env_of_list l =
-  List.fold_left (fun (i, m) l -> (i + 1, IMap.add i l m)) (0, IMap.empty) l
+let find_unsat_certificate (l : Polynomial.cstr list) =
+  let l', env = make_env l in
+  let vr = fresh_var l in
+  match solve false l' (Restricted.make vr) IMap.empty with
+  | Inr c -> Some (make_certificate env c)
+  | Inl _ -> None
 
+open Polynomial
 open ProofFormat
 
-let make_farkas_certificate (env : WithProof.t IMap.t) vm v =
+let make_farkas_certificate (env : PrfEnv.t) v =
   Vect.fold
-    (fun acc x n ->
-      add_proof acc
-        begin
-          try
-            let x', b = IMap.find x vm in
-            mul_cst_proof (if b then n else Q.neg n) (snd (IMap.find x' env))
-          with Not_found ->
-            (* This is an introduced hypothesis *)
-            mul_cst_proof n (snd (IMap.find x env))
-        end)
+    (fun acc x n -> add_proof acc (mul_cst_proof n (snd (PrfEnv.find x env))))
     Zero v
 
-let make_farkas_proof (env : WithProof.t IMap.t) vm v =
+let make_farkas_proof (env : PrfEnv.t) v =
   Vect.fold
     (fun wp x n ->
-      WithProof.addition wp
-        begin
-          try
-            let x', b = IMap.find x vm in
-            let n = if b then n else Q.neg n in
-            let prf = IMap.find x' env in
-            WithProof.mult (Vect.cst n) prf
-          with Not_found ->
-            let prf = IMap.find x env in
-            WithProof.mult (Vect.cst n) prf
-        end)
+      WithProof.addition wp (WithProof.mult (Vect.cst n) (PrfEnv.find x env)))
     WithProof.zero v
 
 let frac_num n = n -/ Q.floor n
@@ -532,7 +607,13 @@ type ('a, 'b) hitkind =
   (* Yes, we have a positive result *)
   | Keep of 'b
 
-let cut env rmin sol vm (rst : Restricted.t) tbl (x, v) =
+let violation sol vect =
+  let sol = Vect.set 0 Q.one sol in
+  let c = Vect.get 0 vect in
+  if Q.zero =/ c then Vect.dotproduct sol vect
+  else Q.abs (Vect.dotproduct sol vect // c)
+
+let cut env rmin sol (rst : Restricted.t) tbl (x, v) =
   let n, r = Vect.decomp_cst v in
   let fn = frac_num (Q.abs n) in
   if fn =/ Q.zero then Forget (* The solution is integral *)
@@ -580,7 +661,7 @@ let cut env rmin sol vm (rst : Restricted.t) tbl (x, v) =
     in
     let lcut =
       ( fst ccoeff
-      , make_farkas_proof env vm (Vect.normalise (cut_vector (snd ccoeff))) )
+      , make_farkas_proof env (Vect.normalise (cut_vector (snd ccoeff))) )
     in
     let check_cutting_plane (p, c) =
       match WithProof.cutting_plane c with
@@ -592,7 +673,9 @@ let cut env rmin sol vm (rst : Restricted.t) tbl (x, v) =
       | Some (v, prf) ->
         if debug then (
           Printf.printf "%s: This is a cutting plane for %a:" p LinPoly.pp_var x;
-          Printf.printf " %a\n" WithProof.output (v, prf) );
+          Printf.printf "(viol %f) %a\n"
+            (Q.to_float (violation sol (fst v)))
+            WithProof.output (v, prf) );
         Some (x, (v, prf))
     in
     match check_cutting_plane lcut with
@@ -621,30 +704,40 @@ let merge_best lt oldr newr =
   | Forget, Keep v -> Keep v
   | Keep v, Keep v' -> Keep v'
 
-let find_cut nb env u sol vm rst tbl =
+(*let size_vect v =
+  let abs z = if Z.compare z Z.zero < 0 then Z.neg z else z in
+  Vect.fold
+    (fun acc _ q -> Z.add (abs (Q.num q)) (Z.add (Q.den q) acc))
+    Z.zero v
+ *)
+
+let find_cut nb env u sol rst tbl =
   if nb = 0 then
     IMap.fold
-      (fun x v acc -> merge_result_old acc (cut env u sol vm rst tbl) (x, v))
+      (fun x v acc -> merge_result_old acc (cut env u sol rst tbl) (x, v))
       tbl Forget
   else
-    let lt (_, (_, p1)) (_, (_, p2)) =
+    let lt (_, ((v1, _), p1)) (_, ((v2, _), p2)) =
+      (*violation sol v1 >/ violation sol v2*)
       ProofFormat.pr_size p1 </ ProofFormat.pr_size p2
     in
     IMap.fold
-      (fun x v acc -> merge_best lt acc (cut env u sol vm rst tbl (x, v)))
+      (fun x v acc -> merge_best lt acc (cut env u sol rst tbl (x, v)))
       tbl Forget
 
 let var_of_vect v = fst (fst (Vect.decomp_fst v))
 
-let eliminate_variable (bounded, vr, env, tbl) x =
+let eliminate_variable (bounded, env, tbl) x =
   if debug then
     Printf.printf "Eliminating variable %a from tableau\n%a\n" LinPoly.pp_var x
       output_tableau tbl;
   (* We identify the new variables with the constraint. *)
-  LinPoly.MonT.reserve vr;
-  let z = LinPoly.var (vr + 1) in
+  let vr = LinPoly.MonT.get_fresh () in
+  let vr1 = LinPoly.MonT.get_fresh () in
+  let vr2 = LinPoly.MonT.get_fresh () in
+  let z = LinPoly.var vr1 in
   let zv = var_of_vect z in
-  let t = LinPoly.var (vr + 2) in
+  let t = LinPoly.var vr2 in
   let tv = var_of_vect t in
   (* x = z - t *)
   let xdef = Vect.add z (Vect.uminus t) in
@@ -653,9 +746,9 @@ let eliminate_variable (bounded, vr, env, tbl) x =
   let tp = ((t, Ge), Def tv) in
   (* Pivot the current tableau using xdef *)
   let tbl = IMap.map (fun v -> Vect.subst x xdef v) tbl in
-  (* Pivot the environment *)
+  (* Pivot the proof environment *)
   let env =
-    IMap.map
+    PrfEnv.map
       (fun lp ->
         let (v, o), p = lp in
         let ai = Vect.get x v in
@@ -664,67 +757,73 @@ let eliminate_variable (bounded, vr, env, tbl) x =
       env
   in
   (* Add the variables to the environment *)
-  let env = IMap.add vr xp (IMap.add zv zp (IMap.add tv tp env)) in
+  let env =
+    PrfEnv.set_prf vr xp (PrfEnv.set_prf zv zp (PrfEnv.set_prf tv tp env))
+  in
   (* Remember the mapping *)
   let bounded = IMap.add x (vr, zv, tv) bounded in
   if debug then (
     Printf.printf "Tableau without\n %a\n" output_tableau tbl;
     Printf.printf "Environment\n %a\n" output_env env );
-  (bounded, vr + 3, env, tbl)
+  (bounded, env, tbl)
 
 let integer_solver lp =
-  let l, _ = List.split lp in
-  let vr0 = 3 * LinPoly.MonT.get_fresh () in
-  let vr, vm, l' = eliminate_equalities vr0 l in
-  let _, env = env_of_list (List.map WithProof.of_cstr lp) in
   let insert_row vr v rst tbl =
     match push_real true vr v rst tbl with
-    | Sat (t', x) -> Inl (Restricted.restrict vr rst, t', x)
+    | Sat (t', x) ->
+      (*pp_smt_goal stdout tbl vr v t';*)
+      Inl (Restricted.restrict vr rst, t', x)
     | Unsat c -> Inr c
   in
+  let vr0 = LinPoly.MonT.get_fresh () in
+  (* Initialise the proof environment mapping variables of the simplex to their proof. *)
+  let l', env =
+    PrfEnv.of_list [] PrfEnv.empty (List.rev_map WithProof.of_cstr lp)
+  in
   let nb = ref 0 in
-  let rec isolve env cr vr res =
+  let rec isolve env cr res =
     incr nb;
     match res with
     | Inr c ->
-      Some (Step (vr, make_farkas_certificate env vm (Vect.normalise c), Done))
+      Some
+        (Step
+           ( LinPoly.MonT.get_fresh ()
+           , make_farkas_certificate env (Vect.normalise c)
+           , Done ))
     | Inl (rst, tbl, x) -> (
       if debug then begin
         Printf.fprintf stdout "Looking for a cut\n";
         Printf.fprintf stdout "Restricted %a ...\n" Restricted.pp rst;
         Printf.fprintf stdout "Current Tableau\n%a\n" output_tableau tbl;
         flush stdout
-        (*           Printf.fprintf stdout "Bounding box\n%a\n" output_box (bounding_box (vr+1) rst tbl l')*)
       end;
       let sol = find_full_solution rst tbl in
-      match find_cut (!nb mod 2) env cr (*x*) sol vm rst tbl with
+      match find_cut (!nb mod 2) env cr (*x*) sol rst tbl with
       | Forget ->
         None (* There is no hope, there should be an integer solution *)
-      | Hit (cr, ((v, op), cut)) ->
+      | Hit (cr, ((v, op), cut)) -> (
+        let vr = LinPoly.MonT.get_fresh () in
         if op = Eq then
           (* This is a contradiction *)
           Some (Step (vr, CutPrf cut, Done))
-        else (
-          LinPoly.MonT.reserve vr;
+        else
           let res = insert_row vr v (Restricted.set_exc vr rst) tbl in
-          let prf =
-            isolve (IMap.add vr ((v, op), Def vr) env) (Some cr) (vr + 1) res
-          in
+          let prf = isolve (IMap.add vr ((v, op), Def vr) env) (Some cr) res in
           match prf with
           | None -> None
           | Some p -> Some (Step (vr, CutPrf cut, p)) )
       | Keep (x, v) -> (
         if debug then
           Printf.fprintf stdout "Remove %a from Tableau\n" LinPoly.pp_var x;
-        let bounded, vr, env, tbl =
+        let bounded, env, tbl =
           Vect.fold
             (fun acc x n ->
               if x <> 0 && not (Restricted.is_restricted x rst) then
                 eliminate_variable acc x
               else acc)
-            (IMap.empty, vr, env, tbl) v
+            (IMap.empty, env, tbl) v
         in
-        let prf = isolve env cr vr (Inl (rst, tbl, None)) in
+        let prf = isolve env cr (Inl (rst, tbl, None)) in
         match prf with
         | None -> None
         | Some pf ->
@@ -734,7 +833,7 @@ let integer_solver lp =
                bounded pf) ) )
   in
   let res = solve true l' (Restricted.make vr0) IMap.empty in
-  isolve env None vr res
+  isolve env None res
 
 let integer_solver lp =
   nb_pivot := 0;
diff --git a/plugins/micromega/vect.ml b/plugins/micromega/vect.ml
index 4df32f2ba4..fe1d721b89 100644
--- a/plugins/micromega/vect.ml
+++ b/plugins/micromega/vect.ml
@@ -57,12 +57,17 @@ let pp_var_num pp_var o {var = v; coe = n} =
   else Printf.fprintf o "%s*%a" (Q.to_string n) pp_var v
 
 let pp_var_num_smt pp_var o {var = v; coe = n} =
-  if Int.equal v 0 then
-    if Q.zero =/ n then () else Printf.fprintf o "%s" (Q.to_string n)
+  let pp_num o q =
+    let nn = Q.num n in
+    let dn = Q.den n in
+    if Z.equal dn Z.one then output_string o (Z.to_string nn)
+    else Printf.fprintf o "(/ %s %s)" (Z.to_string nn) (Z.to_string dn)
+  in
+  if Int.equal v 0 then if Q.zero =/ n then () else pp_num o n
   else if Q.one =/ n then pp_var o v
   else if Q.minus_one =/ n then Printf.fprintf o "(- %a)" pp_var v
   else if Q.zero =/ n then ()
-  else Printf.fprintf o "(* %s %a)" (Q.to_string n) pp_var v
+  else Printf.fprintf o "(* %a %a)" pp_num n pp_var v
 
 let rec pp_gen pp_var o v =
   match v with
diff --git a/plugins/micromega/vect.mli b/plugins/micromega/vect.mli
index 9db6c075f8..b4742430fa 100644
--- a/plugins/micromega/vect.mli
+++ b/plugins/micromega/vect.mli
@@ -56,8 +56,8 @@ val get_cst : t -> Q.t
 (** [decomp_cst v] returns the pair (c,a1.x1+...+an.xn) *)
 val decomp_cst : t -> Q.t * t
 
-(** [decomp_cst v] returns the pair (ai, ai+1.xi+...+an.xn) *)
-val decomp_at : int -> t -> Q.t * t
+(** [decomp_at xi v] returns the pair (ai, ai+1.xi+...+an.xn) *)
+val decomp_at : var -> t -> Q.t * t
 
 val decomp_fst : t -> (var * Q.t) * t