formal_ineqs/taylor/m_taylor.hl

   1 (* =========================================================== *)\r
   2 (* Formal taylor intervals                                     *)\r
   3 (* Author: Alexey Solovyev                                     *)\r
   4 (* Date: 2012-10-27                                            *)\r
   5 (* =========================================================== *)\r
   6 \r
   7 needs "arith/more_float.hl";;\r
   8 needs "arith/float_atn.hl";;\r
   9 needs "arith/eval_interval.hl";;\r
  10 needs "list/list_conversions.hl";;\r
  11 needs "list/list_float.hl";;\r
  12 needs "list/more_list.hl";;\r
  13 needs "misc/vars.hl";;\r
  14 \r
  15 needs "lib/ssreflect/ssreflect.hl";;\r
  16 needs "lib/ssreflect/sections.hl";;\r
  17 needs "taylor/theory/taylor_interval-compiled.hl";;\r
  18 needs "taylor/theory/multivariate_taylor-compiled.hl";;\r
  19 \r
  20 \r
  21 module M_taylor = struct\r
  22 \r
  23 open Arith_misc;;\r
  24 open Arith_float;;\r
  25 open More_float;;\r
  26 open Float_theory;;\r
  27 open Eval_interval;;\r
  28 open List_conversions;;\r
  29 open List_float;;\r
  30 open More_list;;\r
  31 open Interval_arith;;\r
  32 open Misc_vars;;\r
  33 \r
  34 let MY_RULE = UNDISCH_ALL o PURE_REWRITE_RULE[GSYM IMP_IMP] o SPEC_ALL;;\r
  35 let MY_RULE_NUM = UNDISCH_ALL o Arith_nat.NUMERALS_TO_NUM o PURE_REWRITE_RULE[GSYM IMP_IMP] o SPEC_ALL;;\r
  36 let MY_RULE_FLOAT = UNDISCH_ALL o Arith_nat.NUMERALS_TO_NUM o \r
  37   PURE_REWRITE_RULE[FLOAT_OF_NUM; min_exp_def; GSYM IMP_IMP] o SPEC_ALL;;\r
  38 \r
  39 \r
  40 \r
  41 let max_dim = 8;;\r
  42 \r
  43 let inst_first_type_var ty th =\r
  44   let ty_vars = type_vars_in_term (concl th) in\r
  45     if ty_vars = [] then\r
  46       failwith "inst_first_type: no type variables in the theorem"\r
  47     else\r
  48       INST_TYPE [ty, hd ty_vars] th;;\r
  49 \r
  50 \r
  51 let float0 = mk_float 0 0 and\r
  52     interval0 = mk_float_interval_small_num 0;;\r
  53 \r
  54 \r
  55 let has_size_array = Array.init (max_dim + 1) \r
  56   (fun i -> match i with\r
  57      | 0 -> TRUTH\r
  58      | 1 -> HAS_SIZE_1\r
  59      | _ -> define_finite_type i);;\r
  60 \r
  61 let dimindex_array = Array.init (max_dim + 1) \r
  62   (fun i -> if i < 1 then TRUTH else MATCH_MP DIMINDEX_UNIQUE has_size_array.(i));;\r
  63 \r
  64 let n_type_array = Array.init (max_dim + 1)\r
  65   (fun i -> if i < 1 then bool_ty else \r
  66      let dimindex_th = dimindex_array.(i) in\r
  67        (hd o snd o dest_type o snd o dest_const o rand o lhand o concl) dimindex_th);;\r
  68 \r
  69 let n_vector_type_array = Array.init (max_dim + 1)\r
  70   (fun i -> if i < 1 then bool_ty else mk_type ("cart", [real_ty; n_type_array.(i)]));;\r
  71 \r
  72 \r
  73 let x_var_names = Array.init (max_dim + 1) (fun i -> "x"^(string_of_int i)) and\r
  74     y_var_names = Array.init (max_dim + 1) (fun i -> "y"^(string_of_int i)) and\r
  75     z_var_names = Array.init (max_dim + 1) (fun i -> "z"^(string_of_int i)) and\r
  76     w_var_names = Array.init (max_dim + 1) (fun i -> "w"^(string_of_int i));;\r
  77 \r
  78 let x_vars_array = Array.init (max_dim + 1) (fun i -> mk_var(x_var_names.(i), real_ty)) and\r
  79     y_vars_array = Array.init (max_dim + 1) (fun i -> mk_var(y_var_names.(i), real_ty)) and\r
  80     z_vars_array = Array.init (max_dim + 1) (fun i -> mk_var(z_var_names.(i), real_ty)) and\r
  81     w_vars_array = Array.init (max_dim + 1) (fun i -> mk_var(w_var_names.(i), real_ty));;\r
  82 \r
  83 let df_vars_array = Array.init (max_dim + 1) (fun i -> mk_var ("df"^(string_of_int i), real_pair_ty));;\r
  84 let dd_vars_array = Array.init (max_dim + 1) (fun i ->\r
  85    Array.init (max_dim + 1) (fun j -> mk_var ("dd"^(string_of_int i)^(string_of_int j), real_pair_ty)));;\r
  86 \r
  87 let dest_vector = dest_list o rand;;\r
  88 \r
  89 let mk_vector list_tm =\r
  90   let n = (length o dest_list) list_tm in\r
  91   let ty = (hd o snd o dest_type o type_of) list_tm in\r
  92   let vec = mk_const ("vector", [ty, aty; n_type_array.(n), nty]) in\r
  93     mk_comb (vec, list_tm);;\r
  94 \r
  95 let mk_vector_list list =\r
  96     mk_vector (mk_list (list, type_of (hd list)));;\r
  97 \r
  98 let el_thms_array =\r
  99   let el_tm = `EL : num->(A)list->A` in\r
 100   let gen0 n =\r
 101     let e_list = mk_list (map (fun i -> mk_var ("e"^(string_of_int i), aty)) (1--n), aty) in\r
 102     let el0_th = REWRITE_CONV[EL; HD] (mk_binop el_tm `0` e_list) in\r
 103       Array.create n el0_th in\r
 104   let array = Array.init (max_dim + 1) gen0 in\r
 105   let gen_i n i =\r
 106     let e_list = (rand o lhand o concl) array.(n).(i) in\r
 107     let prev_thm = array.(n - 1).(i - 1) in\r
 108     let i_tm = mk_small_numeral i in\r
 109     let prev_i = num_CONV i_tm in\r
 110     let el_th = REWRITE_CONV[prev_i; EL; HD; TL; prev_thm] (mk_binop el_tm i_tm e_list) in\r
 111       array.(n).(i) <- el_th in\r
 112   let _ = map (fun n -> map (fun i -> gen_i n i) (1--(n - 1))) (2--max_dim) in\r
 113     array;;\r
 114 \r
 115 \r
 116 let VECTOR_COMPONENT = prove(`!l i. i IN 1..dimindex (:N) ==>\r
 117                                (vector l:A^N)$i = EL (i - 1) l`,\r
 118 REWRITE_TAC[IN_NUMSEG] THEN REPEAT GEN_TAC THEN DISCH_TAC THEN REWRITE_TAC[vector] THEN\r
 119   MATCH_MP_TAC LAMBDA_BETA THEN ASM_REWRITE_TAC[]);;\r
 120 \r
 121 let gen_comp_thm n i =\r
 122   let i_tm = mk_small_numeral i and\r
 123       x_list = mk_list (map (fun i -> mk_var("x"^(string_of_int i), aty)) (1--n), aty) in\r
 124   let th0 = (ISPECL [x_list; i_tm] o inst_first_type_var (n_type_array.(n))) VECTOR_COMPONENT in\r
 125   let th1 = (CONV_RULE NUM_REDUCE_CONV o REWRITE_RULE[IN_NUMSEG; dimindex_array.(n)]) th0 in\r
 126     REWRITE_RULE[el_thms_array.(n).(i - 1)] th1;;\r
 127 \r
 128 let comp_thms_array = Array.init (max_dim + 1)\r
 129   (fun n -> Array.init (n + 1)\r
 130        (fun i -> if i < 1 or n < 1 then TRUTH else gen_comp_thm n i));;\r
 131 \r
 132 \r
 133 (************************************)\r
 134 (* m_cell_domain *)\r
 135 \r
 136 let ALL2_ALL_ZIP = prove(`!(P:A->B->bool) l1 l2. LENGTH l1 = LENGTH l2 ==> \r
 137     (ALL2 P l1 l2 <=> ALL (\p. P (FST p) (SND p)) (ZIP l1 l2))`,\r
 138   GEN_TAC THEN LIST_INDUCT_TAC THENL\r
 139     [\r
 140       GEN_TAC THEN REWRITE_TAC[LENGTH; EQ_SYM_EQ; LENGTH_EQ_NIL] THEN \r
 141         DISCH_THEN (fun th -> REWRITE_TAC[th]) THEN\r
 142         REWRITE_TAC[ZIP; ALL2; ALL];\r
 143       ALL_TAC\r
 144     ] THEN\r
 145 \r
 146     LIST_INDUCT_TAC THEN REWRITE_TAC[LENGTH] THENL [ARITH_TAC; ALL_TAC] THEN\r
 147     REWRITE_TAC[eqSS] THEN DISCH_TAC THEN\r
 148     REWRITE_TAC[ALL2; ZIP; ALL] THEN\r
 149     FIRST_X_ASSUM (new_rewrite [] []) THEN ASM_REWRITE_TAC[]);;\r
 150 \r
 151 let EL_ZIP = prove(`!(l1:(A)list) (l2:(B)list) i. LENGTH l1 = LENGTH l2 /\ i < LENGTH l1 ==> \r
 152     EL i (ZIP l1 l2) = (EL i l1, EL i l2)`,\r
 153   LIST_INDUCT_TAC THEN LIST_INDUCT_TAC THEN REWRITE_TAC[ZIP; LENGTH] THEN TRY ARITH_TAC THEN\r
 154     case THEN REWRITE_TAC[EL; HD; TL] THEN GEN_TAC THEN\r
 155     REWRITE_TAC[eqSS; ARITH_RULE `SUC n < SUC x <=> n < x`] THEN STRIP_TAC THEN\r
 156     FIRST_X_ASSUM MATCH_MP_TAC THEN ASM_REWRITE_TAC[]);;\r
 157 \r
 158 let LENGTH_ZIP = prove(`!l1 l2. LENGTH l1 = LENGTH l2 ==> LENGTH (ZIP l1 l2) = LENGTH l1`,\r
 159   LIST_INDUCT_TAC THEN LIST_INDUCT_TAC THEN REWRITE_TAC[ZIP; LENGTH] THEN TRY ARITH_TAC THEN\r
 160     REWRITE_TAC[eqSS] THEN DISCH_TAC THEN FIRST_X_ASSUM MATCH_MP_TAC THEN ASM_REWRITE_TAC[]);;\r
 161 \r
 162     \r
 163 \r
 164 let test_domain_xi = new_definition\r
 165   `test_domain_xi xz yw <=> FST xz <= FST yw /\ FST yw <= SND xz /\ \r
 166   FST yw - FST xz <= SND yw /\ SND xz - FST yw <= SND yw`;;\r
 167 \r
 168 \r
 169 let MK_CELL_DOMAIN = prove(`!xz (yw:(real#real)list) x z y w.\r
 170     LENGTH x = dimindex (:N) /\ LENGTH z = dimindex (:N) /\\r
 171     LENGTH y = dimindex (:N) /\ LENGTH w = dimindex (:N) /\\r
 172     ZIP y w = yw /\ ZIP x z = xz /\\r
 173     ALL2 test_domain_xi xz yw ==>\r
 174     m_cell_domain (vector x, vector z:real^N) (vector y) (vector w)`,\r
 175   REPEAT GEN_TAC THEN STRIP_TAC THEN POP_ASSUM MP_TAC THEN\r
 176     SUBGOAL_THEN `LENGTH (xz:(real#real)list) = dimindex (:N) /\ LENGTH (yw:(real#real)list) = dimindex (:N)` ASSUME_TAC THENL\r
 177     [\r
 178       EXPAND_TAC "yw" THEN EXPAND_TAC "xz" THEN\r
 179         REPEAT (new_rewrite [] [] LENGTH_ZIP) THEN ASM_REWRITE_TAC[];\r
 180       ALL_TAC\r
 181     ] THEN\r
 182     rewrite [] [] ALL2_ALL_ZIP THEN ASM_REWRITE_TAC[m_cell_domain; GSYM ALL_EL] THEN DISCH_TAC THEN\r
 183     REWRITE_TAC[m_cell_domain] THEN GEN_TAC THEN DISCH_TAC THEN\r
 184     REPEAT (new_rewrite [] [] VECTOR_COMPONENT) THEN ASM_REWRITE_TAC[] THEN\r
 185     ABBREV_TAC `j = i - 1` THEN\r
 186     SUBGOAL_THEN `j < dimindex (:N)` ASSUME_TAC THENL\r
 187     [\r
 188       POP_ASSUM MP_TAC THEN POP_ASSUM MP_TAC THEN REWRITE_TAC[IN_NUMSEG] THEN ARITH_TAC;\r
 189       ALL_TAC\r
 190     ] THEN\r
 191     FIRST_X_ASSUM (MP_TAC o SPEC `j:num`) THEN REWRITE_TAC[test_domain_xi] THEN\r
 192     rewrite [] [] LENGTH_ZIP THEN ASM_REWRITE_TAC[] THEN\r
 193     rewrite [] [] EL_ZIP THEN ASM_REWRITE_TAC[] THEN\r
 194     EXPAND_TAC "xz" THEN EXPAND_TAC "yw" THEN\r
 195     REPEAT (new_rewrite [] [] EL_ZIP) THEN ASM_REWRITE_TAC[] THEN\r
 196     ARITH_TAC);;\r
 197 \r
 198 \r
 199 \r
 200 (* array of theorems *)\r
 201 let mk_m_domain_array =\r
 202   let mk_m_domain n =\r
 203     let dimindex_th = dimindex_array.(n) in\r
 204     let n_ty = (hd o snd o dest_type o snd o dest_const o rand o lhand o concl) dimindex_th in\r
 205     let nty = `:N` in\r
 206       (UNDISCH_ALL o REWRITE_RULE[float0_eq] o DISCH_ALL o RULE o \r
 207          REWRITE_RULE[dimindex_th] o INST_TYPE[n_ty, nty]) MK_CELL_DOMAIN in\r
 208     Array.init (max_dim + 1) (fun i -> if i < 1 then TRUTH else mk_m_domain i);;\r
 209 \r
 210 \r
 211 let TEST_DOMAIN_XI' = (EQT_INTRO o RULE o prove)(`xz = (x,z) /\ yw = (y,w) /\\r
 212     x <= y /\ y <= z /\ y - x <= w1 /\ z - y <= w2 /\ w1 <= w /\ w2 <= w ==> test_domain_xi xz yw`,\r
 213   SIMP_TAC[test_domain_xi] THEN REAL_ARITH_TAC);;\r
 214 \r
 215 \r
 216 let eval_test_domain_xi pp test_domain_tm =\r
 217   let ltm, yw = dest_comb test_domain_tm in\r
 218   let xz = rand ltm in\r
 219   let x, z = dest_pair xz and\r
 220       y, w = dest_pair yw in\r
 221   let (<=) = (fun t1 t2 -> EQT_ELIM (float_le t1 t2)) and\r
 222       (-) = float_sub_hi pp in\r
 223   let x_le_y = x <= y and\r
 224       y_le_z = y <= z and\r
 225       yx_le_w1 = y - x and\r
 226       zy_le_w2 = z - y in\r
 227   let w1 = (rand o concl) yx_le_w1 and\r
 228       w2 = (rand o concl) zy_le_w2 in\r
 229   let w1_le_w = w1 <= w and\r
 230       w2_le_w = w2 <= w in\r
 231     (MY_PROVE_HYP (REFL xz) o MY_PROVE_HYP (REFL yw) o\r
 232       MY_PROVE_HYP x_le_y o MY_PROVE_HYP y_le_z o \r
 233        MY_PROVE_HYP yx_le_w1 o MY_PROVE_HYP zy_le_w2 o\r
 234        MY_PROVE_HYP w1_le_w o MY_PROVE_HYP w2_le_w o\r
 235        INST[x, x_var_real; y, y_var_real; z, z_var_real; w, w_var_real;\r
 236             w1, w1_var_real; w2, w2_var_real; \r
 237             xz, xz_pair_var; yw, yw_pair_var]) TEST_DOMAIN_XI';;\r
 238     \r
 239 \r
 240 (* mk_m_center_domain *)\r
 241 let mk_m_center_domain n pp x_list_tm z_list_tm =\r
 242   let x_list = dest_list x_list_tm and\r
 243       z_list = dest_list z_list_tm in\r
 244   let y_list =\r
 245     let ( * ) = (fun t1 t2 -> (rand o concl) (float_mul_eq t1 t2)) and\r
 246         (+) = (fun t1 t2 -> (rand o concl) (float_add_hi pp t1 t2)) in\r
 247       map2 (fun x y -> if x = y then x else float_inv2 * (x + y)) x_list z_list in\r
 248 \r
 249   let w_list =\r
 250     let (-) = (fun t1 t2 -> (rand o concl) (float_sub_hi pp t1 t2)) and\r
 251         max = (fun t1 t2 -> (rand o concl) (float_max t1 t2)) in\r
 252     let w1 = map2 (-) y_list x_list and\r
 253         w2 = map2 (-) z_list y_list in\r
 254       map2 max w1 w2 in\r
 255 \r
 256   let y_list_tm = mk_list (y_list, real_ty) and\r
 257       w_list_tm = mk_list (w_list, real_ty) in\r
 258 \r
 259   let yw_zip_th = eval_zip y_list_tm w_list_tm and\r
 260       xz_zip_th = eval_zip x_list_tm z_list_tm in\r
 261 \r
 262   let yw_list_tm = (rand o concl) yw_zip_th and\r
 263       xz_list_tm = (rand o concl) xz_zip_th in\r
 264 \r
 265   let len_x_th = eval_length x_list_tm and\r
 266       len_z_th = eval_length z_list_tm and\r
 267       len_y_th = eval_length y_list_tm and\r
 268       len_w_th = eval_length w_list_tm in\r
 269   let th0 = (MY_PROVE_HYP len_x_th o MY_PROVE_HYP len_z_th o\r
 270                MY_PROVE_HYP len_y_th o MY_PROVE_HYP len_w_th o\r
 271                MY_PROVE_HYP yw_zip_th o MY_PROVE_HYP xz_zip_th o\r
 272                INST[x_list_tm, x_var_real_list; z_list_tm, z_var_real_list;\r
 273                     y_list_tm, y_var_real_list; w_list_tm, w_var_real_list;\r
 274                     yw_list_tm, yw_var; xz_list_tm, xz_var]) mk_m_domain_array.(n) in\r
 275   let all_th = (EQT_ELIM o all2_conv_univ (eval_test_domain_xi pp) o hd o hyp) th0 in\r
 276     MY_PROVE_HYP all_th th0;;\r
 277 \r
 278 \r
 279 \r
 280 (***********************)\r
 281 \r
 282 let MK_M_TAYLOR_INTERVAL' = (RULE o MATCH_MP iffRL o SPEC_ALL) m_taylor_interval;;\r
 283 \r
 284 \r
 285 let get_types_and_vars n =\r
 286   let ty = n_type_array.(n) and\r
 287       xty = n_vector_type_array.(n) in\r
 288   let x_var = mk_var ("x", xty) and\r
 289       f_var = mk_var ("f", mk_fun_ty xty real_ty) and\r
 290       y_var = mk_var ("y", xty) and\r
 291       w_var = mk_var ("w", xty) and\r
 292       domain_var = mk_var ("domain", mk_type ("prod", [xty; xty])) in\r
 293     ty, xty, x_var, f_var, y_var, w_var, domain_var;;\r
 294 \r
 295 \r
 296 \r
 297 let dest_m_cell_domain domain_tm =\r
 298   let lhs, w_tm = dest_comb domain_tm in\r
 299   let lhs2, y_tm = dest_comb lhs in\r
 300     rand lhs2, y_tm, w_tm;;\r
 301 \r
 302 \r
 303 (**************************************************)\r
 304 \r
 305 (* Given a variable of the type `:real^N`, returns the number N *)\r
 306 let get_dim = int_of_string o fst o dest_type o hd o tl o snd o dest_type o type_of;;\r
 307 \r
 308 \r
 309 (**********************)\r
 310 (* eval_m_taylor_poly *)\r
 311 \r
 312 let partial_pow = prove(`!i f n (y:real^N). lift o f differentiable at y ==>\r
 313                           partial i (\x. f x pow n) y = &n * f y pow (n - 1) * partial i f y`,\r
 314   REPEAT STRIP_TAC THEN\r
 315     SUBGOAL_THEN `(\x:real^N. f x pow n) = (\t. t pow n) o f` (fun th -> REWRITE_TAC[th]) THENL\r
 316     [\r
 317       ONCE_REWRITE_TAC[GSYM eq_ext] THEN REWRITE_TAC[o_THM];\r
 318       ALL_TAC\r
 319     ] THEN\r
 320     new_rewrite [] [] partial_uni_compose THENL\r
 321     [\r
 322       ASM_REWRITE_TAC[] THEN\r
 323         new_rewrite [] [] REAL_DIFFERENTIABLE_POW_ATREAL THEN REWRITE_TAC[REAL_DIFFERENTIABLE_ID];\r
 324       ALL_TAC\r
 325     ] THEN\r
 326     new_rewrite [] [] derivative_pow THEN REWRITE_TAC[REAL_DIFFERENTIABLE_ID; derivative_x] THEN\r
 327     REAL_ARITH_TAC);;\r
 328 \r
 329 \r
 330 let nth_diff2_pow = prove(`!n y. nth_diff_strong 2 (\x. x pow n) y`,\r
 331   REWRITE_TAC[nth_diff_strong2_eq] THEN REPEAT GEN_TAC THEN\r
 332     EXISTS_TAC `(:real)` THEN REWRITE_TAC[REAL_OPEN_UNIV; IN_UNIV] THEN GEN_TAC THEN\r
 333     new_rewrite [] [] REAL_DIFFERENTIABLE_POW_ATREAL THEN REWRITE_TAC[REAL_DIFFERENTIABLE_ID] THEN\r
 334     MATCH_MP_TAC differentiable_local THEN\r
 335     EXISTS_TAC `\x. &n * x pow (n - 1)` THEN EXISTS_TAC `(:real)` THEN\r
 336     REWRITE_TAC[REAL_OPEN_UNIV; IN_UNIV] THEN\r
 337     new_rewrite [] [] REAL_DIFFERENTIABLE_MUL_ATREAL THEN REWRITE_TAC[REAL_DIFFERENTIABLE_CONST] THENL\r
 338     [\r
 339       new_rewrite [] [] REAL_DIFFERENTIABLE_POW_ATREAL THEN REWRITE_TAC[REAL_DIFFERENTIABLE_ID];\r
 340       ALL_TAC\r
 341     ] THEN\r
 342     GEN_TAC THEN new_rewrite [] [] derivative_pow THEN REWRITE_TAC[REAL_DIFFERENTIABLE_ID] THEN\r
 343     REWRITE_TAC[derivative_x; REAL_MUL_RID]);;\r
 344 \r
 345 \r
 346 let diff2c_pow = prove(`!f n (x:real^N). diff2c f x ==> diff2c (\x. f x pow n) x`,\r
 347   REPEAT STRIP_TAC THEN\r
 348     SUBGOAL_THEN `(\x:real^N. f x pow n) = (\t. t pow n) o f` (fun th -> REWRITE_TAC[th]) THENL\r
 349     [\r
 350       ONCE_REWRITE_TAC[GSYM eq_ext] THEN REWRITE_TAC[o_THM];\r
 351       ALL_TAC\r
 352     ] THEN\r
 353     apply_tac diff2c_uni_compose THEN ASM_REWRITE_TAC[nth_diff2_pow] THEN\r
 354     REWRITE_TAC[nth_derivative2] THEN\r
 355     SUBGOAL_THEN `!n. derivative (\t. t pow n) = (\t. &n * t pow (n - 1))` ASSUME_TAC THENL\r
 356     [\r
 357       GEN_TAC THEN REWRITE_TAC[FUN_EQ_THM] THEN GEN_TAC THEN\r
 358         new_rewrite [] [] derivative_pow THEN REWRITE_TAC[REAL_DIFFERENTIABLE_ID] THEN\r
 359         REWRITE_TAC[derivative_x; REAL_MUL_RID];\r
 360       ALL_TAC\r
 361     ] THEN\r
 362     ASM_REWRITE_TAC[] THEN\r
 363     SUBGOAL_THEN `!n. derivative (\t. &n * t pow (n - 1)) = (\t. &n * derivative (\t. t pow (n - 1)) t)` ASSUME_TAC THENL\r
 364     [\r
 365       GEN_TAC THEN REWRITE_TAC[FUN_EQ_THM] THEN GEN_TAC THEN\r
 366         new_rewrite [] [] derivative_scale THEN REWRITE_TAC[] THEN\r
 367         MATCH_MP_TAC REAL_DIFFERENTIABLE_POW_ATREAL THEN REWRITE_TAC[REAL_DIFFERENTIABLE_ID];\r
 368       ALL_TAC\r
 369     ] THEN\r
 370     ASM_REWRITE_TAC[] THEN\r
 371     REPEAT (MATCH_MP_TAC REAL_CONTINUOUS_LMUL) THEN\r
 372     MATCH_MP_TAC REAL_CONTINUOUS_POW THEN\r
 373     REWRITE_TAC[REAL_CONTINUOUS_AT_ID]);;\r
 374 \r
 375 \r
 376 let diff2c_domain_pow = prove(`!f n domain. diff2c_domain domain f ==> \r
 377                                diff2c_domain domain (\x. f x pow n)`,\r
 378   REWRITE_TAC[diff2c_domain] THEN REPEAT STRIP_TAC THEN ASM_SIMP_TAC[diff2c_pow]);;\r
 379 \r
 380 \r
 381 let diff2c_domain_tm = `diff2c_domain domain`;;\r
 382 \r
 383 let rec gen_diff2c_domain_poly poly_tm =\r
 384   let x_var, expr = dest_abs poly_tm in\r
 385   let n = (int_of_string o fst o dest_type o hd o tl o snd o dest_type o type_of) x_var in\r
 386   let diff2c_tm = mk_icomb (diff2c_domain_tm, poly_tm) in\r
 387     if frees expr = [] then\r
 388       (* const *)\r
 389       (SPEC_ALL o ISPEC expr o inst_first_type_var (n_type_array.(n))) diff2c_domain_const\r
 390     else\r
 391       let lhs, r_tm = dest_comb expr in\r
 392         if lhs = neg_op_real then\r
 393           (* -- *)\r
 394           let r_th = gen_diff2c_domain_poly (mk_abs (x_var, r_tm)) in\r
 395             prove(diff2c_tm, MATCH_MP_TAC diff2c_domain_neg THEN REWRITE_TAC[r_th])\r
 396         else\r
 397           let op, l_tm = dest_comb lhs in\r
 398           let name = (fst o dest_const) op in\r
 399             if name = "$" then\r
 400               (* x$k *)\r
 401               let dim_th = dimindex_array.(n) in\r
 402                 prove(diff2c_tm, MATCH_MP_TAC diff2c_domain_x THEN\r
 403                         REWRITE_TAC[IN_NUMSEG; dim_th] THEN ARITH_TAC)\r
 404             else\r
 405               let l_th = gen_diff2c_domain_poly (mk_abs (x_var, l_tm)) in\r
 406                 if name = "real_pow" then\r
 407                   (* f pow n *)\r
 408                   prove(diff2c_tm, MATCH_MP_TAC diff2c_domain_pow THEN REWRITE_TAC[l_th])\r
 409                 else\r
 410                   let r_th = gen_diff2c_domain_poly (mk_abs (x_var, r_tm)) in\r
 411                     prove(diff2c_tm,\r
 412                           MAP_FIRST apply_tac [diff2c_domain_add; diff2c_domain_sub; diff2c_domain_mul] THEN\r
 413                             REWRITE_TAC[l_th; r_th]);;\r
 414 \r
 415 \r
 416 let gen_diff2c_poly =\r
 417   let th_imp = prove(`!f. (!domain. diff2c_domain domain f) ==> !x:real^N. diff2c f x`,\r
 418                      REWRITE_TAC[diff2c_domain] THEN REPEAT STRIP_TAC THEN FIRST_X_ASSUM MATCH_MP_TAC THEN\r
 419                        EXISTS_TAC `(x:real^N, x:real^N)` THEN\r
 420                        REWRITE_TAC[INTERVAL_SING; IN_SING]) in\r
 421     fun poly_tm ->\r
 422       (MATCH_MP th_imp o GEN_ALL o gen_diff2c_domain_poly) poly_tm;;\r
 423 \r
 424 \r
 425 let gen_diff_poly =\r
 426   let th_imp = prove(`!f. (!domain. diff2c_domain domain f) ==> !x:real^N. lift o f differentiable at x`,\r
 427                      REWRITE_TAC[diff2c_domain; diff2c; diff2] THEN REPEAT STRIP_TAC THEN\r
 428                        FIRST_X_ASSUM (MP_TAC o SPECL [`x:real^N, x:real^N`; `x:real^N`]) THEN\r
 429                        REWRITE_TAC[INTERVAL_SING; IN_SING] THEN case THEN REPEAT STRIP_TAC THEN\r
 430                        FIRST_X_ASSUM (MP_TAC o SPEC `x:real^N`) THEN ASM_SIMP_TAC[]) in\r
 431     fun poly_tm ->\r
 432       (MATCH_MP th_imp o GEN_ALL o gen_diff2c_domain_poly) poly_tm;;\r
 433 \r
 434 \r
 435 let in_tm = `IN`;;\r
 436 \r
 437 let add_to_hash tbl max_size key value =\r
 438   let _ = if Hashtbl.length tbl >= max_size then Hashtbl.clear tbl else () in\r
 439     Hashtbl.add tbl key value;;\r
 440 \r
 441 (* Formally computes partial derivatives of a polynomial *)\r
 442 let gen_partial_poly =\r
 443   let max_hash = 1000 in\r
 444   let hash = Hashtbl.create max_hash in\r
 445     fun i poly_tm ->\r
 446       let key = (i, poly_tm) in \r
 447         try Hashtbl.find hash (i, poly_tm)\r
 448         with Not_found ->\r
 449           let i_tm = mk_small_numeral i in\r
 450           let rec gen_rec poly_tm =\r
 451             let x_var, expr = dest_abs poly_tm in\r
 452             let n = (int_of_string o fst o dest_type o hd o tl o snd o dest_type o type_of) x_var in\r
 453               if frees expr = [] then\r
 454                 (* const *)\r
 455                 (SPECL [i_tm; expr] o inst_first_type_var (n_type_array.(n))) partial_const\r
 456               else\r
 457                 let lhs, r_tm = dest_comb expr in\r
 458                   if lhs = neg_op_real then\r
 459                     (* -- *)\r
 460                     let r_poly = mk_abs (x_var, r_tm) in\r
 461                     let r_diff = (SPEC_ALL o gen_diff_poly) r_poly and\r
 462                         r_partial = gen_rec r_poly in\r
 463                     let th0 = SPEC i_tm (MATCH_MP partial_neg r_diff) in \r
 464                       REWRITE_RULE[r_partial] th0\r
 465                   else\r
 466                     let op, l_tm = dest_comb lhs in\r
 467                     let name = (fst o dest_const) op in\r
 468                       if name = "$" then\r
 469                         (* comp *)\r
 470                         let dim_th = dimindex_array.(n) in\r
 471                         let dim_tm = (lhand o concl) dim_th in\r
 472                         let i_eq_k = NUM_EQ_CONV (mk_eq (i_tm, r_tm)) in\r
 473                         let int_tm = mk_binop `..` `1` dim_tm in\r
 474                         let k_in_dim = prove(mk_comb (mk_icomb(in_tm, r_tm), int_tm),\r
 475                                              REWRITE_TAC[IN_NUMSEG; dim_th] THEN ARITH_TAC) in\r
 476                           (REWRITE_RULE[i_eq_k] o MATCH_MP (SPECL [r_tm; i_tm] partial_x)) k_in_dim\r
 477                       else\r
 478                         let l_poly = mk_abs (x_var, l_tm) in\r
 479                         let l_partial = gen_rec l_poly in\r
 480                         let l_diff = (SPEC_ALL o gen_diff_poly) l_poly in\r
 481                           if name = "real_pow" then\r
 482                             (* f pow n *)\r
 483                             let th0 = SPECL [i_tm; r_tm] (MATCH_MP partial_pow l_diff) in\r
 484                               REWRITE_RULE[l_partial] th0\r
 485                           else\r
 486                             let r_poly = mk_abs (x_var, r_tm) in\r
 487                             let r_partial = gen_rec r_poly in\r
 488                             let r_diff = (SPEC_ALL o gen_diff_poly) r_poly in\r
 489                             let imp_th = assoc op [add_op_real, partial_add; \r
 490                                                    sub_op_real, partial_sub;\r
 491                                                    mul_op_real, partial_mul] in\r
 492                             let th0 = SPEC i_tm (MATCH_MP (MATCH_MP imp_th l_diff) r_diff) in\r
 493                               REWRITE_RULE[l_partial; r_partial] th0 in\r
 494             \r
 495           let th1 = gen_rec poly_tm in\r
 496           let th2 = ((NUM_REDUCE_CONV THENC REWRITE_CONV[DECIMAL] THENC REAL_POLY_CONV) o rand o concl) th1 in\r
 497           let th3 = (REWRITE_RULE[ETA_AX] o ONCE_REWRITE_RULE[eq_ext] o GEN_ALL) (TRANS th1 th2) in\r
 498           let _ = add_to_hash hash max_hash key th3 in\r
 499             th3;;\r
 500 \r
 501 \r
 502 let gen_partial2_poly i j poly_tm =\r
 503   let partial_j = gen_partial_poly j poly_tm in\r
 504   let partial_ij = gen_partial_poly i (rand (concl partial_j)) in\r
 505   let pi = (rator o lhand o concl) partial_ij in\r
 506     REWRITE_RULE[GSYM partial2] (TRANS (AP_TERM pi partial_j) partial_ij);;\r
 507 \r
 508 \r
 509 (********************************************)\r
 510 \r
 511 let eval_diff2_poly diff2_domain_th =\r
 512   fun xx zz ->\r
 513     let domain_tm = mk_pair (xx, zz) in\r
 514       INST[domain_tm, mk_var ("domain", type_of domain_tm)] diff2_domain_th;;\r
 515 \r
 516 \r
 517 (*****************************)\r
 518 (* m_lin_approx *)\r
 519 \r
 520 let CONST_INTERVAL' = RULE CONST_INTERVAL;;\r
 521 \r
 522 let dest_lin_approx approx_tm =\r
 523   let lhs, df_bounds = dest_comb approx_tm in\r
 524   let lhs2, f_bounds = dest_comb lhs in\r
 525   let lhs3, x_tm = dest_comb lhs2 in\r
 526   let f_tm = rand lhs3 in\r
 527     f_tm, x_tm, f_bounds, df_bounds;;\r
 528 \r
 529 let gen_lin_approx_eq_thm n =\r
 530   let ty = n_type_array.(n) in\r
 531   let df_vars = map (fun i -> df_vars_array.(i)) (1--n) in\r
 532   let df_bounds_list = mk_list (df_vars, real_pair_ty) in\r
 533   let th0 = (SPECL[f_bounds_var; df_bounds_list] o inst_first_type_var ty) m_lin_approx in\r
 534   let th1 = (CONV_RULE NUM_REDUCE_CONV o REWRITE_RULE[all_n]) th0 in\r
 535     th1;;\r
 536 \r
 537 \r
 538 let gen_lin_approx_poly_thm poly_tm diff_th partials =\r
 539   let x_var, _ = dest_abs poly_tm in\r
 540   let n = get_dim x_var in\r
 541   let lin_eq = (REWRITE_RULE partials o SPECL [poly_tm]) (gen_lin_approx_eq_thm n) in\r
 542   let x_vec = mk_vector_list (map (fun i -> x_vars_array.(i)) (1--n)) in\r
 543   let th1 = (REWRITE_RULE (Array.to_list comp_thms_array.(n)) o SPEC x_vec o REWRITE_RULE[diff_th]) lin_eq in\r
 544     th1;;\r
 545 \r
 546 \r
 547 let gen_lin_approx_poly_thm0 poly_tm =\r
 548   let x_var, _ = dest_abs poly_tm in\r
 549   let n = get_dim x_var in\r
 550   let partials = map (fun i -> gen_partial_poly i poly_tm) (1--n) in \r
 551   let diff_th = gen_diff_poly poly_tm in\r
 552     gen_lin_approx_poly_thm poly_tm diff_th partials;;\r
 553 \r
 554 \r
 555 let eval_lin_approx pp0 lin_approx_th =\r
 556   let poly_tm, _, _, _ = (dest_lin_approx o lhand o concl) lin_approx_th in\r
 557   let x_var, _ = dest_abs poly_tm in\r
 558   let n = get_dim x_var in\r
 559   let th0 = lin_approx_th in\r
 560   let th1 = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o MATCH_MP iffRL) th0 in\r
 561   let build_eval int_hyp =\r
 562     let expr, b_var = dest_binary "interval_arith" int_hyp in\r
 563       (eval_constants pp0 o build_interval_fun) expr, b_var in\r
 564   let int_fs = map build_eval (hyp th1) in\r
 565     \r
 566   let rec split_rules i_list =\r
 567     match i_list with\r
 568       | [] -> ([], [])\r
 569       | ((i_fun, var_tm) :: es) -> \r
 570           let th_list, i_list' = split_rules es in\r
 571             match i_fun with\r
 572               | Int_const th -> (var_tm, th) :: th_list, i_list'\r
 573               | Int_var v -> (var_tm, INST[v, x_var_real] CONST_INTERVAL') :: th_list, i_list'\r
 574               | _ -> th_list, (var_tm, i_fun) :: i_list' in\r
 575 \r
 576   let const_th_list, i_list0 = split_rules int_fs in\r
 577   let th2 = itlist (fun (var_tm, th) th0 -> \r
 578                       let b_tm = rand (concl th) in\r
 579                         (MY_PROVE_HYP th o INST[b_tm, var_tm]) th0) const_th_list th1 in\r
 580   let v_list, i_list' = unzip i_list0 in\r
 581   let i_list = find_and_replace_all i_list' [] in\r
 582     fun pp vector_tm ->\r
 583       let x_vals = dest_vector vector_tm in\r
 584       if length x_vals <> n then failwith (sprintf "Wrong vector size; expected size: %d" n)\r
 585       else\r
 586         let x_ints = map mk_const_interval x_vals in\r
 587         let vars = map (fun i -> x_vars_array.(i)) (1--n) in\r
 588         let th3 = INST (zip x_vals vars) th2 in\r
 589         let i_vals = eval_interval_fun_list pp i_list (zip vars x_ints) in\r
 590           itlist2 (fun var_tm th th0 ->\r
 591                      let b_tm = rand (concl th) in\r
 592                        (MY_PROVE_HYP th o INST[b_tm, var_tm]) th0) v_list i_vals th3;;\r
 593 \r
 594 \r
 595 let eval_lin_approx_poly0 pp0 poly_tm =\r
 596   eval_lin_approx pp0 (gen_lin_approx_poly_thm0 poly_tm);;\r
 597 \r
 598 \r
 599 \r
 600 (*************************************)\r
 601 \r
 602 (* 1 <= i /\ i <= n <=> i = 1 \/ i = 2 \/ ... \/ i = n *)\r
 603 let i_int_array =\r
 604   let i_tm = `i:num` in\r
 605   let i_th0 = prove(`1 <= i /\ i <= SUC n <=> (1 <= i /\ i <= n) \/ i = SUC n`, ARITH_TAC) in\r
 606   let th1 = prove(`1 <= i /\ i <= 1 <=> i = 1`, ARITH_TAC) in\r
 607   let array = Array.create (max_dim + 1) th1 in\r
 608   let prove_next n =\r
 609     let n_tm = mk_small_numeral n in\r
 610     let prev_n = num_CONV n_tm in\r
 611     let tm = mk_conj (`1 <= i`, mk_binop le_op_num i_tm n_tm) in\r
 612     let th = REWRITE_CONV[prev_n; i_th0; array.(n - 1)] tm in\r
 613       array.(n) <- REWRITE_RULE[SYM prev_n; GSYM DISJ_ASSOC] th in\r
 614   let _ = map prove_next (2--max_dim) in\r
 615     array;;\r
 616 \r
 617 \r
 618 (* (!i. 1 <= i /\ i <= n ==> P i) <=> P 1 /\ P 2 /\ ... /\ P n *)\r
 619 let gen_in_interval =\r
 620   let th0 = prove(`(!i:num. (i = k \/ Q i) ==> P i) <=> (P k /\ (!i. Q i ==> P i))`, MESON_TAC[]) in\r
 621   let th1 = prove(`(!i:num. (i = k ==> P i)) <=> P k`, MESON_TAC[]) in\r
 622     fun n ->\r
 623       let n_tm = mk_small_numeral n and\r
 624           i_tm = `i:num` in\r
 625       let lhs1 = mk_conj (`1 <= i`, mk_binop le_op_num i_tm n_tm) in\r
 626       let lhs = mk_forall (i_tm, mk_imp (lhs1, `(P:num->bool) i`)) in\r
 627         REWRITE_CONV[i_int_array.(n); th0; th1] lhs;;\r
 628 \r
 629 \r
 630 let gen_second_bounded_eq_thm n =\r
 631   let ty, _, x_var, _, _, _, domain_var = get_types_and_vars n in\r
 632   let dd_vars = map (fun i -> map (fun j -> dd_vars_array.(j).(i)) (1--i)) (1--n) in\r
 633   let dd_bounds_list = mk_list (map (fun l -> mk_list (l, real_pair_ty)) dd_vars, real_pair_list_ty) in\r
 634   let th0 = (SPECL[domain_var; dd_bounds_list] o inst_first_type_var ty) second_bounded in\r
 635   let th1 = (CONV_RULE NUM_REDUCE_CONV o REWRITE_RULE[all_n]) th0 in\r
 636     th1;;\r
 637 \r
 638 \r
 639 let gen_second_bounded_poly_thm poly_tm partials2 =\r
 640   let x_var, _ = dest_abs poly_tm in\r
 641   let n = get_dim x_var in\r
 642   let partials2' = List.flatten partials2 in\r
 643   let second_th = (REWRITE_RULE partials2' o SPECL [poly_tm]) (gen_second_bounded_eq_thm n) in\r
 644     second_th;;\r
 645 \r
 646 \r
 647 let gen_second_bounded_poly_thm0 poly_tm =\r
 648   let x_var, _ = dest_abs poly_tm in\r
 649   let n = get_dim x_var in\r
 650   let partials = map (fun i -> gen_partial_poly i poly_tm) (1--n) in\r
 651   let get_partial i eq_th = \r
 652     let partial_i = gen_partial_poly i (rand (concl eq_th)) in\r
 653     let pi = (rator o lhand o concl) partial_i in\r
 654       REWRITE_RULE[GSYM partial2] (TRANS (AP_TERM pi eq_th) partial_i) in\r
 655   let partials2 = map (fun th, i -> map (fun j -> get_partial j th) (1--i)) (zip partials (1--n)) in\r
 656     gen_second_bounded_poly_thm poly_tm partials2;;\r
 657 \r
 658 (*  let eq_th = TAUT `(P ==> Q /\ R) <=> ((P ==> Q) /\ (P ==> R))` in\r
 659     REWRITE_RULE[eq_th; FORALL_AND_THM; GSYM m_bounded_on_int] second_th;;*)\r
 660 \r
 661 \r
 662 (* eval_second_bounded *)\r
 663 let eval_second_bounded pp0 second_bounded_th =\r
 664   let poly_tm = (lhand o rator o lhand o concl) second_bounded_th in\r
 665   let th0 = second_bounded_th in\r
 666   let n = (get_dim o fst o dest_abs) poly_tm in\r
 667   let x_vector = mk_vector_list (map (fun i -> x_vars_array.(i)) (1--n)) and\r
 668       z_vector = mk_vector_list (map (fun i -> z_vars_array.(i)) (1--n)) in\r
 669   let _, _, _, _, _, _, domain_var = get_types_and_vars n in\r
 670 \r
 671   let th1 = INST[mk_pair (x_vector, z_vector), domain_var] th0 in\r
 672   let th2 = REWRITE_RULE[IN_INTERVAL; dimindex_array.(n)] th1 in\r
 673   let th3 = REWRITE_RULE[gen_in_interval n; GSYM interval_arith] th2 in\r
 674   let th4 = (REWRITE_RULE[CONJ_ACI] o REWRITE_RULE (Array.to_list comp_thms_array.(n))) th3 in\r
 675   let final_th0 = (UNDISCH_ALL o MATCH_MP iffRL) th4 in\r
 676 \r
 677   let x_var, h_tm = (dest_forall o hd o hyp) final_th0 in\r
 678   let _, h2 = dest_imp h_tm in\r
 679   let concl_ints = striplist dest_conj h2 in\r
 680 \r
 681   let i_funs = map (fun int -> \r
 682                       let expr, var = dest_interval_arith int in\r
 683                         (eval_constants pp0 o build_interval_fun) expr, var) concl_ints in\r
 684 \r
 685   let rec split_rules i_list =\r
 686     match i_list with\r
 687       | [] -> ([], [])\r
 688       | ((i_fun, var_tm) :: es) -> \r
 689           let th_list, i_list' = split_rules es in\r
 690             match i_fun with\r
 691               | Int_const th -> (var_tm, th) :: th_list, i_list'\r
 692 (*            | Int_var v -> (var_tm, INST[v, x_var_real] CONST_INTERVAL') :: th_list, i_list' *)\r
 693               | _ -> th_list, (var_tm, i_fun) :: i_list' in\r
 694 \r
 695   let const_th_list, i_list0 = split_rules i_funs in\r
 696   let th5 = itlist (fun (var_tm, th) th0 -> \r
 697                       let b_tm = rand (concl th) in\r
 698                         (REWRITE_RULE[th] o INST[b_tm, var_tm]) th0) const_th_list (SYM th4) in\r
 699   let final_th = REWRITE_RULE[GSYM IMP_IMP] th5 in\r
 700   let v_list, i_list' = unzip i_list0 in\r
 701   let i_list = find_and_replace_all i_list' [] in\r
 702 \r
 703     fun pp x_vector_tm z_vector_tm ->\r
 704       let x_vals = dest_vector x_vector_tm and\r
 705           z_vals = dest_vector z_vector_tm in\r
 706         if length x_vals <> n or length z_vals <> n then \r
 707           failwith (sprintf "Wrong vector size; expected size: %d" n)\r
 708         else\r
 709           let x_vars = map (fun i -> x_vars_array.(i)) (1--n) and\r
 710               z_vars = map (fun i -> z_vars_array.(i)) (1--n) in\r
 711 \r
 712           let inst_th = (INST (zip x_vals x_vars) o INST (zip z_vals z_vars)) final_th in\r
 713             if (not o is_eq) (concl inst_th) then inst_th \r
 714             else\r
 715               let x_var, lhs = (dest_forall o lhand o concl) inst_th in\r
 716               let hs = (butlast o striplist dest_imp) lhs in\r
 717               let vars = map (rand o rator) hs in\r
 718               let int_vars = zip vars (map ASSUME hs) in\r
 719 \r
 720               let dd_ints = eval_interval_fun_list pp i_list int_vars in\r
 721               let inst_dd = map2 (fun var th -> (rand o concl) th, var) v_list dd_ints in\r
 722               let inst_th2 = INST inst_dd inst_th in\r
 723 \r
 724               let conj_th = end_itlist CONJ dd_ints in\r
 725               let lhs_th = GEN x_var (itlist DISCH hs conj_th) in\r
 726                 EQ_MP inst_th2 lhs_th;;\r
 727 \r
 728 \r
 729 \r
 730 let eval_second_bounded_poly0 pp0 poly_tm =\r
 731   eval_second_bounded pp0 (gen_second_bounded_poly_thm0 poly_tm);;\r
 732 \r
 733 \r
 734 \r
 735 (*************************************)\r
 736 (* eval_m_taylor *)\r
 737 \r
 738 let eval_m_taylor pp0 diff2c_th lin_th second_th =\r
 739   let poly_tm = (rand o concl) diff2c_th in\r
 740   let n = (get_dim o fst o dest_abs) poly_tm in\r
 741   let eval_lin = eval_lin_approx pp0 lin_th and\r
 742       eval_second = eval_second_bounded pp0 second_th in\r
 743 \r
 744   let ty, _, x_var, f_var, y_var, w_var, domain_var = get_types_and_vars n in\r
 745   let th0 = (SPEC_ALL o inst_first_type_var ty) m_taylor_interval in\r
 746   let th1 = INST[poly_tm, f_var] th0 in\r
 747   let th2 = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o MATCH_MP iffRL o REWRITE_RULE[diff2c_th]) th1 in\r
 748 \r
 749     fun p_lin p_second domain_th ->\r
 750       let domain_tm, y_tm, w_tm = dest_m_cell_domain (concl domain_th) in\r
 751       let x_tm, z_tm = dest_pair domain_tm in\r
 752 \r
 753       let lin_th = eval_lin p_lin y_tm and\r
 754           second_th = eval_second p_second x_tm z_tm in\r
 755 \r
 756       let _, _, f_bounds, df_bounds_list = dest_lin_approx (concl lin_th) in\r
 757       let dd_bounds_list = (rand o concl) second_th in\r
 758       let df_var = mk_var ("d_bounds_list", type_of df_bounds_list) and\r
 759           dd_var = mk_var ("dd_bounds_list", type_of dd_bounds_list) in\r
 760 \r
 761         (MY_PROVE_HYP domain_th o MY_PROVE_HYP lin_th o MY_PROVE_HYP second_th o\r
 762            INST[domain_tm, domain_var; y_tm, y_var; w_tm, w_var;\r
 763                 f_bounds, f_bounds_var; df_bounds_list, df_var; dd_bounds_list, dd_var]) th2;;\r
 764 \r
 765 \r
 766 let eval_m_taylor_poly0 pp0 poly_tm =\r
 767   let diff2_th = gen_diff2c_domain_poly poly_tm in\r
 768   let x_var, _ = dest_abs poly_tm in\r
 769   let n = get_dim x_var in\r
 770   let partials = map (fun i -> gen_partial_poly i poly_tm) (1--n) in\r
 771   let get_partial i eq_th = \r
 772     let partial_i = gen_partial_poly i (rand (concl eq_th)) in\r
 773     let pi = (rator o lhand o concl) partial_i in\r
 774       REWRITE_RULE[GSYM partial2] (TRANS (AP_TERM pi eq_th) partial_i) in\r
 775   let partials2 = map2 (fun th i -> map (fun j -> get_partial j th) (1--i)) partials (1--n) in\r
 776   let second_th = gen_second_bounded_poly_thm poly_tm partials2 in\r
 777   let diff_th = gen_diff_poly poly_tm in\r
 778   let lin_th = gen_lin_approx_poly_thm poly_tm diff_th partials in\r
 779     eval_m_taylor pp0 diff2_th lin_th second_th;;\r
 780 \r
 781 \r
 782 \r
 783 (******************************************)\r
 784 (* mk_eval_function *)\r
 785 \r
 786 let mk_eval_function_eq pp0 eq_th =\r
 787   let expr_tm = (rand o concl) eq_th in\r
 788   let tm0 = `!x:real^N. x IN interval [domain] ==> interval_arith (f x) f_bounds` in\r
 789   let n = (get_dim o fst o dest_abs) expr_tm in\r
 790   let x_vector = mk_vector_list (map (fun i -> x_vars_array.(i)) (1--n)) and\r
 791       z_vector = mk_vector_list (map (fun i -> z_vars_array.(i)) (1--n)) in\r
 792   let ty, _, _, _, _, _, domain_var = get_types_and_vars n and\r
 793       f_var = mk_var ("f", type_of expr_tm) in\r
 794   let th1 = (REWRITE_CONV[IN_INTERVAL] o subst[mk_pair(x_vector,z_vector), domain_var] o inst[ty, nty]) tm0 in\r
 795   let th2 = REWRITE_RULE [dimindex_array.(n)] th1 in\r
 796   let th3 = REWRITE_RULE [gen_in_interval n; GSYM interval_arith] th2 in\r
 797   let th4 = (REWRITE_RULE[GSYM IMP_IMP; CONJ_ACI] o REWRITE_RULE (Array.to_list comp_thms_array.(n))) th3 in\r
 798   let final_th0 = (CONV_RULE ((RAND_CONV o ONCE_DEPTH_CONV) BETA_CONV) o INST[expr_tm, f_var]) th4 in\r
 799   let x_var, h_tm = (dest_forall o rand o concl) final_th0 in\r
 800   let f_tm = (fst o dest_interval_arith o last o striplist dest_imp) h_tm in\r
 801   let i_fun = (eval_constants pp0 o build_interval_fun) f_tm in\r
 802   let i_list = find_and_replace_all [i_fun] [] in\r
 803   let final_th = (PURE_REWRITE_RULE[SYM eq_th] o SYM) final_th0 in\r
 804     fun pp x_vector_tm z_vector_tm ->\r
 805       let x_vals = dest_vector x_vector_tm and\r
 806           z_vals = dest_vector z_vector_tm in\r
 807         if length x_vals <> n or length z_vals <> n then \r
 808           failwith (sprintf "Wrong vector size; expected size: %d" n)\r
 809         else\r
 810           let x_vars = map (fun i -> x_vars_array.(i)) (1--n) and\r
 811               z_vars = map (fun i -> z_vars_array.(i)) (1--n) in\r
 812 \r
 813           let inst_th = (INST (zip x_vals x_vars) o INST (zip z_vals z_vars)) final_th in\r
 814           let x_var, lhs = (dest_forall o lhand o concl) inst_th in\r
 815           let hs = (butlast o striplist dest_imp) lhs in\r
 816           let vars = map (rand o rator) hs in\r
 817           let int_vars = zip vars (map ASSUME hs) in\r
 818           let eval_th = hd (eval_interval_fun_list pp i_list int_vars) in\r
 819           let f_bounds = (rand o concl) eval_th in\r
 820           let inst_th2 = INST[f_bounds, f_bounds_var] inst_th in\r
 821           let lhs_th = GEN x_var (itlist DISCH hs eval_th) in\r
 822             EQ_MP inst_th2 lhs_th;;\r
 823 \r
 824 \r
 825 let mk_eval_function pp0 expr_tm = mk_eval_function_eq pp0 (REFL expr_tm);;\r
 826 \r
 827 \r
 828 \r
 829 (********************************)\r
 830 (* m_taylor_error *)\r
 831 \r
 832 (* Sum of the list elements *)\r
 833 let ITLIST2_EQ_SUM = prove(`!(f:A->B->real) l1 l2. LENGTH l1 <= LENGTH l2 ==>\r
 834                                ITLIST2 (\x y z. f x y + z) l1 l2 (&0) =\r
 835                                sum (1..(LENGTH l1)) (\i. f (EL (i - 1) l1) (EL (i - 1) l2))`,\r
 836    GEN_TAC THEN\r
 837      LIST_INDUCT_TAC THEN LIST_INDUCT_TAC THEN REWRITE_TAC[LENGTH; ITLIST2_DEF] THEN TRY ARITH_TAC THENL\r
 838      [\r
 839        REWRITE_TAC[SUM_CLAUSES_NUMSEG; ARITH];\r
 840        REWRITE_TAC[SUM_CLAUSES_NUMSEG; ARITH];\r
 841        ALL_TAC\r
 842      ] THEN\r
 843      REWRITE_TAC[leqSS] THEN DISCH_TAC THEN\r
 844      FIRST_X_ASSUM (MP_TAC o SPEC `t':(B)list`) THEN ASM_REWRITE_TAC[] THEN\r
 845      DISCH_THEN (fun th -> REWRITE_TAC[TL; th]) THEN\r
 846      REWRITE_TAC[GSYM add1n] THEN\r
 847      new_rewrite [] [] SUM_ADD_SPLIT THEN REWRITE_TAC[ARITH] THEN\r
 848      REWRITE_TAC[TWO; add1n; SUM_SING_NUMSEG; subnn; EL; HD] THEN\r
 849      REWRITE_TAC[GSYM addn1; SUM_OFFSET; REAL_EQ_ADD_LCANCEL] THEN\r
 850      MATCH_MP_TAC SUM_EQ THEN move ["i"] THEN REWRITE_TAC[IN_NUMSEG] THEN DISCH_TAC THEN\r
 851      ASM_SIMP_TAC[ARITH_RULE `1 <= i ==> (i + 1) - 1 = SUC (i - 1)`; EL; TL]);;\r
 852 \r
 853 \r
 854 let interval_arith_abs_le = prove(`!x int y. interval_arith x int ==> iabs int <= y ==> abs x <= y`,\r
 855    GEN_TAC THEN case THEN REWRITE_TAC[interval_arith; IABS'] THEN REAL_ARITH_TAC);;\r
 856 \r
 857 \r
 858 let ALL_N_ALL2 = prove(`!P (l:(A)list) i0.\r
 859                          (all_n i0 l P <=> if l = [] then T else ALL2 P (l_seq i0 ((i0 + LENGTH l) - 1)) l)`,\r
 860   GEN_TAC THEN LIST_INDUCT_TAC THEN GEN_TAC THEN REWRITE_TAC[all_n; NOT_CONS_NIL] THEN\r
 861     new_rewrite [] [] L_SEQ_CONS THEN REWRITE_TAC[LENGTH; ALL2] THEN TRY ARITH_TAC THEN\r
 862     FIRST_X_ASSUM (new_rewrite [] []) THEN TRY ARITH_TAC THEN\r
 863     REWRITE_TAC[addSn; addnS; addn1] THEN\r
 864     SPEC_TAC (`t:(A)list`, `t:(A)list`) THEN case THEN SIMP_TAC[NOT_CONS_NIL] THEN\r
 865     REWRITE_TAC[LENGTH; addn0] THEN\r
 866     MP_TAC (SPECL [`SUC i0`; `SUC i0 - 1`] L_SEQ_NIL) THEN\r
 867     REWRITE_TAC[ARITH_RULE `SUC i0 - 1 < SUC i0`] THEN DISCH_THEN (fun th -> REWRITE_TAC[th; ALL2]));;\r
 868     \r
 869 \r
 870 let ALL_N_EL = prove(`!P (l:(A)list) i0. all_n i0 l P <=> (!i. i < LENGTH l ==> P (i0 + i) (EL i l))`,\r
 871    REPEAT GEN_TAC THEN REWRITE_TAC[ALL_N_ALL2] THEN\r
 872      SPEC_TAC (`l:(A)list`, `l:(A)list`) THEN case THEN SIMP_TAC[NOT_CONS_NIL; LENGTH; ltn0] THEN\r
 873      REPEAT GEN_TAC THEN\r
 874      new_rewrite [] [] ALL2_ALL_ZIP THENL\r
 875      [\r
 876        REWRITE_TAC[LENGTH_L_SEQ; LENGTH] THEN ARITH_TAC;\r
 877        ALL_TAC\r
 878      ] THEN\r
 879      REWRITE_TAC[GSYM ALL_EL] THEN\r
 880      new_rewrite [] [] LENGTH_ZIP THENL\r
 881      [\r
 882        REWRITE_TAC[LENGTH_L_SEQ; LENGTH] THEN ARITH_TAC;\r
 883        ALL_TAC\r
 884      ] THEN\r
 885      REWRITE_TAC[LENGTH_L_SEQ; ARITH_RULE `((i0 + SUC a) - 1 + 1) - i0 = SUC a`] THEN\r
 886      EQ_TAC THEN REPEAT STRIP_TAC THENL\r
 887      [\r
 888        FIRST_X_ASSUM (MP_TAC o SPEC `i:num`) THEN ASM_REWRITE_TAC[] THEN\r
 889          new_rewrite [] [] EL_ZIP THENL\r
 890          [\r
 891            REWRITE_TAC[LENGTH_L_SEQ; LENGTH] THEN ASM_ARITH_TAC;\r
 892            ALL_TAC\r
 893          ] THEN\r
 894          REWRITE_TAC[] THEN new_rewrite [] [] EL_L_SEQ THEN ASM_ARITH_TAC;\r
 895        ALL_TAC\r
 896      ] THEN\r
 897      new_rewrite [] [] EL_ZIP THENL\r
 898      [\r
 899        REWRITE_TAC[LENGTH_L_SEQ; LENGTH] THEN ASM_ARITH_TAC;\r
 900        ALL_TAC\r
 901      ] THEN\r
 902      REWRITE_TAC[] THEN new_rewrite [] [] EL_L_SEQ THEN TRY ASM_ARITH_TAC THEN\r
 903      FIRST_X_ASSUM MATCH_MP_TAC THEN ASM_REWRITE_TAC[]);;\r
 904      \r
 905 \r
 906 let M_TAYLOR_ERROR_ITLIST2 = prove(`!f domain y w dd_bounds_list error. \r
 907              m_cell_domain domain y (vector w) ==>\r
 908              diff2c_domain domain f ==>\r
 909              second_bounded (f:real^N->real) domain dd_bounds_list ==>\r
 910              LENGTH w = dimindex (:N) ==>\r
 911                  LENGTH dd_bounds_list = dimindex (:N) ==>\r
 912                  all_n 1 dd_bounds_list (\i list. LENGTH list = i) ==>\r
 913                  ITLIST2 (\list x z. x * (x * iabs (LAST list)\r
 914                                           + &2 * ITLIST2 (\a b c. b * iabs a + c) (BUTLAST list) w (&0)) + z)\r
 915                  dd_bounds_list w (&0) <= error ==>\r
 916                  m_taylor_error f domain (vector w) error`,\r
 917    REPEAT GEN_TAC THEN REWRITE_TAC[second_bounded] THEN\r
 918      set_tac "s" `ITLIST2 _1 _2 _3 _4` THEN\r
 919      move ["domain"; "d2f"; "second"; "lw"; "ldd"; "ldd_all"; "s_le"] THEN\r
 920      ASM_SIMP_TAC[m_taylor_error_eq] THEN move ["x"; "x_in"] THEN\r
 921      SUBGOAL_THEN `!i. i IN 1..dimindex (:N) ==> &0 <= EL (i - 1) w` (LABEL_TAC "w_ge0") THENL\r
 922      [\r
 923        GEN_TAC THEN DISCH_TAC THEN REMOVE_THEN "domain" MP_TAC THEN new_rewrite [] [] pair_eq THEN\r
 924          REWRITE_TAC[m_cell_domain] THEN\r
 925          DISCH_THEN (MP_TAC o SPEC `i:num`) THEN ASM_REWRITE_TAC[] THEN\r
 926          ASM_SIMP_TAC[VECTOR_COMPONENT] THEN REAL_ARITH_TAC;\r
 927        ALL_TAC\r
 928      ] THEN\r
 929      MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `s:real` THEN ASM_REWRITE_TAC[] THEN\r
 930      EXPAND_TAC "s" THEN\r
 931      new_rewrite [] [] ITLIST2_EQ_SUM THEN ASM_REWRITE_TAC[LE_REFL] THEN\r
 932      MATCH_MP_TAC SUM_LE THEN REWRITE_TAC[FINITE_NUMSEG] THEN\r
 933      move ["i"; "i_in"] THEN ASM_SIMP_TAC[VECTOR_COMPONENT] THEN\r
 934      USE_THEN "i_in" (ASSUME_TAC o REWRITE_RULE[IN_NUMSEG]) THEN\r
 935      MATCH_MP_TAC REAL_LE_LMUL THEN ASM_SIMP_TAC[] THEN\r
 936      SUBGOAL_THEN `LENGTH (EL (i - 1) dd_bounds_list:(real#real)list) = i` (LABEL_TAC "len_i") THENL\r
 937      [\r
 938        REMOVE_THEN "ldd_all" MP_TAC THEN\r
 939          REWRITE_TAC[ALL_N_EL] THEN DISCH_THEN (MP_TAC o SPEC `i - 1`) THEN\r
 940          ANTS_TAC THENL\r
 941          [\r
 942            ASM_REWRITE_TAC[] THEN POP_ASSUM MP_TAC THEN ARITH_TAC;\r
 943            ALL_TAC\r
 944          ] THEN\r
 945          DISCH_THEN (fun th -> REWRITE_TAC[th]) THEN POP_ASSUM MP_TAC THEN ARITH_TAC;\r
 946        ALL_TAC\r
 947      ] THEN\r
 948      new_rewrite [] [] ITLIST2_EQ_SUM THEN ASM_REWRITE_TAC[] THENL\r
 949      [\r
 950        REWRITE_TAC[LENGTH_BUTLAST] THEN POP_ASSUM MP_TAC THEN POP_ASSUM MP_TAC THEN ARITH_TAC;\r
 951        ALL_TAC\r
 952      ] THEN\r
 953      MATCH_MP_TAC REAL_LE_ADD2 THEN CONJ_TAC THENL\r
 954      [\r
 955        new_rewrite [] [] LAST_EL THENL\r
 956          [\r
 957            ASM_REWRITE_TAC[GSYM LENGTH_EQ_NIL] THEN \r
 958              REMOVE_THEN "i_in" MP_TAC THEN REWRITE_TAC[IN_NUMSEG] THEN ARITH_TAC;\r
 959            ALL_TAC\r
 960          ] THEN\r
 961          MATCH_MP_TAC REAL_LE_LMUL THEN ASM_SIMP_TAC[] THEN\r
 962          FIRST_X_ASSUM (MP_TAC o SPEC `x:real^N`) THEN ASM_REWRITE_TAC[ALL_N_EL] THEN\r
 963          DISCH_THEN (MP_TAC o SPEC `i - 1`) THEN\r
 964          ANTS_TAC THENL\r
 965          [\r
 966            REMOVE_THEN "i_in" MP_TAC THEN REWRITE_TAC[IN_NUMSEG] THEN ARITH_TAC;\r
 967            ALL_TAC\r
 968          ] THEN\r
 969          DISCH_THEN (MP_TAC o SPEC `i - 1`) THEN ANTS_TAC THENL\r
 970          [\r
 971            UNDISCH_TAC `1 <= i /\ i <= dimindex (:N)` THEN ASM_REWRITE_TAC[] THEN ARITH_TAC;\r
 972            ALL_TAC\r
 973          ] THEN\r
 974          SUBGOAL_THEN `1 + i - 1 = i /\ 1 + i - 1 = i` (fun th -> REWRITE_TAC[th]) THENL\r
 975          [\r
 976            REPLICATE_TAC 3 (POP_ASSUM MP_TAC) THEN ARITH_TAC;\r
 977            ALL_TAC\r
 978          ] THEN\r
 979          REWRITE_TAC[partial2] THEN\r
 980          DISCH_THEN (MP_TAC o MATCH_MP interval_arith_abs_le) THEN\r
 981          DISCH_THEN MATCH_MP_TAC THEN REWRITE_TAC[REAL_LE_REFL];\r
 982        ALL_TAC\r
 983      ] THEN\r
 984      MATCH_MP_TAC REAL_LE_LMUL THEN CONJ_TAC THENL [ REAL_ARITH_TAC; ALL_TAC ] THEN\r
 985      ASM_REWRITE_TAC[LENGTH_BUTLAST] THEN\r
 986      MATCH_MP_TAC SUM_LE THEN REWRITE_TAC[FINITE_NUMSEG] THEN\r
 987      move ["j"; "j_in"] THEN\r
 988      SUBGOAL_THEN `j IN 1..dimindex (:N)` (LABEL_TAC "j_in2") THENL\r
 989      [\r
 990        REMOVE_THEN "i_in" MP_TAC THEN REMOVE_THEN "j_in" MP_TAC THEN\r
 991          REWRITE_TAC[IN_NUMSEG] THEN ARITH_TAC;\r
 992        ALL_TAC\r
 993      ] THEN\r
 994      ASM_SIMP_TAC[VECTOR_COMPONENT] THEN\r
 995      USE_THEN "j_in" (ASSUME_TAC o REWRITE_RULE[IN_NUMSEG]) THEN\r
 996      MATCH_MP_TAC REAL_LE_LMUL THEN ASM_SIMP_TAC[] THEN\r
 997      FIRST_X_ASSUM (MP_TAC o SPEC `x:real^N`) THEN ASM_REWRITE_TAC[ALL_N_EL] THEN\r
 998      DISCH_THEN (MP_TAC o SPEC `i - 1`) THEN ANTS_TAC THENL\r
 999      [\r
1000        REMOVE_THEN "i_in" MP_TAC THEN REWRITE_TAC[IN_NUMSEG] THEN ARITH_TAC;\r
1001        ALL_TAC\r
1002      ] THEN\r
1003      DISCH_THEN (MP_TAC o SPEC `j - 1`) THEN ANTS_TAC THENL\r
1004      [\r
1005        ASM_REWRITE_TAC[] THEN POP_ASSUM MP_TAC THEN ARITH_TAC;\r
1006        ALL_TAC\r
1007      ] THEN\r
1008      SUBGOAL_THEN `1 + j - 1 = j /\ 1 + i - 1 = i` (fun th -> REWRITE_TAC[th]) THENL\r
1009      [\r
1010        REPLICATE_TAC 3 (POP_ASSUM MP_TAC) THEN ARITH_TAC;\r
1011        ALL_TAC\r
1012      ] THEN\r
1013      REWRITE_TAC[partial2] THEN\r
1014      DISCH_THEN (MP_TAC o MATCH_MP interval_arith_abs_le) THEN\r
1015      DISCH_THEN MATCH_MP_TAC THEN\r
1016      new_rewrite [] [] EL_BUTLAST THENL\r
1017      [\r
1018        ASM_REWRITE_TAC[] THEN REMOVE_THEN "j_in" MP_TAC THEN REMOVE_THEN "i_in" MP_TAC THEN\r
1019          REWRITE_TAC[IN_NUMSEG] THEN ARITH_TAC;\r
1020        ALL_TAC\r
1021      ] THEN\r
1022      REWRITE_TAC[REAL_LE_REFL]);;\r
1023 \r
1024 \r
1025 let M_TAYLOR_ERROR_ITLIST2_ALT = prove(`!f domain y w f_lo f_hi d_bounds_list dd_bounds_list error. \r
1026          m_taylor_interval f domain y (vector w:real^N) (f_lo, f_hi) d_bounds_list dd_bounds_list ==>\r
1027          LENGTH w = dimindex (:N) ==>\r
1028              LENGTH dd_bounds_list = dimindex (:N) ==>\r
1029              all_n 1 dd_bounds_list (\i list. LENGTH list = i) ==>\r
1030              ITLIST2 (\list x z. x * (x * iabs (LAST list)\r
1031                                       + &2 * ITLIST2 (\a b c. b * iabs a + c) (BUTLAST list) w (&0)) + z)\r
1032              dd_bounds_list w (&0) <= error ==>\r
1033              m_taylor_error f domain (vector w) error`,\r
1034    REWRITE_TAC[m_taylor_interval] THEN REPEAT STRIP_TAC THEN\r
1035      MP_TAC (SPEC_ALL M_TAYLOR_ERROR_ITLIST2) THEN ASM_REWRITE_TAC[]);;\r
1036 \r
1037 \r
1038 \r
1039 (****************************)\r
1040 \r
1041 let M_TAYLOR_INTERVAL' = MY_RULE m_taylor_interval;;\r
1042 \r
1043 let dest_m_taylor m_taylor_tm =\r
1044   let ltm1, dd_bounds_list = dest_comb m_taylor_tm in\r
1045   let ltm2, d_bounds_list = dest_comb ltm1 in\r
1046   let ltm3, f_bounds = dest_comb ltm2 in\r
1047   let ltm4, w = dest_comb ltm3 in\r
1048   let ltm5, y = dest_comb ltm4 in\r
1049   let ltm6, domain = dest_comb ltm5 in\r
1050     rand ltm6, domain, y, w, f_bounds, d_bounds_list, dd_bounds_list;;\r
1051 \r
1052 let dest_m_taylor_thms n =\r
1053   let ty, xty, x_var, f_var, y_var, w_var, domain_var = get_types_and_vars n in\r
1054     fun m_taylor_th ->\r
1055       let f, domain, y, w, f_bounds, d_bounds_list, dd_bounds_list = dest_m_taylor (concl m_taylor_th) in\r
1056       let th0 = (INST[f, f_var; domain, domain_var; y, y_var; w, w_var; f_bounds, f_bounds_var;\r
1057                       d_bounds_list, d_bounds_list_var; dd_bounds_list, dd_bounds_list_var] o\r
1058                    inst_first_type_var ty) M_TAYLOR_INTERVAL' in\r
1059       let th1 = EQ_MP th0 m_taylor_th in\r
1060       let [domain_th; d2_th; lin_th; second_th] = CONJUNCTS th1 in\r
1061         domain_th, d2_th, lin_th, second_th;;\r
1062   \r
1063 \r
1064 \r
1065 (**********************)\r
1066 \r
1067 (* bound *)\r
1068 let M_TAYLOR_BOUND' = \r
1069   prove(`m_taylor_interval f domain y (vector w:real^N) (f_lo, f_hi) d_bounds_list dd_bounds_list /\\r
1070           m_taylor_error f domain (vector w) error /\\r
1071           ITLIST2 (\a b c. b * iabs a + c) d_bounds_list w (&0) <= b /\\r
1072           b + inv(&2) * error <= a /\\r
1073           lo <= f_lo - a /\\r
1074           f_hi + a <= hi /\\r
1075           LENGTH w = dimindex (:N) /\ LENGTH d_bounds_list = dimindex (:N) ==>\r
1076               (!x. x IN interval [domain] ==> interval_arith (f x) (lo, hi))`,\r
1077    REWRITE_TAC[GSYM m_bounded_on_int; m_taylor_interval; m_lin_approx; ALL_N_EL] THEN\r
1078      set_tac "s" `ITLIST2 _1 _2 _3 _4` THEN STRIP_TAC THEN\r
1079      SUBGOAL_THEN `diff2_domain domain (f:real^N->real)` ASSUME_TAC THENL\r
1080      [\r
1081        UNDISCH_TAC `diff2c_domain domain (f:real^N->real)` THEN\r
1082          SIMP_TAC[diff2_domain; diff2c_domain; diff2c];\r
1083        ALL_TAC\r
1084      ] THEN\r
1085      apply_tac m_taylor_bounds THEN\r
1086      MAP_EVERY EXISTS_TAC [`y:real^N`; `vector w:real^N`; `error:real`; `f_lo:real`; `f_hi:real`; `a:real`] THEN\r
1087      ASM_REWRITE_TAC[] THEN\r
1088      MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b + inv (&2) * error` THEN ASM_REWRITE_TAC[] THEN\r
1089      REWRITE_TAC[real_div; REAL_MUL_AC] THEN\r
1090      MATCH_MP_TAC REAL_LE_ADD2 THEN REWRITE_TAC[REAL_LE_REFL] THEN\r
1091      MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `s:real` THEN ASM_REWRITE_TAC[] THEN\r
1092      EXPAND_TAC "s" THEN new_rewrite [] [] ITLIST2_EQ_SUM THEN ASM_REWRITE_TAC[LE_REFL] THEN\r
1093      MATCH_MP_TAC SUM_LE THEN REWRITE_TAC[FINITE_NUMSEG] THEN REPEAT STRIP_TAC THEN\r
1094      MATCH_MP_TAC REAL_LE_MUL2 THEN ASM_SIMP_TAC[VECTOR_COMPONENT; REAL_LE_REFL; REAL_ABS_POS] THEN\r
1095      CONJ_TAC THENL\r
1096      [\r
1097        UNDISCH_TAC `m_cell_domain domain (y:real^N) (vector w)` THEN\r
1098          new_rewrite [] [] pair_eq THEN REWRITE_TAC[m_cell_domain] THEN\r
1099          DISCH_THEN (MP_TAC o SPEC `x:num`) THEN ASM_REWRITE_TAC[] THEN\r
1100          ASM_SIMP_TAC[VECTOR_COMPONENT] THEN REAL_ARITH_TAC;\r
1101        ALL_TAC\r
1102      ] THEN\r
1103      FIRST_X_ASSUM (MP_TAC o SPEC `x - 1`) THEN\r
1104      ANTS_TAC THENL\r
1105      [\r
1106        POP_ASSUM MP_TAC THEN ASM_REWRITE_TAC[IN_NUMSEG] THEN ARITH_TAC;\r
1107        ALL_TAC\r
1108      ] THEN\r
1109      DISCH_THEN (MP_TAC o MATCH_MP interval_arith_abs_le) THEN\r
1110      SUBGOAL_THEN `1 + x - 1 = x` (fun th -> REWRITE_TAC[th]) THENL\r
1111      [\r
1112        POP_ASSUM MP_TAC THEN REWRITE_TAC[IN_NUMSEG] THEN ARITH_TAC;\r
1113        ALL_TAC\r
1114      ] THEN\r
1115      DISCH_THEN MATCH_MP_TAC THEN REWRITE_TAC[REAL_LE_REFL]);;\r
1116 \r
1117 \r
1118 (* upper *)\r
1119 let M_TAYLOR_UPPER_BOUND' = prove(`m_taylor_interval f domain y (vector w) (f_lo, f_hi) \r
1120                                     d_bounds_list dd_bounds_list /\\r
1121                                     m_taylor_error f domain (vector w:real^N) error /\\r
1122                                     ITLIST2 (\a b c. b * iabs a + c) d_bounds_list w (&0) <= b /\\r
1123                                     b + inv(&2) * error <= a /\\r
1124                                     f_hi + a <= hi /\\r
1125                                     LENGTH w = dimindex (:N) /\ LENGTH d_bounds_list = dimindex (:N) ==>\r
1126                                         (!p. p IN interval [domain] ==> f p <= hi)`,\r
1127    STRIP_TAC THEN\r
1128      MP_TAC (INST[`f_lo - a:real`, `lo:real`] M_TAYLOR_BOUND') THEN\r
1129      ASM_SIMP_TAC[interval_arith; REAL_LE_REFL]);;\r
1130 \r
1131 \r
1132 (* lower *)\r
1133 let M_TAYLOR_LOWER_BOUND' = prove(`m_taylor_interval f domain y (vector w:real^N) (f_lo, f_hi) \r
1134                                     d_bounds_list dd_bounds_list /\\r
1135                                     m_taylor_error f domain (vector w) error /\\r
1136                                     ITLIST2 (\a b c. b * iabs a + c) d_bounds_list w (&0) <= b /\\r
1137                                     b + inv(&2) * error <= a /\\r
1138                                     lo <= f_lo - a /\\r
1139                                     LENGTH w = dimindex (:N) /\ LENGTH d_bounds_list = dimindex (:N) ==>\r
1140                                     (!p. p IN interval [domain] ==> lo <= f p)`,\r
1141    STRIP_TAC THEN\r
1142      MP_TAC (INST[`f_hi + a:real`, `hi:real`] M_TAYLOR_BOUND') THEN\r
1143      ASM_SIMP_TAC[interval_arith; REAL_LE_REFL]);;\r
1144 \r
1145 \r
1146 \r
1147 (* arrays *)\r
1148 let gen_taylor_bound_th bound_th n =\r
1149   let th0 = (DISCH_ALL o MY_RULE o REWRITE_RULE[MY_RULE M_TAYLOR_ERROR_ITLIST2_ALT]) bound_th in\r
1150   let ns = 1--n in\r
1151   let mk_list_hd l = mk_list (l, type_of (hd l)) in\r
1152   let w_list = mk_list_hd (map (fun i -> w_vars_array.(i)) ns) in\r
1153   let d_bounds_list = mk_list_hd (map (fun i -> df_vars_array.(i)) ns) in\r
1154   let dd_bounds_list = mk_list_hd (map (fun i -> mk_list_hd (map (fun j -> dd_vars_array.(i).(j)) (1--i))) ns) in\r
1155   let th1 = (INST[w_list, w_var_real_list; d_bounds_list, d_bounds_list_var; \r
1156                   dd_bounds_list, dd_bounds_list_var] o INST_TYPE[n_type_array.(n), nty]) th0 in\r
1157   let th2 = REWRITE_RULE[LAST; NOT_CONS_NIL; BUTLAST; all_n; ITLIST2_DEF; LENGTH; ARITH; dimindex_array.(n)] th1 in\r
1158   let th3 = REWRITE_RULE[HD; TL; REAL_MUL_RZERO; REAL_ADD_RID; GSYM error_mul_f2] th2 in\r
1159     (MY_RULE o REWRITE_RULE[float_inv2_th; SYM float2_eq]) th3;;\r
1160 \r
1161 \r
1162 let m_taylor_upper_array = Array.init (max_dim + 1)\r
1163   (fun i -> if i < 1 then TRUTH else gen_taylor_bound_th M_TAYLOR_UPPER_BOUND' i);;\r
1164 \r
1165 let m_taylor_lower_array = Array.init (max_dim + 1)\r
1166   (fun i -> if i < 1 then TRUTH else gen_taylor_bound_th M_TAYLOR_LOWER_BOUND' i);;\r
1167 \r
1168 let m_taylor_bound_array = Array.init (max_dim + 1)\r
1169   (fun i -> if i < 1 then TRUTH else gen_taylor_bound_th M_TAYLOR_BOUND' i);;\r
1170 \r
1171 \r
1172 (***************************)\r
1173 (* eval_m_taylor_bounds0 *)\r
1174 \r
1175 let eval_m_taylor_bounds0 mode n pp m_taylor_th =\r
1176   let bound_th = \r
1177     if mode = "upper" then \r
1178       m_taylor_upper_array.(n)\r
1179     else if mode = "bound" then\r
1180       m_taylor_bound_array.(n)\r
1181     else\r
1182       m_taylor_lower_array.(n) in\r
1183 \r
1184   let f_tm, domain_tm, y_tm, w_tm, f_bounds, d_bounds_list, dd_bounds_list = dest_m_taylor (concl m_taylor_th) in\r
1185   let f_lo, f_hi = dest_pair f_bounds and\r
1186       ws = dest_list (rand w_tm) and\r
1187       dfs = dest_list d_bounds_list and\r
1188       dds = map dest_list (dest_list dd_bounds_list) in\r
1189   let ns = 1--n in\r
1190   let df_vars = map (fun i -> df_vars_array.(i)) ns and\r
1191       dd_vars = map (fun i -> map (fun j -> dd_vars_array.(i).(j)) (1--i)) ns and\r
1192       w_vars = map (fun i -> w_vars_array.(i)) ns and\r
1193       y_var = mk_var ("y", type_of y_tm) and\r
1194       f_var = mk_var ("f", type_of f_tm) and\r
1195       domain_var = mk_var ("domain", type_of domain_tm) in\r
1196     \r
1197   (* sum of first partials *)\r
1198   let d_th =\r
1199     let mul_wd = map2 (error_mul_f2_le_conv2 pp) ws dfs in\r
1200       end_itlist (add_ineq_hi pp) mul_wd in\r
1201     \r
1202   let b_tm = (rand o concl) d_th in\r
1203     \r
1204   (* sum of second partials *)\r
1205   let dd_th =\r
1206     let ( * ), ( + ) = mul_ineq_pos_const_hi pp, add_ineq_hi pp in\r
1207     let mul_wdd = map2 (fun list i -> my_map2 (error_mul_f2_le_conv2 pp) ws list) dds ns in\r
1208     let sums1 = map (end_itlist ( + ) o butlast) (tl mul_wdd) in\r
1209     let sums2 = (hd o hd) mul_wdd :: map2 (fun list th1 -> last list + two_float * th1) (tl mul_wdd) sums1 in\r
1210     let sums = map2 ( * ) ws sums2 in\r
1211       end_itlist ( + ) sums in\r
1212 \r
1213   let error_tm = (rand o concl) dd_th in\r
1214 \r
1215   (* additional inequalities *)\r
1216   let ineq1_th = \r
1217     let ( * ), ( + ) = float_mul_hi pp, add_ineq_hi pp in\r
1218       mk_refl_ineq b_tm + float_inv2 * error_tm in\r
1219 \r
1220   let a_tm = (rand o concl) ineq1_th in\r
1221 \r
1222   let prove_ineq2, bounds_inst = \r
1223     if mode = "upper" then\r
1224       let ineq2 = float_add_hi pp f_hi a_tm in\r
1225         MY_PROVE_HYP ineq2, [(rand o concl) ineq2, hi_var_real]\r
1226     else if mode = "bound" then\r
1227       let ineq2 = float_add_hi pp f_hi a_tm in\r
1228       let ineq3 = float_sub_lo pp f_lo a_tm in\r
1229         MY_PROVE_HYP ineq2 o MY_PROVE_HYP ineq3, \r
1230     [(rand o concl) ineq2, hi_var_real; (lhand o concl) ineq3, lo_var_real]\r
1231     else\r
1232       let ineq2 = float_sub_lo pp f_lo a_tm in\r
1233         MY_PROVE_HYP ineq2, [(lhand o concl) ineq2, lo_var_real] in\r
1234 \r
1235   (* final step *)\r
1236   let inst_list =\r
1237     let inst1 = zip dfs df_vars in\r
1238     let inst2 = zip (List.flatten dds) (List.flatten dd_vars) in\r
1239     let inst3 = zip ws w_vars in\r
1240       inst1 @ inst2 @ inst3 in\r
1241     \r
1242     (MY_PROVE_HYP m_taylor_th o MY_PROVE_HYP dd_th o MY_PROVE_HYP d_th o \r
1243        MY_PROVE_HYP ineq1_th o prove_ineq2 o\r
1244        INST ([f_hi, f_hi_var; f_lo, f_lo_var; error_tm, error_var; \r
1245               a_tm, a_var_real; b_tm, b_var_real;\r
1246               y_tm, y_var; domain_tm, domain_var;\r
1247               f_tm, f_var] @ bounds_inst @ inst_list)) bound_th;;\r
1248 \r
1249 \r
1250 (* upper *)\r
1251 let eval_m_taylor_upper_bound = eval_m_taylor_bounds0 "upper";;\r
1252 \r
1253 (* lower *)\r
1254 let eval_m_taylor_lower_bound = eval_m_taylor_bounds0 "lower";;\r
1255 \r
1256 (* bound *)\r
1257 let eval_m_taylor_bound = eval_m_taylor_bounds0 "bound";;\r
1258 \r
1259 \r
1260 \r
1261 \r
1262 (******************************)\r
1263 (* taylor_upper_partial_bound *)\r
1264 (* taylor_lower_partial_bound *)\r
1265 \r
1266 \r
1267 (* bound *)\r
1268 let M_TAYLOR_PARTIAL_BOUND' = \r
1269   prove(`m_taylor_interval f domain (y:real^N) (vector w) f_bounds d_bounds_list dd_bounds_list /\\r
1270           i IN 1..dimindex (:N) /\\r
1271           EL (i - 1) d_bounds_list = (df_lo, df_hi) /\\r
1272           (!x. x IN interval [domain] ==> all_n 1 dd_list \r
1273              (\j int. interval_arith (if j <= i then partial2 j i f x else partial2 i j f x) int)) /\\r
1274           LENGTH dd_list = dimindex (:N) /\\r
1275           LENGTH d_bounds_list = dimindex (:N) /\ \r
1276           LENGTH w = dimindex (:N) /\\r
1277           ITLIST2 (\a b c. b * iabs a + c) dd_list w (&0) <= error /\\r
1278           df_hi + error <= hi ==>\r
1279           lo <= df_lo - error ==>\r
1280           (!x. x IN interval [domain] ==> interval_arith (partial i f x) (lo, hi))`,\r
1281    REWRITE_TAC[m_taylor_interval; m_lin_approx; ALL_N_EL; GSYM m_bounded_on_int] THEN\r
1282      set_tac "s" `ITLIST2 _1 _2 _3 _4` THEN REPEAT STRIP_TAC THEN\r
1283      SUBGOAL_THEN `1 <= i /\ i <= dimindex (:N)` (LABEL_TAC "i_in") THENL\r
1284      [\r
1285        ASM_REWRITE_TAC[GSYM IN_NUMSEG];\r
1286        ALL_TAC\r
1287      ] THEN\r
1288      SUBGOAL_THEN `diff2_domain domain (f:real^N->real)` ASSUME_TAC THENL\r
1289      [\r
1290        UNDISCH_TAC `diff2c_domain domain (f:real^N->real)` THEN\r
1291          SIMP_TAC[diff2_domain; diff2c_domain; diff2c];\r
1292        ALL_TAC\r
1293      ] THEN\r
1294      REWRITE_TAC[ETA_AX] THEN apply_tac m_taylor_partial_bounds THEN\r
1295      MAP_EVERY EXISTS_TAC [`y:real^N`; `vector w:real^N`; `error:real`; `df_lo:real`; `df_hi:real`] THEN\r
1296      ASM_REWRITE_TAC[] THEN\r
1297      FIRST_X_ASSUM (MP_TAC o SPEC `i - 1`) THEN\r
1298      ASM_SIMP_TAC[ARITH_RULE `1 <= i /\ i <= n ==> i - 1 < n`] THEN\r
1299      ASM_SIMP_TAC[ARITH_RULE `1 <= i ==> 1 + i - 1 = i`; interval_arith] THEN\r
1300      DISCH_THEN (fun th -> ALL_TAC) THEN\r
1301      REWRITE_TAC[m_taylor_partial_error] THEN REPEAT STRIP_TAC THEN\r
1302      MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `s:real` THEN ASM_REWRITE_TAC[] THEN\r
1303      EXPAND_TAC "s" THEN new_rewrite [] [] ITLIST2_EQ_SUM THEN ASM_REWRITE_TAC[LE_REFL] THEN\r
1304      MATCH_MP_TAC SUM_LE THEN REWRITE_TAC[FINITE_NUMSEG] THEN\r
1305      move ["j"; "j_in"] THEN\r
1306      ASM_SIMP_TAC[VECTOR_COMPONENT] THEN\r
1307      MATCH_MP_TAC REAL_LE_LMUL THEN CONJ_TAC THENL\r
1308      [\r
1309        UNDISCH_TAC `m_cell_domain domain (y:real^N) (vector w)` THEN\r
1310          new_rewrite [] [] pair_eq THEN REWRITE_TAC[m_cell_domain] THEN\r
1311          DISCH_THEN (MP_TAC o SPEC `j:num`) THEN ASM_REWRITE_TAC[] THEN\r
1312          ASM_SIMP_TAC[VECTOR_COMPONENT] THEN REAL_ARITH_TAC;\r
1313        ALL_TAC\r
1314      ] THEN\r
1315      FIRST_X_ASSUM (MP_TAC o SPEC `x:real^N`) THEN ASM_REWRITE_TAC[] THEN\r
1316      DISCH_THEN (MP_TAC o SPEC `j - 1`) THEN\r
1317      POP_ASSUM MP_TAC THEN REWRITE_TAC[IN_NUMSEG] THEN DISCH_TAC THEN\r
1318      ASM_SIMP_TAC[ARITH_RULE `1 <= j /\ j <= n ==> j - 1 < n`] THEN\r
1319      ASM_SIMP_TAC[ARITH_RULE `!i. 1 <= i ==> 1 + i - 1 = i`; GSYM partial2] THEN\r
1320      DISCH_THEN (MP_TAC o MATCH_MP interval_arith_abs_le) THEN\r
1321      COND_CASES_TAC THEN TRY (DISCH_THEN MATCH_MP_TAC THEN REWRITE_TAC[REAL_LE_REFL]) THEN\r
1322      new_rewrite [] [] mixed_second_partials THENL\r
1323      [\r
1324        UNDISCH_TAC `diff2c_domain domain (f:real^N->real)` THEN\r
1325          ASM_SIMP_TAC[diff2c_domain];\r
1326        ALL_TAC\r
1327      ] THEN\r
1328      DISCH_THEN MATCH_MP_TAC THEN REWRITE_TAC[REAL_LE_REFL]);;\r
1329 \r
1330 \r
1331 (* upper *)\r
1332 let M_TAYLOR_PARTIAL_UPPER' = \r
1333   prove(`m_taylor_interval f domain (y:real^N) (vector w) f_bounds d_bounds_list dd_bounds_list /\\r
1334           i IN 1..dimindex (:N) /\\r
1335           EL (i - 1) d_bounds_list = (df_lo, df_hi) /\\r
1336           (!x. x IN interval [domain] ==> all_n 1 dd_list \r
1337              (\j int. interval_arith (if j <= i then partial2 j i f x else partial2 i j f x) int)) /\\r
1338           LENGTH dd_list = dimindex (:N) /\\r
1339           LENGTH d_bounds_list = dimindex (:N) /\ \r
1340           LENGTH w = dimindex (:N) /\\r
1341           ITLIST2 (\a b c. b * iabs a + c) dd_list w (&0) <= error /\\r
1342           df_hi + error <= hi ==>\r
1343           (!x. x IN interval [domain] ==> partial i f x <= hi)`,\r
1344    REPEAT STRIP_TAC THEN\r
1345      MP_TAC (INST[`df_lo - error`, `lo:real`] M_TAYLOR_PARTIAL_BOUND') THEN ASM_REWRITE_TAC[REAL_LE_REFL] THEN\r
1346      ASM_SIMP_TAC[interval_arith]);;\r
1347 \r
1348 \r
1349 (* lower *)\r
1350 let M_TAYLOR_PARTIAL_LOWER' = \r
1351   prove(`m_taylor_interval f domain (y:real^N) (vector w) f_bounds d_bounds_list dd_bounds_list /\\r
1352           i IN 1..dimindex (:N) /\\r
1353           EL (i - 1) d_bounds_list = (df_lo, df_hi) /\\r
1354           (!x. x IN interval [domain] ==> all_n 1 dd_list \r
1355              (\j int. interval_arith (if j <= i then partial2 j i f x else partial2 i j f x) int)) /\\r
1356           LENGTH dd_list = dimindex (:N) /\\r
1357           LENGTH d_bounds_list = dimindex (:N) /\ \r
1358           LENGTH w = dimindex (:N) /\\r
1359           ITLIST2 (\a b c. b * iabs a + c) dd_list w (&0) <= error /\\r
1360           lo <= df_lo - error ==>\r
1361           (!x. x IN interval [domain] ==> lo <= partial i f x)`,\r
1362    REPEAT STRIP_TAC THEN\r
1363      MP_TAC (INST[`df_hi + error`, `hi:real`] M_TAYLOR_PARTIAL_BOUND') THEN ASM_REWRITE_TAC[REAL_LE_REFL] THEN\r
1364      ASM_SIMP_TAC[interval_arith]);;\r
1365 \r
1366 \r
1367 \r
1368 \r
1369 (* arrays *)\r
1370 let gen_taylor_partial_bound_th =\r
1371   let imp_and_eq = TAUT `((P ==> Q) /\ (P ==> R)) <=> (P ==> Q /\ R)` in\r
1372   let mk_list_hd l = mk_list (l, type_of (hd l)) in\r
1373   let dd_list_var = `dd_list : (real#real)list` in\r
1374     fun bound_th n i ->\r
1375       let ns = 1--n in\r
1376       let i_tm = mk_small_numeral i in\r
1377       let w_list = mk_list_hd (map (fun i -> w_vars_array.(i)) ns) in\r
1378       let d_bounds_list = mk_list_hd (map (fun i -> df_vars_array.(i)) ns) in\r
1379       let dd_bounds_list = mk_list_hd\r
1380         (map (fun i -> mk_list_hd (map (fun j -> dd_vars_array.(i).(j)) (1--i))) ns) in\r
1381       let dd_list = mk_list_hd \r
1382         (map (fun j -> if j <= i then dd_vars_array.(i).(j) else dd_vars_array.(j).(i)) ns) in\r
1383 \r
1384       let th1 = (INST[w_list, w_var_real_list; d_bounds_list, d_bounds_list_var; dd_list, dd_list_var; i_tm, `i:num`;\r
1385                       dd_bounds_list, dd_bounds_list_var] o INST_TYPE[n_type_array.(n), nty]) bound_th in\r
1386       let th2 = REWRITE_RULE[REAL_ADD_RID; HD; TL; ITLIST2_DEF; LENGTH; ARITH; dimindex_array.(n)] th1 in\r
1387       let th3 = REWRITE_RULE[IN_NUMSEG; ARITH; el_thms_array.(n).(i - 1)] th2 in\r
1388       let th4 = (REWRITE_RULE[] o INST[`df_lo:real, df_hi:real`, df_vars_array.(i)]) th3 in\r
1389       let th5 = (MY_RULE o REWRITE_RULE[GSYM imp_and_eq; GSYM AND_FORALL_THM; all_n; ARITH]) th4 in\r
1390 \r
1391       let m_taylor_hyp = find (can dest_m_taylor) (hyp th5) in\r
1392       let t_th0 = (REWRITE_RULE[ARITH; all_n; second_bounded; m_taylor_interval] o ASSUME) m_taylor_hyp in\r
1393       let t_th1 = REWRITE_RULE[GSYM imp_and_eq; GSYM AND_FORALL_THM] t_th0 in\r
1394         (MY_RULE_NUM o REWRITE_RULE[GSYM error_mul_f2; t_th1] o DISCH_ALL) th5;;\r
1395 \r
1396 \r
1397 (* The (n, i)-th element is the theorem |- i IN 1..dimindex (:n) *)\r
1398 let i_in_array = Array.init (max_dim + 1)\r
1399   (fun i -> Array.init (i + 1)\r
1400      (fun j ->\r
1401         if j < 1 then TRUTH else\r
1402           let j_tm = mk_small_numeral j in\r
1403           let tm0 = `j IN 1..dimindex (:N)` in\r
1404           let tm1 = (subst [j_tm, `j:num`] o inst [n_type_array.(i), nty]) tm0 in\r
1405             prove(tm1, REWRITE_TAC[dimindex_array.(i); IN_NUMSEG] THEN ARITH_TAC)));;\r
1406 \r
1407 let m_taylor_partial_upper_array, m_taylor_partial_lower_array, m_taylor_partial_bound_array =\r
1408   let gen_array bound_th = Array.init (max_dim + 1)\r
1409     (fun i -> Array.init (i + 1)\r
1410        (fun j -> if j < 1 then TRUTH else gen_taylor_partial_bound_th bound_th i j)) in\r
1411     gen_array M_TAYLOR_PARTIAL_UPPER', gen_array M_TAYLOR_PARTIAL_LOWER', gen_array M_TAYLOR_PARTIAL_BOUND';;\r
1412 \r
1413 \r
1414 (***************************)\r
1415 \r
1416 \r
1417 let eval_m_taylor_partial_bounds0 mode n pp i m_taylor_th =\r
1418   let bound_th = \r
1419     if mode = "upper" then \r
1420       m_taylor_partial_upper_array.(n).(i)\r
1421     else if mode = "bound" then\r
1422       m_taylor_partial_bound_array.(n).(i)\r
1423     else\r
1424       m_taylor_partial_lower_array.(n).(i) in\r
1425 \r
1426   let f_tm, domain_tm, y_tm, w_tm, f_bounds, d_bounds_list, dd_bounds_list = dest_m_taylor (concl m_taylor_th) in\r
1427   let ws = dest_list (rand w_tm) and\r
1428       dfs = dest_list d_bounds_list and\r
1429       dds = map dest_list (dest_list dd_bounds_list) in\r
1430   let ns = 1--n in\r
1431   let df_lo, df_hi = dest_pair (List.nth dfs (i - 1)) and\r
1432       df_vars = map (fun i -> df_vars_array.(i)) ns and\r
1433       dd_vars = map (fun i -> map (fun j -> dd_vars_array.(i).(j)) (1--i)) ns and\r
1434       w_vars = map (fun i -> w_vars_array.(i)) ns and\r
1435       y_var = mk_var ("y", type_of y_tm) and\r
1436       f_var = mk_var ("f", type_of f_tm) and\r
1437       domain_var = mk_var ("domain", type_of domain_tm) in\r
1438 \r
1439   (* sum of second partials *)\r
1440   let dd_list = map \r
1441     (fun j -> if j <= i then List.nth (List.nth dds (i-1)) (j-1) else List.nth (List.nth dds (j-1)) (i-1)) ns in\r
1442 \r
1443   let dd_th =\r
1444     let mul_dd = map2 (error_mul_f2_le_conv2 pp) ws dd_list in\r
1445       end_itlist (add_ineq_hi pp) mul_dd in\r
1446 \r
1447   let error_tm = (rand o concl) dd_th in\r
1448 \r
1449   (* additional inequalities *)\r
1450   let prove_ineq, bounds_inst = \r
1451     if mode = "upper" then\r
1452       let ineq2 = float_add_hi pp df_hi error_tm in\r
1453         MY_PROVE_HYP ineq2, [(rand o concl) ineq2, hi_var_real]\r
1454     else if mode = "bound" then\r
1455       let ineq2 = float_add_hi pp df_hi error_tm in\r
1456       let ineq3 = float_sub_lo pp df_lo error_tm in\r
1457         MY_PROVE_HYP ineq2 o MY_PROVE_HYP ineq3, \r
1458     [(rand o concl) ineq2, hi_var_real; (lhand o concl) ineq3, lo_var_real]\r
1459     else\r
1460       let ineq2 = float_sub_lo pp df_lo error_tm in\r
1461         MY_PROVE_HYP ineq2, [(lhand o concl) ineq2, lo_var_real] in\r
1462 \r
1463   (* final step *)\r
1464   let inst_list =\r
1465     let inst1 = zip dfs df_vars in\r
1466     let inst2 = zip (List.flatten dds) (List.flatten dd_vars) in\r
1467     let inst3 = zip ws w_vars in\r
1468       inst1 @ inst2 @ inst3 in\r
1469     \r
1470     (MY_PROVE_HYP m_taylor_th o MY_PROVE_HYP dd_th o prove_ineq o\r
1471        INST ([df_hi, df_hi_var; df_lo, df_lo_var; error_tm, error_var; \r
1472               y_tm, y_var; domain_tm, domain_var; f_bounds, f_bounds_var;\r
1473               f_tm, f_var] @ bounds_inst @ inst_list)) bound_th;;\r
1474 \r
1475 \r
1476 (* upper *)\r
1477 let eval_m_taylor_partial_upper = eval_m_taylor_partial_bounds0 "upper";;\r
1478 \r
1479 (* lower *)\r
1480 let eval_m_taylor_partial_lower = eval_m_taylor_partial_bounds0 "lower";;\r
1481 \r
1482 (* bound *)\r
1483 let eval_m_taylor_partial_bound = eval_m_taylor_partial_bounds0 "bound";;\r
1484 \r
1485 \r
1486 \r
1487 end;;\r