formal_lp/old/formal_interval/interval_1d/recurse.hl

   1 (* port of recurse.cc *)
   2
   3 (*
   4 This is the code that verifies a disjunct of nonlinear inequalities.
   5 The are given as a list   (tf:tfunction list).  If tf = [f1;....;fk], then
   6 the list represents the inequality (f1 < 0 \/ f2 < 0 .... fk < 0).
   7
   8 The end user should only need to define a cell option,
   9 and then call recursive_verifier, which recursively bisects the domain
  10 until a partition of the domain is found on which verifier_cell gives
  11 a pass on each piece of the partition.
  12
  13 *)
  14
  15 flyspeck_needs "../formal_lp/formal_interval/interval_1d/types.hl";;
  16 flyspeck_needs "../formal_lp/formal_interval/interval_1d/report.hl";;
  17 flyspeck_needs "../formal_lp/formal_interval/interval_1d/interval.hl";;
  18 flyspeck_needs "../formal_lp/formal_interval/interval_1d/univariate.hl";;
  19 flyspeck_needs "../formal_lp/formal_interval/interval_1d/line_interval.hl";;
  20 flyspeck_needs "../formal_lp/formal_interval/interval_1d/taylor.hl";;
  21
  22 module Recurse = struct
  23
  24 open Interval;;
  25 open Univariate;;
  26 open Line_interval;;
  27 open Taylor;;
  28
  29 type cellOption = {
  30   allow_sharp : bool;
  31   mutable iteration_count : int;
  32   iteration_limit : int;
  33   recursion_depth : int;
  34 };;
  35
  36 (* cell verification is complex, and we use exceptions to
  37     exit as soon as the status has been determined.   *)
  38
  39 type cell_status =
  40   | Cell_pass
  41   | Cell_counterexample
  42   | Cell_inconclusive of (float * float * tfunction);;
  43
  44 type result_tree =
  45         | Result_false
  46         | Result_pass of (float * float)
  47         | Result_glue of result_tree * result_tree;;
  48
  49
  50 exception Return of  cell_status;;
  51
  52
  53 let return c = raise (Return c);;
  54
  55
  56 (* error checking and reporting functions *)
  57
  58 let boolify _ = true;;
  59
  60 let string_of_domain x =
  61   Printf.sprintf "{%f}" x;;
  62
  63 let string3 (x,z,s) =  (string_of_domain x ^"\n"^ string_of_domain z ^ "\n" ^ s);;
  64
  65 let report_current = boolify o Report.report_timed o string3;;
  66
  67 let report_error = boolify o Report.report_error o string3;;
  68
  69 let report_fatal = boolify o Report.report_fatal o string3;;
  70
  71
  72 (* let t = [0.1;0.2;0.3;0.4;0.5;0.6] in report_error (t,t,"ok");; *)
  73
  74 let periodic_count =
  75   let end_count = ref 0 in
  76     fun () ->
  77       let _ = end_count := !end_count + 1 in
  78         (0 = ( !end_count mod 80000));;
  79
  80 let check_limit opt depth =
  81   let _ = opt.iteration_count <- opt.iteration_count + 1 in
  82    ( opt.iteration_count < opt.iteration_limit or opt.iteration_limit = 0 ) &&
  83      (depth < opt.recursion_depth);;
  84
  85 let sgn x = if (x.lo > 0.0) then 1 else if (x.hi < 0.0) then -1 else 0;;
  86
  87
  88 (* a bit tricky because unstable exceptions are common early on,
  89     and evaluations are very expensive.
  90    We don't want a single unstable disjunct to ruin everything.
  91
  92    When boxes are small, then we throw away unstable disjuncts and continue w/o them.
  93    When the boxes are still big, we return inconclusive to force a bisection.
  94
  95    Starting with this procedure, we can throw an exception to return early,
  96    as soon as we are able to determine the cell_status.  All these exceptions
  97    get caught at the last line of verify_cell.
  98 *)
  99
 100 let rec set_targets (x,z,tf) =
 101     try (
 102                 let target = evalf tf x z in
 103                 let _ =  upper_bound target >= 0.0 or return Cell_pass in
 104                         [target, tf]
 105     )
 106     with Unstable -> (
 107 (*              let _ = (2.0 *. maxwidth <= opt.width_cutoff) or return (Cell_inconclusive (x, z, x0, z0, tf)) in []*)
 108       return (Cell_inconclusive (x, z, tf))
 109     );
 110 ;;
 111
 112
 113 let rec delete_false acc tis =
 114   if (tis=[]) then List.rev acc
 115   else if (lower_bound (fst (hd tis)) > 0.0) then delete_false acc (tl tis)
 116   else delete_false (hd tis::acc) (tl tis);;
 117
 118
 119
 120 (* loop as long as monotonicity keeps making progress.  *)
 121
 122 let rec going_strong(x,z,tf,opt) =
 123   let (y,w) = center_form (x,z) in
 124   let maxwidth = w in
 125   let tis = set_targets(x,z,tf) in
 126   let epsilon_width = 1.0e-8 in
 127   let _ = (maxwidth >= epsilon_width) or return
 128         (if (opt.allow_sharp) then (Report.inc_corner_count(); Cell_pass)
 129                                    else Cell_counterexample) in
 130   let tis = delete_false [] tis in
 131   let _ = (List.length tis > 0) or return Cell_counterexample in
 132         (x, z, maxwidth, tis);;
 133
 134 (*
 135 This procedure is mostly guided by heuristics that don't require formal
 136 verification. In particular, no justification is required for tossing out inequalities
 137 (since they appear as disjuncts, we can choose which one to prove).
 138
 139 Formal verification is required whenever a Cell_passes is issued,
 140 and whenever the domain (x,z) is restricted.
 141
 142 The record (x0,z0) of the current outer boundary must be restricted to (x,z)
 143 whenever an inequality is tossed out.
 144 *)
 145
 146 let rec verify_cell (x,z,tf,opt) =
 147   try (
 148   let _ = not(periodic_count ()) or report_current (x,z,"periodic report") in
 149   let (x,z,maxwidth,tis) = going_strong(x,z,tf,opt) in
 150     Cell_inconclusive (x,z,hd (map snd tis));
 151   )
 152   with Return c -> c;;
 153
 154 let rec recursive_verifier (depth,x,z,tf,opt) =
 155   let _ = check_limit opt depth or report_fatal(x,z,Printf.sprintf "check_limit: depth %d" depth) in
 156     match verify_cell(x,z,tf,opt)  with
 157       | Cell_counterexample -> Result_false
 158       | Cell_pass -> Result_pass(x,z)
 159       | Cell_inconclusive(x,z,tf) ->
 160           (let ( ++ ), ( / ) = up(); upadd, updiv in
 161            let yj = (x ++ z) / 2.0 in
 162            let r1 = recursive_verifier(depth + 1, x, yj, tf, opt) in
 163              if r1 = Result_false then Result_false else
 164                let r2 = recursive_verifier(depth + 1, yj, z, tf, opt) in
 165                  if r2 = Result_false then Result_false else
 166                    Result_glue (r1, r2)
 167           );;
 168
 169
 170  end;;