summaryrefslogtreecommitdiff
path: root/plugins/omega/omega.ml
diff options
context:
space:
mode:
Diffstat (limited to 'plugins/omega/omega.ml')
-rw-r--r--plugins/omega/omega.ml716
1 files changed, 716 insertions, 0 deletions
diff --git a/plugins/omega/omega.ml b/plugins/omega/omega.ml
new file mode 100644
index 00000000..11ab9c03
--- /dev/null
+++ b/plugins/omega/omega.ml
@@ -0,0 +1,716 @@
+(************************************************************************)
+(* v * The Coq Proof Assistant / The Coq Development Team *)
+(* <O___,, * CNRS-Ecole Polytechnique-INRIA Futurs-Universite Paris Sud *)
+(* \VV/ **************************************************************)
+(* // * This file is distributed under the terms of the *)
+(* * GNU Lesser General Public License Version 2.1 *)
+(************************************************************************)
+(**************************************************************************)
+(* *)
+(* Omega: a solver of quantifier-free problems in Presburger Arithmetic *)
+(* *)
+(* Pierre Crégut (CNET, Lannion, France) *)
+(* *)
+(* 13/10/2002 : modified to cope with an external numbering of equations *)
+(* and hypothesis. Its use for Omega is not more complex and it makes *)
+(* things much simpler for the reflexive version where we should limit *)
+(* the number of source of numbering. *)
+(**************************************************************************)
+
+open Names
+
+module type INT = sig
+ type bigint
+ val less_than : bigint -> bigint -> bool
+ val add : bigint -> bigint -> bigint
+ val sub : bigint -> bigint -> bigint
+ val mult : bigint -> bigint -> bigint
+ val euclid : bigint -> bigint -> bigint * bigint
+ val neg : bigint -> bigint
+ val zero : bigint
+ val one : bigint
+ val to_string : bigint -> string
+end
+
+let debug = ref false
+
+module MakeOmegaSolver (Int:INT) = struct
+
+type bigint = Int.bigint
+let (<?) = Int.less_than
+let (<=?) x y = Int.less_than x y or x = y
+let (>?) x y = Int.less_than y x
+let (>=?) x y = Int.less_than y x or x = y
+let (=?) = (=)
+let (+) = Int.add
+let (-) = Int.sub
+let ( * ) = Int.mult
+let (/) x y = fst (Int.euclid x y)
+let (mod) x y = snd (Int.euclid x y)
+let zero = Int.zero
+let one = Int.one
+let two = one + one
+let negone = Int.neg one
+let abs x = if Int.less_than x zero then Int.neg x else x
+let string_of_bigint = Int.to_string
+let neg = Int.neg
+
+(* To ensure that polymorphic (<) is not used mistakenly on big integers *)
+(* Warning: do not use (=) either on big int *)
+let (<) = ((<) : int -> int -> bool)
+let (>) = ((>) : int -> int -> bool)
+let (<=) = ((<=) : int -> int -> bool)
+let (>=) = ((>=) : int -> int -> bool)
+
+let pp i = print_int i; print_newline (); flush stdout
+
+let push v l = l := v :: !l
+
+let rec pgcd x y = if y =? zero then x else pgcd y (x mod y)
+
+let pgcd_l = function
+ | [] -> failwith "pgcd_l"
+ | x :: l -> List.fold_left pgcd x l
+
+let floor_div a b =
+ match a >=? zero , b >? zero with
+ | true,true -> a / b
+ | false,false -> a / b
+ | true, false -> (a-one) / b - one
+ | false,true -> (a+one) / b - one
+
+type coeff = {c: bigint ; v: int}
+
+type linear = coeff list
+
+type eqn_kind = EQUA | INEQ | DISE
+
+type afine = {
+ (* a number uniquely identifying the equation *)
+ id: int ;
+ (* a boolean true for an eq, false for an ineq (Sigma a_i x_i >= 0) *)
+ kind: eqn_kind;
+ (* the variables and their coefficient *)
+ body: coeff list;
+ (* a constant *)
+ constant: bigint }
+
+type state_action = {
+ st_new_eq : afine;
+ st_def : afine;
+ st_orig : afine;
+ st_coef : bigint;
+ st_var : int }
+
+type action =
+ | DIVIDE_AND_APPROX of afine * afine * bigint * bigint
+ | NOT_EXACT_DIVIDE of afine * bigint
+ | FORGET_C of int
+ | EXACT_DIVIDE of afine * bigint
+ | SUM of int * (bigint * afine) * (bigint * afine)
+ | STATE of state_action
+ | HYP of afine
+ | FORGET of int * int
+ | FORGET_I of int * int
+ | CONTRADICTION of afine * afine
+ | NEGATE_CONTRADICT of afine * afine * bool
+ | MERGE_EQ of int * afine * int
+ | CONSTANT_NOT_NUL of int * bigint
+ | CONSTANT_NUL of int
+ | CONSTANT_NEG of int * bigint
+ | SPLIT_INEQ of afine * (int * action list) * (int * action list)
+ | WEAKEN of int * bigint
+
+exception UNSOLVABLE
+
+exception NO_CONTRADICTION
+
+let display_eq print_var (l,e) =
+ let _ =
+ List.fold_left
+ (fun not_first f ->
+ print_string
+ (if f.c <? zero then "- " else if not_first then "+ " else "");
+ let c = abs f.c in
+ if c =? one then
+ Printf.printf "%s " (print_var f.v)
+ else
+ Printf.printf "%s %s " (string_of_bigint c) (print_var f.v);
+ true)
+ false l
+ in
+ if e >? zero then
+ Printf.printf "+ %s " (string_of_bigint e)
+ else if e <? zero then
+ Printf.printf "- %s " (string_of_bigint (abs e))
+
+let rec trace_length l =
+ let action_length accu = function
+ | SPLIT_INEQ (_,(_,l1),(_,l2)) ->
+ accu + one + trace_length l1 + trace_length l2
+ | _ -> accu + one in
+ List.fold_left action_length zero l
+
+let operator_of_eq = function
+ | EQUA -> "=" | DISE -> "!=" | INEQ -> ">="
+
+let kind_of = function
+ | EQUA -> "equation" | DISE -> "disequation" | INEQ -> "inequation"
+
+let display_system print_var l =
+ List.iter
+ (fun { kind=b; body=e; constant=c; id=id} ->
+ Printf.printf "E%d: " id;
+ display_eq print_var (e,c);
+ Printf.printf "%s 0\n" (operator_of_eq b))
+ l;
+ print_string "------------------------\n\n"
+
+let display_inequations print_var l =
+ List.iter (fun e -> display_eq print_var e;print_string ">= 0\n") l;
+ print_string "------------------------\n\n"
+
+let sbi = string_of_bigint
+
+let rec display_action print_var = function
+ | act :: l -> begin match act with
+ | DIVIDE_AND_APPROX (e1,e2,k,d) ->
+ Printf.printf
+ "Inequation E%d is divided by %s and the constant coefficient is \
+ rounded by substracting %s.\n" e1.id (sbi k) (sbi d)
+ | NOT_EXACT_DIVIDE (e,k) ->
+ Printf.printf
+ "Constant in equation E%d is not divisible by the pgcd \
+ %s of its other coefficients.\n" e.id (sbi k)
+ | EXACT_DIVIDE (e,k) ->
+ Printf.printf
+ "Equation E%d is divided by the pgcd \
+ %s of its coefficients.\n" e.id (sbi k)
+ | WEAKEN (e,k) ->
+ Printf.printf
+ "To ensure a solution in the dark shadow \
+ the equation E%d is weakened by %s.\n" e (sbi k)
+ | SUM (e,(c1,e1),(c2,e2)) ->
+ Printf.printf
+ "We state %s E%d = %s %s E%d + %s %s E%d.\n"
+ (kind_of e1.kind) e (sbi c1) (kind_of e1.kind) e1.id (sbi c2)
+ (kind_of e2.kind) e2.id
+ | STATE { st_new_eq = e } ->
+ Printf.printf "We define a new equation E%d: " e.id;
+ display_eq print_var (e.body,e.constant);
+ print_string (operator_of_eq e.kind); print_string " 0"
+ | HYP e ->
+ Printf.printf "We define E%d: " e.id;
+ display_eq print_var (e.body,e.constant);
+ print_string (operator_of_eq e.kind); print_string " 0\n"
+ | FORGET_C e -> Printf.printf "E%d is trivially satisfiable.\n" e
+ | FORGET (e1,e2) -> Printf.printf "E%d subsumes E%d.\n" e1 e2
+ | FORGET_I (e1,e2) -> Printf.printf "E%d subsumes E%d.\n" e1 e2
+ | MERGE_EQ (e,e1,e2) ->
+ Printf.printf "E%d and E%d can be merged into E%d.\n" e1.id e2 e
+ | CONTRADICTION (e1,e2) ->
+ Printf.printf
+ "Equations E%d and E%d imply a contradiction on their \
+ constant factors.\n" e1.id e2.id
+ | NEGATE_CONTRADICT(e1,e2,b) ->
+ Printf.printf
+ "Equations E%d and E%d state that their body is at the same time
+ equal and different\n" e1.id e2.id
+ | CONSTANT_NOT_NUL (e,k) ->
+ Printf.printf "Equation E%d states %s = 0.\n" e (sbi k)
+ | CONSTANT_NEG(e,k) ->
+ Printf.printf "Equation E%d states %s >= 0.\n" e (sbi k)
+ | CONSTANT_NUL e ->
+ Printf.printf "Inequation E%d states 0 != 0.\n" e
+ | SPLIT_INEQ (e,(e1,l1),(e2,l2)) ->
+ Printf.printf "Equation E%d is split in E%d and E%d\n\n" e.id e1 e2;
+ display_action print_var l1;
+ print_newline ();
+ display_action print_var l2;
+ print_newline ()
+ end; display_action print_var l
+ | [] ->
+ flush stdout
+
+let default_print_var v = Printf.sprintf "X%d" v (* For debugging *)
+
+(*""*)
+let add_event, history, clear_history =
+ let accu = ref [] in
+ (fun (v:action) -> if !debug then display_action default_print_var [v]; push v accu),
+ (fun () -> !accu),
+ (fun () -> accu := [])
+
+let nf_linear = Sort.list (fun x y -> x.v > y.v)
+
+let nf ((b : bool),(e,(x : int))) = (b,(nf_linear e,x))
+
+let map_eq_linear f =
+ let rec loop = function
+ | x :: l -> let c = f x.c in if c=?zero then loop l else {v=x.v; c=c} :: loop l
+ | [] -> []
+ in
+ loop
+
+let map_eq_afine f e =
+ { id = e.id; kind = e.kind; body = map_eq_linear f e.body;
+ constant = f e.constant }
+
+let negate_eq = map_eq_afine (fun x -> neg x)
+
+let rec sum p0 p1 = match (p0,p1) with
+ | ([], l) -> l | (l, []) -> l
+ | (((x1::l1) as l1'), ((x2::l2) as l2')) ->
+ if x1.v = x2.v then
+ let c = x1.c + x2.c in
+ if c =? zero then sum l1 l2 else {v=x1.v;c=c} :: sum l1 l2
+ else if x1.v > x2.v then
+ x1 :: sum l1 l2'
+ else
+ x2 :: sum l1' l2
+
+let sum_afine new_eq_id eq1 eq2 =
+ { kind = eq1.kind; id = new_eq_id ();
+ body = sum eq1.body eq2.body; constant = eq1.constant + eq2.constant }
+
+exception FACTOR1
+
+let rec chop_factor_1 = function
+ | x :: l ->
+ if abs x.c =? one then x,l else let (c',l') = chop_factor_1 l in (c',x::l')
+ | [] -> raise FACTOR1
+
+exception CHOPVAR
+
+let rec chop_var v = function
+ | f :: l -> if f.v = v then f,l else let (f',l') = chop_var v l in (f',f::l')
+ | [] -> raise CHOPVAR
+
+let normalize ({id=id; kind=eq_flag; body=e; constant =x} as eq) =
+ if e = [] then begin
+ match eq_flag with
+ | EQUA ->
+ if x =? zero then [] else begin
+ add_event (CONSTANT_NOT_NUL(id,x)); raise UNSOLVABLE
+ end
+ | DISE ->
+ if x <> zero then [] else begin
+ add_event (CONSTANT_NUL id); raise UNSOLVABLE
+ end
+ | INEQ ->
+ if x >=? zero then [] else begin
+ add_event (CONSTANT_NEG(id,x)); raise UNSOLVABLE
+ end
+ end else
+ let gcd = pgcd_l (List.map (fun f -> abs f.c) e) in
+ if eq_flag=EQUA & x mod gcd <> zero then begin
+ add_event (NOT_EXACT_DIVIDE (eq,gcd)); raise UNSOLVABLE
+ end else if eq_flag=DISE & x mod gcd <> zero then begin
+ add_event (FORGET_C eq.id); []
+ end else if gcd <> one then begin
+ let c = floor_div x gcd in
+ let d = x - c * gcd in
+ let new_eq = {id=id; kind=eq_flag; constant=c;
+ body=map_eq_linear (fun c -> c / gcd) e} in
+ add_event (if eq_flag=EQUA or eq_flag = DISE then EXACT_DIVIDE(eq,gcd)
+ else DIVIDE_AND_APPROX(eq,new_eq,gcd,d));
+ [new_eq]
+ end else [eq]
+
+let eliminate_with_in new_eq_id {v=v;c=c_unite} eq2
+ ({body=e1; constant=c1} as eq1) =
+ try
+ let (f,_) = chop_var v e1 in
+ let coeff = if c_unite=?one then neg f.c else if c_unite=? negone then f.c
+ else failwith "eliminate_with_in" in
+ let res = sum_afine new_eq_id eq1 (map_eq_afine (fun c -> c * coeff) eq2) in
+ add_event (SUM (res.id,(one,eq1),(coeff,eq2))); res
+ with CHOPVAR -> eq1
+
+let omega_mod a b = a - b * floor_div (two * a + b) (two * b)
+let banerjee_step (new_eq_id,new_var_id,print_var) original l1 l2 =
+ let e = original.body in
+ let sigma = new_var_id () in
+ let smallest,var =
+ try
+ List.fold_left (fun (v,p) c -> if v >? (abs c.c) then abs c.c,c.v else (v,p))
+ (abs (List.hd e).c, (List.hd e).v) (List.tl e)
+ with Failure "tl" -> display_system print_var [original] ; failwith "TL" in
+ let m = smallest + one in
+ let new_eq =
+ { constant = omega_mod original.constant m;
+ body = {c= neg m;v=sigma} ::
+ map_eq_linear (fun a -> omega_mod a m) original.body;
+ id = new_eq_id (); kind = EQUA } in
+ let definition =
+ { constant = neg (floor_div (two * original.constant + m) (two * m));
+ body = map_eq_linear (fun a -> neg (floor_div (two * a + m) (two * m)))
+ original.body;
+ id = new_eq_id (); kind = EQUA } in
+ add_event (STATE {st_new_eq = new_eq; st_def = definition;
+ st_orig = original; st_coef = m; st_var = sigma});
+ let new_eq = List.hd (normalize new_eq) in
+ let eliminated_var, def = chop_var var new_eq.body in
+ let other_equations =
+ Util.list_map_append
+ (fun e ->
+ normalize (eliminate_with_in new_eq_id eliminated_var new_eq e)) l1 in
+ let inequations =
+ Util.list_map_append
+ (fun e ->
+ normalize (eliminate_with_in new_eq_id eliminated_var new_eq e)) l2 in
+ let original' = eliminate_with_in new_eq_id eliminated_var new_eq original in
+ let mod_original = map_eq_afine (fun c -> c / m) original' in
+ add_event (EXACT_DIVIDE (original',m));
+ List.hd (normalize mod_original),other_equations,inequations
+
+let rec eliminate_one_equation ((new_eq_id,new_var_id,print_var) as new_ids) (e,other,ineqs) =
+ if !debug then display_system print_var (e::other);
+ try
+ let v,def = chop_factor_1 e.body in
+ (Util.list_map_append
+ (fun e' -> normalize (eliminate_with_in new_eq_id v e e')) other,
+ Util.list_map_append
+ (fun e' -> normalize (eliminate_with_in new_eq_id v e e')) ineqs)
+ with FACTOR1 ->
+ eliminate_one_equation new_ids (banerjee_step new_ids e other ineqs)
+
+let rec banerjee ((_,_,print_var) as new_ids) (sys_eq,sys_ineq) =
+ let rec fst_eq_1 = function
+ (eq::l) ->
+ if List.exists (fun x -> abs x.c =? one) eq.body then eq,l
+ else let (eq',l') = fst_eq_1 l in (eq',eq::l')
+ | [] -> raise Not_found in
+ match sys_eq with
+ [] -> if !debug then display_system print_var sys_ineq; sys_ineq
+ | (e1::rest) ->
+ let eq,other = try fst_eq_1 sys_eq with Not_found -> (e1,rest) in
+ if eq.body = [] then
+ if eq.constant =? zero then begin
+ add_event (FORGET_C eq.id); banerjee new_ids (other,sys_ineq)
+ end else begin
+ add_event (CONSTANT_NOT_NUL(eq.id,eq.constant)); raise UNSOLVABLE
+ end
+ else
+ banerjee new_ids
+ (eliminate_one_equation new_ids (eq,other,sys_ineq))
+
+type kind = INVERTED | NORMAL
+
+let redundancy_elimination new_eq_id system =
+ let normal = function
+ ({body=f::_} as e) when f.c <? zero -> negate_eq e, INVERTED
+ | e -> e,NORMAL in
+ let table = Hashtbl.create 7 in
+ List.iter
+ (fun e ->
+ let ({body=ne} as nx) ,kind = normal e in
+ if ne = [] then
+ if nx.constant <? zero then begin
+ add_event (CONSTANT_NEG(nx.id,nx.constant)); raise UNSOLVABLE
+ end else add_event (FORGET_C nx.id)
+ else
+ try
+ let (optnormal,optinvert) = Hashtbl.find table ne in
+ let final =
+ if kind = NORMAL then begin
+ match optnormal with
+ Some v ->
+ let kept =
+ if v.constant <? nx.constant
+ then begin add_event (FORGET (v.id,nx.id));v end
+ else begin add_event (FORGET (nx.id,v.id));nx end in
+ (Some(kept),optinvert)
+ | None -> Some nx,optinvert
+ end else begin
+ match optinvert with
+ Some v ->
+ let _kept =
+ if v.constant >? nx.constant
+ then begin add_event (FORGET_I (v.id,nx.id));v end
+ else begin add_event (FORGET_I (nx.id,v.id));nx end in
+ (optnormal,Some(if v.constant >? nx.constant then v else nx))
+ | None -> optnormal,Some nx
+ end in
+ begin match final with
+ (Some high, Some low) ->
+ if high.constant <? low.constant then begin
+ add_event(CONTRADICTION (high,negate_eq low));
+ raise UNSOLVABLE
+ end
+ | _ -> () end;
+ Hashtbl.remove table ne;
+ Hashtbl.add table ne final
+ with Not_found ->
+ Hashtbl.add table ne
+ (if kind = NORMAL then (Some nx,None) else (None,Some nx)))
+ system;
+ let accu_eq = ref [] in
+ let accu_ineq = ref [] in
+ Hashtbl.iter
+ (fun p0 p1 -> match (p0,p1) with
+ | (e, (Some x, Some y)) when x.constant =? y.constant ->
+ let id=new_eq_id () in
+ add_event (MERGE_EQ(id,x,y.id));
+ push {id=id; kind=EQUA; body=x.body; constant=x.constant} accu_eq
+ | (e, (optnorm,optinvert)) ->
+ begin match optnorm with
+ Some x -> push x accu_ineq | _ -> () end;
+ begin match optinvert with
+ Some x -> push (negate_eq x) accu_ineq | _ -> () end)
+ table;
+ !accu_eq,!accu_ineq
+
+exception SOLVED_SYSTEM
+
+let select_variable system =
+ let table = Hashtbl.create 7 in
+ let push v c=
+ try let r = Hashtbl.find table v in r := max !r (abs c)
+ with Not_found -> Hashtbl.add table v (ref (abs c)) in
+ List.iter (fun {body=l} -> List.iter (fun f -> push f.v f.c) l) system;
+ let vmin,cmin = ref (-1), ref zero in
+ let var_cpt = ref 0 in
+ Hashtbl.iter
+ (fun v ({contents = c}) ->
+ incr var_cpt;
+ if c <? !cmin or !vmin = (-1) then begin vmin := v; cmin := c end)
+ table;
+ if !var_cpt < 1 then raise SOLVED_SYSTEM;
+ !vmin
+
+let classify v system =
+ List.fold_left
+ (fun (not_occ,below,over) eq ->
+ try let f,eq' = chop_var v eq.body in
+ if f.c >=? zero then (not_occ,((f.c,eq) :: below),over)
+ else (not_occ,below,((neg f.c,eq) :: over))
+ with CHOPVAR -> (eq::not_occ,below,over))
+ ([],[],[]) system
+
+let product new_eq_id dark_shadow low high =
+ List.fold_left
+ (fun accu (a,eq1) ->
+ List.fold_left
+ (fun accu (b,eq2) ->
+ let eq =
+ sum_afine new_eq_id (map_eq_afine (fun c -> c * b) eq1)
+ (map_eq_afine (fun c -> c * a) eq2) in
+ add_event(SUM(eq.id,(b,eq1),(a,eq2)));
+ match normalize eq with
+ | [eq] ->
+ let final_eq =
+ if dark_shadow then
+ let delta = (a - one) * (b - one) in
+ add_event(WEAKEN(eq.id,delta));
+ {id = eq.id; kind=INEQ; body = eq.body;
+ constant = eq.constant - delta}
+ else eq
+ in final_eq :: accu
+ | (e::_) -> failwith "Product dardk"
+ | [] -> accu)
+ accu high)
+ [] low
+
+let fourier_motzkin (new_eq_id,_,print_var) dark_shadow system =
+ let v = select_variable system in
+ let (ineq_out, ineq_low,ineq_high) = classify v system in
+ let expanded = ineq_out @ product new_eq_id dark_shadow ineq_low ineq_high in
+ if !debug then display_system print_var expanded; expanded
+
+let simplify ((new_eq_id,new_var_id,print_var) as new_ids) dark_shadow system =
+ if List.exists (fun e -> e.kind = DISE) system then
+ failwith "disequation in simplify";
+ clear_history ();
+ List.iter (fun e -> add_event (HYP e)) system;
+ let system = Util.list_map_append normalize system in
+ let eqs,ineqs = List.partition (fun e -> e.kind=EQUA) system in
+ let simp_eq,simp_ineq = redundancy_elimination new_eq_id ineqs in
+ let system = (eqs @ simp_eq,simp_ineq) in
+ let rec loop1a system =
+ let sys_ineq = banerjee new_ids system in
+ loop1b sys_ineq
+ and loop1b sys_ineq =
+ let simp_eq,simp_ineq = redundancy_elimination new_eq_id sys_ineq in
+ if simp_eq = [] then simp_ineq else loop1a (simp_eq,simp_ineq)
+ in
+ let rec loop2 system =
+ try
+ let expanded = fourier_motzkin new_ids dark_shadow system in
+ loop2 (loop1b expanded)
+ with SOLVED_SYSTEM ->
+ if !debug then display_system print_var system; system
+ in
+ loop2 (loop1a system)
+
+let rec depend relie_on accu = function
+ | act :: l ->
+ begin match act with
+ | DIVIDE_AND_APPROX (e,_,_,_) ->
+ if List.mem e.id relie_on then depend relie_on (act::accu) l
+ else depend relie_on accu l
+ | EXACT_DIVIDE (e,_) ->
+ if List.mem e.id relie_on then depend relie_on (act::accu) l
+ else depend relie_on accu l
+ | WEAKEN (e,_) ->
+ if List.mem e relie_on then depend relie_on (act::accu) l
+ else depend relie_on accu l
+ | SUM (e,(_,e1),(_,e2)) ->
+ if List.mem e relie_on then
+ depend (e1.id::e2.id::relie_on) (act::accu) l
+ else
+ depend relie_on accu l
+ | STATE {st_new_eq=e;st_orig=o} ->
+ if List.mem e.id relie_on then depend (o.id::relie_on) (act::accu) l
+ else depend relie_on accu l
+ | HYP e ->
+ if List.mem e.id relie_on then depend relie_on (act::accu) l
+ else depend relie_on accu l
+ | FORGET_C _ -> depend relie_on accu l
+ | FORGET _ -> depend relie_on accu l
+ | FORGET_I _ -> depend relie_on accu l
+ | MERGE_EQ (e,e1,e2) ->
+ if List.mem e relie_on then
+ depend (e1.id::e2::relie_on) (act::accu) l
+ else
+ depend relie_on accu l
+ | NOT_EXACT_DIVIDE (e,_) -> depend (e.id::relie_on) (act::accu) l
+ | CONTRADICTION (e1,e2) ->
+ depend (e1.id::e2.id::relie_on) (act::accu) l
+ | CONSTANT_NOT_NUL (e,_) -> depend (e::relie_on) (act::accu) l
+ | CONSTANT_NEG (e,_) -> depend (e::relie_on) (act::accu) l
+ | CONSTANT_NUL e -> depend (e::relie_on) (act::accu) l
+ | NEGATE_CONTRADICT (e1,e2,_) ->
+ depend (e1.id::e2.id::relie_on) (act::accu) l
+ | SPLIT_INEQ _ -> failwith "depend"
+ end
+ | [] -> relie_on, accu
+
+(*
+let depend relie_on accu trace =
+ Printf.printf "Longueur de la trace initiale : %d\n"
+ (trace_length trace + trace_length accu);
+ let rel',trace' = depend relie_on accu trace in
+ Printf.printf "Longueur de la trace simplifiée : %d\n" (trace_length trace');
+ rel',trace'
+*)
+
+let solve (new_eq_id,new_eq_var,print_var) system =
+ try let _ = simplify new_eq_id false system in failwith "no contradiction"
+ with UNSOLVABLE -> display_action print_var (snd (depend [] [] (history ())))
+
+let negation (eqs,ineqs) =
+ let diseq,_ = List.partition (fun e -> e.kind = DISE) ineqs in
+ let normal = function
+ | ({body=f::_} as e) when f.c <? zero -> negate_eq e, INVERTED
+ | e -> e,NORMAL in
+ let table = Hashtbl.create 7 in
+ List.iter (fun e ->
+ let {body=ne;constant=c} ,kind = normal e in
+ Hashtbl.add table (ne,c) (kind,e)) diseq;
+ List.iter (fun e ->
+ assert (e.kind = EQUA);
+ let {body=ne;constant=c},kind = normal e in
+ try
+ let (kind',e') = Hashtbl.find table (ne,c) in
+ add_event (NEGATE_CONTRADICT (e,e',kind=kind'));
+ raise UNSOLVABLE
+ with Not_found -> ()) eqs
+
+exception FULL_SOLUTION of action list * int list
+
+let simplify_strong ((new_eq_id,new_var_id,print_var) as new_ids) system =
+ clear_history ();
+ List.iter (fun e -> add_event (HYP e)) system;
+ (* Initial simplification phase *)
+ let rec loop1a system =
+ negation system;
+ let sys_ineq = banerjee new_ids system in
+ loop1b sys_ineq
+ and loop1b sys_ineq =
+ let dise,ine = List.partition (fun e -> e.kind = DISE) sys_ineq in
+ let simp_eq,simp_ineq = redundancy_elimination new_eq_id ine in
+ if simp_eq = [] then dise @ simp_ineq
+ else loop1a (simp_eq,dise @ simp_ineq)
+ in
+ let rec loop2 system =
+ try
+ let expanded = fourier_motzkin new_ids false system in
+ loop2 (loop1b expanded)
+ with SOLVED_SYSTEM -> if !debug then display_system print_var system; system
+ in
+ let rec explode_diseq = function
+ | (de::diseq,ineqs,expl_map) ->
+ let id1 = new_eq_id ()
+ and id2 = new_eq_id () in
+ let e1 =
+ {id = id1; kind=INEQ; body = de.body; constant = de.constant -one} in
+ let e2 =
+ {id = id2; kind=INEQ; body = map_eq_linear neg de.body;
+ constant = neg de.constant - one} in
+ let new_sys =
+ List.map (fun (what,sys) -> ((de.id,id1,true)::what, e1::sys))
+ ineqs @
+ List.map (fun (what,sys) -> ((de.id,id2,false)::what,e2::sys))
+ ineqs
+ in
+ explode_diseq (diseq,new_sys,(de.id,(de,id1,id2))::expl_map)
+ | ([],ineqs,expl_map) -> ineqs,expl_map
+ in
+ try
+ let system = Util.list_map_append normalize system in
+ let eqs,ineqs = List.partition (fun e -> e.kind=EQUA) system in
+ let dise,ine = List.partition (fun e -> e.kind = DISE) ineqs in
+ let simp_eq,simp_ineq = redundancy_elimination new_eq_id ine in
+ let system = (eqs @ simp_eq,simp_ineq @ dise) in
+ let system' = loop1a system in
+ let diseq,ineq = List.partition (fun e -> e.kind = DISE) system' in
+ let first_segment = history () in
+ let sys_exploded,explode_map = explode_diseq (diseq,[[],ineq],[]) in
+ let all_solutions =
+ List.map
+ (fun (decomp,sys) ->
+ clear_history ();
+ try let _ = loop2 sys in raise NO_CONTRADICTION
+ with UNSOLVABLE ->
+ let relie_on,path = depend [] [] (history ()) in
+ let dc,_ = List.partition (fun (_,id,_) -> List.mem id relie_on) decomp in
+ let red = List.map (fun (x,_,_) -> x) dc in
+ (red,relie_on,decomp,path))
+ sys_exploded
+ in
+ let max_count sys =
+ let tbl = Hashtbl.create 7 in
+ let augment x =
+ try incr (Hashtbl.find tbl x)
+ with Not_found -> Hashtbl.add tbl x (ref 1) in
+ let eq = ref (-1) and c = ref 0 in
+ List.iter (function
+ | ([],r_on,_,path) -> raise (FULL_SOLUTION (path,r_on))
+ | (l,_,_,_) -> List.iter augment l) sys;
+ Hashtbl.iter (fun x v -> if !v > !c then begin eq := x; c := !v end) tbl;
+ !eq
+ in
+ let rec solve systems =
+ try
+ let id = max_count systems in
+ let rec sign = function
+ | ((id',_,b)::l) -> if id=id' then b else sign l
+ | [] -> failwith "solve" in
+ let s1,s2 =
+ List.partition (fun (_,_,decomp,_) -> sign decomp) systems in
+ let s1' =
+ List.map (fun (dep,ro,dc,pa) -> (Util.list_except id dep,ro,dc,pa)) s1 in
+ let s2' =
+ List.map (fun (dep,ro,dc,pa) -> (Util.list_except id dep,ro,dc,pa)) s2 in
+ let (r1,relie1) = solve s1'
+ and (r2,relie2) = solve s2' in
+ let (eq,id1,id2) = List.assoc id explode_map in
+ [SPLIT_INEQ(eq,(id1,r1),(id2, r2))], eq.id :: Util.list_union relie1 relie2
+ with FULL_SOLUTION (x0,x1) -> (x0,x1)
+ in
+ let act,relie_on = solve all_solutions in
+ snd(depend relie_on act first_segment)
+ with UNSOLVABLE -> snd (depend [] [] (history ()))
+
+end
generated by cgit v1.2.3 (git 2.39.1) at 2024-05-14 09:44:07 +0000