Prove soundness of approxPoly, add to realZeroLib

checkerScript.sml

@@ -349,6 +349,150 @@ Proof
   Induct_on ‘hints’ >> gs[getExpHint_def, CaseEq"hint"]
 QED
 
+Definition approxPolySideCond_def:
+  approxPolySideCond approxSteps =
+    (0 < approxSteps ∧ EVEN approxSteps ∧ EVEN (approxSteps DIV 2))
+End
+
+Theorem approxPoly_soundness:
+  ∀ transc iv hs approxSteps p err.
+    approxPolySideCond approxSteps ∧
+    approxPoly transc iv hs approxSteps = SOME (p, err) ⇒
+    ∀ x.
+      FST iv ≤ x ∧ x ≤ SND iv ⇒
+      abs (interp transc x - evalPoly p x) ≤ err
+Proof
+  simp[approxPoly_def, approxPolySideCond_def, CaseEq"transc"]
+  >> rpt strip_tac
+  >> ‘(tr = "exp" ∧
+       ((iv = (0, 1 * inv 2) ∧ getExpHint hs = NONE) ∨
+        ∃ n. getExpHint hs = SOME n ∧ iv = (0,&n * inv 2))) ∨
+      tr = "cos" ∨
+      tr = "sin"’
+    by (every_case_tac >> gs[getExpHint_SOME_MEM])
+  (* exp function, 0 to 1/2 *)
+  >- (
+    gs[interp_def, getFun_def]
+    >> qspecl_then [‘x’, ‘approxSteps’] strip_assume_tac MCLAURIN_EXP_LE
+    >> pop_assum $ rewrite_tac o single
+    >> rpt VAR_EQ_TAC
+    >> rewrite_tac[exp_sum_to_poly]
+    >> qmatch_goalsub_abbrev_tac ‘abs (exp_taylor + taylor_rem - exp_taylor) ≤ _’
+    >> ‘exp_taylor + taylor_rem - exp_taylor = taylor_rem’ by real_tac
+    >> pop_assum $ rewrite_tac o single
+    >> unabbrev_all_tac
+    >> ‘exp_err_small approxSteps = inv (&FACT approxSteps * 2 pow (approxSteps - 1))’ by EVAL_TAC
+    >> qspecl_then [‘approxSteps’, ‘x’,‘t’] mp_tac exp_remainder_bounded_small
+    >> impl_tac >> gs[]
+    >> real_tac)
+  (* exp function, 0 to 1 *)
+  >- (
+    gs[interp_def, getFun_def]
+    >> ‘1 ≠ inv 2’
+      by (once_rewrite_tac [GSYM REAL_INV1]
+          >> CCONTR_TAC
+          >> pop_assum $ mp_tac o SIMP_RULE std_ss []
+          >> rewrite_tac[REAL_INV_INJ] >> real_tac)
+    >> rpt VAR_EQ_TAC
+    >> rewrite_tac[GSYM poly_compat, eval_simps]
+    >> pop_assum $ rewrite_tac o single
+    >> rewrite_tac[exp_sum_to_poly]
+    >> qspecl_then [‘x’, ‘approxSteps’] strip_assume_tac MCLAURIN_EXP_LE
+    >> pop_assum $ rewrite_tac o single
+    >> qmatch_goalsub_abbrev_tac ‘abs (exp_taylor + taylor_rem - exp_taylor) ≤ _’
+    >> ‘exp_taylor + taylor_rem - exp_taylor = taylor_rem’ by real_tac
+    >> pop_assum $ rewrite_tac o single
+    >> unabbrev_all_tac
+    >> ‘exp_err_big n approxSteps = 2 pow n * &n pow approxSteps * inv (&FACT approxSteps * 2 pow approxSteps)’
+      by (rewrite_tac[] >> EVAL_TAC)
+    >> pop_assum $ rewrite_tac o single
+    >> qspecl_then [‘approxSteps’, ‘n’, ‘x’,‘t’] mp_tac exp_remainder_bounded_big
+    >> impl_tac
+    >- (rpt conj_tac >> gs[] >> real_tac)
+    >> rewrite_tac[])
+  (* cos function *)
+  >- (
+    gs[interp_def, getFun_def] >> rpt VAR_EQ_TAC
+    >> qspecl_then [‘x’, ‘approxSteps’] strip_assume_tac MCLAURIN_COS_LE
+    >> gs[]
+    >> pop_assum $ rewrite_tac o single
+    >> gs[cos_sum_to_poly]
+    >> qmatch_goalsub_abbrev_tac ‘abs (cos_taylor + taylor_rem - cos_taylor) ≤ _’
+    >> ‘cos_taylor + taylor_rem - cos_taylor = taylor_rem’ by real_tac
+    >> pop_assum $ rewrite_tac o single
+    >> unabbrev_all_tac
+    >> ‘(x pow approxSteps) * cos t * inv (&FACT approxSteps) =
+        (cos t * ((x pow approxSteps) * inv (&FACT approxSteps)))’
+      by real_tac
+    >> ‘-(x pow approxSteps) * cos t * inv (&FACT approxSteps) =
+        -(cos t * ((x pow approxSteps) * inv (&FACT approxSteps)))’
+      by real_tac
+    >> rewrite_tac []
+    >> ntac 2 $ pop_assum $ rewrite_tac o single
+    >> rewrite_tac [GSYM REAL_MUL_ASSOC]
+    >> qmatch_goalsub_abbrev_tac ‘abs (cos _ * err_cos_concr)’
+    >> irule REAL_LE_TRANS
+    >> qexists_tac ‘ 1 * abs err_cos_concr’ >> conj_tac
+    >- (rewrite_tac[ABS_MUL] >> irule REAL_LE_RMUL_IMP >> unabbrev_all_tac >> gs[COS_BOUND, ABS_POS])
+    >> rewrite_tac[REAL_MUL_LID]
+    >> ‘abs err_cos_concr = err_cos_concr’
+      by (unabbrev_all_tac
+          >> rewrite_tac[ABS_REFL]
+          >> irule REAL_LE_MUL >> conj_tac
+          >- (irule REAL_LE_INV >> gs[REAL_POS])
+          >> irule REAL_LE_MUL >> conj_tac
+          >> gs[REAL_POW_GE0])
+    >> pop_assum $ rewrite_tac o single
+    >> unabbrev_all_tac
+    >> rewrite_tac [cos_err_def]
+    >> imp_res_tac EVEN_ODD_EXISTS >> gs[POW_MINUS1]
+    >> irule REAL_LE_LMUL_IMP >> gs[GSYM POW_ABS]
+    >> irule REAL_LE_TRANS
+    >> qexists_tac ‘abs (x pow (2 * m'))’ >> gs[ABS_LE, GSYM POW_ABS]
+    >> irule POW_LE >> gs[ABS_POS]
+    >> irule RealSimpsTheory.maxAbs >> gs[])
+  (* sin *)
+  >> gs[interp_def, getFun_def] >> rpt VAR_EQ_TAC
+  >> qspecl_then [‘x’, ‘approxSteps’] strip_assume_tac MCLAURIN_SIN_LE
+  >> gs[]
+  >> pop_assum $ rewrite_tac o single
+  >> gs[sin_sum_to_poly]
+  >> qmatch_goalsub_abbrev_tac ‘abs (sin_taylor + taylor_rem - sin_taylor) ≤ _’
+  >> ‘sin_taylor + taylor_rem - sin_taylor = taylor_rem’ by real_tac
+  >> pop_assum $ rewrite_tac o single
+  >> unabbrev_all_tac
+  >> ‘inv (&FACT approxSteps) * sin t  * x pow approxSteps * -1 pow (approxSteps DIV 2) =
+      (sin t * ((x pow approxSteps) * inv (&FACT approxSteps) * -1 pow (approxSteps DIV 2)))’
+    by real_tac
+  >> ‘-(x pow approxSteps) * inv (&FACT approxSteps) * sin t =
+      -(sin t * ((x pow approxSteps) * inv (&FACT approxSteps)))’
+    by real_tac
+  >> rewrite_tac []
+  >> ntac 2 $ pop_assum $ rewrite_tac o single
+  >> rewrite_tac[GSYM REAL_MUL_ASSOC]
+  >> qmatch_goalsub_abbrev_tac ‘_ * err_sin_concr’
+  >> rewrite_tac [ABS_NEG, Once ABS_MUL]
+  >> irule REAL_LE_TRANS
+  >> qexists_tac ‘ 1 * abs err_sin_concr’ >> conj_tac
+  >- (irule REAL_LE_RMUL_IMP >> unabbrev_all_tac >> gs[SIN_BOUND, ABS_POS])
+  >> rewrite_tac [REAL_MUL_LID, sin_err_def, ABS_MUL]
+  >> ‘abs err_sin_concr = err_sin_concr’
+    by (unabbrev_all_tac
+        >> rewrite_tac[ABS_REFL]
+        >> irule REAL_LE_MUL >> conj_tac
+        >> gs[REAL_POW_GE0]
+        >> irule REAL_LE_MUL >> gs[REAL_POS, REAL_POW_GE0])
+  >> pop_assum $ rewrite_tac o single
+  >> unabbrev_all_tac
+  >> rewrite_tac [sin_err_def]
+  >> imp_res_tac EVEN_ODD_EXISTS >> gs[POW_MINUS1]
+  >> irule REAL_LE_LMUL_IMP >> gs[GSYM POW_ABS]
+  >> irule REAL_LE_TRANS
+  >> qexists_tac ‘abs (x pow (2 * m'))’ >> gs[ABS_LE, GSYM POW_ABS]
+  >> irule POW_LE >> gs[ABS_POS]
+  >> irule RealSimpsTheory.maxAbs >> gs[]
+QED
+
 (*
 Theorem checker_soundness:
   ∀ cert approxSteps.

realZeroLib.sml

 structure realZeroLib =
 struct
   open RealArith realTheory realLib realSyntax polyTheory;
-  open realPolyTheory checkerTheory sturmComputeTheory;
+  open realPolyTheory realPolyProofsTheory checkerTheory sturmComputeTheory;
   open bossLib preamble;
 
     exception ZeroLibErr of string;
@@ -146,6 +146,7 @@ struct
             (let val (instr, outstr) = Unix.streamsOf(Unix.execute("/usr/bin/which", ["sollya"]))
             in TextIO.inputAll(instr) end)
             handle SysErr _ => (print "Could not get path for Sollya\n"; "");
+      val sollyaPath = explode sollyaPath |> List.rev |> (fn ls => if hd ls = #"\n" then tl ls else ls) |> List.rev |> implode
       val (instr, outStr) =
           (Unix.streamsOf(Unix.execute(sollyaPath, ["--warnonstderr", "/tmp/sollya_input.sollya"])))
           handle SysErr _ => (print ("Could not run Sollya at "^sollyaPath ^ "\n"); raise ZeroLibErr "")
@@ -275,6 +276,7 @@ struct
           (let val (instr, outstr) = Unix.streamsOf(Unix.execute("/usr/bin/which", ["sollya"]))
           in TextIO.inputAll(instr) end)
           handle SysErr _ => (print "Could not get path for Sollya\n"; "");
+    val sollyaPath = explode sollyaPath |> List.rev |> tl |> List.rev |> implode
     val (instr, outStr) =
         (Unix.streamsOf(Unix.execute(sollyaPath, ["--warnonstderr", "/tmp/sollya_input.sollya"])))
         handle SysErr _ => (print ("Could not run Sollya at "^sollyaPath ^ "\n"); raise ZeroLibErr "")
@@ -283,15 +285,35 @@ struct
     print res
   end;
 
+val REAL_ABS_TRIANGLE_PRE = Q.prove (‘! (lb:real) ub f g h err1 err2.
+  (! x. lb <= x /\ x <= ub ==> abs (f x - g x) <= err1) ==>
+  (! x. lb <= x /\ x <= ub ==> abs (g x - h x) <= err2) ==>
+  (! x. lb <= x /\ x <= ub ==> abs (f x - h x) <= err1 + err2)’,
+  rpt strip_tac
+  >> ‘f x - h x = (f x - g x) + (g x - h x)’ by real_tac
+  >> pop_assum $ rewrite_tac o single
+  >> irule REAL_LE_TRANS >> qexists_tac ‘abs (f x - g x) + abs (g x - h x)’
+  >> gs[REAL_ABS_TRIANGLE] >> irule REAL_LE_ADD2
+  >> gs[]);
+
   fun validateCert (defTh:thm) =
   let
      val (transc, poly, eps, iv, hints) = destCert (defTh |> concl |> rhs)
+     val iv_FST = EVAL “FST ^iv”
+     val iv_SND = EVAL “SND ^iv”
+     val approxSideThm = Parse.Term ‘approxPolySideCond 16’ |> EVAL |> SIMP_RULE std_ss [EQ_CLAUSES]
      val approxThm = Parse.Term ‘approxPoly ^transc ^iv  ^hints 16’ |> EVAL
+     val approxSoundThm = MATCH_MP (MATCH_MP (REWRITE_RULE [GSYM AND_IMP_INTRO] approxPoly_soundness) approxSideThm) approxThm
+                          |> REWRITE_RULE [iv_FST, iv_SND] |> REWRITE_RULE [AND_IMP_INTRO]
      val (transp, err) =  approxThm |> concl |> rhs |> optionSyntax.dest_some |> pairSyntax.dest_pair
-     val errorp = Parse.Term ‘^transp -p ^poly’ |> EVAL |> rhs o concl
+     val errorpThm = Parse.Term ‘^transp -p ^poly’ |> EVAL
+     val errorp = errorpThm |> rhs o concl
+     val polyErrThm = (testSollya();
+      REAL_ZERO_CONV (Parse.Term ‘! x. FST (^iv) <= x /\ x <= SND (^iv) ==> evalPoly ^errorp x <= ^eps - ^err’))
   in
-    (testSollya();
-    REAL_ZERO_CONV (Parse.Term ‘! x. FST (^iv) <= x /\ x <= SND (^iv) ==> evalPoly ^errorp x <= ^eps - ^err’))
+    MATCH_MP
+      (MATCH_MP REAL_ABS_TRIANGLE_PRE approxSoundThm)
+      (REWRITE_RULE [GSYM errorpThm, eval_simps, GSYM poly_compat] polyErrThm)
   end;
 
   (** Some tests **)