1 (* Recursive verification of inequalities using the basic interval arithmetic only *)
3 flyspeck_needs "../formal_lp/formal_interval/interval_m/recurse.hl";;
5 module Recurse0 = struct
17 let rec has_monotone0 opt tf domain0 x z x0 z0 is found = match is with
18 | [] -> (x,z,x0,z0,List.rev found)
19 | j::js when (mth x j >= mth z j) ->
20 has_monotone0 opt tf domain0 x z x0 z0 js found
22 let df_int = try evalf0 tf (j + 1) (fst domain0) (snd domain0) with Unstable -> mk_interval (-1.0,1.0) in
23 let allpos_df0 = df_int.lo >= opt.eps in
24 let allneg_df0 = df_int.hi < ~-.(opt.eps) in
27 {variable = j + 1; decr_flag = false; df0_flag = allpos_df0; ti_flag = false} in
28 if opt.mono_pass && mth z j < mth z0 j then return (Cell_pass_mono ([], status))
30 let setj u = table (fun i -> (if i=j then mth z j else mth u i)) in
31 has_monotone0 opt tf domain0 (setj x) (setj z)
32 (setj x0) (setj z0) js (status :: found)
33 else if allneg_df0 then
35 {variable = j + 1; decr_flag = true; df0_flag = allneg_df0; ti_flag = false} in
36 if opt.mono_pass && mth x j > mth x0 j then return (Cell_pass_mono ([], status))
38 let setj u = table (fun i -> (if i=j then mth x j else mth u i)) in
39 has_monotone0 opt tf domain0 (setj x) (setj z)
40 (setj x0) (setj z0) js (status :: found)
41 else has_monotone0 opt tf domain0 x z x0 z0 js found;;
43 (* loop as long as monotonicity keeps making progress. *)
45 let rec going_strong0(x,z,x0,z0,tf,opt,mono) =
46 let (y,w) = center_form (x,z) in
47 let maxwidth = maxl w in
48 let target0 = try evalf0 tf 0 x z with Unstable -> return (Cell_inconclusive (mono,x,z,x0,z0)) in
49 let _ = target0.hi >= ~-.(opt.eps) or return (Cell_pass (mono, true)) in
50 let epsilon_width = 1.0e-8 in
51 let _ = (maxwidth >= epsilon_width) or return Cell_counterexample in
52 let (x,z,x0,z0,strong) =
53 if (opt.allow_derivatives) then
55 has_monotone0 opt tf (x,z) x z x0 z0 iter8 []
56 with Return (Cell_pass_mono (_, status)) -> return (Cell_pass_mono (mono, status))
57 else (x,z,x0,z0,[]) in
58 if (strong <> []) then
59 going_strong0(x,z,x0,z0,tf,opt,mono @ [strong])
61 (x,z,x0,z0,maxwidth,mono);;
65 This procedure is mostly guided by heuristics that don't require formal
66 verification. In particular, no justification is required for tossing out inequalities
67 (since they appear as disjuncts, we can choose which one to prove).
69 Formal verification is required whenever a Cell_passes is issued,
70 and whenever the domain (x,z) is restricted.
72 The record (x0,z0) of the current outer boundary must be restricted to (x,z)
73 whenever an inequality is tossed out.
76 let rec verify_cell0 (x,z,x0,z0,tf,opt) =
78 let _ = not(periodic_count ()) or report_current (x,z,"periodic report") in
79 let (x,z,x0,z0,maxwidth,mono) = going_strong0(x,z,x0,z0,tf,opt,[]) in
80 if opt.convex_flag then
81 let ti = try evalf tf x z with Unstable -> return (Cell_inconclusive (mono,x,z,x0,z0)) in
82 Cell_inconclusive_ti (mono,ti,x,z,x0,z0)
84 Cell_inconclusive (mono,x,z,x0,z0)
88 let rec recursive_verifier0 (depth,x,z,x0,z0,tf,opt) =
89 let _ = check_limit opt depth or report_fatal(x,z,Printf.sprintf "depth %d" depth) in
90 let split_and_verify j x z x0 z0 convex_flag =
91 let ( ++ ), ( / ) = up(); upadd, updiv in
92 let yj = (mth x j ++ mth z j) / 2.0 in
93 let delta b v = table (fun i-> if (i = j && b) then yj else mth v i) in
96 x, table (fun i -> if i = j then mth x i else mth z i)
98 delta false x, delta true z in
101 table (fun i -> if i = j then mth z i else mth x i), z
103 delta true x, delta false z in
104 let r1 = recursive_verifier0(depth+1,x1,z1,x0,z0,tf,opt) in
106 | Result_false t -> Result_false t
108 (let r2 = recursive_verifier0(depth+1,x2,z2,x0,z0,tf,opt) in
110 | Result_false t -> Result_false t
111 | _ -> Result_glue (j, convex_flag, r1, r2)) in
113 let add_mono mono r1 =
114 itlist (fun m r -> Result_mono (m, r)) mono r1 in
117 match verify_cell0(x,z,x0,z0,tf,opt) with
118 | Cell_counterexample -> Result_false (x,z)
119 | Cell_pass (mono, f0_flag) -> add_mono mono (Result_pass (f0_flag,x,z))
120 | Cell_pass_mono (mono, status) -> add_mono mono (Result_pass_mono status)
121 | Cell_inconclusive_ti(mono,ti,x,z,x0,z0) ->
122 let dds = map (fun i -> mth (mth ti.dd i) i, i) iter8 in
123 let convex_dds = filter (fun dd, i -> dd.lo >= opt.eps && mth x i < mth z i) dds in
124 let convex_i = map snd convex_dds in
125 let w2 = List.map2 upsub z x in
126 let convex_flag, ws, ws_i =
127 if convex_dds = [] then
130 true, map (mth w2) convex_i, convex_i in
131 let maxwidth2 = maxl ws in
132 let j_wide = try( find (fun i -> mth w2 i = maxwidth2) ws_i) with
133 | _ -> failwith "recursive_verifier find" in
134 add_mono mono (split_and_verify j_wide x z x0 z0 convex_flag)
136 | Cell_inconclusive(mono,x,z,x0,z0) ->
137 let w2 = List.map2 upsub z x in
138 let maxwidth2 = maxl w2 in
139 let j_wide = try( find (fun i -> mth w2 i = maxwidth2) iter8) with
140 | _ -> failwith "recursive_verifier find" in
141 add_mono mono (split_and_verify j_wide x z x0 z0 false);;