Extend definition of sturm sequences to prove a first soundness theorem about num of zeroes

checkerScript.sml

@@ -3,8 +3,8 @@
   Dandelion
 **)
 open realTheory realLib RealArith stringTheory;
-open renameTheory realPolyTheory transcLangTheory sturmComputeTheory
-     checkerDefsTheory pointCheckerTheory;
+open renameTheory realPolyTheory transcLangTheory sturmComputeTheory sturmTheory
+     drangTheory checkerDefsTheory pointCheckerTheory;
 open preamble;
 
 val _ = new_theory "checker";
@@ -30,18 +30,20 @@ End
 **)
 (* TODO: Is this sufficient for a check? *)
 Definition checkZeroes_def:
-  checkZeroes numZeroes zeroes =
-  if LENGTH zeroes ≠ numZeroes
-  then Invalid "Incorrect number of zeroes found"
-  else Valid
+  checkZeroes deriv1 deriv2 iv sseq zeroes =
+  let numZeroes =
+      (variation (MAP (λp. poly p (FST iv)) (deriv1::deriv2::sseq)) -
+       variation (MAP (λp. poly p (SND iv)) (deriv1::deriv2::sseq)))
+  in
+    if (LENGTH zeroes ≠ numZeroes ∨ poly deriv1 (FST iv) = 0 ∨ poly deriv1 (SND iv) = 0)
+    then Invalid "Incorrect number of zeroes found"
+    else Valid
 End
 
 (**
   Checks that the value at the extrema of the error polynomial are less than or
   equal to the actual error encoded in the certificate **)
 (* TODO: Is this the correct way of checking the zeroes?? Check Harrison paper *)
-
-
 (** In the second part of the definition,
     the first point line verifies that the zero interval (z1, z2) is inside
     the large interval vi
@@ -58,17 +60,17 @@ End
     Tolerable bound for size of the inteval (z1, z2) i.e. e in the paper
  **)
 Definition validateZeroesLeqErr_def:
-  validateZeroesLeqErr diffPoly deriv_poly iv [] eps =
+  validateZeroesLeqErr diffPoly iv [] eps =
     (if abs(poly diffPoly (FST iv)) ≤ eps ∧ abs(poly diffPoly (SND iv)) ≤ eps
     then Valid
     else Invalid "Outer points")
   ∧
-  validateZeroesLeqErr diffPoly deriv_poly iv ((z1,z2)::zeroes) eps =
+  validateZeroesLeqErr diffPoly iv ((z1,z2)::zeroes) eps =
   ( if ((FST iv) < z1 ∧ z2 < (SND iv) ∧ z1 ≤ z2 ∧
-        ((poly deriv_poly z1) * (poly deriv_poly z2) < &0) ∧
+        ((poly diffPoly z1) * (poly diffPoly z2) < &0) ∧
         abs(poly diffPoly z1) ≤ eps ∧ abs(poly diffPoly z2) ≤ eps)
-        then validateZeroesLeqErr diffPoly deriv_poly iv zeroes eps
-     else Invalid "Zero of derivate not an extremal point or extrema too big")
+        then validateZeroesLeqErr diffPoly iv zeroes eps
+     else Invalid "Zero of derivativ not an extremal point or extrema too big")
 End
 
 (**
@@ -81,20 +83,78 @@ Definition checker_def:
     | Valid =>
       let (transp, err) = approx cert.transc cert.iv;
           errorp = transp -p cert.poly;
-          deriv = derive errorp
+          deriv1 = diff errorp;
+          deriv2 = diff deriv1
       in
-        case sturm_seq errorp deriv of
+        case sturm_seq deriv1 deriv2 of
         NONE => Invalid "Could not compute sturm sequence"
         | SOME sseq =>
-          let numZeroes =
-              (variation (MAP (λp. poly p (FST cert.iv)) sseq) -
-               variation (MAP (λp. poly p (SND cert.iv)) sseq));
-          in
-            case checkZeroes numZeroes cert.zeroes of
+            case checkZeroes deriv1 deriv2 cert.iv sseq cert.zeroes of
             | Valid =>
-                validateZeroesLeqErr errorp deriv
-                                     cert.iv cert.zeroes (cert.eps - err)
+                validateZeroesLeqErr deriv1 cert.iv cert.zeroes (cert.eps - err)
             | Invalid s => Invalid s
 End
 
+Definition getMaxWidth_def:
+  getMaxWidth [] = 0 ∧
+  getMaxWidth ((u,v)::xs) = max (abs(u-v)) (getMaxWidth xs)
+End
+
+Theorem getMaxWidth_is_max:
+  EVERY (λ (u,v). abs(u-v) ≤ getMaxWidth l)
+Proof
+  cheat
+QED
+
+Theorem checkZeroes_sound:
+  ∀ sseq deriv1 deriv2 iv zeroes.
+    sturm_seq deriv1 (diff deriv1) = SOME sseq ∧
+    checkZeroes deriv1 (diff deriv1) iv sseq zeroes = Valid ∧
+    FST iv ≤ SND iv ⇒
+    {x | FST iv ≤ x ∧ x ≤ SND iv ∧ (poly deriv1 x = &0)} HAS_SIZE
+                                                         LENGTH zeroes
+Proof
+  rpt gen_tac
+  >> rewrite_tac[checkZeroes_def]
+  >> CONV_TAC $ DEPTH_CONV let_CONV
+  >> cond_cases_tac >- gs[]
+  >> rpt $ pop_assum $ mp_tac o SIMP_RULE std_ss []
+  >> rpt strip_tac
+  >> qpat_x_assum ‘LENGTH zeroes = _’ $ rewrite_tac o single
+  >> imp_res_tac sturm_seq_equiv
+  >> irule STURM_THEOREM >> gs[]
+QED
+
+Theorem validateZeroesLeqErr_EVERY:
+  validateZeroesLeqErr deriv iv zeroes err = Valid ⇒
+  EVERY (λ (u,v). FST iv ≤ u ∧ v ≤ SND iv ∧ abs (poly errorp u) ≤ ub) l
+Proof
+  cheat
+QED
+
+Theorem validateZeroesLeqErr_sound:
+  ∀ deriv errorp iv zeroes err.
+    deriv = diff errorp ∧
+    {x | FST iv ≤ x ∧ x ≤ SND iv ∧ (poly deriv x = &0)} HAS_SIZE LENGTH zeroes ∧
+    validateZeroesLeqErr deriv iv zeroes err = Valid ⇒
+    ∀ x. FST iv ≤ x ∧ x ≤ SND iv ⇒ poly errorp x ≤ err
+Proof
+  rpt strip_tac
+  >> ‘∀ x. FST iv ≤ x ∧ x ≤ SND iv ⇒ ((λ x. poly errorp x) diffl poly deriv x) x’
+     by (rpt strip_tac >> gs[polyTheory.POLY_DIFF])
+  >> drule (GEN_ALL BOUND_THEOREM_INEXACT |> SPEC_ALL |> GEN “e:real”)
+  >> gs[validateZeroesLeqErr_def]
+  >> cheat
+QED
+
+Theorem checker_soundness:
+  ∀ cert.
+    checker cert = Valid ⇒
+    ∀ x. FST(cert.iv) ≤ x ∧ x ≤ SND (cert.iv) ⇒
+         (interp cert.transc x - cert.poly) ≤ cert.err
+Proof
+  cheat
+QED
+
+
 val _ = export_theory();

sturmComputeScript.sml

@@ -25,8 +25,15 @@ End
 
 Definition sturm_seq_def:
   sturm_seq (p:poly) (q:poly) : poly list option =
-  if (deg p = 1 ∨ deg p = 0) then NONE
-  else sturm_seq_aux (deg p - 1) p q
+  if zerop q ∨ deg p ≤ 1 then NONE
+  else
+    case sturm_seq_aux (deg p - 1) p q of
+    | NONE => NONE
+    | SOME sseq =>
+        case LAST sseq of
+        | [ x ] => if x ≠ 0 then SOME sseq
+                   else NONE
+        | _ => NONE
 End
 
 Theorem poly_neg_evals:
@@ -316,4 +323,55 @@ Proof
   >> pop_assum kall_tac >> last_x_assum kall_tac >> real_tac
 QED
 
+Theorem reduce_nonzero:
+  reduce p ≠ [] ⇒
+  ~ (EVERY (λ c. c = 0) p)
+Proof
+  Induct_on ‘p’ >> gs[reduce_def]
+  >> cond_cases_tac >> gs[]
+QED
+
+Theorem reduce_not_zero:
+  reduce p ≠ [] ⇒
+  ∃ x. evalPoly (reduce p) x ≠ 0
+Proof
+  gs[reduce_preserving] >> gs[poly_compat] >> rpt strip_tac
+  >> CCONTR_TAC >> gs[]
+  >> ‘poly p = poly []’ by (gs[FUN_EQ_THM, poly_def])
+  >> gs[POLY_ZERO]
+  >> imp_res_tac reduce_nonzero
+QED
+
+Theorem sturm_seq_aux_length:
+  ∀ n p q sseq.
+    sturm_seq_aux n p q = SOME sseq ⇒
+    LENGTH sseq = n
+Proof
+  Induct_on ‘n’ >> gs[sturm_seq_aux_def]
+  >> rpt strip_tac >> gs[CaseEq"option"]
+  >> rpt VAR_EQ_TAC >> res_tac >> gs[]
+QED
+
+Theorem sturm_seq_equiv:
+  sturm_seq p q = SOME sseq ⇒
+  poly q ≠ poly [] ∧ sseq ≠ [] ∧
+  (∃ d. d ≠ 0 ∧ LAST sseq = [d]) ∧
+  STURM p q sseq
+Proof
+  gs[sturm_seq_def, CaseEq"option", CaseEq"list"]
+  >> rpt $ disch_then strip_assume_tac
+  >> imp_res_tac sturm_equiv >> gs[]
+  >> gs[zerop_def]
+  >> rpt conj_tac
+  >- (
+    gs[FUN_EQ_THM, GSYM poly_compat, Once $ GSYM reduce_preserving,
+       evalPoly_def]
+    >> imp_res_tac reduce_not_zero
+    >> fsrw_tac [SATISFY_ss] [reduce_idempotent])
+  >> CCONTR_TAC >> ‘LENGTH sseq = 0’ by gs[]
+  >> ‘0 < deg p - 1’ by gs[]
+  >> ‘0 < LENGTH sseq’ suffices_by gs[]
+  >> imp_res_tac sturm_seq_aux_length >> gs[]
+QED
+
 val _ = export_theory();