← Back to branch summary

~ubuntu-branches/ubuntu/wily/coq-doc/wily

« back to all changes in this revision

Viewing changes to contrib/micromega/csdpcert.ml

  • Committer: Bazaar Package Importer
  • Author(s): Stéphane Glondu, Stéphane Glondu, Samuel Mimram
  • Date: 2010-01-07 22:50:39 UTC
  • mfrom: (1.2.2 upstream)
  • Revision ID: james.westby@ubuntu.com-20100107225039-n3cq82589u0qt0s2
Tags: 8.2pl1-1
http://bugs.debian.org/563669
http://bugs.debian.org/543545
[ Stéphane Glondu ]
* New upstream release (Closes: #563669)
  - remove patches
* Packaging overhaul:
  - use git, advertise it in Vcs-* fields of debian/control
  - use debhelper 7 and dh with override
  - use source format 3.0 (quilt)
* debian/control:
  - set Maintainer to d-o-m, set Uploaders to Sam and myself
  - add Homepage field
  - bump Standards-Version to 3.8.3
* Register PDF documentation into doc-base
* Add debian/watch
* Update debian/copyright

[ Samuel Mimram ]
* Change coq-doc's description to mention that it provides documentation in
  pdf format, not postscript, closes: #543545.

Show diffs side-by-side

added added

removed removed

Lines of Context:
contrib/micromega/csdpcert.ml
 
1
(************************************************************************)
 
2
(*  v      *   The Coq Proof Assistant  /  The Coq Development Team     *)
 
3
(* <O___,, * CNRS-Ecole Polytechnique-INRIA Futurs-Universite Paris Sud *)
 
4
(*   \VV/  **************************************************************)
 
5
(*    //   *      This file is distributed under the terms of the       *)
 
6
(*         *       GNU Lesser General Public License Version 2.1        *)
 
7
(************************************************************************)
 
8
(*                                                                      *)
 
9
(* Micromega: A reflexive tactic using the Positivstellensatz           *)
 
10
(*                                                                      *)
 
11
(*  Frédéric Besson (Irisa/Inria) 2006-2008                             *)
 
12
(*                                                                      *)
 
13
(************************************************************************)
 
14
 
 
15
open Big_int
 
16
open Num
 
17
open Sos
 
18
 
 
19
module Mc = Micromega
 
20
module Ml2C = Mutils.CamlToCoq
 
21
module C2Ml = Mutils.CoqToCaml
 
22
 
 
23
let debug = false
 
24
 
 
25
module M =
 
26
struct
 
27
 open Mc
 
28
 
 
29
 let rec expr_to_term = function
 
30
  |  PEc z ->  Const  (C2Ml.q_to_num z)
 
31
  |  PEX v ->  Var ("x"^(string_of_int (C2Ml.index v)))
 
32
  |  PEmul(p1,p2) -> 
 
33
      let p1 = expr_to_term p1 in
 
34
      let p2 = expr_to_term p2 in
 
35
      let res = Mul(p1,p2) in      res
 
36
 
 
37
  | PEadd(p1,p2) -> Add(expr_to_term p1, expr_to_term p2)
 
38
  | PEsub(p1,p2) -> Sub(expr_to_term p1, expr_to_term p2)
 
39
  | PEpow(p,n)   -> Pow(expr_to_term p , C2Ml.n n)
 
40
  |  PEopp p ->  Opp (expr_to_term p)
 
41
 
 
42
      
 
43
 
 
44
 
 
45
(* let term_to_expr e =
 
46
  let e' = term_to_expr e in
 
47
   if debug 
 
48
   then Printf.printf "term_to_expr : %s - %s\n"  
 
49
    (string_of_poly (poly_of_term  e)) 
 
50
    (string_of_poly (poly_of_term (expr_to_term e')));
 
51
   e' *)
 
52
 
 
53
end 
 
54
open M      
 
55
 
 
56
open List
 
57
open Mutils 
 
58
 
 
59
 
 
60
 
 
61
 
 
62
let rec canonical_sum_to_string = function s -> failwith "not implemented"
 
63
 
 
64
let print_canonical_sum m = Format.print_string (canonical_sum_to_string m)
 
65
 
 
66
let print_list_term l = 
 
67
 print_string "print_list_term\n";
 
68
 List.iter (fun (Mc.Pair(e,k)) -> Printf.printf "q: %s %s ;" 
 
69
  (string_of_poly (poly_of_term (expr_to_term e))) 
 
70
  (match k with 
 
71
    Mc.Equal -> "= " 
 
72
   | Mc.Strict -> "> " 
 
73
   | Mc.NonStrict -> ">= " 
 
74
   | _ -> failwith "not_implemented")) l ;
 
75
 print_string "\n"
 
76
 
 
77
 
 
78
let partition_expr l = 
 
79
 let rec f i = function
 
80
  | [] -> ([],[],[])
 
81
  | Mc.Pair(e,k)::l -> 
 
82
     let (eq,ge,neq) = f (i+1) l in
 
83
      match k with 
 
84
       | Mc.Equal -> ((e,i)::eq,ge,neq)
 
85
       | Mc.NonStrict -> (eq,(e,Axiom_le i)::ge,neq)
 
86
       | Mc.Strict    -> (* e > 0 == e >= 0 /\ e <> 0 *) 
 
87
          (eq, (e,Axiom_lt i)::ge,(e,Axiom_lt i)::neq)
 
88
       | Mc.NonEqual -> (eq,ge,(e,Axiom_eq i)::neq) 
 
89
          (* Not quite sure -- Coq interface has changed *)
 
90
 in f 0 l
 
91
 
 
92
 
 
93
let rec sets_of_list l =
 
94
 match l with
 
95
  | [] -> [[]]
 
96
  | e::l -> let s = sets_of_list l in
 
97
             s@(List.map (fun s0 -> e::s0) s) 
 
98
 
 
99
(* The exploration is probably not complete - for simple cases, it works... *)
 
100
let real_nonlinear_prover d l =
 
101
 try 
 
102
  let (eq,ge,neq) = partition_expr l in
 
103
 
 
104
  let rec elim_const  = function
 
105
    [] ->  []
 
106
   | (x,y)::l -> let p = poly_of_term (expr_to_term x) in
 
107
                  if poly_isconst p 
 
108
                  then elim_const l 
 
109
                  else (p,y)::(elim_const l) in
 
110
 
 
111
  let eq = elim_const eq in
 
112
  let peq = List.map  fst eq in
 
113
   
 
114
  let pge = List.map 
 
115
   (fun (e,psatz) -> poly_of_term (expr_to_term e),psatz) ge in
 
116
   
 
117
  let monoids = List.map (fun m ->  (List.fold_right (fun (p,kd) y -> 
 
118
   let p = poly_of_term (expr_to_term p) in
 
119
    match kd with
 
120
     | Axiom_lt i -> poly_mul p y
 
121
     | Axiom_eq i -> poly_mul (poly_pow p 2) y
 
122
     |   _        -> failwith "monoids") m (poly_const (Int 1)) , map  snd m))
 
123
   (sets_of_list neq) in
 
124
 
 
125
  let (cert_ideal, cert_cone,monoid) = deepen_until d (fun d -> 
 
126
   list_try_find (fun m -> let (ci,cc) = 
 
127
    real_positivnullstellensatz_general false d peq pge (poly_neg (fst m) ) in
 
128
                      (ci,cc,snd m)) monoids) 0 in
 
129
   
 
130
  let proofs_ideal = map2 (fun q i -> Eqmul(term_of_poly q,Axiom_eq i))  
 
131
   cert_ideal (List.map snd eq) in
 
132
 
 
133
  let proofs_cone = map term_of_sos cert_cone in
 
134
   
 
135
  let proof_ne =  
 
136
   let (neq , lt) = List.partition 
 
137
    (function  Axiom_eq _ -> true | _ -> false ) monoid in
 
138
   let sq = match 
 
139
     (List.map (function Axiom_eq i -> i | _ -> failwith "error") neq) 
 
140
    with
 
141
    | []  -> Rational_lt (Int 1)
 
142
    | l   -> Monoid l in
 
143
    List.fold_right (fun x y -> Product(x,y)) lt sq in
 
144
 
 
145
  let proof = list_fold_right_elements 
 
146
   (fun s t -> Sum(s,t)) (proof_ne :: proofs_ideal @ proofs_cone) in
 
147
   Some proof
 
148
 with 
 
149
  | Sos.TooDeep -> None
 
150
 
 
151
 
 
152
(* This is somewhat buggy, over Z, strict inequality vanish... *)
 
153
let pure_sos  l =
 
154
 (* If there is no strict inequality, 
 
155
    I should nonetheless be able to try something - over Z  > is equivalent to -1  >= *)
 
156
 try 
 
157
  let l = List.combine l (interval 0 (length l -1)) in
 
158
  let (lt,i) =  try (List.find (fun (x,_) -> snd' x =  Mc.Strict) l) 
 
159
   with Not_found -> List.hd l in
 
160
  let plt = poly_neg (poly_of_term (expr_to_term (fst' lt))) in
 
161
  let (n,polys) = sumofsquares plt in (* n * (ci * pi^2) *)
 
162
  let pos = Product (Rational_lt n, 
 
163
                    List.fold_right (fun (c,p) rst -> Sum (Product (Rational_lt c, Square
 
164
                     (term_of_poly p)), rst)) 
 
165
                     polys (Rational_lt (Int 0))) in
 
166
  let proof = Sum(Axiom_lt i, pos) in
 
167
(*  let s,proof' = scale_certificate proof in
 
168
  let cert  = snd (cert_of_pos proof') in *)
 
169
   Some proof
 
170
 with
 
171
  | Not_found -> (* This is no strict inequality *) None
 
172
  |  x        -> None
 
173
 
 
174
 
 
175
type micromega_polys = (Micromega.q Mc.pExpr, Mc.op1) Micromega.prod list
 
176
type csdp_certificate = Sos.positivstellensatz option
 
177
type provername = string * int option
 
178
 
 
179
let main () =
 
180
  if Array.length Sys.argv <> 3 then
 
181
    (Printf.printf "Usage: csdpcert inputfile outputfile\n"; exit 1);
 
182
  let input_file = Sys.argv.(1) in
 
183
  let output_file = Sys.argv.(2) in
 
184
  let inch = open_in input_file in
 
185
  let (prover,poly) = (input_value inch : provername * micromega_polys) in
 
186
  close_in inch;
 
187
  let cert =
 
188
    match prover with
 
189
    | "real_nonlinear_prover", Some d -> real_nonlinear_prover d poly
 
190
    | "pure_sos", None -> pure_sos poly
 
191
    | prover, _ -> (Printf.printf "unknown prover: %s\n" prover; exit 1) in
 
192
  let outch = open_out output_file in
 
193
  output_value outch (cert:csdp_certificate);
 
194
  close_out outch;
 
195
  exit 0;;
 
196
 
 
197
let _ = main () in ()