needs "tame/ArcProperties.hl";;
needs "../formal_lp/arith/float_atn.hl";;


let REAL_LE_POW_2 = 
prove(`!x. &0 <= x pow 2`,
REWRITE_TAC[REAL_POW_2; REAL_LE_SQUARE]);;
let REAL_CONTINUOUS_ON_EQ_REAL_CONTINUOUS_AT = 
prove(`!f s. real_open s ==> (f real_continuous_on s <=> (!x. x IN s ==> f real_continuous atreal x))`,
REWRITE_TAC[REAL_OPEN; REAL_CONTINUOUS_ON; REAL_CONTINUOUS_CONTINUOUS_ATREAL] THEN REPEAT STRIP_TAC THEN MP_TAC (ISPECL[`lift o f o drop`; `IMAGE lift s`] CONTINUOUS_ON_EQ_CONTINUOUS_AT) THEN ASM_REWRITE_TAC[] THEN DISCH_THEN (fun th -> REWRITE_TAC[th; IN_IMAGE_LIFT_DROP]) THEN EQ_TAC THEN REPEAT STRIP_TAC THENL [ FIRST_X_ASSUM MATCH_MP_TAC THEN ASM_REWRITE_TAC[LIFT_DROP]; FIRST_X_ASSUM (MP_TAC o SPEC `drop x`) THEN ASM_REWRITE_TAC[LIFT_DROP] ]);;
let ALL_EL = 
prove (`!P l. (!i. i < LENGTH l ==> P (EL i l)) <=> ALL P l`,
REWRITE_TAC[GSYM ALL_MEM; MEM_EXISTS_EL] THEN MESON_TAC[]);;
needs "../formal_lp/formal_interval/theory/taylor_interval.hl";; open Arith_misc;; open Float_theory;; open Interval_arith;; open Arith_float;; open Float_atn;; let RULE = UNDISCH_ALL o Arith_nat.NUMERALS_TO_NUM o REWRITE_RULE[FLOAT_OF_NUM; min_exp_def; GSYM IMP_IMP] o SPEC_ALL;; let x_var_real = `x:real` and y_var_real = `y:real` and w_var_real = `w:real` and f_var_fun = `f : real->real` and g_var_fun = `g : real->real` and f1_var_fun = `f1 : real->real` and f2_var_fun = `f2 : real->real` and int_var = `int : real#real` and f_bounds_var = `f_bounds : real#real` and df_bounds_var = `df_bounds : real#real` and dd_bounds_var = `dd_bounds : real#real` and x_lo_var = `x_lo : real` and x_hi_var = `x_hi : real` and dd_var_real = `dd : real` and df_lo_var = `df_lo : real` and df_hi_var = `df_hi : real` and df_var_real = `df : real` and f_lo_var = `f_lo : real` and f_hi_var = `f_hi : real` and s_var_bool = `s:bool` and n_var_num = `n:num` and e_var_num = `e:num`;; let add_op_real = `(+):real->real->real` and mul_op_real = `( * ):real->real->real` and sub_op_real = `(-):real->real->real` and div_op_real = `(/):real->real->real` and inv_op_real = `inv` and neg_op_real = `--` and eq_op_real = `(=):real->real->bool` and lt_op_real = `(<):real->real->bool`;; let w1_var_real = `w1 : real` and w2_var_real = `w2 : real` and t_var_real = `t : real` and g_bounds_var = `g_bounds : real#real` and dg_bounds_var = `dg_bounds : real#real` and bounds_var = `bounds : real#real` and d_bounds_var = `d_bounds : real#real` and x0_var_real = `x0 : real` and z0_var_real = `z0 : real` and w0_var_real = `w0 : real`;; (*************************************) (* More float *) (* Converts a float term to the corresponding rational number *) let num_of_float_tm tm = let s, n_tm, e_tm = dest_float tm in let b = Num.num_of_int Arith_options.base in let m = Num.num_of_int Arith_options.min_exp in let ( * ), (^), (-), (!) = ( */ ), ( **/ ), (-/), Arith_nat.raw_dest_hash in let r = !n_tm * (b ^ (!e_tm - m)) in if s = "T" then minus_num r else r;; let float_of_float_tm tm = float_of_num (num_of_float_tm tm);; let mk_float n e = let n, s = if n < 0 then -n, `T` else n, `F` in let n_tm = rand (Arith_nat.mk_small_numeral_array n) in let e_tm = rand (Arith_nat.mk_small_numeral_array e) in mk_comb(mk_comb(mk_comb (`float`, s), n_tm), e_tm);; let float0_eq = FLOAT_TO_NUM_CONV (mk_float 0 Arith_options.min_exp);; let float1_eq = FLOAT_TO_NUM_CONV (mk_float 1 Arith_options.min_exp);; let float2_eq = FLOAT_TO_NUM_CONV (mk_float 2 Arith_options.min_exp);; let float3_eq = FLOAT_TO_NUM_CONV (mk_float 3 Arith_options.min_exp);; let float4_eq = FLOAT_TO_NUM_CONV (mk_float 4 Arith_options.min_exp);; let num1_eq = (SYM o REWRITE_RULE[Arith_hash.NUM_THM] o Arith_nat.NUMERAL_TO_NUM_CONV) `1`;; let num2_eq = (SYM o REWRITE_RULE[Arith_hash.NUM_THM] o Arith_nat.NUMERAL_TO_NUM_CONV) `2`;; let num3_eq = (SYM o REWRITE_RULE[Arith_hash.NUM_THM] o Arith_nat.NUMERAL_TO_NUM_CONV) `3`;; (*********************) let interval_tm = `interval_arith`;; let mk_interval tm bounds = mk_comb(mk_comb(interval_tm, tm), bounds);; (* const interval *) let CONST_INTERVAL' = SPEC_ALL CONST_INTERVAL;; let mk_const_interval tm = INST[tm, x_var_real] CONST_INTERVAL';; let float_F_const = `float F`;; let mk_float_small n = let n_tm0 = mk_small_numeral n in let n_th = Arith_nat.NUMERAL_TO_NUM_CONV n_tm0 in let n_tm = rand(rand(concl n_th)) in mk_comb(mk_comb(float_F_const, n_tm), min_exp_num_const);; let one_float = mk_float_small 1;; let two_float = mk_float_small 2;; let one_interval = mk_float_interval_small_num 1;; let two_interval = mk_float_interval_small_num 2;; let neg_two_interval = float_interval_neg two_interval;; (***********************************) (* float_eq0 *) let FLOAT_EQ_0' = (GEN_REWRITE_RULE (RAND_CONV o RAND_CONV) [NUMERAL] o SPEC_ALL) FLOAT_EQ_0;; let float_eq0 f_tm = let lhs, e_tm = dest_comb f_tm in let lhs2, n_tm = dest_comb lhs in let th0 = INST[rand lhs2, s_var_bool; n_tm, n_var_num; e_tm, e_var_num] FLOAT_EQ_0' in let eq_th = Arith_nat.raw_eq0_hash_conv n_tm in TRANS th0 eq_th;; (***********************************) (* float_interval_scale *) let float_interval_scale pp c_tm th = let c_th = mk_const_interval c_tm in float_interval_mul pp c_th th;; (**********************************) (* float_pos *) let FLOAT_F_POS' = SPEC_ALL FLOAT_F_POS;; (* Returns &0 <= float F n e *) let float_pos tm = let _, n_tm, e_tm = dest_float tm in INST[n_tm, n_var_num; e_tm, e_var_num] FLOAT_F_POS';; (************************************) (* float_iabs *) let FLOAT_NEG_F' = (RULE) FLOAT_NEG_T;; let FLOAT_NEG_T' = (RULE) FLOAT_NEG_F;; let float_neg tm = let sign, n_tm, e_tm = dest_float tm in if sign = "T" then INST[n_tm, n_var_num; e_tm, e_var_num] FLOAT_NEG_T' else INST[n_tm, n_var_num; e_tm, e_var_num] FLOAT_NEG_F';; let IABS' = RULE iabs;; let float_iabs int_tm = let lo_tm, hi_tm = dest_pair int_tm in let neg_lo_th = float_neg lo_tm in let max_th = SYM (float_max hi_tm ((rand o rator o concl) neg_lo_th)) in let lhs, rhs = dest_comb (concl max_th) in let th0 = SYM (EQ_MP (AP_TERM lhs (AP_TERM (rator rhs) neg_lo_th)) max_th) in let th1 = INST[lo_tm, x_lo_var; hi_tm, x_hi_var] IABS' in TRANS th1 th0;;
let FLOAT_IABS_FF = 
prove(`iabs (float F n1 e1, float F n2 e2) = float F n2 e2`,
REWRITE_TAC[iabs] THEN MP_TAC (SPECL [`n1:num`; `e1:num`] FLOAT_F_POS) THEN MP_TAC (SPECL [`n2:num`; `e2:num`] FLOAT_F_POS) THEN REAL_ARITH_TAC);;
let FLOAT_IABS_TT = 
prove(`iabs (float T n1 e1, float T n2 e2) = float F n1 e1`,
REWRITE_TAC[iabs; GSYM FLOAT_NEG_F] THEN MP_TAC (SPECL [`n1:num`; `e1:num`] FLOAT_F_POS) THEN MP_TAC (SPECL [`n2:num`; `e2:num`] FLOAT_T_NEG) THEN REAL_ARITH_TAC);;
(* let FLOAT_IABS_FT = prove(`iabs (float F n1 e1, float T n2 e2) *) (****************************) (* interval_not_zero *) let INTERVAL_NOT_ZERO1' = (UNDISCH_ALL o prove) (`(&0 < lo <=> T) ==> interval_not_zero (lo, hi)`, SIMP_TAC[interval_not_zero]);; let INTERVAL_NOT_ZERO2' = (UNDISCH_ALL o prove) (`(hi < &0 <=> T) ==> interval_not_zero (lo, hi)`, SIMP_TAC[interval_not_zero]);; let lo_var_real = `lo:real` and hi_var_real = `hi:real`;; let check_interval_not_zero int_tm = let lo, hi = dest_pair int_tm in let inst = INST[lo, lo_var_real; hi, hi_var_real] in let s1, _, _ = dest_float lo in if s1 = "F" then let gt_th = float_gt0 lo in if (fst o dest_const o rand o concl) gt_th <> "T" then failwith "check_interval_not_zero: &0 < lo <=> F" else (MY_PROVE_HYP gt_th o inst) INTERVAL_NOT_ZERO1' else let lt_th = float_lt0 hi in if (fst o dest_const o rand o concl) lt_th <> "T" then failwith "check_interval_not_zero: hi < &0 <=> F" else (MY_PROVE_HYP lt_th o inst) INTERVAL_NOT_ZERO2';; (* let int1 = mk_pair (mk_float 3 50, mk_float 4 50);; let int2 = mk_pair (mk_float (-4) 50, mk_float (-3) 50);; check_interval_not_zero int1;; check_interval_not_zero int2;; *) (*************************************) (* interval_pos *) let INTERVAL_POS' = (UNDISCH_ALL o prove) (`(&0 < lo <=> T) ==> interval_pos (lo, hi:real)`, SIMP_TAC[interval_pos]);; let check_interval_pos int_tm = let lo, hi = dest_pair int_tm in let gt_th = float_gt0 lo in if (fst o dest_const o rand o concl) gt_th <> "T" then failwith "check_interval_pos: &0 < lo <=> F" else (MY_PROVE_HYP gt_th o INST[lo, lo_var_real; hi, hi_var_real]) INTERVAL_POS';; (* check_interval_pos int1;; check_interval_pos int2;; *) (************************************) (* check_interval_iabs *) (* proves iabs int < rhs *) let check_interval_iabs int_tm rhs_tm = let iabs_eq = float_iabs int_tm in let lt_th = float_lt (rand (concl iabs_eq)) rhs_tm in if (fst o dest_const o rand o concl) lt_th <> "T" then failwith "check_interval_iabs: iabs < rhs <=> F" else let th0 = AP_THM (AP_TERM lt_op_real iabs_eq) rhs_tm in TRANS th0 lt_th;; (* check_interval_iabs int2 (mk_float 41 49);; *) (****************************) (* inv *)
let INV_EQ_DIV_LEMMA = 
prove(`&1 / x = inv x`,
REWRITE_TAC[real_div; REAL_MUL_LID]);;
let float_interval_inv pp th = let x_tm = (rand o rator o concl) th in let div_th = INST[x_tm, x_var_real] INV_EQ_DIV_LEMMA in let th0 = float_interval_div pp one_interval th in let lhs, rhs = dest_comb (concl th0) in let lhs2, rhs2 = dest_comb lhs in EQ_MP (AP_THM (AP_TERM lhs2 div_th) rhs) th0;; (*****************************************) (* bounded_on_int *) let norm_derivative d_th eq_th = let lhs, rhs = (dest_eq o concl) d_th in let lhs2, rhs2 = dest_comb lhs in let th0 = AP_THM (AP_TERM (rator lhs2) eq_th) rhs2 in TRANS (SYM th0) d_th;; let norm_diff d_th eq_th = let lhs, rhs = (dest_comb o concl) d_th in let th0 = AP_THM (AP_TERM (rator lhs) eq_th) rhs in EQ_MP th0 d_th;; let norm_interval int_th eq_th = let lhs, rhs = (dest_comb o concl) int_th in let th0 = AP_THM (AP_TERM (rator lhs) eq_th) rhs in EQ_MP (SYM th0) int_th;; let norm_second_derivative th eq_th = let lhs, dd_bounds = dest_comb (concl th) in let lhs2, int_tm = dest_comb lhs in let th0 = AP_THM (AP_THM (AP_TERM (rator lhs2) eq_th) int_tm) dd_bounds in EQ_MP th0 th;; let norm_lin_approx th eq_th = let lhs, df_bounds = dest_comb (concl th) in let lhs2, f_bounds = dest_comb lhs in let lhs3, x_tm = dest_comb lhs2 in let th0 = AP_THM (AP_THM (AP_THM (AP_TERM (rator lhs3) eq_th) x_tm) f_bounds) df_bounds in EQ_MP th0 th;; let BOUNDED_ON_INT = (UNDISCH_ALL o prove)(`(!x. interval_arith x int ==> interval_arith (f x) f_bounds) ==> bounded_on_int f int f_bounds`, REWRITE_TAC[bounded_on_int; interval_arith; IN_REAL_INTERVAL]);; let BOUNDED_ON_INT_DEST = (UNDISCH_ALL o prove)(`bounded_on_int f int f_bounds ==> (!x. interval_arith x int ==> interval_arith (f x) f_bounds)`, REWRITE_TAC[bounded_on_int; interval_arith; IN_REAL_INTERVAL]);; let mk_bounded_on_int th = let int_tm = (rand o hd o hyp) th in let lhs, f_bounds_tm = dest_comb (concl th) in let lhs2, rhs2 = dest_comb lhs in let f_tm = mk_abs (x_var_real, rhs2) in let b_th0 = (SYM o BETA_CONV) (mk_comb (f_tm , x_var_real)) in let b_th1 = AP_THM (AP_TERM lhs2 b_th0) f_bounds_tm in let th2 = EQ_MP b_th1 th in let th3 = DISCH_ALL th2 in let th4 = GEN x_var_real th3 in let th_int = INST[int_tm, int_var; f_bounds_tm, f_bounds_var; f_tm, f_var_fun] BOUNDED_ON_INT in MY_PROVE_HYP th4 th_int;; let dest_bounded_on_int th = let lhs, f_bounds = dest_comb (concl th) in let lhs2, int_tm = dest_comb lhs in let f_tm = rand lhs2 in let th0 = INST[f_tm, f_var_fun; int_tm, int_var; f_bounds, f_bounds_var] BOUNDED_ON_INT_DEST in let th1 = UNDISCH_ALL (SPEC x_var_real (MY_PROVE_HYP th th0)) in if is_abs f_tm then let f_tm = (rand o rator o concl) th1 in let eq_th = BETA_CONV f_tm in norm_interval th1 (SYM eq_th) else th1;; let dest_bounded_on_int_raw th = let lhs, f_bounds = dest_comb (concl th) in let lhs2, int_tm = dest_comb lhs in let f_tm = rand lhs2 in let th0 = INST[f_tm, f_var_fun; int_tm, int_var; f_bounds, f_bounds_var] BOUNDED_ON_INT_DEST in UNDISCH_ALL (SPEC x_var_real (MY_PROVE_HYP th th0));; (* let pp = 5;; let th0 = (RULE o ASSUME) `interval_arith x (&2, &3)`;; let th = float_interval_atn pp th0;; let b_th = mk_bounded_on_int th;; let d_th = dest_bounded_on_int b_th;; *) (***********************************) (* bounded_on_int arithmetic *) let bounded_on_int_scale pp c_tm th = let i_th = dest_bounded_on_int th in let th0 = float_interval_scale pp c_tm i_th in mk_bounded_on_int th0;; let bounded_on_int_mul_int pp int_th th = let i_th = dest_bounded_on_int th in let th0 = float_interval_mul pp int_th i_th in mk_bounded_on_int th0;; let bounded_on_int_neg th1 = let i_th = dest_bounded_on_int th1 in let th0 = float_interval_neg i_th in mk_bounded_on_int th0;; let bounded_on_int_add pp th1 th2 = let i_th1, i_th2 = dest_bounded_on_int th1, dest_bounded_on_int th2 in let th0 = float_interval_add pp i_th1 i_th2 in mk_bounded_on_int th0;; let bounded_on_int_sub pp th1 th2 = let i_th1, i_th2 = dest_bounded_on_int th1, dest_bounded_on_int th2 in let th0 = float_interval_sub pp i_th1 i_th2 in mk_bounded_on_int th0;; let bounded_on_int_mul pp th1 th2 = let i_th1, i_th2 = dest_bounded_on_int th1, dest_bounded_on_int th2 in let th0 = float_interval_mul pp i_th1 i_th2 in mk_bounded_on_int th0;; let bounded_on_int_mul_raw pp th1 th2 = let i_th1, i_th2 = dest_bounded_on_int_raw th1, dest_bounded_on_int_raw th2 in let th0 = float_interval_mul pp i_th1 i_th2 in mk_bounded_on_int th0;; let bounded_on_int_div pp th1 th2 = let i_th1, i_th2 = dest_bounded_on_int th1, dest_bounded_on_int th2 in let th0 = float_interval_div pp i_th1 i_th2 in mk_bounded_on_int th0;; (* let b_sub = bounded_on_int_sub pp b_th b_th;; bounded_on_int_mul pp b_sub b_th;; *) (************************************) let ADD_INEQ_HI = (RULE o REAL_ARITH) `x1 <= y1 /\ x2 <= y2 /\ y1 + y2 <= y ==> x1 + x2 <= y`;; let ADD_INEQ_LO = (RULE o REAL_ARITH) `x1 <= y1 /\ x2 <= y2 /\ x <= x1 + x2 ==> x <= y1 + y2`;; let SUB_INEQ_HI = (RULE o REAL_ARITH) `x1 <= y1 /\ y2 <= x2 /\ y1 - y2 <= y ==> x1 - x2 <= y`;; let SUB_INEQ_LO = (RULE o REAL_ARITH) `x1 <= y1 /\ y2 <= x2 /\ x <= x1 - x2 ==> x <= y1 - y2`;; let MUL_INEQ_HI = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove) (`&0 <= x1 /\ &0 <= x2 /\ x1 <= y1 /\ x2 <= y2 /\ y1 * y2 <= y ==> x1 * x2 <= y`, DISCH_TAC THEN MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `y1 * y2` THEN ASM_REWRITE_TAC[] THEN MATCH_MP_TAC REAL_LE_MUL2 THEN ASM_REWRITE_TAC[]);; let MUL_INEQ_POS_CONST_HI = (UNDISCH_ALL o prove) (`(&0 <= x <=> T) ==> y1 <= y2 ==> x * y2 <= z ==> x * y1 <= z`, REPEAT STRIP_TAC THEN MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `x * y2` THEN ASM_REWRITE_TAC[] THEN MATCH_MP_TAC REAL_LE_LMUL THEN ASM_REWRITE_TAC[]);; let REAL_LE_REFL' = RULE REAL_LE_REFL;; let mk_refl_ineq tm = INST[tm, x_var_real] REAL_LE_REFL';; let dest_le_op ineq = let lhs, y_tm = dest_comb ineq in (rand lhs, y_tm);; let y1_var_real = `y1:real` and x1_var_real = `x1:real` and y2_var_real = `y2:real` and x2_var_real = `x2:real` and z_var_real = `z:real`;; let mul_ineq_pos_const_hi pp c_tm ineq = let y1_tm, y2_tm = dest_le_op (concl ineq) in let ge0_th = float_ge0 c_tm in let mul_hi_th = float_mul_hi pp c_tm y2_tm in let z_tm = (rand o concl) mul_hi_th in (MY_PROVE_HYP ge0_th o MY_PROVE_HYP ineq o MY_PROVE_HYP mul_hi_th o INST[c_tm, x_var_real; y1_tm, y1_var_real; y2_tm, y2_var_real; z_tm, z_var_real]) MUL_INEQ_POS_CONST_HI;; let mul_ineq_hi pp ineq1 ineq2 = let x1_tm, y1_tm = dest_le_op (concl ineq1) in let x2_tm, y2_tm = dest_le_op (concl ineq2) in let x1_pos, x2_pos = float_pos x1_tm, float_pos x2_tm in let rhs_mul = float_mul_hi pp y1_tm y2_tm in let y_tm = (rand o concl) rhs_mul in let th0 = INST[x1_tm, x1_var_real; y1_tm, y1_var_real; x2_tm, x2_var_real; y2_tm, y2_var_real; y_tm, y_var_real] MUL_INEQ_HI in (MY_PROVE_HYP x1_pos o MY_PROVE_HYP x2_pos o MY_PROVE_HYP ineq1 o MY_PROVE_HYP ineq2 o MY_PROVE_HYP rhs_mul) th0;; let sub_ineq_hi pp ineq1 ineq2 = let x1_tm, y1_tm = dest_le_op (concl ineq1) in let y2_tm, x2_tm = dest_le_op (concl ineq2) in let rhs_sub = float_sub_hi pp y1_tm y2_tm in let y_tm = (rand o concl) rhs_sub in let th0 = INST[x1_tm, x1_var_real; y1_tm, y1_var_real; x2_tm, x2_var_real; y2_tm, y2_var_real; y_tm, y_var_real] SUB_INEQ_HI in MY_PROVE_HYP ineq1 (MY_PROVE_HYP ineq2 (MY_PROVE_HYP rhs_sub th0));; let sub_ineq_lo pp ineq1 ineq2 = let x1_tm, y1_tm = dest_le_op (concl ineq1) in let y2_tm, x2_tm = dest_le_op (concl ineq2) in let lhs_sub = float_sub_lo pp x1_tm x2_tm in let x_tm = (lhand o concl) lhs_sub in let th0 = INST[x1_tm, x1_var_real; y1_tm, y1_var_real; x2_tm, x2_var_real; y2_tm, y2_var_real; x_tm, x_var_real] SUB_INEQ_LO in MY_PROVE_HYP ineq1 (MY_PROVE_HYP ineq2 (MY_PROVE_HYP lhs_sub th0));; let add_ineq_hi pp ineq1 ineq2 = let x1_tm, y1_tm = dest_le_op (concl ineq1) in let x2_tm, y2_tm = dest_le_op (concl ineq2) in let rhs_sum = float_add_hi pp y1_tm y2_tm in let y_tm = (rand o concl) rhs_sum in let th0 = INST[x1_tm, x1_var_real; y1_tm, y1_var_real; x2_tm, x2_var_real; y2_tm, y2_var_real; y_tm, y_var_real] ADD_INEQ_HI in MY_PROVE_HYP ineq1 (MY_PROVE_HYP ineq2 (MY_PROVE_HYP rhs_sum th0));; let add_ineq_lo pp ineq1 ineq2 = let x1_tm, y1_tm = dest_le_op (concl ineq1) in let x2_tm, y2_tm = dest_le_op (concl ineq2) in let lhs_sum = float_add_lo pp x1_tm x2_tm in let x_tm = (lhand o concl) lhs_sum in let th0 = INST[x1_tm, x1_var_real; y1_tm, y1_var_real; x2_tm, x2_var_real; y2_tm, y2_var_real; x_tm, x_var_real] ADD_INEQ_LO in MY_PROVE_HYP ineq1 (MY_PROVE_HYP ineq2 (MY_PROVE_HYP lhs_sum th0));;