summaryrefslogtreecommitdiff
path: root/plugins/nsatz/nsatz.ml
diff options
context:
space:
mode:
Diffstat (limited to 'plugins/nsatz/nsatz.ml')
-rw-r--r--plugins/nsatz/nsatz.ml288
1 files changed, 104 insertions, 184 deletions
diff --git a/plugins/nsatz/nsatz.ml b/plugins/nsatz/nsatz.ml
index 36bce780..81b44ffa 100644
--- a/plugins/nsatz/nsatz.ml
+++ b/plugins/nsatz/nsatz.ml
@@ -1,14 +1,16 @@
(************************************************************************)
-(* v * The Coq Proof Assistant / The Coq Development Team *)
-(* <O___,, * INRIA - CNRS - LIX - LRI - PPS - Copyright 1999-2016 *)
+(* * The Coq Proof Assistant / The Coq Development Team *)
+(* v * INRIA, CNRS and contributors - Copyright 1999-2018 *)
+(* <O___,, * (see CREDITS file for the list of authors) *)
(* \VV/ **************************************************************)
-(* // * This file is distributed under the terms of the *)
-(* * GNU Lesser General Public License Version 2.1 *)
+(* // * This file is distributed under the terms of the *)
+(* * GNU Lesser General Public License Version 2.1 *)
+(* * (see LICENSE file for the text of the license) *)
(************************************************************************)
open CErrors
open Util
-open Term
+open Constr
open Tactics
open Coqlib
@@ -22,7 +24,6 @@ open Utile
let num_0 = Int 0
and num_1 = Int 1
and num_2 = Int 2
-and num_10 = Int 10
let numdom r =
let r' = Ratio.normalize_ratio (ratio_of_num r) in
@@ -35,7 +36,6 @@ module BigInt = struct
type t = big_int
let of_int = big_int_of_int
let coef0 = of_int 0
- let coef1 = of_int 1
let of_num = Num.big_int_of_num
let to_num = Num.num_of_big_int
let equal = eq_big_int
@@ -49,7 +49,6 @@ module BigInt = struct
let div = div_big_int
let modulo = mod_big_int
let to_string = string_of_big_int
- let to_int x = int_of_big_int x
let hash x =
try (int_of_big_int x)
with Failure _ -> 1
@@ -61,15 +60,6 @@ module BigInt = struct
then a
else if lt a b then pgcd b a else pgcd b (modulo a b)
-
- (* signe du pgcd = signe(a)*signe(b) si non nuls. *)
- let pgcd2 a b =
- if equal a coef0 then b
- else if equal b coef0 then a
- else let c = pgcd (abs a) (abs b) in
- if ((lt coef0 a)&&(lt b coef0))
- ||((lt coef0 b)&&(lt a coef0))
- then opp c else c
end
(*
@@ -146,10 +136,10 @@ let mul = function
| (Const n,q) when eq_num n num_1 -> q
| (p,q) -> Mul(p,q)
-let unconstr = mkRel 1
+let gen_constant msg path s = Universes.constr_of_global @@
+ coq_reference msg path s
-let tpexpr =
- lazy (gen_constant "CC" ["setoid_ring";"Ring_polynom"] "PExpr")
+let tpexpr = lazy (gen_constant "CC" ["setoid_ring";"Ring_polynom"] "PExpr")
let ttconst = lazy (gen_constant "CC" ["setoid_ring";"Ring_polynom"] "PEc")
let ttvar = lazy (gen_constant "CC" ["setoid_ring";"Ring_polynom"] "PEX")
let ttadd = lazy (gen_constant "CC" ["setoid_ring";"Ring_polynom"] "PEadd")
@@ -216,76 +206,60 @@ else
mkt_app ttpow [Lazy.force tz; mkt_term t1; mkt_n (num_of_int n)]
let rec parse_pos p =
- match kind_of_term p with
+ match Constr.kind p with
| App (a,[|p2|]) ->
- if eq_constr a (Lazy.force pxO) then num_2 */ (parse_pos p2)
+ if Constr.equal a (Lazy.force pxO) then num_2 */ (parse_pos p2)
else num_1 +/ (num_2 */ (parse_pos p2))
| _ -> num_1
let parse_z z =
- match kind_of_term z with
+ match Constr.kind z with
| App (a,[|p2|]) ->
- if eq_constr a (Lazy.force zpos) then parse_pos p2 else (num_0 -/ (parse_pos p2))
+ if Constr.equal a (Lazy.force zpos) then parse_pos p2 else (num_0 -/ (parse_pos p2))
| _ -> num_0
let parse_n z =
- match kind_of_term z with
+ match Constr.kind z with
| App (a,[|p2|]) ->
parse_pos p2
| _ -> num_0
let rec parse_term p =
- match kind_of_term p with
+ match Constr.kind p with
| App (a,[|_;p2|]) ->
- if eq_constr a (Lazy.force ttvar) then Var (string_of_num (parse_pos p2))
- else if eq_constr a (Lazy.force ttconst) then Const (parse_z p2)
- else if eq_constr a (Lazy.force ttopp) then Opp (parse_term p2)
+ if Constr.equal a (Lazy.force ttvar) then Var (string_of_num (parse_pos p2))
+ else if Constr.equal a (Lazy.force ttconst) then Const (parse_z p2)
+ else if Constr.equal a (Lazy.force ttopp) then Opp (parse_term p2)
else Zero
| App (a,[|_;p2;p3|]) ->
- if eq_constr a (Lazy.force ttadd) then Add (parse_term p2, parse_term p3)
- else if eq_constr a (Lazy.force ttsub) then Sub (parse_term p2, parse_term p3)
- else if eq_constr a (Lazy.force ttmul) then Mul (parse_term p2, parse_term p3)
- else if eq_constr a (Lazy.force ttpow) then
+ if Constr.equal a (Lazy.force ttadd) then Add (parse_term p2, parse_term p3)
+ else if Constr.equal a (Lazy.force ttsub) then Sub (parse_term p2, parse_term p3)
+ else if Constr.equal a (Lazy.force ttmul) then Mul (parse_term p2, parse_term p3)
+ else if Constr.equal a (Lazy.force ttpow) then
Pow (parse_term p2, int_of_num (parse_n p3))
else Zero
| _ -> Zero
let rec parse_request lp =
- match kind_of_term lp with
+ match Constr.kind lp with
| App (_,[|_|]) -> []
| App (_,[|_;p;lp1|]) ->
(parse_term p)::(parse_request lp1)
|_-> assert false
-let nvars = ref 0
-
-let set_nvars_term t =
- let rec aux t =
+let set_nvars_term nvars t =
+ let rec aux t nvars =
match t with
- | Zero -> ()
- | Const r -> ()
+ | Zero -> nvars
+ | Const r -> nvars
| Var v -> let n = int_of_string v in
- nvars:= max (!nvars) n
- | Opp t1 -> aux t1
- | Add (t1,t2) -> aux t1; aux t2
- | Sub (t1,t2) -> aux t1; aux t2
- | Mul (t1,t2) -> aux t1; aux t2
- | Pow (t1,n) -> aux t1
- in aux t
-
-let string_of_term p =
- let rec aux p =
- match p with
- | Zero -> "0"
- | Const r -> string_of_num r
- | Var v -> "x"^v
- | Opp t1 -> "(-"^(aux t1)^")"
- | Add (t1,t2) -> "("^(aux t1)^"+"^(aux t2)^")"
- | Sub (t1,t2) -> "("^(aux t1)^"-"^(aux t2)^")"
- | Mul (t1,t2) -> "("^(aux t1)^"*"^(aux t2)^")"
- | Pow (t1,n) -> (aux t1)^"^"^(string_of_int n)
- in aux p
-
+ max nvars n
+ | Opp t1 -> aux t1 nvars
+ | Add (t1,t2) -> aux t2 (aux t1 nvars)
+ | Sub (t1,t2) -> aux t2 (aux t1 nvars)
+ | Mul (t1,t2) -> aux t2 (aux t1 nvars)
+ | Pow (t1,n) -> aux t1 nvars
+ in aux t nvars
(***********************************************************************
Coefficients: recursive polynomials
@@ -301,8 +275,8 @@ open PIdeal
varaibles <=np are in the coefficients
*)
-let term_pol_sparse np t=
- let d = !nvars in
+let term_pol_sparse nvars np t=
+ let d = nvars in
let rec aux t =
(* info ("conversion de: "^(string_of_term t)^"\n");*)
let res =
@@ -336,14 +310,8 @@ let polrec_to_term p =
match p with
|Poly.Pint n -> const (Coef.to_num n)
|Poly.Prec (v,coefs) ->
- let res = ref Zero in
- Array.iteri
- (fun i c ->
- res:=add(!res, mul(aux c,
- pow (Var (string_of_int v),
- i))))
- coefs;
- !res
+ let fold i c res = add (res, mul (aux c, pow (Var (string_of_int v), i))) in
+ Array.fold_right_i fold coefs Zero
in aux p
(* approximation of the Horner form used in the tactic ring *)
@@ -355,9 +323,11 @@ let pol_sparse_to_term n2 p =
match p with
[] -> const (num_of_string "0")
| (a,m)::p1 ->
+ let m = Ideal.Monomial.repr m in
let n = (Array.length m)-1 in
let (i0,e0) =
List.fold_left (fun (r,d) (a,m) ->
+ let m = Ideal.Monomial.repr m in
let i0= ref 0 in
for k=1 to n do
if m.(k)>0
@@ -374,45 +344,28 @@ let pol_sparse_to_term n2 p =
p in
if Int.equal i0 0
then
- let mp = ref (polrec_to_term a) in
- if List.is_empty p1
- then !mp
- else add(!mp,aux p1)
- else (
- let p1=ref [] in
- let p2=ref [] in
- List.iter
- (fun (a,m) ->
- if m.(i0)>=e0
- then (m.(i0)<-m.(i0)-e0;
- p1:=(a,m)::(!p1))
- else p2:=(a,m)::(!p2))
- p;
+ let mp = polrec_to_term a in
+ if List.is_empty p1 then mp else add (mp, aux p1)
+ else
+ let fold (p1, p2) (a, m) =
+ if (Ideal.Monomial.repr m).(i0) >= e0 then begin
+ let m0 = Array.copy (Ideal.Monomial.repr m) in
+ let () = m0.(i0) <- m0.(i0) - e0 in
+ let m0 = Ideal.Monomial.make m0 in
+ ((a, m0) :: p1, p2)
+ end else
+ (p1, (a, m) :: p2)
+ in
+ let (p1, p2) = List.fold_left fold ([], []) p in
let vm =
if Int.equal e0 1
then Var (string_of_int (i0))
else pow (Var (string_of_int (i0)),e0) in
- add(mul(vm, aux (List.rev (!p1))), aux (List.rev (!p2))))
+ add (mul(vm, aux (List.rev p1)), aux (List.rev p2))
in (*info "-> pol_sparse_to_term\n";*)
aux p
-let remove_list_tail l i =
- let rec aux l i =
- if List.is_empty l
- then []
- else if i<0
- then l
- else if Int.equal i 0
- then List.tl l
- else
- match l with
- |(a::l1) ->
- a::(aux l1 (i-1))
- |_ -> assert false
- in
- List.rev (aux (List.rev l) i)
-
(*
lq = [cn+m+1 n+m ...cn+m+1 1]
lci=[[cn+1 n,...,cn1 1]
@@ -422,49 +375,35 @@ let remove_list_tail l i =
removes intermediate polynomials not useful to compute the last one.
*)
-let remove_zeros zero lci =
- let n = List.length (List.hd lci) in
- let m=List.length lci in
+let remove_zeros lci =
+ let m = List.length lci in
let u = Array.make m false in
let rec utiles k =
- if k>=m
- then ()
- else (
- u.(k)<-true;
+ (** TODO: Find a more reasonable implementation of this traversal. *)
+ if k >= m || u.(k) then ()
+ else
+ let () = u.(k) <- true in
let lc = List.nth lci k in
- for i=0 to List.length lc - 1 do
- if not (zero (List.nth lc i))
- then utiles (i+k+1);
- done)
- in utiles 0;
- let lr = ref [] in
- for i=0 to m-1 do
- if u.(i)
- then lr:=(List.nth lci i)::(!lr)
- done;
- let lr=List.rev !lr in
- let lr = List.map
- (fun lc ->
- let lcr=ref lc in
- for i=0 to m-1 do
- if not u.(i)
- then lcr:=remove_list_tail !lcr (m-i+(n-m))
- done;
- !lcr)
- lr in
- info ("useless spolynomials: "
- ^string_of_int (m-List.length lr)^"\n");
- info ("useful spolynomials: "
- ^string_of_int (List.length lr)^"\n");
+ let iter i c = if not (PIdeal.equal c zeroP) then utiles (i + k + 1) in
+ List.iteri iter lc
+ in
+ let () = utiles 0 in
+ let filter i l =
+ let f j l = if m <= i + j + 1 then true else u.(i + j + 1) in
+ if u.(i) then Some (List.filteri f l)
+ else None
+ in
+ let lr = CList.map_filter_i filter lci in
+ info (fun () -> Printf.sprintf "useless spolynomials: %i" (m-List.length lr));
+ info (fun () -> Printf.sprintf "useful spolynomials: %i " (List.length lr));
lr
-let theoremedeszeros lpol p =
+let theoremedeszeros metadata nvars lpol p =
let t1 = Unix.gettimeofday() in
- let m = !nvars in
- let (lp0,p,cert) = in_ideal m lpol p in
- let lpc = List.rev !poldepcontent in
- info ("time: "^Format.sprintf "@[%10.3f@]s\n" (Unix.gettimeofday ()-.t1));
- (cert,lp0,p,lpc)
+ let m = nvars in
+ let cert = in_ideal metadata m lpol p in
+ info (fun () -> Printf.sprintf "time: @[%10.3f@]s" (Unix.gettimeofday ()-.t1));
+ cert
open Ideal
@@ -474,7 +413,7 @@ open Ideal
that has the same size than lp and where true indicates an
element that has been removed
*)
-let rec clean_pol lp =
+let clean_pol lp =
let t = Hashpol.create 12 in
let find p = try Hashpol.find t p
with
@@ -507,51 +446,33 @@ let expand_pol lb lp =
in List.rev (aux lb (List.rev lp))
let theoremedeszeros_termes lp =
- nvars:=0;(* mise a jour par term_pol_sparse *)
- List.iter set_nvars_term lp;
+ let nvars = List.fold_left set_nvars_term 0 lp in
match lp with
| Const (Int sugarparam)::Const (Int nparam)::lp ->
((match sugarparam with
- |0 -> info "computation without sugar\n";
+ |0 -> sinfo "computation without sugar";
lexico:=false;
- sugar_flag := false;
- divide_rem_with_critical_pair := false
- |1 -> info "computation with sugar\n";
+ |1 -> sinfo "computation with sugar";
lexico:=false;
- sugar_flag := true;
- divide_rem_with_critical_pair := false
- |2 -> info "ordre lexico computation without sugar\n";
+ |2 -> sinfo "ordre lexico computation without sugar";
lexico:=true;
- sugar_flag := false;
- divide_rem_with_critical_pair := false
- |3 -> info "ordre lexico computation with sugar\n";
+ |3 -> sinfo "ordre lexico computation with sugar";
lexico:=true;
- sugar_flag := true;
- divide_rem_with_critical_pair := false
- |4 -> info "computation without sugar, division by pairs\n";
+ |4 -> sinfo "computation without sugar, division by pairs";
lexico:=false;
- sugar_flag := false;
- divide_rem_with_critical_pair := true
- |5 -> info "computation with sugar, division by pairs\n";
+ |5 -> sinfo "computation with sugar, division by pairs";
lexico:=false;
- sugar_flag := true;
- divide_rem_with_critical_pair := true
- |6 -> info "ordre lexico computation without sugar, division by pairs\n";
+ |6 -> sinfo "ordre lexico computation without sugar, division by pairs";
lexico:=true;
- sugar_flag := false;
- divide_rem_with_critical_pair := true
- |7 -> info "ordre lexico computation with sugar, division by pairs\n";
+ |7 -> sinfo "ordre lexico computation with sugar, division by pairs";
lexico:=true;
- sugar_flag := true;
- divide_rem_with_critical_pair := true
- | _ -> error "nsatz: bad parameter"
+ | _ -> user_err Pp.(str "nsatz: bad parameter")
);
- let m= !nvars in
- let lvar=ref [] in
- for i=m downto 1 do lvar:=["x"^(string_of_int i)^""]@(!lvar); done;
- lvar:=["a";"b";"c";"d";"e";"f";"g";"h";"i";"j";"k";"l";"m";"n";"o";"p";"q";"r";"s";"t";"u";"v";"w";"x";"y";"z"] @ (!lvar); (* pour macaulay *)
- name_var:=!lvar;
- let lp = List.map (term_pol_sparse nparam) lp in
+ let lvar = List.init nvars (fun i -> Printf.sprintf "x%i" (i + 1)) in
+ let lvar = ["a";"b";"c";"d";"e";"f";"g";"h";"i";"j";"k";"l";"m";"n";"o";"p";"q";"r";"s";"t";"u";"v";"w";"x";"y";"z"] @ lvar in
+ (* pour macaulay *)
+ let metadata = { name_var = lvar } in
+ let lp = List.map (term_pol_sparse nvars nparam) lp in
match lp with
| [] -> assert false
| p::lp1 ->
@@ -561,16 +482,16 @@ let theoremedeszeros_termes lp =
lb is kept in order to fix the certificate in the post-processing
*)
let lpol, lb = clean_pol lpol in
- let (cert,lp0,p,_lct) = theoremedeszeros lpol p in
- info "cert ok\n";
+ let cert = theoremedeszeros metadata nvars lpol p in
+ sinfo "cert ok";
let lc = cert.last_comb::List.rev cert.gb_comb in
- match remove_zeros (fun x -> equal x zeroP) lc with
+ match remove_zeros lc with
| [] -> assert false
| (lq::lci) ->
(* post-processing : we apply the correction for the last line *)
let lq = expand_pol lb lq in
(* lci commence par les nouveaux polynomes *)
- let m = !nvars in
+ let m = nvars in
let c = pol_sparse_to_term m (polconst m cert.coef) in
let r = Pow(Zero,cert.power) in
let lci = List.rev lci in
@@ -578,8 +499,8 @@ let theoremedeszeros_termes lp =
let lci = List.map (expand_pol lb) lci in
let lci = List.map (List.map (pol_sparse_to_term m)) lci in
let lq = List.map (pol_sparse_to_term m) lq in
- info ("number of parameters: "^string_of_int nparam^"\n");
- info "term computed\n";
+ info (fun () -> Printf.sprintf "number of parameters: %i" nparam);
+ sinfo "term computed";
(c,r,lci,lq)
)
|_ -> assert false
@@ -619,19 +540,18 @@ let nsatz lpol =
mkt_app lcons [tlp ();ltterm;r])
res
(mkt_app lnil [tlp ()]) in
- info "term computed\n";
+ sinfo "term computed";
res
let return_term t =
let a =
- mkApp(gen_constant "CC" ["Init";"Logic"] "refl_equal",[|tllp ();t|]) in
+ mkApp(gen_constant "CC" ["Init";"Logic"] "eq_refl",[|tllp ();t|]) in
+ let a = EConstr.of_constr a in
generalize [a]
let nsatz_compute t =
let lpol =
try nsatz t
with Ideal.NotInIdeal ->
- error "nsatz cannot solve this problem" in
+ user_err Pp.(str "nsatz cannot solve this problem") in
return_term lpol
-
-