LP_arith/lp_arith.ml

   1
   2 (* small LP-based prover, to convert the HOL-terms to a coefficient
   3    matrix and back it uses the code of REAL_LINEAR_PROVER in the HOL
   4    Light distribution *)
   5
   6
   7 let cddwrapper = "cdd_cert";;
   8
   9 (* in lin_of_hol one can replace the call to linear_add to a call to lin_add *)
  10
  11 let lin_of_hol =
  12   let one_tm = `&1:real`
  13   and zero_tm = `&0:real`
  14   and add_tm = `(+):real->real->real`
  15   and mul_tm = `(*):real->real->real`
  16   and lin_add = combine (+/) (fun x -> x =/ num_0) in
  17   let rec lin_of_hol tm =
  18     if tm = zero_tm then undefined
  19     else if not (is_comb tm) then (tm |=> Int 1)
  20     else if is_ratconst tm then (one_tm |=> rat_of_term tm) else
  21       let lop,r = dest_comb tm in
  22         if not (is_comb lop) then (tm |=> Int 1) else
  23           let op,l = dest_comb lop in
  24             if op = add_tm then lin_add (lin_of_hol l) (lin_of_hol r)
  25             else if op = mul_tm & is_ratconst l then (r |=> rat_of_term l)
  26             else (tm |=> Int 1) in
  27     lin_of_hol;;
  28
  29 let words s =
  30   let stre = Stream.of_string s in
  31   let is_empty st = match Stream.peek st with
  32       None -> true
  33     | Some _ -> false in
  34   let rec sb acc st =
  35     if is_empty st
  36     then [acc]
  37     else
  38       let t = Stream.next st in
  39         if t = ' ' then acc :: (sb "" st) else sb (acc ^ Char.escaped t) st
  40   in filter (fun x -> x <> "") (sb "" stre);;
  41
  42 let cdd ins =
  43   let outfn = Filename.temp_file "cdd" ".res"
  44   and infn = Filename.temp_file "cdd" ".ine" in
  45   let s = "cat " ^ infn ^ "| " ^ cddwrapper ^ " 2> /dev/null > " ^ outfn in
  46   let inch = open_out infn in
  47     output_string inch ins;
  48     close_out inch;
  49     if Sys.command s <> 0 then failwith "cdd" else
  50       let fd = Pervasives.open_in outfn in let data = input_line fd in close_in fd; Sys.remove infn; Sys.remove outfn; data;;
  51
  52 let rec take n l =
  53   match l with
  54       x :: xs -> if n = 0 then [] else x :: (take (n-1) xs)
  55     | [] -> [];;
  56
  57 let rec drop n l =
  58   match l with
  59       x :: xs -> if n = 0 then l else (drop (n-1) xs)
  60     | [] -> [];;
  61
  62 let lp_prover (eq,le,lt) =
  63   let one_tm = `&1:real` in
  64   let vars = (subtract (itlist (union o dom) (eq@le@lt) []) [one_tm]) in
  65   let neq = length eq
  66   and nle = length le
  67   and nlt = length lt
  68   and nr = length (eq@le@lt) in
  69   let get_row v = map (fun x -> applyd x (fun _ -> num_0) v) (eq@le@lt) in
  70   let rec rep n e = if n = 0 then [] else e :: (rep (n-1) e) in
  71   let one_at n =
  72     map (fun i -> (rep i (num_0))@[num_1]@(rep (n-i-1) (num_0))) (0--(n-1)) in
  73   let main_rows = map ((fun l -> num_0::l) o get_row) vars
  74   and lt_row = [minus_num num_1] @ (rep (length eq) num_0) @
  75     (rep (length le) num_0) @ (rep (length lt) num_1)
  76   and pos_rows = map (fun l -> (rep (length eq + 1 ) num_0) @ l)
  77     (one_at (length (le@lt)))
  78   and bvec = (num_0 :: (get_row one_tm)) in
  79   let mat = main_rows@[lt_row]@pos_rows in
  80   let string_of_row = (String.concat " ") o (map string_of_num) in
  81   let cddlp = (String.concat "\n"
  82                  ["H-representation";
  83                   "linearity "^(string_of_int (length main_rows))^" "^
  84                     (String.concat " " (map string_of_int (1--(length main_rows))));
  85                   "begin";
  86                   String.concat " " [string_of_int (length mat);string_of_int (nr+1);"rational"];
  87                   String.concat "\n" (map string_of_row mat);
  88                   "end";
  89                   String.concat " " ["minimize";string_of_row bvec]]) in
  90   let outp = (cdd cddlp) in
  91   let res = (* print_string cddlp; print_newline();   *)(* print_string outp; print_newline(); *)
  92     if outp = "No Contradiction" then failwith "No contradiction" else map Num.num_of_string (words outp) in
  93   let (req,rle,rlt) = (take neq res,
  94                        take nle (drop neq res),
  95                        take nlt (drop (nle+neq) res)) in
  96   let peq = map2 (fun r e -> if (r =/ num_0) then [] else [Eqmul (term_of_rat r, Axiom_eq e)]) req (0--(neq-1))
  97   and ple = map2 (fun r e -> if (r =/ num_0) then [] else [Product (Rational_lt r,Axiom_le e)]) rle (0--(nle-1))
  98   and plt = map2 (fun r e -> if (r =/ num_0) then [] else [Product (Rational_lt r,Axiom_lt e)]) rlt (0--(nlt-1)) in
  99   let pp = List.flatten (peq@ple@plt) in
 100   let refu = itlist (fun acc x -> Sum (acc,x)) (tl pp) (hd pp) in
 101     (*     print_string outp; *)
 102     (*     print_newline(); *)
 103     refu;;
 104
 105 let LP_PROVER =
 106   let is_alien tm =
 107     match tm with
 108         Comb(Const("real_of_num",_),n) when not(is_numeral n) -> true
 109       | _ -> false in
 110   let n_tm = `n:num` in
 111   let pth = REWRITE_RULE[GSYM real_ge] (SPEC n_tm REAL_POS) in
 112     fun translator (eq,le,lt) ->
 113       let eq_pols = map (lin_of_hol o lhand o concl) eq
 114       and le_pols = map (lin_of_hol o lhand o concl) le
 115       and lt_pols = map (lin_of_hol o lhand o concl) lt in
 116       let aliens =  filter is_alien
 117         (itlist (union o dom) (eq_pols @ le_pols @ lt_pols) []) in
 118       let le_pols' = le_pols @ map (fun v -> (v |=> Int 1)) aliens in
 119       let proof = lp_prover(eq_pols,le_pols',lt_pols) in
 120       let le' = le @ map (fun a -> INST [rand a,n_tm] pth) aliens in
 121         translator (eq,le',lt) proof;;
 122
 123 let LP_ARITH =
 124   let init = GEN_REWRITE_CONV ONCE_DEPTH_CONV [DECIMAL]
 125   and pure = GEN_REAL_ARITH LP_PROVER in
 126     fun tm -> let th = init tm in EQ_MP (SYM th) (pure(rand(concl th)));;
 127
 128 let LP_ARITH_TAC = CONV_TAC LP_ARITH;;