diff options
Diffstat (limited to 'contrib')
| -rw-r--r-- | contrib/fourier/Fourier.v | 27 | ||||
| -rw-r--r-- | contrib/fourier/Fourier_util.v | 261 | ||||
| -rw-r--r-- | contrib/fourier/fourier.ml | 202 | ||||
| -rw-r--r-- | contrib/fourier/fourierR.ml | 543 |
4 files changed, 1033 insertions, 0 deletions
diff --git a/contrib/fourier/Fourier.v b/contrib/fourier/Fourier.v new file mode 100644 index 0000000000..67d734281d --- /dev/null +++ b/contrib/fourier/Fourier.v @@ -0,0 +1,27 @@ +(***********************************************************************) +(* v * The Coq Proof Assistant / The Coq Development Team *) +(* <O___,, * INRIA-Rocquencourt & LRI-CNRS-Orsay *) +(* \VV/ *************************************************************) +(* // * This file is distributed under the terms of the *) +(* * GNU Lesser General Public License Version 2.1 *) +(***********************************************************************) + +(* $Id$ *) + +(* "Fourier's method to solve linear inequations/equations systems.".*) + +Declare ML Module "quote". +Declare ML Module "ring". +Declare ML Module "fourier". +Declare ML Module "fourierR". + +Require Export Fourier_util. + +Grammar tactic simple_tactic:ast:= + fourier + ["Fourier" constrarg_list($arg)] -> + [(Fourier ($LIST $arg))]. + +Tactic Definition FourierEq := + Apply Rge_ge_eq ; Fourier. + diff --git a/contrib/fourier/Fourier_util.v b/contrib/fourier/Fourier_util.v new file mode 100644 index 0000000000..91134f022f --- /dev/null +++ b/contrib/fourier/Fourier_util.v @@ -0,0 +1,261 @@ +(***********************************************************************) +(* v * The Coq Proof Assistant / The Coq Development Team *) +(* <O___,, * INRIA-Rocquencourt & LRI-CNRS-Orsay *) +(* \VV/ *************************************************************) +(* // * This file is distributed under the terms of the *) +(* * GNU Lesser General Public License Version 2.1 *) +(***********************************************************************) + +(* $Id$ *) + +Require Export Rbase. +Comments "Lemmas used by the tactic Fourier". + +Lemma Rfourier_lt: + (x1, y1, a : R) (Rlt x1 y1) -> (Rlt R0 a) -> (Rlt (Rmult a x1) (Rmult a y1)). +Intros x1 y1 a H H0; Try Assumption. +Apply Rlt_monotony; Auto with real. +Qed. + +Lemma Rfourier_le: + (x1, y1, a : R) (Rle x1 y1) -> (Rlt R0 a) -> (Rle (Rmult a x1) (Rmult a y1)). +Red. +Intros. +Case H; Auto with real. +Qed. + +Lemma Rfourier_lt_lt: + (x1, y1, x2, y2, a : R) + (Rlt x1 y1) -> + (Rlt x2 y2) -> + (Rlt R0 a) -> (Rlt (Rplus x1 (Rmult a x2)) (Rplus y1 (Rmult a y2))). +Intros x1 y1 x2 y2 a H H0 H1; Try Assumption. +Apply Rplus_lt. +Try Exact H. +Apply Rfourier_lt. +Try Exact H0. +Try Exact H1. +Qed. + +Lemma Rfourier_lt_le: + (x1, y1, x2, y2, a : R) + (Rlt x1 y1) -> + (Rle x2 y2) -> + (Rlt R0 a) -> (Rlt (Rplus x1 (Rmult a x2)) (Rplus y1 (Rmult a y2))). +Intros x1 y1 x2 y2 a H H0 H1; Try Assumption. +Case H0; Intros. +Apply Rplus_lt. +Try Exact H. +Apply Rfourier_lt; Auto with real. +Rewrite H2. +Rewrite (Rplus_sym y1 (Rmult a y2)). +Rewrite (Rplus_sym x1 (Rmult a y2)). +Apply Rlt_compatibility. +Try Exact H. +Qed. + +Lemma Rfourier_le_lt: + (x1, y1, x2, y2, a : R) + (Rle x1 y1) -> + (Rlt x2 y2) -> + (Rlt R0 a) -> (Rlt (Rplus x1 (Rmult a x2)) (Rplus y1 (Rmult a y2))). +Intros x1 y1 x2 y2 a H H0 H1; Try Assumption. +Case H; Intros. +Apply Rfourier_lt_le; Auto with real. +Rewrite H2. +Apply Rlt_compatibility. +Apply Rfourier_lt; Auto with real. +Qed. + +Lemma Rfourier_le_le: + (x1, y1, x2, y2, a : R) + (Rle x1 y1) -> + (Rle x2 y2) -> + (Rlt R0 a) -> (Rle (Rplus x1 (Rmult a x2)) (Rplus y1 (Rmult a y2))). +Intros x1 y1 x2 y2 a H H0 H1; Try Assumption. +Case H0; Intros. +Red. +Left; Try Assumption. +Apply Rfourier_le_lt; Auto with real. +Rewrite H2. +Case H; Intros. +Red. +Left; Try Assumption. +Rewrite (Rplus_sym x1 (Rmult a y2)). +Rewrite (Rplus_sym y1 (Rmult a y2)). +Apply Rlt_compatibility. +Try Exact H3. +Rewrite H3. +Red. +Right; Try Assumption. +Auto with real. +Qed. + +Lemma Rlt_zero_pos_plus1: (x : R) (Rlt R0 x) -> (Rlt R0 (Rplus R1 x)). +Intros x H; Try Assumption. +Rewrite Rplus_sym. +Apply Rlt_r_plus_R1. +Red; Auto with real. +Qed. + +Lemma Rlt_mult_inv_pos: + (x, y : R) (Rlt R0 x) -> (Rlt R0 y) -> (Rlt R0 (Rmult x (Rinv y))). +Intros x y H H0; Try Assumption. +Replace R0 with (Rmult x R0). +Apply Rlt_monotony; Auto with real. +Ring. +Qed. + +Lemma Rlt_zero_1: (Rlt R0 R1). +Exact Rlt_R0_R1. +Qed. + +Lemma Rle_zero_pos_plus1: (x : R) (Rle R0 x) -> (Rle R0 (Rplus R1 x)). +Intros x H; Try Assumption. +Case H; Intros. +Red. +Left; Try Assumption. +Apply Rlt_zero_pos_plus1; Auto with real. +Rewrite <- H0. +Replace (Rplus R1 R0) with R1. +Red; Left. +Exact Rlt_zero_1. +Ring. +Qed. + +Lemma Rle_mult_inv_pos: + (x, y : R) (Rle R0 x) -> (Rlt R0 y) -> (Rle R0 (Rmult x (Rinv y))). +Intros x y H H0; Try Assumption. +Case H; Intros. +Red; Left. +Apply Rlt_mult_inv_pos; Auto with real. +Rewrite <- H1. +Red; Right; Ring. +Qed. + +Lemma Rle_zero_1: (Rle R0 R1). +Red; Left. +Exact Rlt_zero_1. +Qed. + +Lemma Rle_not_lt: + (n, d : R) (Rle R0 (Rmult n (Rinv d))) -> ~ (Rlt R0 (Rmult (Ropp n) (Rinv d))). +Intros n d H; Red; Intros H0; Try Exact H0. +Generalize (Rle_not R0 (Rmult n (Rinv d))). +Intros H1; Elim H1; Try Assumption. +Replace (Rmult n (Rinv d)) with (Ropp (Ropp (Rmult n (Rinv d)))). +Replace R0 with (Ropp (Ropp R0)). +Replace (Ropp (Rmult n (Rinv d))) with (Rmult (Ropp n) (Rinv d)). +Replace (Ropp R0) with R0. +Red. +Try Exact H0. +Apply Rgt_Ropp. +Replace (Ropp (Rmult n (Rinv d))) with (Rmult (Ropp n) (Rinv d)). +Replace (Ropp R0) with R0. +Red. +Try Exact H0. +Ring. +Ring. +Ring. +Ring. +Ring. +Ring. +Qed. + +Lemma Rnot_lt0: (x : R) ~ (Rlt R0 (Rmult R0 x)). +Intros x; Try Assumption. +Replace (Rmult R0 x) with R0. +Apply Rlt_antirefl. +Ring. +Qed. + +Lemma Rlt_not_le: + (n, d : R) (Rlt R0 (Rmult n (Rinv d))) -> ~ (Rle R0 (Rmult (Ropp n) (Rinv d))). +Intros n d H; Try Assumption. +Apply Rle_not. +Replace R0 with (Ropp R0). +Replace (Rmult (Ropp n) (Rinv d)) with (Ropp (Rmult n (Rinv d))). +Apply Rlt_Ropp. +Try Exact H. +Ring. +Ring. +Qed. + +Lemma Rnot_lt_lt: (x, y : R) ~ (Rlt R0 (Rminus y x)) -> ~ (Rlt x y). +Unfold not; Intros. +Apply H. +Apply Rlt_anti_compatibility with x. +Replace (Rplus x R0) with x. +Replace (Rplus x (Rminus y x)) with y. +Try Exact H0. +Ring. +Ring. +Qed. + +Lemma Rnot_le_le: (x, y : R) ~ (Rle R0 (Rminus y x)) -> ~ (Rle x y). +Unfold not; Intros. +Apply H. +Case H0; Intros. +Left. +Apply Rlt_anti_compatibility with x. +Replace (Rplus x R0) with x. +Replace (Rplus x (Rminus y x)) with y. +Try Exact H1. +Ring. +Ring. +Right. +Rewrite H1; Ring. +Qed. + +Lemma Rfourier_gt_to_lt: (x, y : R) (Rgt y x) -> (Rlt x y). +Unfold Rgt; Auto with *. +Qed. + +Lemma Rfourier_ge_to_le: (x, y : R) (Rge y x) -> (Rle x y). +Unfold Rge Rle; Intuition. +Qed. + +Lemma Rfourier_eqLR_to_le: (x, y : R) x == y -> (Rle x y). +Intros. +Rewrite H; Red. +Right; Auto with *. +Qed. + +Lemma Rfourier_eqRL_to_le: (x, y : R) y == x -> (Rle x y). +Intros. +Rewrite H; Red. +Right; Auto with *. +Qed. + +Lemma Rfourier_not_ge_lt: (x, y : R) ((Rge x y) -> False) -> (Rlt x y). +Intros. +Unfold Rge in H. +Elim (total_order x y); Intros. +Try Exact H0. +Intuition. +Qed. + +Lemma Rfourier_not_gt_le: (x, y : R) ((Rgt x y) -> False) -> (Rle x y). +Intros. +Red. +Elim (total_order x y); Intros. +Intuition. +Intuition. +Qed. + +Lemma Rfourier_not_le_gt: (x, y : R) ((Rle x y) -> False) -> (Rgt x y). +Unfold Rle. +Intros. +Red. +Elim (total_order x y); Intros. +Intuition. +Intuition. +Qed. + +Lemma Rfourier_not_lt_ge: (x, y : R) ((Rlt x y) -> False) -> (Rge x y). +Intros. +Red. +Elim (total_order x y); Intros. +Intuition. +Intuition. +Qed. diff --git a/contrib/fourier/fourier.ml b/contrib/fourier/fourier.ml new file mode 100644 index 0000000000..ea124ed17e --- /dev/null +++ b/contrib/fourier/fourier.ml @@ -0,0 +1,202 @@ +(***********************************************************************) +(* v * The Coq Proof Assistant / The Coq Development Team *) +(* <O___,, * INRIA-Rocquencourt & LRI-CNRS-Orsay *) +(* \VV/ *************************************************************) +(* // * This file is distributed under the terms of the *) +(* * GNU Lesser General Public License Version 2.1 *) +(***********************************************************************) + +(* $Id$ *) + +(* Méthode d'élimination de Fourier *) +(* Référence: +Auteur(s) : Fourier, Jean-Baptiste-Joseph + +Titre(s) : Oeuvres de Fourier [Document électronique]. Tome second. Mémoires publiés dans divers recueils / publ. par les soins de M. Gaston Darboux,... + +Publication : Numérisation BnF de l'édition de Paris : Gauthier-Villars, 1890 + +Pages: 326-327 + +http://gallica.bnf.fr/ +*) + +(* Un peu de calcul sur les rationnels... +Les opérations rendent des rationnels normalisés, +i.e. le numérateur et le dénominateur sont premiers entre eux. +*) +type rational = {num:int; + den:int} +;; +let print_rational x = + print_int x.num; + print_string "/"; + print_int x.den +;; + +let rec pgcd x y = if y = 0 then x else pgcd y (x mod y);; + + +let r0 = {num=0;den=1};; +let r1 = {num=1;den=1};; + +let rnorm x = let x = (if x.den<0 then {num=(-x.num);den=(-x.den)} else x) in + if x.num=0 then r0 + else (let d=pgcd x.num x.den in + {num=(x.num)/d;den=(x.den)/d});; + +let rop x = rnorm {num=(-x.num);den=x.den};; + +let rplus x y = rnorm {num=x.num*y.den + y.num*x.den;den=x.den*y.den};; + +let rminus x y = rnorm {num=x.num*y.den - y.num*x.den;den=x.den*y.den};; + +let rmult x y = rnorm {num=x.num*y.num;den=x.den*y.den};; + +let rinv x = rnorm {num=x.den;den=x.num};; + +let rinf x y = x.num*y.den < y.num*x.den;; +let rinfeq x y = x.num*y.den <= y.num*x.den;; + +(* {coef;hist;strict}, où coef=[c1; ...; cn; d], représente l'inéquation +c1x1+...+cnxn < d si strict=true, <= sinon, +hist donnant les coefficients (positifs) d'une combinaison linéaire qui permet d'obtenir l'inéquation à partir de celles du départ. +*) + +type ineq = {coef:rational list; + hist:rational list; + strict:bool};; + +let pop x l = l:=x::(!l);; + +(* sépare la liste d'inéquations s selon que leur premier coefficient est +négatif, nul ou positif. *) +let partitionne s = + let lpos=ref [] in + let lneg=ref [] in + let lnul=ref [] in + List.iter (fun ie -> match ie.coef with + [] -> raise (Failure "empty ineq") + |(c::r) -> if rinf c r0 + then pop ie lneg + else if rinf r0 c then pop ie lpos + else pop ie lnul) + s; + [!lneg;!lnul;!lpos] +;; +(* initialise les histoires d'une liste d'inéquations données par leurs listes de coefficients et leurs strictitudes (!): +(add_hist [(equation 1, s1);...;(équation n, sn)]) += +[{équation 1, [1;0;...;0], s1}; + {équation 2, [0;1;...;0], s2}; + ... + {équation n, [0;0;...;1], sn}] +*) +let add_hist le = + let n = List.length le in + let i=ref 0 in + List.map (fun (ie,s) -> + let h =ref [] in + for k=1 to (n-(!i)-1) do pop r0 h; done; + pop r1 h; + for k=1 to !i do pop r0 h; done; + i:=!i+1; + {coef=ie;hist=(!h);strict=s}) + le +;; +(* additionne deux inéquations *) +let ie_add ie1 ie2 = {coef=List.map2 rplus ie1.coef ie2.coef; + hist=List.map2 rplus ie1.hist ie2.hist; + strict=ie1.strict || ie2.strict} +;; +(* multiplication d'une inéquation par un rationnel (positif) *) +let ie_emult a ie = {coef=List.map (fun x -> rmult a x) ie.coef; + hist=List.map (fun x -> rmult a x) ie.hist; + strict= ie.strict} +;; +(* on enlève le premier coefficient *) +let ie_tl ie = {coef=List.tl ie.coef;hist=ie.hist;strict=ie.strict} +;; +(* le premier coefficient: "tête" de l'inéquation *) +let hd_coef ie = List.hd ie.coef +;; + +(* calcule toutes les combinaisons entre inéquations de tête négative et inéquations de tête positive qui annulent le premier coefficient. +*) +let deduce_add lneg lpos = + let res=ref [] in + List.iter (fun i1 -> + List.iter (fun i2 -> + let a = rop (hd_coef i1) in + let b = hd_coef i2 in + pop (ie_tl (ie_add (ie_emult b i1) + (ie_emult a i2))) res) + lpos) + lneg; + !res +;; +(* élimination de la première variable à partir d'une liste d'inéquations: +opération qu'on itère dans l'algorithme de Fourier. +*) +let deduce1 s = + match (partitionne s) with + [lneg;lnul;lpos] -> + let lnew = deduce_add lneg lpos in + (List.map ie_tl lnul)@lnew + |_->assert false +;; +(* algorithme de Fourier: on élimine successivement toutes les variables. +*) +let deduce lie = + let n = List.length (fst (List.hd lie)) in + let lie=ref (add_hist lie) in + for i=1 to n-1 do + lie:= deduce1 !lie; + done; + !lie +;; + +(* donne [] si le système a des solutions, +sinon donne [c,s,lc] +où lc est la combinaison linéaire des inéquations de départ +qui donne 0 < c si s=true + ou 0 <= c sinon +cette inéquation étant absurde. +*) +let unsolvable lie = + let lr = deduce lie in + let res = ref [] in + (try (List.iter (fun e -> + match e with + {coef=[c];hist=lc;strict=s} -> + if (rinf c r0 && (not s)) || (rinfeq c r0 && s) + then (res := [c,s,lc]; + raise (Failure "contradiction found")) + |_->assert false) + lr) + with _ -> ()); + !res +;; + +(* Exemples: + +let test1=[[r1;r1;r0],true;[rop r1;r1;r1],false;[r0;rop r1;rop r1],false];; +deduce test1;; +unsolvable test1;; + +let test2=[ +[r1;r1;r0;r0;r0],false; +[r0;r1;r1;r0;r0],false; +[r0;r0;r1;r1;r0],false; +[r0;r0;r0;r1;r1],false; +[r1;r0;r0;r0;r1],false; +[rop r1;rop r1;r0;r0;r0],false; +[r0;rop r1;rop r1;r0;r0],false; +[r0;r0;rop r1;rop r1;r0],false; +[r0;r0;r0;rop r1;rop r1],false; +[rop r1;r0;r0;r0;rop r1],false +];; +deduce test2;; +unsolvable test2;; + +*)
\ No newline at end of file diff --git a/contrib/fourier/fourierR.ml b/contrib/fourier/fourierR.ml new file mode 100644 index 0000000000..c85ead5362 --- /dev/null +++ b/contrib/fourier/fourierR.ml @@ -0,0 +1,543 @@ +(***********************************************************************) +(* v * The Coq Proof Assistant / The Coq Development Team *) +(* <O___,, * INRIA-Rocquencourt & LRI-CNRS-Orsay *) +(* \VV/ *************************************************************) +(* // * This file is distributed under the terms of the *) +(* * GNU Lesser General Public License Version 2.1 *) +(***********************************************************************) + +(* $Id$ *) + + + +(* La tactique Fourier ne fonctionne de manière sûre que si les coefficients +des inéquations et équations sont entiers. En attendant la tactique Field. +*) + +open Term +open Tactics +open Clenv +open Names +open Tacmach +open Fourier + +(****************************************************************************** +Opérations sur les combinaisons linéaires affines. +La partie homogène d'une combinaison linéaire est en fait une table de hash +qui donne le coefficient d'un terme du calcul des constructions, +qui est zéro si le terme n'y est pas. +*) + +type flin = {fhom:(constr , rational)Hashtbl.t; + fcste:rational};; + +let flin_zero () = {fhom=Hashtbl.create 50;fcste=r0};; + +let flin_coef f x = try (Hashtbl.find f.fhom x) with _-> r0;; + +let flin_add f x c = + let cx = flin_coef f x in + Hashtbl.remove f.fhom x; + Hashtbl.add f.fhom x (rplus cx c); + f +;; +let flin_add_cste f c = + {fhom=f.fhom; + fcste=rplus f.fcste c} +;; + +let flin_one () = flin_add_cste (flin_zero()) r1;; + +let flin_plus f1 f2 = + let f3 = flin_zero() in + Hashtbl.iter (fun x c -> let _=flin_add f3 x c in ()) f1.fhom; + Hashtbl.iter (fun x c -> let _=flin_add f3 x c in ()) f2.fhom; + flin_add_cste (flin_add_cste f3 f1.fcste) f2.fcste; +;; + +let flin_minus f1 f2 = + let f3 = flin_zero() in + Hashtbl.iter (fun x c -> let _=flin_add f3 x c in ()) f1.fhom; + Hashtbl.iter (fun x c -> let _=flin_add f3 x (rop c) in ()) f2.fhom; + flin_add_cste (flin_add_cste f3 f1.fcste) (rop f2.fcste); +;; +let flin_emult a f = + let f2 = flin_zero() in + Hashtbl.iter (fun x c -> let _=flin_add f2 x (rmult a c) in ()) f.fhom; + flin_add_cste f2 (rmult a f.fcste); +;; + +(*****************************************************************************) +let parse_ast = Pcoq.parse_string Pcoq.Constr.constr;; +let parse s = Astterm.interp_constr Evd.empty (Global.env()) (parse_ast s);; +let pf_parse_constr gl s = + Astterm.interp_constr Evd.empty (pf_env gl) (parse_ast s);; + +let rec string_of_constr c = + match kind_of_term c with + IsCast (c,t) -> string_of_constr c + |IsConst (c,l) -> string_of_path c + |IsVar(c) -> string_of_id c + | _ -> "not_of_constant" +;; + +let rec rational_of_constr c = + match kind_of_term c with + | IsCast (c,t) -> (rational_of_constr c) + | IsApp (c,args) -> + (match kind_of_term c with + IsConst (c,l) -> + (match (string_of_path c) with + "Coq.Reals.Rdefinitions.Ropp" -> + rop (rational_of_constr args.(0)) + |"Coq.Reals.Rdefinitions.Rinv" -> + rinv (rational_of_constr args.(0)) + |"Coq.Reals.Rdefinitions.Rmult" -> + rmult (rational_of_constr args.(0)) + (rational_of_constr args.(1)) + |"Coq.Reals.Rdefinitions.Rplus" -> + rplus (rational_of_constr args.(0)) + (rational_of_constr args.(1)) + |"Coq.Reals.Rdefinitions.Rminus" -> + rminus (rational_of_constr args.(0)) + (rational_of_constr args.(1)) + | _ -> failwith "not a rational") + | _ -> failwith "not a rational") + | IsConst (c,l) -> + (match (string_of_path c) with + "Coq.Reals.Rdefinitions.R1" -> r1 + |"Coq.Reals.Rdefinitions.R0" -> r0 + | _ -> failwith "not a rational") + | _ -> failwith "not a rational" +;; + +let rec flin_of_constr c = + try( + match kind_of_term c with + | IsCast (c,t) -> (flin_of_constr c) + | IsApp (c,args) -> + (match kind_of_term c with + IsConst (c,l) -> + (match (string_of_path c) with + "Coq.Reals.Rdefinitions.Ropp" -> + flin_emult (rop r1) (flin_of_constr args.(0)) + |"Coq.Reals.Rdefinitions.Rplus"-> + flin_plus (flin_of_constr args.(0)) + (flin_of_constr args.(1)) + |"Coq.Reals.Rdefinitions.Rminus"-> + flin_minus (flin_of_constr args.(0)) + (flin_of_constr args.(1)) + |"Coq.Reals.Rdefinitions.Rmult"-> + (try (let a=(rational_of_constr args.(0)) in + try (let b = (rational_of_constr args.(1)) in + (flin_add_cste (flin_zero()) (rmult a b))) + with _-> (flin_add (flin_zero()) + args.(1) + a)) + with _-> (flin_add (flin_zero()) + args.(0) + (rational_of_constr args.(1)))) + |_->assert false) + |_ -> assert false) + | IsConst (c,l) -> + (match (string_of_path c) with + "Coq.Reals.Rdefinitions.R1" -> flin_one () + |"Coq.Reals.Rdefinitions.R0" -> flin_zero () + |_-> assert false) + |_-> assert false) + with _ -> flin_add (flin_zero()) + c + r1 +;; + +let flin_to_alist f = + let res=ref [] in + Hashtbl.iter (fun x c -> res:=(c,x)::(!res)) f; + !res +;; + +(* Représentation des hypothèses qui sont des inéquations ou des équations. +*) +type hineq={hname:constr; (* le nom de l'hypothèse *) + htype:string; (* Rlt, Rgt, Rle, Rge, eqTLR ou eqTRL *) + hleft:constr; + hright:constr; + hflin:flin; + hstrict:bool} +;; + +(* Transforme une hypothese h:t en inéquation flin<0 ou flin<=0 +*) +let ineq1_of_constr (h,t) = + match (kind_of_term t) with + IsApp (f,args) -> + let t1= args.(0) in + let t2= args.(1) in + (match kind_of_term f with + IsConst (c,l) -> + (match (string_of_path c) with + "Coq.Reals.Rdefinitions.Rlt" -> [{hname=h; + htype="Rlt"; + hleft=t1; + hright=t2; + hflin= flin_minus (flin_of_constr t1) + (flin_of_constr t2); + hstrict=true}] + |"Coq.Reals.Rdefinitions.Rgt" -> [{hname=h; + htype="Rgt"; + hleft=t2; + hright=t1; + hflin= flin_minus (flin_of_constr t2) + (flin_of_constr t1); + hstrict=true}] + |"Coq.Reals.Rdefinitions.Rle" -> [{hname=h; + htype="Rle"; + hleft=t1; + hright=t2; + hflin= flin_minus (flin_of_constr t1) + (flin_of_constr t2); + hstrict=false}] + |"Coq.Reals.Rdefinitions.Rge" -> [{hname=h; + htype="Rge"; + hleft=t2; + hright=t1; + hflin= flin_minus (flin_of_constr t2) + (flin_of_constr t1); + hstrict=false}] + |_->assert false) + | IsMutInd ((sp,i),l) -> + (match (string_of_path sp) with + "Coq.Init.Logic_Type.eqT" -> let t0= args.(0) in + let t1= args.(1) in + let t2= args.(2) in + (match (kind_of_term t0) with + IsConst (c,l) -> + (match (string_of_path c) with + "Coq.Reals.Rdefinitions.R"-> + [{hname=h; + htype="eqTLR"; + hleft=t1; + hright=t2; + hflin= flin_minus (flin_of_constr t1) + (flin_of_constr t2); + hstrict=false}; + {hname=h; + htype="eqTRL"; + hleft=t2; + hright=t1; + hflin= flin_minus (flin_of_constr t2) + (flin_of_constr t1); + hstrict=false}] + |_-> assert false) + |_-> assert false) + |_-> assert false) + |_-> assert false) + |_-> assert false +;; + +(* Applique la méthode de Fourier à une liste d'hypothèses (type hineq) +*) + +let fourier_lineq lineq1 = + let nvar=ref (-1) in + let hvar=Hashtbl.create 50 in (* la table des variables des inéquations *) + List.iter (fun f -> + Hashtbl.iter (fun x c -> + try (Hashtbl.find hvar x;()) + with _-> nvar:=(!nvar)+1; + Hashtbl.add hvar x (!nvar)) + f.hflin.fhom) + lineq1; + let sys= List.map (fun h-> + let v=Array.create ((!nvar)+1) r0 in + Hashtbl.iter (fun x c -> v.(Hashtbl.find hvar x)<-c) + h.hflin.fhom; + ((Array.to_list v)@[rop h.hflin.fcste],h.hstrict)) + lineq1 in + unsolvable sys +;; + +(****************************************************************************** +Construction de la preuve en cas de succès de la méthode de Fourier, +i.e. on obtient une contradiction. +*) +let is_int x = (x.den)=1 +;; + +(* fraction = couple (num,den) *) +let rec rational_to_fraction x= (x.num,x.den) +;; + +(* traduction -3 -> (Ropp (Rplus R1 (Rplus R1 R1))) +*) +let int_to_real n = + let nn=abs n in + let s=ref "" in + if nn=0 + then s:="R0" + else (s:="R1"; + for i=1 to (nn-1) do s:="(Rplus R1 "^(!s)^")"; done;); + if n<0 then s:="(Ropp "^(!s)^")"; + !s +;; +(* -1/2 -> (Rmult (Ropp R1) (Rinv (Rplus R1 R1))) +*) +let rational_to_real x = + let (n,d)=rational_to_fraction x in + ("(Rmult "^(int_to_real n)^" (Rinv "^(int_to_real d)^"))") +;; + +(* preuve que 0<n*1/d +*) +let tac_zero_inf_pos gl (n,d) = + let cste = pf_parse_constr gl in + let tacn=ref (apply (cste "Rlt_zero_1")) in + let tacd=ref (apply (cste "Rlt_zero_1")) in + for i=1 to n-1 do + tacn:=(tclTHEN (apply (cste "Rlt_zero_pos_plus1")) !tacn); done; + for i=1 to d-1 do + tacd:=(tclTHEN (apply (cste "Rlt_zero_pos_plus1")) !tacd); done; + (tclTHENS (apply (cste "Rlt_mult_inv_pos")) [!tacn;!tacd]) +;; + +(* preuve que 0<=n*1/d +*) +let tac_zero_infeq_pos gl (n,d)= + let cste = pf_parse_constr gl in + let tacn=ref (if n=0 + then (apply (cste "Rle_zero_zero")) + else (apply (cste "Rle_zero_1"))) in + let tacd=ref (apply (cste "Rlt_zero_1")) in + for i=1 to n-1 do + tacn:=(tclTHEN (apply (cste "Rle_zero_pos_plus1")) !tacn); done; + for i=1 to d-1 do + tacd:=(tclTHEN (apply (cste "Rlt_zero_pos_plus1")) !tacd); done; + (tclTHENS (apply (cste "Rle_mult_inv_pos")) [!tacn;!tacd]) +;; + +(* preuve que 0<(-n)*(1/d) => False +*) +let tac_zero_inf_false gl (n,d) = + let cste = pf_parse_constr gl in + if n=0 then (apply (cste "Rnot_lt0")) + else + (tclTHEN (apply (cste "Rle_not_lt")) + (tac_zero_infeq_pos gl (-n,d))) +;; + +(* preuve que 0<=(-n)*(1/d) => False +*) +let tac_zero_infeq_false gl (n,d) = + let cste = pf_parse_constr gl in + (tclTHEN (apply (cste "Rlt_not_le")) + (tac_zero_inf_pos gl (-n,d))) +;; + +let create_meta () = mkMeta(new_meta());; + +let my_cut c gl= + let concl = pf_concl gl in + apply_type (mkProd(Anonymous,c,concl)) [create_meta()] gl +;; +let exact = exact_check;; + +let tac_use h = match h.htype with + "Rlt" -> exact h.hname + |"Rle" -> exact h.hname + |"Rgt" -> (tclTHEN (apply (parse "Rfourier_gt_to_lt")) + (exact h.hname)) + |"Rge" -> (tclTHEN (apply (parse "Rfourier_ge_to_le")) + (exact h.hname)) + |"eqTLR" -> (tclTHEN (apply (parse "Rfourier_eqLR_to_le")) + (exact h.hname)) + |"eqTRL" -> (tclTHEN (apply (parse "Rfourier_eqRL_to_le")) + (exact h.hname)) + |_->assert false +;; + +let is_ineq (h,t) = + match (kind_of_term t) with + IsApp (f,args) -> + (match (string_of_constr f) with + "Coq.Reals.Rdefinitions.Rlt" -> true + |"Coq.Reals.Rdefinitions.Rgt" -> true + |"Coq.Reals.Rdefinitions.Rle" -> true + |"Coq.Reals.Rdefinitions.Rge" -> true + |"Coq.Init.Logic_Type.eqT" -> (match (string_of_constr args.(0)) with + "Coq.Reals.Rdefinitions.R"->true + |_->false) + |_->false) + |_->false +;; + +let list_of_sign s = List.map (fun (x,_,z)->(x,z)) s;; + +let mkAppL a = + let l = Array.to_list a in + mkApp(List.hd l, Array.of_list (List.tl l)) +;; + +(* Résolution d'inéquations linéaires dans R *) +let rec fourier gl= + let parse = pf_parse_constr gl in + let goal = strip_outer_cast (pf_concl gl) in + let fhyp=id_of_string "new_hyp_for_fourier" in + (* si le but est une inéquation, on introduit son contraire, + et le but à prouver devient False *) + try (let tac = + match (kind_of_term goal) with + IsApp (f,args) -> + (match (string_of_constr f) with + "Coq.Reals.Rdefinitions.Rlt" -> + (tclTHEN + (tclTHEN (apply (parse "Rfourier_not_ge_lt")) + (intro_using fhyp)) + fourier) + |"Coq.Reals.Rdefinitions.Rle" -> + (tclTHEN + (tclTHEN (apply (parse "Rfourier_not_gt_le")) + (intro_using fhyp)) + fourier) + |"Coq.Reals.Rdefinitions.Rgt" -> + (tclTHEN + (tclTHEN (apply (parse "Rfourier_not_le_gt")) + (intro_using fhyp)) + fourier) + |"Coq.Reals.Rdefinitions.Rge" -> + (tclTHEN + (tclTHEN (apply (parse "Rfourier_not_lt_ge")) + (intro_using fhyp)) + fourier) + |_->assert false) + |_->assert false + in tac gl) + with _ -> + (* les hypothèses *) + let hyps = List.map (fun (h,t)-> (mkVar h,(body_of_type t))) + (list_of_sign (pf_hyps gl)) in + let lineq =ref [] in + List.iter (fun h -> try (lineq:=(ineq1_of_constr h)@(!lineq)) + with _-> ()) + hyps; + (* lineq = les inéquations découlant des hypothèses *) + let res=fourier_lineq (!lineq) in + let tac=ref tclIDTAC in + if res=[] + then (print_string "Tactic Fourier fails, perhaps because it cannot prove some equality on rationnal numbers with the tactic Ring."; + flush stdout) + (* l'algorithme de Fourier a réussi: on va en tirer une preuve Coq *) + else (match res with + [(cres,sres,lc)]-> + (* lc=coefficients multiplicateurs des inéquations + qui donnent 0<cres ou 0<=cres selon sres *) + (*print_string "Fourier's method can prove the goal...";flush stdout;*) + let lutil=ref [] in + List.iter + (fun (h,c) -> + if c<>r0 + then (lutil:=(h,c)::(!lutil)(*; + print_rational(c);print_string " "*))) + (List.combine (!lineq) lc); + (* on construit la combinaison linéaire des inéquation *) + (match (!lutil) with + (h1,c1)::lutil -> + let s=ref (h1.hstrict) in + let t1=ref (mkAppL [|parse "Rmult"; + parse (rational_to_real c1); + h1.hleft|]) in + let t2=ref (mkAppL [|parse "Rmult"; + parse (rational_to_real c1); + h1.hright|]) in + List.iter (fun (h,c) -> + s:=(!s)||(h.hstrict); + t1:=(mkAppL [|parse "Rplus"; + !t1; + mkAppL [|parse "Rmult"; + parse (rational_to_real c); + h.hleft|] |]); + t2:=(mkAppL [|parse "Rplus"; + !t2; + mkAppL [|parse "Rmult"; + parse (rational_to_real c); + h.hright|] |])) + lutil; + let ineq=mkAppL [|parse (if (!s) then "Rlt" else "Rle"); + !t1; + !t2 |] in + let tc=parse (rational_to_real cres) in + (* puis sa preuve *) + let tac1=ref (if h1.hstrict + then (tclTHENS (apply (parse "Rfourier_lt")) + [tac_use h1; + tac_zero_inf_pos gl + (rational_to_fraction c1)]) + else (tclTHENS (apply (parse "Rfourier_le")) + [tac_use h1; + tac_zero_inf_pos gl + (rational_to_fraction c1)])) in + s:=h1.hstrict; + List.iter (fun (h,c)-> + (if (!s) + then (if h.hstrict + then tac1:=(tclTHENS (apply (parse "Rfourier_lt_lt")) + [!tac1;tac_use h; + tac_zero_inf_pos gl + (rational_to_fraction c)]) + else tac1:=(tclTHENS (apply (parse "Rfourier_lt_le")) + [!tac1;tac_use h; + tac_zero_inf_pos gl + (rational_to_fraction c)])) + else (if h.hstrict + then tac1:=(tclTHENS (apply (parse "Rfourier_le_lt")) + [!tac1;tac_use h; + tac_zero_inf_pos gl + (rational_to_fraction c)]) + else tac1:=(tclTHENS (apply (parse "Rfourier_le_le")) + [!tac1;tac_use h; + tac_zero_inf_pos gl + (rational_to_fraction c)]))); + s:=(!s)||(h.hstrict)) + lutil; + let tac2= if sres + then tac_zero_inf_false gl (rational_to_fraction cres) + else tac_zero_infeq_false gl (rational_to_fraction cres) + in + tac:=(tclTHENS (my_cut ineq) + [tclTHEN (change_in_concl + (mkAppL [| parse "not"; ineq|] + )) + (tclTHEN (apply (if sres then parse "Rnot_lt_lt" + else parse "Rnot_le_le")) + (tclTHENS (Equality.replace + (mkAppL [|parse "Rminus";!t2;!t1|] + ) + tc) + [tac2; + (tclTHENS (Equality.replace (parse "(Rinv R1)") + (parse "R1")) +(* en attendant Field, ça peut aider Ring de remplacer 1/1 par 1 ... *) + + [(Ring.polynom []); + (tclTHEN (apply (parse "sym_eqT")) + (apply (parse "Rinv_R1")))] + + ) + ])); + !tac1]); + tac:=(tclTHENS (cut (parse "False")) + [tclTHEN intro contradiction; + !tac]) + |_-> assert false) |_-> assert false + ); +(* ((tclTHEN !tac (tclFAIL 1 (* 1 au hasard... *))) gl) *) + ((tclABSTRACT None !tac) gl) + +;; + +let fourier_tac x gl = + fourier gl +;; + +let v_fourier = add_tactic "Fourier" fourier_tac + + |