src/HOL/Decision_Procs/Approximation.thy
author hoelzl
Mon, 06 Dec 2010 19:54:48 +0100
changeset 41024 ba961a606c67
parent 41022 81d337539d57
child 41413 64cd30d6b0b8
permissions -rw-r--r--
move coercions to appropriate places
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
     1
(* Author:     Johannes Hoelzl, TU Muenchen
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
     2
   Coercions removed by Dmitriy Traytel *)
30122
1c912a9d8200 standard headers;
wenzelm
parents: 29823
diff changeset
     3
30886
dda08b76fa99 updated official title of contribution by Johannes Hoelzl;
wenzelm
parents: 30549
diff changeset
     4
header {* Prove Real Valued Inequalities by Computation *}
30122
1c912a9d8200 standard headers;
wenzelm
parents: 29823
diff changeset
     5
40892
6f7292b94652 proper theory name (cf. e84f82418e09);
wenzelm
parents: 40881
diff changeset
     6
theory Approximation
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
     7
imports Complex_Main Float Reflection "~~/src/HOL/Decision_Procs/Dense_Linear_Order" Efficient_Nat
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
     8
begin
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
     9
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    10
section "Horner Scheme"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    11
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    12
subsection {* Define auxiliary helper @{text horner} function *}
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    13
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    14
primrec horner :: "(nat \<Rightarrow> nat) \<Rightarrow> (nat \<Rightarrow> nat \<Rightarrow> nat) \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> real \<Rightarrow> real" where
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    15
"horner F G 0 i k x       = 0" |
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    16
"horner F G (Suc n) i k x = 1 / k - x * horner F G n (F i) (G i k) x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    17
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    18
lemma horner_schema': fixes x :: real  and a :: "nat \<Rightarrow> real"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    19
  shows "a 0 - x * (\<Sum> i=0..<n. (-1)^i * a (Suc i) * x^i) = (\<Sum> i=0..<Suc n. (-1)^i * a i * x^i)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    20
proof -
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    21
  have shift_pow: "\<And>i. - (x * ((-1)^i * a (Suc i) * x ^ i)) = (-1)^(Suc i) * a (Suc i) * x ^ (Suc i)" by auto
37887
2ae085b07f2f diff_minus subsumes diff_def
haftmann
parents: 37677
diff changeset
    22
  show ?thesis unfolding setsum_right_distrib shift_pow diff_minus setsum_negf[symmetric] setsum_head_upt_Suc[OF zero_less_Suc]
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    23
    setsum_reindex[OF inj_Suc, unfolded comp_def, symmetric, of "\<lambda> n. (-1)^n  *a n * x^n"] by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    24
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    25
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    26
lemma horner_schema: fixes f :: "nat \<Rightarrow> nat" and G :: "nat \<Rightarrow> nat \<Rightarrow> nat" and F :: "nat \<Rightarrow> nat"
30971
7fbebf75b3ef funpow and relpow with shared "^^" syntax
haftmann
parents: 30968
diff changeset
    27
  assumes f_Suc: "\<And>n. f (Suc n) = G ((F ^^ n) s) (f n)"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    28
  shows "horner F G n ((F ^^ j') s) (f j') x = (\<Sum> j = 0..< n. -1 ^ j * (1 / (f (j' + j))) * x ^ j)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    29
proof (induct n arbitrary: i k j')
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    30
  case (Suc n)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    31
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    32
  show ?case unfolding horner.simps Suc[where j'="Suc j'", unfolded funpow.simps comp_def f_Suc]
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    33
    using horner_schema'[of "\<lambda> j. 1 / (f (j' + j))"] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    34
qed auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    35
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    36
lemma horner_bounds':
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    37
  fixes lb :: "nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> float \<Rightarrow> float" and ub :: "nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> float \<Rightarrow> float"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    38
  assumes "0 \<le> real x" and f_Suc: "\<And>n. f (Suc n) = G ((F ^^ n) s) (f n)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    39
  and lb_0: "\<And> i k x. lb 0 i k x = 0"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    40
  and lb_Suc: "\<And> n i k x. lb (Suc n) i k x = lapprox_rat prec 1 k - x * (ub n (F i) (G i k) x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    41
  and ub_0: "\<And> i k x. ub 0 i k x = 0"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    42
  and ub_Suc: "\<And> n i k x. ub (Suc n) i k x = rapprox_rat prec 1 k - x * (lb n (F i) (G i k) x)"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    43
  shows "(lb n ((F ^^ j') s) (f j') x) \<le> horner F G n ((F ^^ j') s) (f j') x \<and>
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    44
         horner F G n ((F ^^ j') s) (f j') x \<le> (ub n ((F ^^ j') s) (f j') x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    45
  (is "?lb n j' \<le> ?horner n j' \<and> ?horner n j' \<le> ?ub n j'")
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    46
proof (induct n arbitrary: j')
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    47
  case 0 thus ?case unfolding lb_0 ub_0 horner.simps by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    48
next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    49
  case (Suc n)
37887
2ae085b07f2f diff_minus subsumes diff_def
haftmann
parents: 37677
diff changeset
    50
  have "?lb (Suc n) j' \<le> ?horner (Suc n) j'" unfolding lb_Suc ub_Suc horner.simps real_of_float_sub diff_minus
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    51
  proof (rule add_mono)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    52
    show "(lapprox_rat prec 1 (f j')) \<le> 1 / (f j')" using lapprox_rat[of prec 1  "f j'"] by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    53
    from Suc[where j'="Suc j'", unfolded funpow.simps comp_def f_Suc, THEN conjunct2] `0 \<le> real x`
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    54
    show "- real (x * ub n (F ((F ^^ j') s)) (G ((F ^^ j') s) (f j')) x) \<le>
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    55
          - (x * horner F G n (F ((F ^^ j') s)) (G ((F ^^ j') s) (f j')) x)"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    56
      unfolding real_of_float_mult neg_le_iff_le by (rule mult_left_mono)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    57
  qed
37887
2ae085b07f2f diff_minus subsumes diff_def
haftmann
parents: 37677
diff changeset
    58
  moreover have "?horner (Suc n) j' \<le> ?ub (Suc n) j'" unfolding ub_Suc ub_Suc horner.simps real_of_float_sub diff_minus
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    59
  proof (rule add_mono)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    60
    show "1 / (f j') \<le> (rapprox_rat prec 1 (f j'))" using rapprox_rat[of 1 "f j'" prec] by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    61
    from Suc[where j'="Suc j'", unfolded funpow.simps comp_def f_Suc, THEN conjunct1] `0 \<le> real x`
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    62
    show "- (x * horner F G n (F ((F ^^ j') s)) (G ((F ^^ j') s) (f j')) x) \<le>
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    63
          - real (x * lb n (F ((F ^^ j') s)) (G ((F ^^ j') s) (f j')) x)"
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    64
      unfolding real_of_float_mult neg_le_iff_le by (rule mult_left_mono)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    65
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    66
  ultimately show ?case by blast
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    67
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    68
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    69
subsection "Theorems for floating point functions implementing the horner scheme"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    70
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    71
text {*
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    72
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    73
Here @{term_type "f :: nat \<Rightarrow> nat"} is the sequence defining the Taylor series, the coefficients are
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    74
all alternating and reciprocs. We use @{term G} and @{term F} to describe the computation of @{term f}.
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    75
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    76
*}
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    77
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    78
lemma horner_bounds: fixes F :: "nat \<Rightarrow> nat" and G :: "nat \<Rightarrow> nat \<Rightarrow> nat"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    79
  assumes "0 \<le> real x" and f_Suc: "\<And>n. f (Suc n) = G ((F ^^ n) s) (f n)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    80
  and lb_0: "\<And> i k x. lb 0 i k x = 0"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    81
  and lb_Suc: "\<And> n i k x. lb (Suc n) i k x = lapprox_rat prec 1 k - x * (ub n (F i) (G i k) x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    82
  and ub_0: "\<And> i k x. ub 0 i k x = 0"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    83
  and ub_Suc: "\<And> n i k x. ub (Suc n) i k x = rapprox_rat prec 1 k - x * (lb n (F i) (G i k) x)"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    84
  shows "(lb n ((F ^^ j') s) (f j') x) \<le> (\<Sum>j=0..<n. -1 ^ j * (1 / (f (j' + j))) * (x ^ j))" (is "?lb") and
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    85
    "(\<Sum>j=0..<n. -1 ^ j * (1 / (f (j' + j))) * (x ^ j)) \<le> (ub n ((F ^^ j') s) (f j') x)" (is "?ub")
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    86
proof -
31809
hoelzl
parents: 31790
diff changeset
    87
  have "?lb  \<and> ?ub"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    88
    using horner_bounds'[where lb=lb, OF `0 \<le> real x` f_Suc lb_0 lb_Suc ub_0 ub_Suc]
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    89
    unfolding horner_schema[where f=f, OF f_Suc] .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    90
  thus "?lb" and "?ub" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    91
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    92
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    93
lemma horner_bounds_nonpos: fixes F :: "nat \<Rightarrow> nat" and G :: "nat \<Rightarrow> nat \<Rightarrow> nat"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
    94
  assumes "real x \<le> 0" and f_Suc: "\<And>n. f (Suc n) = G ((F ^^ n) s) (f n)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    95
  and lb_0: "\<And> i k x. lb 0 i k x = 0"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    96
  and lb_Suc: "\<And> n i k x. lb (Suc n) i k x = lapprox_rat prec 1 k + x * (ub n (F i) (G i k) x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
    97
  and ub_0: "\<And> i k x. ub 0 i k x = 0"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    98
  and ub_Suc: "\<And> n i k x. ub (Suc n) i k x = rapprox_rat prec 1 k + x * (lb n (F i) (G i k) x)"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
    99
  shows "(lb n ((F ^^ j') s) (f j') x) \<le> (\<Sum>j=0..<n. (1 / (f (j' + j))) * real x ^ j)" (is "?lb") and
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   100
    "(\<Sum>j=0..<n. (1 / (f (j' + j))) * real x ^ j) \<le> (ub n ((F ^^ j') s) (f j') x)" (is "?ub")
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   101
proof -
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   102
  { fix x y z :: float have "x - y * z = x + - y * z"
30968
10fef94f40fc adaptions due to rearrangment of power operation
haftmann
parents: 30952
diff changeset
   103
      by (cases x, cases y, cases z, simp add: plus_float.simps minus_float_def uminus_float.simps times_float.simps algebra_simps)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   104
  } note diff_mult_minus = this
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   105
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   106
  { fix x :: float have "- (- x) = x" by (cases x, auto simp add: uminus_float.simps) } note minus_minus = this
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   107
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   108
  have move_minus: "(-x) = -1 * real x" by auto (* coercion "inside" is necessary *)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   109
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   110
  have sum_eq: "(\<Sum>j=0..<n. (1 / (f (j' + j))) * real x ^ j) =
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   111
    (\<Sum>j = 0..<n. -1 ^ j * (1 / (f (j' + j))) * real (- x) ^ j)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   112
  proof (rule setsum_cong, simp)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   113
    fix j assume "j \<in> {0 ..< n}"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   114
    show "1 / (f (j' + j)) * real x ^ j = -1 ^ j * (1 / (f (j' + j))) * real (- x) ^ j"
36778
739a9379e29b avoid using real-specific versions of generic lemmas
huffman
parents: 36531
diff changeset
   115
      unfolding move_minus power_mult_distrib mult_assoc[symmetric]
739a9379e29b avoid using real-specific versions of generic lemmas
huffman
parents: 36531
diff changeset
   116
      unfolding mult_commute unfolding mult_assoc[of "-1 ^ j", symmetric] power_mult_distrib[symmetric]
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   117
      by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   118
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   119
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   120
  have "0 \<le> real (-x)" using assms by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   121
  from horner_bounds[where G=G and F=F and f=f and s=s and prec=prec
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   122
    and lb="\<lambda> n i k x. lb n i k (-x)" and ub="\<lambda> n i k x. ub n i k (-x)", unfolded lb_Suc ub_Suc diff_mult_minus,
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   123
    OF this f_Suc lb_0 refl ub_0 refl]
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   124
  show "?lb" and "?ub" unfolding minus_minus sum_eq
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   125
    by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   126
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   127
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   128
subsection {* Selectors for next even or odd number *}
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   129
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   130
text {*
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   131
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   132
The horner scheme computes alternating series. To get the upper and lower bounds we need to
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   133
guarantee to access a even or odd member. To do this we use @{term get_odd} and @{term get_even}.
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   134
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   135
*}
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   136
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   137
definition get_odd :: "nat \<Rightarrow> nat" where
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   138
  "get_odd n = (if odd n then n else (Suc n))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   139
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   140
definition get_even :: "nat \<Rightarrow> nat" where
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   141
  "get_even n = (if even n then n else (Suc n))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   142
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   143
lemma get_odd[simp]: "odd (get_odd n)" unfolding get_odd_def by (cases "odd n", auto)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   144
lemma get_even[simp]: "even (get_even n)" unfolding get_even_def by (cases "even n", auto)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   145
lemma get_odd_ex: "\<exists> k. Suc k = get_odd n \<and> odd (Suc k)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   146
proof (cases "odd n")
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   147
  case True hence "0 < n" by (rule odd_pos)
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   148
  from gr0_implies_Suc[OF this] obtain k where "Suc k = n" by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   149
  thus ?thesis unfolding get_odd_def if_P[OF True] using True[unfolded `Suc k = n`[symmetric]] by blast
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   150
next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   151
  case False hence "odd (Suc n)" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   152
  thus ?thesis unfolding get_odd_def if_not_P[OF False] by blast
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   153
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   154
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   155
lemma get_even_double: "\<exists>i. get_even n = 2 * i" using get_even[unfolded even_mult_two_ex] .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   156
lemma get_odd_double: "\<exists>i. get_odd n = 2 * i + 1" using get_odd[unfolded odd_Suc_mult_two_ex] by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   157
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   158
section "Power function"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   159
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   160
definition float_power_bnds :: "nat \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float * float" where
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   161
"float_power_bnds n l u = (if odd n \<or> 0 < l then (l ^ n, u ^ n)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   162
                      else if u < 0         then (u ^ n, l ^ n)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   163
                                            else (0, (max (-l) u) ^ n))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   164
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   165
lemma float_power_bnds: fixes x :: real
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   166
  assumes "(l1, u1) = float_power_bnds n l u" and "x \<in> {l .. u}"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   167
  shows "x ^ n \<in> {l1..u1}"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   168
proof (cases "even n")
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   169
  case True
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   170
  show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   171
  proof (cases "0 < l")
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   172
    case True hence "odd n \<or> 0 < l" and "0 \<le> real l" unfolding less_float_def by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   173
    have u1: "u1 = u ^ n" and l1: "l1 = l ^ n" using assms unfolding float_power_bnds_def if_P[OF `odd n \<or> 0 < l`] by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   174
    have "real l ^ n \<le> x ^ n" and "x ^ n \<le> real u ^ n " using `0 \<le> real l` and assms unfolding atLeastAtMost_iff using power_mono[of l x] power_mono[of x u] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   175
    thus ?thesis using assms `0 < l` unfolding atLeastAtMost_iff l1 u1 float_power less_float_def by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   176
  next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   177
    case False hence P: "\<not> (odd n \<or> 0 < l)" using `even n` by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   178
    show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   179
    proof (cases "u < 0")
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   180
      case True hence "0 \<le> - real u" and "- real u \<le> - x" and "0 \<le> - x" and "-x \<le> - real l" using assms unfolding less_float_def by auto
31809
hoelzl
parents: 31790
diff changeset
   181
      hence "real u ^ n \<le> x ^ n" and "x ^ n \<le> real l ^ n" using power_mono[of  "-x" "-real l" n] power_mono[of "-real u" "-x" n]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   182
        unfolding power_minus_even[OF `even n`] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   183
      moreover have u1: "u1 = l ^ n" and l1: "l1 = u ^ n" using assms unfolding float_power_bnds_def if_not_P[OF P] if_P[OF True] by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   184
      ultimately show ?thesis using float_power by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   185
    next
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   186
      case False
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   187
      have "\<bar>x\<bar> \<le> real (max (-l) u)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   188
      proof (cases "-l \<le> u")
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   189
        case True thus ?thesis unfolding max_def if_P[OF True] using assms unfolding le_float_def by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   190
      next
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   191
        case False thus ?thesis unfolding max_def if_not_P[OF False] using assms unfolding le_float_def by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   192
      qed
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   193
      hence x_abs: "\<bar>x\<bar> \<le> \<bar>real (max (-l) u)\<bar>" by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   194
      have u1: "u1 = (max (-l) u) ^ n" and l1: "l1 = 0" using assms unfolding float_power_bnds_def if_not_P[OF P] if_not_P[OF False] by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   195
      show ?thesis unfolding atLeastAtMost_iff l1 u1 float_power using zero_le_even_power[OF `even n`] power_mono_even[OF `even n` x_abs] by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   196
    qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   197
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   198
next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   199
  case False hence "odd n \<or> 0 < l" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   200
  have u1: "u1 = u ^ n" and l1: "l1 = l ^ n" using assms unfolding float_power_bnds_def if_P[OF `odd n \<or> 0 < l`] by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   201
  have "real l ^ n \<le> x ^ n" and "x ^ n \<le> real u ^ n " using assms unfolding atLeastAtMost_iff using power_mono_odd[OF False] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   202
  thus ?thesis unfolding atLeastAtMost_iff l1 u1 float_power less_float_def by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   203
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   204
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   205
lemma bnds_power: "\<forall> (x::real) l u. (l1, u1) = float_power_bnds n l u \<and> x \<in> {l .. u} \<longrightarrow> l1 \<le> x ^ n \<and> x ^ n \<le> u1"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   206
  using float_power_bnds by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   207
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   208
section "Square root"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   209
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   210
text {*
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   211
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   212
The square root computation is implemented as newton iteration. As first first step we use the
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   213
nearest power of two greater than the square root.
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   214
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   215
*}
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   216
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   217
fun sqrt_iteration :: "nat \<Rightarrow> nat \<Rightarrow> float \<Rightarrow> float" where
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   218
"sqrt_iteration prec 0 (Float m e) = Float 1 ((e + bitlen m) div 2 + 1)" |
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   219
"sqrt_iteration prec (Suc m) x = (let y = sqrt_iteration prec m x
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   220
                                  in Float 1 -1 * (y + float_divr prec x y))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   221
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   222
function ub_sqrt lb_sqrt :: "nat \<Rightarrow> float \<Rightarrow> float" where
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   223
"ub_sqrt prec x = (if 0 < x then (sqrt_iteration prec prec x)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   224
              else if x < 0 then - lb_sqrt prec (- x)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   225
                            else 0)" |
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   226
"lb_sqrt prec x = (if 0 < x then (float_divl prec x (sqrt_iteration prec prec x))
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   227
              else if x < 0 then - ub_sqrt prec (- x)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   228
                            else 0)"
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   229
by pat_completeness auto
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   230
termination by (relation "measure (\<lambda> v. let (prec, x) = sum_case id id v in (if x < 0 then 1 else 0))", auto simp add: less_float_def)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   231
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   232
declare lb_sqrt.simps[simp del]
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   233
declare ub_sqrt.simps[simp del]
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   234
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   235
lemma sqrt_ub_pos_pos_1:
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   236
  assumes "sqrt x < b" and "0 < b" and "0 < x"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   237
  shows "sqrt x < (b + x / b)/2"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   238
proof -
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   239
  from assms have "0 < (b - sqrt x) ^ 2 " by simp
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   240
  also have "\<dots> = b ^ 2 - 2 * b * sqrt x + (sqrt x) ^ 2" by algebra
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   241
  also have "\<dots> = b ^ 2 - 2 * b * sqrt x + x" using assms by (simp add: real_sqrt_pow2)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   242
  finally have "0 < b ^ 2 - 2 * b * sqrt x + x" by assumption
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   243
  hence "0 < b / 2 - sqrt x + x / (2 * b)" using assms
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   244
    by (simp add: field_simps power2_eq_square)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   245
  thus ?thesis by (simp add: field_simps)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   246
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   247
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   248
lemma sqrt_iteration_bound: assumes "0 < real x"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   249
  shows "sqrt x < (sqrt_iteration prec n x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   250
proof (induct n)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   251
  case 0
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   252
  show ?case
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   253
  proof (cases x)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   254
    case (Float m e)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   255
    hence "0 < m" using float_pos_m_pos[unfolded less_float_def] assms by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   256
    hence "0 < sqrt m" by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   257
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   258
    have int_nat_bl: "(nat (bitlen m)) = bitlen m" using bitlen_ge0 by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   259
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   260
    have "x = (m / 2^nat (bitlen m)) * pow2 (e + (nat (bitlen m)))"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   261
      unfolding pow2_add pow2_int Float real_of_float_simp by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   262
    also have "\<dots> < 1 * pow2 (e + nat (bitlen m))"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   263
    proof (rule mult_strict_right_mono, auto)
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   264
      show "real m < 2^nat (bitlen m)" using bitlen_bounds[OF `0 < m`, THEN conjunct2]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   265
        unfolding real_of_int_less_iff[of m, symmetric] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   266
    qed
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   267
    finally have "sqrt x < sqrt (pow2 (e + bitlen m))" unfolding int_nat_bl by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   268
    also have "\<dots> \<le> pow2 ((e + bitlen m) div 2 + 1)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   269
    proof -
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   270
      let ?E = "e + bitlen m"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   271
      have E_mod_pow: "pow2 (?E mod 2) < 4"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   272
      proof (cases "?E mod 2 = 1")
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   273
        case True thus ?thesis by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   274
      next
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   275
        case False
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   276
        have "0 \<le> ?E mod 2" by auto
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   277
        have "?E mod 2 < 2" by auto
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   278
        from this[THEN zless_imp_add1_zle]
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   279
        have "?E mod 2 \<le> 0" using False by auto
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   280
        from xt1(5)[OF `0 \<le> ?E mod 2` this]
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   281
        show ?thesis by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   282
      qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   283
      hence "sqrt (pow2 (?E mod 2)) < sqrt (2 * 2)" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   284
      hence E_mod_pow: "sqrt (pow2 (?E mod 2)) < 2" unfolding real_sqrt_abs2 by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   285
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   286
      have E_eq: "pow2 ?E = pow2 (?E div 2 + ?E div 2 + ?E mod 2)" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   287
      have "sqrt (pow2 ?E) = sqrt (pow2 (?E div 2) * pow2 (?E div 2) * pow2 (?E mod 2))"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   288
        unfolding E_eq unfolding pow2_add ..
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   289
      also have "\<dots> = pow2 (?E div 2) * sqrt (pow2 (?E mod 2))"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   290
        unfolding real_sqrt_mult[of _ "pow2 (?E mod 2)"] real_sqrt_abs2 by auto
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   291
      also have "\<dots> < pow2 (?E div 2) * 2"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   292
        by (rule mult_strict_left_mono, auto intro: E_mod_pow)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   293
      also have "\<dots> = pow2 (?E div 2 + 1)" unfolding zadd_commute[of _ 1] pow2_add1 by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   294
      finally show ?thesis by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   295
    qed
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   296
    finally show ?thesis
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   297
      unfolding Float sqrt_iteration.simps real_of_float_simp by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   298
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   299
next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   300
  case (Suc n)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   301
  let ?b = "sqrt_iteration prec n x"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   302
  have "0 < sqrt x" using `0 < real x` by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   303
  also have "\<dots> < real ?b" using Suc .
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   304
  finally have "sqrt x < (?b + x / ?b)/2" using sqrt_ub_pos_pos_1[OF Suc _ `0 < real x`] by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   305
  also have "\<dots> \<le> (?b + (float_divr prec x ?b))/2" by (rule divide_right_mono, auto simp add: float_divr)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   306
  also have "\<dots> = (Float 1 -1) * (?b + (float_divr prec x ?b))" by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   307
  finally show ?case unfolding sqrt_iteration.simps Let_def real_of_float_mult real_of_float_add right_distrib .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   308
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   309
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   310
lemma sqrt_iteration_lower_bound: assumes "0 < real x"
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   311
  shows "0 < real (sqrt_iteration prec n x)" (is "0 < ?sqrt")
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   312
proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   313
  have "0 < sqrt x" using assms by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   314
  also have "\<dots> < ?sqrt" using sqrt_iteration_bound[OF assms] .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   315
  finally show ?thesis .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   316
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   317
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   318
lemma lb_sqrt_lower_bound: assumes "0 \<le> real x"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   319
  shows "0 \<le> real (lb_sqrt prec x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   320
proof (cases "0 < x")
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   321
  case True hence "0 < real x" and "0 \<le> x" using `0 \<le> real x` unfolding less_float_def le_float_def by auto
31809
hoelzl
parents: 31790
diff changeset
   322
  hence "0 < sqrt_iteration prec prec x" unfolding less_float_def using sqrt_iteration_lower_bound by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   323
  hence "0 \<le> real (float_divl prec x (sqrt_iteration prec prec x))" using float_divl_lower_bound[OF `0 \<le> x`] unfolding le_float_def by auto
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   324
  thus ?thesis unfolding lb_sqrt.simps using True by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   325
next
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   326
  case False with `0 \<le> real x` have "real x = 0" unfolding less_float_def by auto
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   327
  thus ?thesis unfolding lb_sqrt.simps less_float_def by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   328
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   329
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   330
lemma bnds_sqrt':
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   331
  shows "sqrt x \<in> {(lb_sqrt prec x) .. (ub_sqrt prec x) }"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   332
proof -
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   333
  { fix x :: float assume "0 < x"
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   334
    hence "0 < real x" and "0 \<le> real x" unfolding less_float_def by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   335
    hence sqrt_gt0: "0 < sqrt x" by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   336
    hence sqrt_ub: "sqrt x < sqrt_iteration prec prec x" using sqrt_iteration_bound by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   337
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   338
    have "(float_divl prec x (sqrt_iteration prec prec x)) \<le>
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   339
          x / (sqrt_iteration prec prec x)" by (rule float_divl)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   340
    also have "\<dots> < x / sqrt x"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   341
      by (rule divide_strict_left_mono[OF sqrt_ub `0 < real x`
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   342
               mult_pos_pos[OF order_less_trans[OF sqrt_gt0 sqrt_ub] sqrt_gt0]])
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   343
    also have "\<dots> = sqrt x"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   344
      unfolding inverse_eq_iff_eq[of _ "sqrt x", symmetric]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   345
                sqrt_divide_self_eq[OF `0 \<le> real x`, symmetric] by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   346
    finally have "lb_sqrt prec x \<le> sqrt x"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   347
      unfolding lb_sqrt.simps if_P[OF `0 < x`] by auto }
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   348
  note lb = this
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   349
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   350
  { fix x :: float assume "0 < x"
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   351
    hence "0 < real x" unfolding less_float_def by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   352
    hence "0 < sqrt x" by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   353
    hence "sqrt x < sqrt_iteration prec prec x"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   354
      using sqrt_iteration_bound by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   355
    hence "sqrt x \<le> ub_sqrt prec x"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   356
      unfolding ub_sqrt.simps if_P[OF `0 < x`] by auto }
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   357
  note ub = this
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   358
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   359
  show ?thesis
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   360
  proof (cases "0 < x")
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   361
    case True with lb ub show ?thesis by auto
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   362
  next case False show ?thesis
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   363
  proof (cases "real x = 0")
31809
hoelzl
parents: 31790
diff changeset
   364
    case True thus ?thesis
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   365
      by (auto simp add: less_float_def lb_sqrt.simps ub_sqrt.simps)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   366
  next
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   367
    case False with `\<not> 0 < x` have "x < 0" and "0 < -x"
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   368
      by (auto simp add: less_float_def)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   369
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   370
    with `\<not> 0 < x`
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   371
    show ?thesis using lb[OF `0 < -x`] ub[OF `0 < -x`]
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   372
      by (auto simp add: real_sqrt_minus lb_sqrt.simps ub_sqrt.simps)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   373
  qed qed
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   374
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   375
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   376
lemma bnds_sqrt: "\<forall> (x::real) lx ux. (l, u) = (lb_sqrt prec lx, ub_sqrt prec ux) \<and> x \<in> {lx .. ux} \<longrightarrow> l \<le> sqrt x \<and> sqrt x \<le> u"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   377
proof ((rule allI) +, rule impI, erule conjE, rule conjI)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   378
  fix x :: real fix lx ux
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   379
  assume "(l, u) = (lb_sqrt prec lx, ub_sqrt prec ux)"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   380
    and x: "x \<in> {lx .. ux}"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   381
  hence l: "l = lb_sqrt prec lx " and u: "u = ub_sqrt prec ux" by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   382
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   383
  have "sqrt lx \<le> sqrt x" using x by auto
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   384
  from order_trans[OF _ this]
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   385
  show "l \<le> sqrt x" unfolding l using bnds_sqrt'[of lx prec] by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   386
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   387
  have "sqrt x \<le> sqrt ux" using x by auto
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   388
  from order_trans[OF this]
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   389
  show "sqrt x \<le> u" unfolding u using bnds_sqrt'[of ux prec] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   390
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   391
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   392
section "Arcus tangens and \<pi>"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   393
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   394
subsection "Compute arcus tangens series"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   395
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   396
text {*
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   397
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   398
As first step we implement the computation of the arcus tangens series. This is only valid in the range
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   399
@{term "{-1 :: real .. 1}"}. This is used to compute \<pi> and then the entire arcus tangens.
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   400
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   401
*}
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   402
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   403
fun ub_arctan_horner :: "nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> float \<Rightarrow> float"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   404
and lb_arctan_horner :: "nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> float \<Rightarrow> float" where
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   405
  "ub_arctan_horner prec 0 k x = 0"
31809
hoelzl
parents: 31790
diff changeset
   406
| "ub_arctan_horner prec (Suc n) k x =
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   407
    (rapprox_rat prec 1 k) - x * (lb_arctan_horner prec n (k + 2) x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   408
| "lb_arctan_horner prec 0 k x = 0"
31809
hoelzl
parents: 31790
diff changeset
   409
| "lb_arctan_horner prec (Suc n) k x =
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   410
    (lapprox_rat prec 1 k) - x * (ub_arctan_horner prec n (k + 2) x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   411
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   412
lemma arctan_0_1_bounds': assumes "0 \<le> real x" "real x \<le> 1" and "even n"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   413
  shows "arctan x \<in> {(x * lb_arctan_horner prec n 1 (x * x)) .. (x * ub_arctan_horner prec (Suc n) 1 (x * x))}"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   414
proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   415
  let "?c i" = "-1^i * (1 / (i * 2 + (1::nat)) * real x ^ (i * 2 + 1))"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   416
  let "?S n" = "\<Sum> i=0..<n. ?c i"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   417
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   418
  have "0 \<le> real (x * x)" by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   419
  from `even n` obtain m where "2 * m = n" unfolding even_mult_two_ex by auto
31809
hoelzl
parents: 31790
diff changeset
   420
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   421
  have "arctan x \<in> { ?S n .. ?S (Suc n) }"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   422
  proof (cases "real x = 0")
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   423
    case False
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   424
    hence "0 < real x" using `0 \<le> real x` by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   425
    hence prem: "0 < 1 / (0 * 2 + (1::nat)) * real x ^ (0 * 2 + 1)" by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   426
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   427
    have "\<bar> real x \<bar> \<le> 1"  using `0 \<le> real x` `real x \<le> 1` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   428
    from mp[OF summable_Leibniz(2)[OF zeroseq_arctan_series[OF this] monoseq_arctan_series[OF this]] prem, THEN spec, of m, unfolded `2 * m = n`]
31790
05c92381363c corrected and unified thm names
nipkow
parents: 31508
diff changeset
   429
    show ?thesis unfolding arctan_series[OF `\<bar> real x \<bar> \<le> 1`] Suc_eq_plus1  .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   430
  qed auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   431
  note arctan_bounds = this[unfolded atLeastAtMost_iff]
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   432
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   433
  have F: "\<And>n. 2 * Suc n + 1 = 2 * n + 1 + 2" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   434
31809
hoelzl
parents: 31790
diff changeset
   435
  note bounds = horner_bounds[where s=1 and f="\<lambda>i. 2 * i + 1" and j'=0
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   436
    and lb="\<lambda>n i k x. lb_arctan_horner prec n k x"
31809
hoelzl
parents: 31790
diff changeset
   437
    and ub="\<lambda>n i k x. ub_arctan_horner prec n k x",
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   438
    OF `0 \<le> real (x*x)` F lb_arctan_horner.simps ub_arctan_horner.simps]
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   439
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   440
  { have "(x * lb_arctan_horner prec n 1 (x*x)) \<le> ?S n"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   441
      using bounds(1) `0 \<le> real x`
36778
739a9379e29b avoid using real-specific versions of generic lemmas
huffman
parents: 36531
diff changeset
   442
      unfolding real_of_float_mult power_add power_one_right mult_assoc[symmetric] setsum_left_distrib[symmetric]
739a9379e29b avoid using real-specific versions of generic lemmas
huffman
parents: 36531
diff changeset
   443
      unfolding mult_commute[where 'a=real] mult_commute[of _ "2::nat"] power_mult power2_eq_square[of "real x"]
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   444
      by (auto intro!: mult_left_mono)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   445
    also have "\<dots> \<le> arctan x" using arctan_bounds ..
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   446
    finally have "(x * lb_arctan_horner prec n 1 (x*x)) \<le> arctan x" . }
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   447
  moreover
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   448
  { have "arctan x \<le> ?S (Suc n)" using arctan_bounds ..
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   449
    also have "\<dots> \<le> (x * ub_arctan_horner prec (Suc n) 1 (x*x))"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   450
      using bounds(2)[of "Suc n"] `0 \<le> real x`
36778
739a9379e29b avoid using real-specific versions of generic lemmas
huffman
parents: 36531
diff changeset
   451
      unfolding real_of_float_mult power_add power_one_right mult_assoc[symmetric] setsum_left_distrib[symmetric]
739a9379e29b avoid using real-specific versions of generic lemmas
huffman
parents: 36531
diff changeset
   452
      unfolding mult_commute[where 'a=real] mult_commute[of _ "2::nat"] power_mult power2_eq_square[of "real x"]
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   453
      by (auto intro!: mult_left_mono)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   454
    finally have "arctan x \<le> (x * ub_arctan_horner prec (Suc n) 1 (x*x))" . }
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   455
  ultimately show ?thesis by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   456
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   457
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   458
lemma arctan_0_1_bounds: assumes "0 \<le> real x" "real x \<le> 1"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   459
  shows "arctan x \<in> {(x * lb_arctan_horner prec (get_even n) 1 (x * x)) .. (x * ub_arctan_horner prec (get_odd n) 1 (x * x))}"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   460
proof (cases "even n")
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   461
  case True
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   462
  obtain n' where "Suc n' = get_odd n" and "odd (Suc n')" using get_odd_ex by auto
31148
7ba7c1f8bc22 Cleaned up Parity a little
nipkow
parents: 31099
diff changeset
   463
  hence "even n'" unfolding even_Suc by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   464
  have "arctan x \<le> x * ub_arctan_horner prec (get_odd n) 1 (x * x)"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   465
    unfolding `Suc n' = get_odd n`[symmetric] using arctan_0_1_bounds'[OF `0 \<le> real x` `real x \<le> 1` `even n'`] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   466
  moreover
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   467
  have "x * lb_arctan_horner prec (get_even n) 1 (x * x) \<le> arctan x"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   468
    unfolding get_even_def if_P[OF True] using arctan_0_1_bounds'[OF `0 \<le> real x` `real x \<le> 1` `even n`] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   469
  ultimately show ?thesis by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   470
next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   471
  case False hence "0 < n" by (rule odd_pos)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   472
  from gr0_implies_Suc[OF this] obtain n' where "n = Suc n'" ..
31148
7ba7c1f8bc22 Cleaned up Parity a little
nipkow
parents: 31099
diff changeset
   473
  from False[unfolded this even_Suc]
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   474
  have "even n'" and "even (Suc (Suc n'))" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   475
  have "get_odd n = Suc n'" unfolding get_odd_def if_P[OF False] using `n = Suc n'` .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   476
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   477
  have "arctan x \<le> x * ub_arctan_horner prec (get_odd n) 1 (x * x)"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   478
    unfolding `get_odd n = Suc n'` using arctan_0_1_bounds'[OF `0 \<le> real x` `real x \<le> 1` `even n'`] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   479
  moreover
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   480
  have "(x * lb_arctan_horner prec (get_even n) 1 (x * x)) \<le> arctan x"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   481
    unfolding get_even_def if_not_P[OF False] unfolding `n = Suc n'` using arctan_0_1_bounds'[OF `0 \<le> real x` `real x \<le> 1` `even (Suc (Suc n'))`] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   482
  ultimately show ?thesis by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   483
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   484
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   485
subsection "Compute \<pi>"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   486
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   487
definition ub_pi :: "nat \<Rightarrow> float" where
31809
hoelzl
parents: 31790
diff changeset
   488
  "ub_pi prec = (let A = rapprox_rat prec 1 5 ;
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   489
                     B = lapprox_rat prec 1 239
31809
hoelzl
parents: 31790
diff changeset
   490
                 in ((Float 1 2) * ((Float 1 2) * A * (ub_arctan_horner prec (get_odd (prec div 4 + 1)) 1 (A * A)) -
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   491
                                                  B * (lb_arctan_horner prec (get_even (prec div 14 + 1)) 1 (B * B)))))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   492
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   493
definition lb_pi :: "nat \<Rightarrow> float" where
31809
hoelzl
parents: 31790
diff changeset
   494
  "lb_pi prec = (let A = lapprox_rat prec 1 5 ;
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   495
                     B = rapprox_rat prec 1 239
31809
hoelzl
parents: 31790
diff changeset
   496
                 in ((Float 1 2) * ((Float 1 2) * A * (lb_arctan_horner prec (get_even (prec div 4 + 1)) 1 (A * A)) -
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   497
                                                  B * (ub_arctan_horner prec (get_odd (prec div 14 + 1)) 1 (B * B)))))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   498
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   499
lemma pi_boundaries: "pi \<in> {(lb_pi n) .. (ub_pi n)}"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   500
proof -
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   501
  have machin_pi: "pi = 4 * (4 * arctan (1 / 5) - arctan (1 / 239))" unfolding machin[symmetric] by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   502
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   503
  { fix prec n :: nat fix k :: int assume "1 < k" hence "0 \<le> k" and "0 < k" and "1 \<le> k" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   504
    let ?k = "rapprox_rat prec 1 k"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   505
    have "1 div k = 0" using div_pos_pos_trivial[OF _ `1 < k`] by auto
31809
hoelzl
parents: 31790
diff changeset
   506
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   507
    have "0 \<le> real ?k" by (rule order_trans[OF _ rapprox_rat], auto simp add: `0 \<le> k`)
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   508
    have "real ?k \<le> 1" unfolding rapprox_rat.simps(2)[OF zero_le_one `0 < k`]
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   509
      by (rule rapprox_posrat_le1, auto simp add: `0 < k` `1 \<le> k`)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   510
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   511
    have "1 / k \<le> ?k" using rapprox_rat[where x=1 and y=k] by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   512
    hence "arctan (1 / k) \<le> arctan ?k" by (rule arctan_monotone')
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   513
    also have "\<dots> \<le> (?k * ub_arctan_horner prec (get_odd n) 1 (?k * ?k))"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   514
      using arctan_0_1_bounds[OF `0 \<le> real ?k` `real ?k \<le> 1`] by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   515
    finally have "arctan (1 / k) \<le> ?k * ub_arctan_horner prec (get_odd n) 1 (?k * ?k)" .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   516
  } note ub_arctan = this
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   517
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   518
  { fix prec n :: nat fix k :: int assume "1 < k" hence "0 \<le> k" and "0 < k" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   519
    let ?k = "lapprox_rat prec 1 k"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   520
    have "1 div k = 0" using div_pos_pos_trivial[OF _ `1 < k`] by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   521
    have "1 / k \<le> 1" using `1 < k` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   522
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   523
    have "\<And>n. 0 \<le> real ?k" using lapprox_rat_bottom[where x=1 and y=k, OF zero_le_one `0 < k`] by (auto simp add: `1 div k = 0`)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   524
    have "\<And>n. real ?k \<le> 1" using lapprox_rat by (rule order_trans, auto simp add: `1 / k \<le> 1`)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   525
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   526
    have "?k \<le> 1 / k" using lapprox_rat[where x=1 and y=k] by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   527
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   528
    have "?k * lb_arctan_horner prec (get_even n) 1 (?k * ?k) \<le> arctan ?k"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   529
      using arctan_0_1_bounds[OF `0 \<le> real ?k` `real ?k \<le> 1`] by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   530
    also have "\<dots> \<le> arctan (1 / k)" using `?k \<le> 1 / k` by (rule arctan_monotone')
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   531
    finally have "?k * lb_arctan_horner prec (get_even n) 1 (?k * ?k) \<le> arctan (1 / k)" .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   532
  } note lb_arctan = this
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   533
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   534
  have "pi \<le> ub_pi n"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   535
    unfolding ub_pi_def machin_pi Let_def real_of_float_mult real_of_float_sub unfolding Float_num
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   536
    using lb_arctan[of 239] ub_arctan[of 5]
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   537
    by (auto intro!: mult_left_mono add_mono simp add: diff_minus simp del: lapprox_rat.simps rapprox_rat.simps)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   538
  moreover
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   539
  have "lb_pi n \<le> pi"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   540
    unfolding lb_pi_def machin_pi Let_def real_of_float_mult real_of_float_sub Float_num
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   541
    using lb_arctan[of 5] ub_arctan[of 239]
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   542
    by (auto intro!: mult_left_mono add_mono simp add: diff_minus simp del: lapprox_rat.simps rapprox_rat.simps)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   543
  ultimately show ?thesis by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   544
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   545
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   546
subsection "Compute arcus tangens in the entire domain"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   547
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   548
function lb_arctan :: "nat \<Rightarrow> float \<Rightarrow> float" and ub_arctan :: "nat \<Rightarrow> float \<Rightarrow> float" where
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   549
  "lb_arctan prec x = (let ub_horner = \<lambda> x. x * ub_arctan_horner prec (get_odd (prec div 4 + 1)) 1 (x * x) ;
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   550
                           lb_horner = \<lambda> x. x * lb_arctan_horner prec (get_even (prec div 4 + 1)) 1 (x * x)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   551
    in (if x < 0          then - ub_arctan prec (-x) else
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   552
        if x \<le> Float 1 -1 then lb_horner x else
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   553
        if x \<le> Float 1 1  then Float 1 1 * lb_horner (float_divl prec x (1 + ub_sqrt prec (1 + x * x)))
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   554
                          else (let inv = float_divr prec 1 x
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   555
                                in if inv > 1 then 0
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   556
                                              else lb_pi prec * Float 1 -1 - ub_horner inv)))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   557
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   558
| "ub_arctan prec x = (let lb_horner = \<lambda> x. x * lb_arctan_horner prec (get_even (prec div 4 + 1)) 1 (x * x) ;
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   559
                           ub_horner = \<lambda> x. x * ub_arctan_horner prec (get_odd (prec div 4 + 1)) 1 (x * x)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   560
    in (if x < 0          then - lb_arctan prec (-x) else
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   561
        if x \<le> Float 1 -1 then ub_horner x else
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   562
        if x \<le> Float 1 1  then let y = float_divr prec x (1 + lb_sqrt prec (1 + x * x))
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   563
                               in if y > 1 then ub_pi prec * Float 1 -1
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   564
                                           else Float 1 1 * ub_horner y
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   565
                          else ub_pi prec * Float 1 -1 - lb_horner (float_divl prec 1 x)))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   566
by pat_completeness auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   567
termination by (relation "measure (\<lambda> v. let (prec, x) = sum_case id id v in (if x < 0 then 1 else 0))", auto simp add: less_float_def)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   568
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   569
declare ub_arctan_horner.simps[simp del]
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   570
declare lb_arctan_horner.simps[simp del]
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   571
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   572
lemma lb_arctan_bound': assumes "0 \<le> real x"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   573
  shows "lb_arctan prec x \<le> arctan x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   574
proof -
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   575
  have "\<not> x < 0" and "0 \<le> x" unfolding less_float_def le_float_def using `0 \<le> real x` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   576
  let "?ub_horner x" = "x * ub_arctan_horner prec (get_odd (prec div 4 + 1)) 1 (x * x)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   577
    and "?lb_horner x" = "x * lb_arctan_horner prec (get_even (prec div 4 + 1)) 1 (x * x)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   578
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   579
  show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   580
  proof (cases "x \<le> Float 1 -1")
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   581
    case True hence "real x \<le> 1" unfolding le_float_def Float_num by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   582
    show ?thesis unfolding lb_arctan.simps Let_def if_not_P[OF `\<not> x < 0`] if_P[OF True]
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   583
      using arctan_0_1_bounds[OF `0 \<le> real x` `real x \<le> 1`] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   584
  next
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   585
    case False hence "0 < real x" unfolding le_float_def Float_num by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   586
    let ?R = "1 + sqrt (1 + real x * real x)"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   587
    let ?fR = "1 + ub_sqrt prec (1 + x * x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   588
    let ?DIV = "float_divl prec x ?fR"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   589
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   590
    have sqr_ge0: "0 \<le> 1 + real x * real x" using sum_power2_ge_zero[of 1 "real x", unfolded numeral_2_eq_2] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   591
    hence divisor_gt0: "0 < ?R" by (auto intro: add_pos_nonneg)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   592
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   593
    have "sqrt (1 + x * x) \<le> ub_sqrt prec (1 + x * x)"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   594
      using bnds_sqrt'[of "1 + x * x"] by auto
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   595
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   596
    hence "?R \<le> ?fR" by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   597
    hence "0 < ?fR" and "0 < real ?fR" unfolding less_float_def using `0 < ?R` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   598
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   599
    have monotone: "(float_divl prec x ?fR) \<le> x / ?R"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   600
    proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   601
      have "?DIV \<le> real x / ?fR" by (rule float_divl)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   602
      also have "\<dots> \<le> x / ?R" by (rule divide_left_mono[OF `?R \<le> ?fR` `0 \<le> real x` mult_pos_pos[OF order_less_le_trans[OF divisor_gt0 `?R \<le> real ?fR`] divisor_gt0]])
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   603
      finally show ?thesis .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   604
    qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   605
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   606
    show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   607
    proof (cases "x \<le> Float 1 1")
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   608
      case True
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   609
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   610
      have "x \<le> sqrt (1 + x * x)" using real_sqrt_sum_squares_ge2[where x=1, unfolded numeral_2_eq_2] by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   611
      also have "\<dots> \<le> (ub_sqrt prec (1 + x * x))"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   612
        using bnds_sqrt'[of "1 + x * x"] by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   613
      finally have "real x \<le> ?fR" by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   614
      moreover have "?DIV \<le> real x / ?fR" by (rule float_divl)
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   615
      ultimately have "real ?DIV \<le> 1" unfolding divide_le_eq_1_pos[OF `0 < real ?fR`, symmetric] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   616
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   617
      have "0 \<le> real ?DIV" using float_divl_lower_bound[OF `0 \<le> x` `0 < ?fR`] unfolding le_float_def by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   618
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   619
      have "(Float 1 1 * ?lb_horner ?DIV) \<le> 2 * arctan (float_divl prec x ?fR)" unfolding real_of_float_mult[of "Float 1 1"] Float_num
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   620
        using arctan_0_1_bounds[OF `0 \<le> real ?DIV` `real ?DIV \<le> 1`] by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   621
      also have "\<dots> \<le> 2 * arctan (x / ?R)"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   622
        using arctan_monotone'[OF monotone] by (auto intro!: mult_left_mono)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   623
      also have "2 * arctan (x / ?R) = arctan x" using arctan_half[symmetric] unfolding numeral_2_eq_2 power_Suc2 power_0 mult_1_left .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   624
      finally show ?thesis unfolding lb_arctan.simps Let_def if_not_P[OF `\<not> x < 0`] if_not_P[OF `\<not> x \<le> Float 1 -1`] if_P[OF True] .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   625
    next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   626
      case False
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   627
      hence "2 < real x" unfolding le_float_def Float_num by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   628
      hence "1 \<le> real x" by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   629
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   630
      let "?invx" = "float_divr prec 1 x"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   631
      have "0 \<le> arctan x" using arctan_monotone'[OF `0 \<le> real x`] using arctan_tan[of 0, unfolded tan_zero] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   632
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   633
      show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   634
      proof (cases "1 < ?invx")
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   635
        case True
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   636
        show ?thesis unfolding lb_arctan.simps Let_def if_not_P[OF `\<not> x < 0`] if_not_P[OF `\<not> x \<le> Float 1 -1`] if_not_P[OF False] if_P[OF True]
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   637
          using `0 \<le> arctan x` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   638
      next
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   639
        case False
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   640
        hence "real ?invx \<le> 1" unfolding less_float_def by auto
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   641
        have "0 \<le> real ?invx" by (rule order_trans[OF _ float_divr], auto simp add: `0 \<le> real x`)
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   642
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   643
        have "1 / x \<noteq> 0" and "0 < 1 / x" using `0 < real x` by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   644
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   645
        have "arctan (1 / x) \<le> arctan ?invx" unfolding real_of_float_1[symmetric] by (rule arctan_monotone', rule float_divr)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   646
        also have "\<dots> \<le> (?ub_horner ?invx)" using arctan_0_1_bounds[OF `0 \<le> real ?invx` `real ?invx \<le> 1`] by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   647
        finally have "pi / 2 - (?ub_horner ?invx) \<le> arctan x"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   648
          using `0 \<le> arctan x` arctan_inverse[OF `1 / x \<noteq> 0`]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   649
          unfolding real_sgn_pos[OF `0 < 1 / real x`] le_diff_eq by auto
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   650
        moreover
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   651
        have "lb_pi prec * Float 1 -1 \<le> pi / 2" unfolding real_of_float_mult Float_num times_divide_eq_right mult_1_left using pi_boundaries by auto
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   652
        ultimately
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   653
        show ?thesis unfolding lb_arctan.simps Let_def if_not_P[OF `\<not> x < 0`] if_not_P[OF `\<not> x \<le> Float 1 -1`] if_not_P[OF `\<not> x \<le> Float 1 1`] if_not_P[OF False]
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   654
          by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   655
      qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   656
    qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   657
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   658
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   659
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   660
lemma ub_arctan_bound': assumes "0 \<le> real x"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   661
  shows "arctan x \<le> ub_arctan prec x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   662
proof -
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   663
  have "\<not> x < 0" and "0 \<le> x" unfolding less_float_def le_float_def using `0 \<le> real x` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   664
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   665
  let "?ub_horner x" = "x * ub_arctan_horner prec (get_odd (prec div 4 + 1)) 1 (x * x)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   666
    and "?lb_horner x" = "x * lb_arctan_horner prec (get_even (prec div 4 + 1)) 1 (x * x)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   667
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   668
  show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   669
  proof (cases "x \<le> Float 1 -1")
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   670
    case True hence "real x \<le> 1" unfolding le_float_def Float_num by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   671
    show ?thesis unfolding ub_arctan.simps Let_def if_not_P[OF `\<not> x < 0`] if_P[OF True]
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   672
      using arctan_0_1_bounds[OF `0 \<le> real x` `real x \<le> 1`] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   673
  next
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   674
    case False hence "0 < real x" unfolding le_float_def Float_num by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   675
    let ?R = "1 + sqrt (1 + real x * real x)"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   676
    let ?fR = "1 + lb_sqrt prec (1 + x * x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   677
    let ?DIV = "float_divr prec x ?fR"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   678
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   679
    have sqr_ge0: "0 \<le> 1 + real x * real x" using sum_power2_ge_zero[of 1 "real x", unfolded numeral_2_eq_2] by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   680
    hence "0 \<le> real (1 + x*x)" by auto
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   681
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   682
    hence divisor_gt0: "0 < ?R" by (auto intro: add_pos_nonneg)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   683
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   684
    have "lb_sqrt prec (1 + x * x) \<le> sqrt (1 + x * x)"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   685
      using bnds_sqrt'[of "1 + x * x"] by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   686
    hence "?fR \<le> ?R" by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   687
    have "0 < real ?fR" unfolding real_of_float_add real_of_float_1 by (rule order_less_le_trans[OF zero_less_one], auto simp add: lb_sqrt_lower_bound[OF `0 \<le> real (1 + x*x)`])
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   688
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   689
    have monotone: "x / ?R \<le> (float_divr prec x ?fR)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   690
    proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   691
      from divide_left_mono[OF `?fR \<le> ?R` `0 \<le> real x` mult_pos_pos[OF divisor_gt0 `0 < real ?fR`]]
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   692
      have "x / ?R \<le> x / ?fR" .
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   693
      also have "\<dots> \<le> ?DIV" by (rule float_divr)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   694
      finally show ?thesis .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   695
    qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   696
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   697
    show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   698
    proof (cases "x \<le> Float 1 1")
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   699
      case True
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   700
      show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   701
      proof (cases "?DIV > 1")
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   702
        case True
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   703
        have "pi / 2 \<le> ub_pi prec * Float 1 -1" unfolding real_of_float_mult Float_num times_divide_eq_right mult_1_left using pi_boundaries by auto
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   704
        from order_less_le_trans[OF arctan_ubound this, THEN less_imp_le]
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   705
        show ?thesis unfolding ub_arctan.simps Let_def if_not_P[OF `\<not> x < 0`] if_not_P[OF `\<not> x \<le> Float 1 -1`] if_P[OF `x \<le> Float 1 1`] if_P[OF True] .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   706
      next
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   707
        case False
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   708
        hence "real ?DIV \<le> 1" unfolding less_float_def by auto
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   709
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   710
        have "0 \<le> x / ?R" using `0 \<le> real x` `0 < ?R` unfolding real_0_le_divide_iff by auto
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   711
        hence "0 \<le> real ?DIV" using monotone by (rule order_trans)
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   712
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   713
        have "arctan x = 2 * arctan (x / ?R)" using arctan_half unfolding numeral_2_eq_2 power_Suc2 power_0 mult_1_left .
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   714
        also have "\<dots> \<le> 2 * arctan (?DIV)"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   715
          using arctan_monotone'[OF monotone] by (auto intro!: mult_left_mono)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   716
        also have "\<dots> \<le> (Float 1 1 * ?ub_horner ?DIV)" unfolding real_of_float_mult[of "Float 1 1"] Float_num
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   717
          using arctan_0_1_bounds[OF `0 \<le> real ?DIV` `real ?DIV \<le> 1`] by auto
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   718
        finally show ?thesis unfolding ub_arctan.simps Let_def if_not_P[OF `\<not> x < 0`] if_not_P[OF `\<not> x \<le> Float 1 -1`] if_P[OF `x \<le> Float 1 1`] if_not_P[OF False] .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   719
      qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   720
    next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   721
      case False
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   722
      hence "2 < real x" unfolding le_float_def Float_num by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   723
      hence "1 \<le> real x" by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   724
      hence "0 < real x" by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   725
      hence "0 < x" unfolding less_float_def by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   726
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   727
      let "?invx" = "float_divl prec 1 x"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   728
      have "0 \<le> arctan x" using arctan_monotone'[OF `0 \<le> real x`] using arctan_tan[of 0, unfolded tan_zero] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   729
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   730
      have "real ?invx \<le> 1" unfolding less_float_def by (rule order_trans[OF float_divl], auto simp add: `1 \<le> real x` divide_le_eq_1_pos[OF `0 < real x`])
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   731
      have "0 \<le> real ?invx" unfolding real_of_float_0[symmetric] by (rule float_divl_lower_bound[unfolded le_float_def], auto simp add: `0 < x`)
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
   732
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   733
      have "1 / x \<noteq> 0" and "0 < 1 / x" using `0 < real x` by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   734
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   735
      have "(?lb_horner ?invx) \<le> arctan (?invx)" using arctan_0_1_bounds[OF `0 \<le> real ?invx` `real ?invx \<le> 1`] by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   736
      also have "\<dots> \<le> arctan (1 / x)" unfolding real_of_float_1[symmetric] by (rule arctan_monotone', rule float_divl)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   737
      finally have "arctan x \<le> pi / 2 - (?lb_horner ?invx)"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   738
        using `0 \<le> arctan x` arctan_inverse[OF `1 / x \<noteq> 0`]
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   739
        unfolding real_sgn_pos[OF `0 < 1 / x`] le_diff_eq by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   740
      moreover
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   741
      have "pi / 2 \<le> ub_pi prec * Float 1 -1" unfolding real_of_float_mult Float_num times_divide_eq_right mult_1_right using pi_boundaries by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   742
      ultimately
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   743
      show ?thesis unfolding ub_arctan.simps Let_def if_not_P[OF `\<not> x < 0`] if_not_P[OF `\<not> x \<le> Float 1 -1`] if_not_P[OF `\<not> x \<le> Float 1 1`] if_not_P[OF False]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   744
        by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   745
    qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   746
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   747
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   748
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   749
lemma arctan_boundaries:
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   750
  "arctan x \<in> {(lb_arctan prec x) .. (ub_arctan prec x)}"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   751
proof (cases "0 \<le> x")
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   752
  case True hence "0 \<le> real x" unfolding le_float_def by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   753
  show ?thesis using ub_arctan_bound'[OF `0 \<le> real x`] lb_arctan_bound'[OF `0 \<le> real x`] unfolding atLeastAtMost_iff by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   754
next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   755
  let ?mx = "-x"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   756
  case False hence "x < 0" and "0 \<le> real ?mx" unfolding le_float_def less_float_def by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   757
  hence bounds: "lb_arctan prec ?mx \<le> arctan ?mx \<and> arctan ?mx \<le> ub_arctan prec ?mx"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   758
    using ub_arctan_bound'[OF `0 \<le> real ?mx`] lb_arctan_bound'[OF `0 \<le> real ?mx`] by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   759
  show ?thesis unfolding real_of_float_minus arctan_minus lb_arctan.simps[where x=x] ub_arctan.simps[where x=x] Let_def if_P[OF `x < 0`]
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   760
    unfolding atLeastAtMost_iff using bounds[unfolded real_of_float_minus arctan_minus] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   761
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   762
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   763
lemma bnds_arctan: "\<forall> (x::real) lx ux. (l, u) = (lb_arctan prec lx, ub_arctan prec ux) \<and> x \<in> {lx .. ux} \<longrightarrow> l \<le> arctan x \<and> arctan x \<le> u"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   764
proof (rule allI, rule allI, rule allI, rule impI)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   765
  fix x :: real fix lx ux
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   766
  assume "(l, u) = (lb_arctan prec lx, ub_arctan prec ux) \<and> x \<in> {lx .. ux}"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   767
  hence l: "lb_arctan prec lx = l " and u: "ub_arctan prec ux = u" and x: "x \<in> {lx .. ux}" by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   768
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   769
  { from arctan_boundaries[of lx prec, unfolded l]
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   770
    have "l \<le> arctan lx" by (auto simp del: lb_arctan.simps)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   771
    also have "\<dots> \<le> arctan x" using x by (auto intro: arctan_monotone')
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   772
    finally have "l \<le> arctan x" .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   773
  } moreover
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   774
  { have "arctan x \<le> arctan ux" using x by (auto intro: arctan_monotone')
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   775
    also have "\<dots> \<le> u" using arctan_boundaries[of ux prec, unfolded u] by (auto simp del: ub_arctan.simps)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   776
    finally have "arctan x \<le> u" .
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   777
  } ultimately show "l \<le> arctan x \<and> arctan x \<le> u" ..
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   778
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   779
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   780
section "Sinus and Cosinus"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   781
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   782
subsection "Compute the cosinus and sinus series"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   783
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   784
fun ub_sin_cos_aux :: "nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> float \<Rightarrow> float"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   785
and lb_sin_cos_aux :: "nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> float \<Rightarrow> float" where
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   786
  "ub_sin_cos_aux prec 0 i k x = 0"
31809
hoelzl
parents: 31790
diff changeset
   787
| "ub_sin_cos_aux prec (Suc n) i k x =
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   788
    (rapprox_rat prec 1 k) - x * (lb_sin_cos_aux prec n (i + 2) (k * i * (i + 1)) x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   789
| "lb_sin_cos_aux prec 0 i k x = 0"
31809
hoelzl
parents: 31790
diff changeset
   790
| "lb_sin_cos_aux prec (Suc n) i k x =
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   791
    (lapprox_rat prec 1 k) - x * (ub_sin_cos_aux prec n (i + 2) (k * i * (i + 1)) x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   792
lemma cos_aux:
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   793
  shows "(lb_sin_cos_aux prec n 1 1 (x * x)) \<le> (\<Sum> i=0..<n. -1^i * (1/real (fact (2 * i))) * x ^(2 * i))" (is "?lb")
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   794
  and "(\<Sum> i=0..<n. -1^i * (1/real (fact (2 * i))) * x^(2 * i)) \<le> (ub_sin_cos_aux prec n 1 1 (x * x))" (is "?ub")
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   795
proof -
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   796
  have "0 \<le> real (x * x)" unfolding real_of_float_mult by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   797
  let "?f n" = "fact (2 * n)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   798
31809
hoelzl
parents: 31790
diff changeset
   799
  { fix n
30971
7fbebf75b3ef funpow and relpow with shared "^^" syntax
haftmann
parents: 30968
diff changeset
   800
    have F: "\<And>m. ((\<lambda>i. i + 2) ^^ n) m = m + 2 * n" by (induct n arbitrary: m, auto)
7fbebf75b3ef funpow and relpow with shared "^^" syntax
haftmann
parents: 30968
diff changeset
   801
    have "?f (Suc n) = ?f n * ((\<lambda>i. i + 2) ^^ n) 1 * (((\<lambda>i. i + 2) ^^ n) 1 + 1)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   802
      unfolding F by auto } note f_eq = this
31809
hoelzl
parents: 31790
diff changeset
   803
hoelzl
parents: 31790
diff changeset
   804
  from horner_bounds[where lb="lb_sin_cos_aux prec" and ub="ub_sin_cos_aux prec" and j'=0,
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   805
    OF `0 \<le> real (x * x)` f_eq lb_sin_cos_aux.simps ub_sin_cos_aux.simps]
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   806
  show "?lb" and "?ub" by (auto simp add: power_mult power2_eq_square[of "real x"])
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   807
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   808
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   809
lemma cos_boundaries: assumes "0 \<le> real x" and "x \<le> pi / 2"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   810
  shows "cos x \<in> {(lb_sin_cos_aux prec (get_even n) 1 1 (x * x)) .. (ub_sin_cos_aux prec (get_odd n) 1 1 (x * x))}"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   811
proof (cases "real x = 0")
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   812
  case False hence "real x \<noteq> 0" by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   813
  hence "0 < x" and "0 < real x" using `0 \<le> real x` unfolding less_float_def by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   814
  have "0 < x * x" using `0 < x` unfolding less_float_def real_of_float_mult real_of_float_0
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   815
    using mult_pos_pos[where a="real x" and b="real x"] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   816
30952
7ab2716dd93b power operation on functions with syntax o^; power operation on relations with syntax ^^
haftmann
parents: 30886
diff changeset
   817
  { fix x n have "(\<Sum> i=0..<n. -1^i * (1/real (fact (2 * i))) * x ^ (2 * i))
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   818
    = (\<Sum> i = 0 ..< 2 * n. (if even(i) then (-1 ^ (i div 2))/(real (fact i)) else 0) * x ^ i)" (is "?sum = ?ifsum")
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   819
  proof -
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   820
    have "?sum = ?sum + (\<Sum> j = 0 ..< n. 0)" by auto
31809
hoelzl
parents: 31790
diff changeset
   821
    also have "\<dots> =
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   822
      (\<Sum> j = 0 ..< n. -1 ^ ((2 * j) div 2) / (real (fact (2 * j))) * x ^(2 * j)) + (\<Sum> j = 0 ..< n. 0)" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   823
    also have "\<dots> = (\<Sum> i = 0 ..< 2 * n. if even i then -1 ^ (i div 2) / (real (fact i)) * x ^ i else 0)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   824
      unfolding sum_split_even_odd ..
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   825
    also have "\<dots> = (\<Sum> i = 0 ..< 2 * n. (if even i then -1 ^ (i div 2) / (real (fact i)) else 0) * x ^ i)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   826
      by (rule setsum_cong2) auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   827
    finally show ?thesis by assumption
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   828
  qed } note morph_to_if_power = this
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   829
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   830
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   831
  { fix n :: nat assume "0 < n"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   832
    hence "0 < 2 * n" by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   833
    obtain t where "0 < t" and "t < real x" and
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   834
      cos_eq: "cos x = (\<Sum> i = 0 ..< 2 * n. (if even(i) then (-1 ^ (i div 2))/(real (fact i)) else 0) * (real x) ^ i)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   835
      + (cos (t + 1/2 * (2 * n) * pi) / real (fact (2*n))) * (real x)^(2*n)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   836
      (is "_ = ?SUM + ?rest / ?fact * ?pow")
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   837
      using Maclaurin_cos_expansion2[OF `0 < real x` `0 < 2 * n`] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   838
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   839
    have "cos t * -1^n = cos t * cos (n * pi) + sin t * sin (n * pi)" by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   840
    also have "\<dots> = cos (t + n * pi)"  using cos_add by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   841
    also have "\<dots> = ?rest" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   842
    finally have "cos t * -1^n = ?rest" .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   843
    moreover
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   844
    have "t \<le> pi / 2" using `t < real x` and `x \<le> pi / 2` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   845
    hence "0 \<le> cos t" using `0 < t` and cos_ge_zero by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   846
    ultimately have even: "even n \<Longrightarrow> 0 \<le> ?rest" and odd: "odd n \<Longrightarrow> 0 \<le> - ?rest " by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   847
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   848
    have "0 < ?fact" by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   849
    have "0 < ?pow" using `0 < real x` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   850
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   851
    {
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   852
      assume "even n"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   853
      have "(lb_sin_cos_aux prec n 1 1 (x * x)) \<le> ?SUM"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   854
        unfolding morph_to_if_power[symmetric] using cos_aux by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   855
      also have "\<dots> \<le> cos x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   856
      proof -
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   857
        from even[OF `even n`] `0 < ?fact` `0 < ?pow`
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   858
        have "0 \<le> (?rest / ?fact) * ?pow" by (metis mult_nonneg_nonneg divide_nonneg_pos less_imp_le)
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   859
        thus ?thesis unfolding cos_eq by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   860
      qed
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   861
      finally have "(lb_sin_cos_aux prec n 1 1 (x * x)) \<le> cos x" .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   862
    } note lb = this
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   863
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   864
    {
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   865
      assume "odd n"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   866
      have "cos x \<le> ?SUM"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   867
      proof -
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   868
        from `0 < ?fact` and `0 < ?pow` and odd[OF `odd n`]
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   869
        have "0 \<le> (- ?rest) / ?fact * ?pow"
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   870
          by (metis mult_nonneg_nonneg divide_nonneg_pos less_imp_le)
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   871
        thus ?thesis unfolding cos_eq by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   872
      qed
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   873
      also have "\<dots> \<le> (ub_sin_cos_aux prec n 1 1 (x * x))"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   874
        unfolding morph_to_if_power[symmetric] using cos_aux by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   875
      finally have "cos x \<le> (ub_sin_cos_aux prec n 1 1 (x * x))" .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   876
    } note ub = this and lb
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   877
  } note ub = this(1) and lb = this(2)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   878
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   879
  have "cos x \<le> (ub_sin_cos_aux prec (get_odd n) 1 1 (x * x))" using ub[OF odd_pos[OF get_odd] get_odd] .
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   880
  moreover have "(lb_sin_cos_aux prec (get_even n) 1 1 (x * x)) \<le> cos x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   881
  proof (cases "0 < get_even n")
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   882
    case True show ?thesis using lb[OF True get_even] .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   883
  next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   884
    case False
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   885
    hence "get_even n = 0" by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   886
    have "- (pi / 2) \<le> x" by (rule order_trans[OF _ `0 < real x`[THEN less_imp_le]], auto)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   887
    with `x \<le> pi / 2`
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   888
    show ?thesis unfolding `get_even n = 0` lb_sin_cos_aux.simps real_of_float_minus real_of_float_0 using cos_ge_zero by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   889
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   890
  ultimately show ?thesis by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   891
next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   892
  case True
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   893
  show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   894
  proof (cases "n = 0")
31809
hoelzl
parents: 31790
diff changeset
   895
    case True
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   896
    thus ?thesis unfolding `n = 0` get_even_def get_odd_def using `real x = 0` lapprox_rat[where x="-1" and y=1] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   897
  next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   898
    case False with not0_implies_Suc obtain m where "n = Suc m" by blast
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   899
    thus ?thesis unfolding `n = Suc m` get_even_def get_odd_def using `real x = 0` rapprox_rat[where x=1 and y=1] lapprox_rat[where x=1 and y=1] by (cases "even (Suc m)", auto)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   900
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   901
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   902
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   903
lemma sin_aux: assumes "0 \<le> real x"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   904
  shows "(x * lb_sin_cos_aux prec n 2 1 (x * x)) \<le> (\<Sum> i=0..<n. -1^i * (1/real (fact (2 * i + 1))) * x^(2 * i + 1))" (is "?lb")
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   905
  and "(\<Sum> i=0..<n. -1^i * (1/real (fact (2 * i + 1))) * x^(2 * i + 1)) \<le> (x * ub_sin_cos_aux prec n 2 1 (x * x))" (is "?ub")
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   906
proof -
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   907
  have "0 \<le> real (x * x)" unfolding real_of_float_mult by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   908
  let "?f n" = "fact (2 * n + 1)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   909
31809
hoelzl
parents: 31790
diff changeset
   910
  { fix n
30971
7fbebf75b3ef funpow and relpow with shared "^^" syntax
haftmann
parents: 30968
diff changeset
   911
    have F: "\<And>m. ((\<lambda>i. i + 2) ^^ n) m = m + 2 * n" by (induct n arbitrary: m, auto)
7fbebf75b3ef funpow and relpow with shared "^^" syntax
haftmann
parents: 30968
diff changeset
   912
    have "?f (Suc n) = ?f n * ((\<lambda>i. i + 2) ^^ n) 2 * (((\<lambda>i. i + 2) ^^ n) 2 + 1)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   913
      unfolding F by auto } note f_eq = this
31809
hoelzl
parents: 31790
diff changeset
   914
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   915
  from horner_bounds[where lb="lb_sin_cos_aux prec" and ub="ub_sin_cos_aux prec" and j'=0,
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   916
    OF `0 \<le> real (x * x)` f_eq lb_sin_cos_aux.simps ub_sin_cos_aux.simps]
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   917
  show "?lb" and "?ub" using `0 \<le> real x` unfolding real_of_float_mult
36778
739a9379e29b avoid using real-specific versions of generic lemmas
huffman
parents: 36531
diff changeset
   918
    unfolding power_add power_one_right mult_assoc[symmetric] setsum_left_distrib[symmetric]
739a9379e29b avoid using real-specific versions of generic lemmas
huffman
parents: 36531
diff changeset
   919
    unfolding mult_commute[where 'a=real]
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   920
    by (auto intro!: mult_left_mono simp add: power_mult power2_eq_square[of "real x"])
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   921
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   922
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   923
lemma sin_boundaries: assumes "0 \<le> real x" and "x \<le> pi / 2"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   924
  shows "sin x \<in> {(x * lb_sin_cos_aux prec (get_even n) 2 1 (x * x)) .. (x * ub_sin_cos_aux prec (get_odd n) 2 1 (x * x))}"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   925
proof (cases "real x = 0")
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   926
  case False hence "real x \<noteq> 0" by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   927
  hence "0 < x" and "0 < real x" using `0 \<le> real x` unfolding less_float_def by auto
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   928
  have "0 < x * x" using `0 < x` unfolding less_float_def real_of_float_mult real_of_float_0
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   929
    using mult_pos_pos[where a="real x" and b="real x"] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   930
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   931
  { fix x n have "(\<Sum> j = 0 ..< n. -1 ^ (((2 * j + 1) - Suc 0) div 2) / (real (fact (2 * j + 1))) * x ^(2 * j + 1))
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   932
    = (\<Sum> i = 0 ..< 2 * n. (if even(i) then 0 else (-1 ^ ((i - Suc 0) div 2))/(real (fact i))) * x ^ i)" (is "?SUM = _")
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   933
    proof -
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   934
      have pow: "!!i. x ^ (2 * i + 1) = x * x ^ (2 * i)" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   935
      have "?SUM = (\<Sum> j = 0 ..< n. 0) + ?SUM" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   936
      also have "\<dots> = (\<Sum> i = 0 ..< 2 * n. if even i then 0 else -1 ^ ((i - Suc 0) div 2) / (real (fact i)) * x ^ i)"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   937
        unfolding sum_split_even_odd ..
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   938
      also have "\<dots> = (\<Sum> i = 0 ..< 2 * n. (if even i then 0 else -1 ^ ((i - Suc 0) div 2) / (real (fact i))) * x ^ i)"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   939
        by (rule setsum_cong2) auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   940
      finally show ?thesis by assumption
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   941
    qed } note setsum_morph = this
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   942
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   943
  { fix n :: nat assume "0 < n"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   944
    hence "0 < 2 * n + 1" by auto
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   945
    obtain t where "0 < t" and "t < real x" and
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   946
      sin_eq: "sin x = (\<Sum> i = 0 ..< 2 * n + 1. (if even(i) then 0 else (-1 ^ ((i - Suc 0) div 2))/(real (fact i))) * (real x) ^ i)
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   947
      + (sin (t + 1/2 * (2 * n + 1) * pi) / real (fact (2*n + 1))) * (real x)^(2*n + 1)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   948
      (is "_ = ?SUM + ?rest / ?fact * ?pow")
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   949
      using Maclaurin_sin_expansion3[OF `0 < 2 * n + 1` `0 < real x`] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   950
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   951
    have "?rest = cos t * -1^n" unfolding sin_add cos_add real_of_nat_add left_distrib right_distrib by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   952
    moreover
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   953
    have "t \<le> pi / 2" using `t < real x` and `x \<le> pi / 2` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   954
    hence "0 \<le> cos t" using `0 < t` and cos_ge_zero by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   955
    ultimately have even: "even n \<Longrightarrow> 0 \<le> ?rest" and odd: "odd n \<Longrightarrow> 0 \<le> - ?rest " by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   956
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   957
    have "0 < ?fact" by (rule real_of_nat_fact_gt_zero)
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   958
    have "0 < ?pow" using `0 < real x` by (rule zero_less_power)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   959
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   960
    {
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   961
      assume "even n"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   962
      have "(x * lb_sin_cos_aux prec n 2 1 (x * x)) \<le>
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   963
            (\<Sum> i = 0 ..< 2 * n. (if even(i) then 0 else (-1 ^ ((i - Suc 0) div 2))/(real (fact i))) * (real x) ^ i)"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   964
        using sin_aux[OF `0 \<le> real x`] unfolding setsum_morph[symmetric] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   965
      also have "\<dots> \<le> ?SUM" by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   966
      also have "\<dots> \<le> sin x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   967
      proof -
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   968
        from even[OF `even n`] `0 < ?fact` `0 < ?pow`
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   969
        have "0 \<le> (?rest / ?fact) * ?pow" by (metis mult_nonneg_nonneg divide_nonneg_pos less_imp_le)
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   970
        thus ?thesis unfolding sin_eq by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   971
      qed
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   972
      finally have "(x * lb_sin_cos_aux prec n 2 1 (x * x)) \<le> sin x" .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   973
    } note lb = this
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   974
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   975
    {
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   976
      assume "odd n"
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   977
      have "sin x \<le> ?SUM"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   978
      proof -
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   979
        from `0 < ?fact` and `0 < ?pow` and odd[OF `odd n`]
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   980
        have "0 \<le> (- ?rest) / ?fact * ?pow"
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   981
          by (metis mult_nonneg_nonneg divide_nonneg_pos less_imp_le)
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   982
        thus ?thesis unfolding sin_eq by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   983
      qed
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
   984
      also have "\<dots> \<le> (\<Sum> i = 0 ..< 2 * n. (if even(i) then 0 else (-1 ^ ((i - Suc 0) div 2))/(real (fact i))) * (real x) ^ i)"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   985
         by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   986
      also have "\<dots> \<le> (x * ub_sin_cos_aux prec n 2 1 (x * x))"
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
   987
        using sin_aux[OF `0 \<le> real x`] unfolding setsum_morph[symmetric] by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   988
      finally have "sin x \<le> (x * ub_sin_cos_aux prec n 2 1 (x * x))" .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   989
    } note ub = this and lb
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   990
  } note ub = this(1) and lb = this(2)
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   991
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   992
  have "sin x \<le> (x * ub_sin_cos_aux prec (get_odd n) 2 1 (x * x))" using ub[OF odd_pos[OF get_odd] get_odd] .
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   993
  moreover have "(x * lb_sin_cos_aux prec (get_even n) 2 1 (x * x)) \<le> sin x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   994
  proof (cases "0 < get_even n")
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   995
    case True show ?thesis using lb[OF True get_even] .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   996
  next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   997
    case False
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
   998
    hence "get_even n = 0" by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
   999
    with `x \<le> pi / 2` `0 \<le> real x`
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
  1000
    show ?thesis unfolding `get_even n = 0` ub_sin_cos_aux.simps real_of_float_minus real_of_float_0 using sin_ge_zero by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1001
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1002
  ultimately show ?thesis by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1003
next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1004
  case True
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1005
  show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1006
  proof (cases "n = 0")
31809
hoelzl
parents: 31790
diff changeset
  1007
    case True
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
  1008
    thus ?thesis unfolding `n = 0` get_even_def get_odd_def using `real x = 0` lapprox_rat[where x="-1" and y=1] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1009
  next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1010
    case False with not0_implies_Suc obtain m where "n = Suc m" by blast
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
  1011
    thus ?thesis unfolding `n = Suc m` get_even_def get_odd_def using `real x = 0` rapprox_rat[where x=1 and y=1] lapprox_rat[where x=1 and y=1] by (cases "even (Suc m)", auto)
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1012
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1013
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1014
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1015
subsection "Compute the cosinus in the entire domain"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1016
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1017
definition lb_cos :: "nat \<Rightarrow> float \<Rightarrow> float" where
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1018
"lb_cos prec x = (let
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1019
    horner = \<lambda> x. lb_sin_cos_aux prec (get_even (prec div 4 + 1)) 1 1 (x * x) ;
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1020
    half = \<lambda> x. if x < 0 then - 1 else Float 1 1 * x * x - 1
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1021
  in if x < Float 1 -1 then horner x
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1022
else if x < 1          then half (horner (x * Float 1 -1))
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1023
                       else half (half (horner (x * Float 1 -2))))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1024
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1025
definition ub_cos :: "nat \<Rightarrow> float \<Rightarrow> float" where
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1026
"ub_cos prec x = (let
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1027
    horner = \<lambda> x. ub_sin_cos_aux prec (get_odd (prec div 4 + 1)) 1 1 (x * x) ;
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1028
    half = \<lambda> x. Float 1 1 * x * x - 1
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1029
  in if x < Float 1 -1 then horner x
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1030
else if x < 1          then half (horner (x * Float 1 -1))
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1031
                       else half (half (horner (x * Float 1 -2))))"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1032
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1033
lemma lb_cos: assumes "0 \<le> real x" and "x \<le> pi"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1034
  shows "cos x \<in> {(lb_cos prec x) .. (ub_cos prec x)}" (is "?cos x \<in> {(?lb x) .. (?ub x) }")
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1035
proof -
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1036
  { fix x :: real
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1037
    have "cos x = cos (x / 2 + x / 2)" by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1038
    also have "\<dots> = cos (x / 2) * cos (x / 2) + sin (x / 2) * sin (x / 2) - sin (x / 2) * sin (x / 2) + cos (x / 2) * cos (x / 2) - 1"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1039
      unfolding cos_add by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1040
    also have "\<dots> = 2 * cos (x / 2) * cos (x / 2) - 1" by algebra
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1041
    finally have "cos x = 2 * cos (x / 2) * cos (x / 2) - 1" .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1042
  } note x_half = this[symmetric]
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1043
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
  1044
  have "\<not> x < 0" using `0 \<le> real x` unfolding less_float_def by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1045
  let "?ub_horner x" = "ub_sin_cos_aux prec (get_odd (prec div 4 + 1)) 1 1 (x * x)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1046
  let "?lb_horner x" = "lb_sin_cos_aux prec (get_even (prec div 4 + 1)) 1 1 (x * x)"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1047
  let "?ub_half x" = "Float 1 1 * x * x - 1"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1048
  let "?lb_half x" = "if x < 0 then - 1 else Float 1 1 * x * x - 1"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1049
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1050
  show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1051
  proof (cases "x < Float 1 -1")
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1052
    case True hence "x \<le> pi / 2" unfolding less_float_def using pi_ge_two by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1053
    show ?thesis unfolding lb_cos_def[where x=x] ub_cos_def[where x=x] if_not_P[OF `\<not> x < 0`] if_P[OF `x < Float 1 -1`] Let_def
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1054
      using cos_boundaries[OF `0 \<le> real x` `x \<le> pi / 2`] .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1055
  next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1056
    case False
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1057
    { fix y x :: float let ?x2 = "(x * Float 1 -1)"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1058
      assume "y \<le> cos ?x2" and "-pi \<le> x" and "x \<le> pi"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
  1059
      hence "- (pi / 2) \<le> ?x2" and "?x2 \<le> pi / 2" using pi_ge_two unfolding real_of_float_mult Float_num by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1060
      hence "0 \<le> cos ?x2" by (rule cos_ge_zero)
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1061
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1062
      have "(?lb_half y) \<le> cos x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1063
      proof (cases "y < 0")
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1064
        case True show ?thesis using cos_ge_minus_one unfolding if_P[OF True] by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1065
      next
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1066
        case False
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1067
        hence "0 \<le> real y" unfolding less_float_def by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1068
        from mult_mono[OF `y \<le> cos ?x2` `y \<le> cos ?x2` `0 \<le> cos ?x2` this]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1069
        have "real y * real y \<le> cos ?x2 * cos ?x2" .
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1070
        hence "2 * real y * real y \<le> 2 * cos ?x2 * cos ?x2" by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1071
        hence "2 * real y * real y - 1 \<le> 2 * cos (x / 2) * cos (x / 2) - 1" unfolding Float_num real_of_float_mult by auto
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1072
        thus ?thesis unfolding if_not_P[OF False] x_half Float_num real_of_float_mult real_of_float_sub by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1073
      qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1074
    } note lb_half = this
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1075
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1076
    { fix y x :: float let ?x2 = "(x * Float 1 -1)"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1077
      assume ub: "cos ?x2 \<le> y" and "- pi \<le> x" and "x \<le> pi"
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
  1078
      hence "- (pi / 2) \<le> ?x2" and "?x2 \<le> pi / 2" using pi_ge_two unfolding real_of_float_mult Float_num by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1079
      hence "0 \<le> cos ?x2" by (rule cos_ge_zero)
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1080
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1081
      have "cos x \<le> (?ub_half y)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1082
      proof -
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1083
        have "0 \<le> real y" using `0 \<le> cos ?x2` ub by (rule order_trans)
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1084
        from mult_mono[OF ub ub this `0 \<le> cos ?x2`]
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1085
        have "cos ?x2 * cos ?x2 \<le> real y * real y" .
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1086
        hence "2 * cos ?x2 * cos ?x2 \<le> 2 * real y * real y" by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1087
        hence "2 * cos (x / 2) * cos (x / 2) - 1 \<le> 2 * real y * real y - 1" unfolding Float_num real_of_float_mult by auto
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1088
        thus ?thesis unfolding x_half real_of_float_mult Float_num real_of_float_sub by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1089
      qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1090
    } note ub_half = this
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1091
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1092
    let ?x2 = "x * Float 1 -1"
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1093
    let ?x4 = "x * Float 1 -1 * Float 1 -1"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1094
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1095
    have "-pi \<le> x" using pi_ge_zero[THEN le_imp_neg_le, unfolded minus_zero] `0 \<le> real x` by (rule order_trans)
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1096
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1097
    show ?thesis
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1098
    proof (cases "x < 1")
31098
73dd67adf90a replaced Ifloat => real_of_float and real, renamed ApproxEq => inequality, uneq => interpret_inequality, uneq' => approx_inequality, Ifloatarith => interpret_floatarith
hoelzl
parents: 30971
diff changeset
  1099
      case True hence "real x \<le> 1" unfolding less_float_def by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1100
      have "0 \<le> real ?x2" and "?x2 \<le> pi / 2" using pi_ge_two `0 \<le> real x` unfolding real_of_float_mult Float_num using assms by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1101
      from cos_boundaries[OF this]
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1102
      have lb: "(?lb_horner ?x2) \<le> ?cos ?x2" and ub: "?cos ?x2 \<le> (?ub_horner ?x2)" by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1103
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1104
      have "(?lb x) \<le> ?cos x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1105
      proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1106
        from lb_half[OF lb `-pi \<le> x` `x \<le> pi`]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1107
        show ?thesis unfolding lb_cos_def[where x=x] Let_def using `\<not> x < 0` `\<not> x < Float 1 -1` `x < 1` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1108
      qed
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1109
      moreover have "?cos x \<le> (?ub x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1110
      proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1111
        from ub_half[OF ub `-pi \<le> x` `x \<le> pi`]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1112
        show ?thesis unfolding ub_cos_def[where x=x] Let_def using `\<not> x < 0` `\<not> x < Float 1 -1` `x < 1` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1113
      qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1114
      ultimately show ?thesis by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1115
    next
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1116
      case False
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1117
      have "0 \<le> real ?x4" and "?x4 \<le> pi / 2" using pi_ge_two `0 \<le> real x` `x \<le> pi` unfolding real_of_float_mult Float_num by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1118
      from cos_boundaries[OF this]
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1119
      have lb: "(?lb_horner ?x4) \<le> ?cos ?x4" and ub: "?cos ?x4 \<le> (?ub_horner ?x4)" by auto
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1120
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1121
      have eq_4: "?x2 * Float 1 -1 = x * Float 1 -2" by (cases x, auto simp add: times_float.simps)
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1122
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1123
      have "(?lb x) \<le> ?cos x"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1124
      proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1125
        have "-pi \<le> ?x2" and "?x2 \<le> pi" unfolding real_of_float_mult Float_num using pi_ge_two `0 \<le> real x` `x \<le> pi` by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1126
        from lb_half[OF lb_half[OF lb this] `-pi \<le> x` `x \<le> pi`, unfolded eq_4]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1127
        show ?thesis unfolding lb_cos_def[where x=x] if_not_P[OF `\<not> x < 0`] if_not_P[OF `\<not> x < Float 1 -1`] if_not_P[OF `\<not> x < 1`] Let_def .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1128
      qed
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1129
      moreover have "?cos x \<le> (?ub x)"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1130
      proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1131
        have "-pi \<le> ?x2" and "?x2 \<le> pi" unfolding real_of_float_mult Float_num using pi_ge_two `0 \<le> real x` ` x \<le> pi` by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1132
        from ub_half[OF ub_half[OF ub this] `-pi \<le> x` `x \<le> pi`, unfolded eq_4]
32960
69916a850301 eliminated hard tabulators, guessing at each author's individual tab-width;
wenzelm
parents: 32920
diff changeset
  1133
        show ?thesis unfolding ub_cos_def[where x=x] if_not_P[OF `\<not> x < 0`] if_not_P[OF `\<not> x < Float 1 -1`] if_not_P[OF `\<not> x < 1`] Let_def .
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1134
      qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1135
      ultimately show ?thesis by auto
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1136
    qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1137
  qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1138
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1139
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1140
lemma lb_cos_minus: assumes "-pi \<le> x" and "real x \<le> 0"
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1141
  shows "cos (real(-x)) \<in> {(lb_cos prec (-x)) .. (ub_cos prec (-x))}"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1142
proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1143
  have "0 \<le> real (-x)" and "(-x) \<le> pi" using `-pi \<le> x` `real x \<le> 0` by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1144
  from lb_cos[OF this] show ?thesis .
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1145
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1146
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1147
definition bnds_cos :: "nat \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float * float" where
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1148
"bnds_cos prec lx ux = (let
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1149
    lpi = round_down prec (lb_pi prec) ;
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1150
    upi = round_up prec (ub_pi prec) ;
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1151
    k = floor_fl (float_divr prec (lx + lpi) (2 * lpi)) ;
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1152
    lx = lx - k * 2 * (if k < 0 then lpi else upi) ;
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1153
    ux = ux - k * 2 * (if k < 0 then upi else lpi)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1154
  in   if - lpi \<le> lx \<and> ux \<le> 0    then (lb_cos prec (-lx), ub_cos prec (-ux))
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1155
  else if 0 \<le> lx \<and> ux \<le> lpi      then (lb_cos prec ux, ub_cos prec lx)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1156
  else if - lpi \<le> lx \<and> ux \<le> lpi  then (min (lb_cos prec (-lx)) (lb_cos prec ux), Float 1 0)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1157
  else if 0 \<le> lx \<and> ux \<le> 2 * lpi  then (Float -1 0, max (ub_cos prec lx) (ub_cos prec (- (ux - 2 * lpi))))
31508
1ea1c70aba00 Better approximation of cos around pi.
hoelzl
parents: 31468
diff changeset
  1158
  else if -2 * lpi \<le> lx \<and> ux \<le> 0 then (Float -1 0, max (ub_cos prec (lx + 2 * lpi)) (ub_cos prec (-ux)))
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1159
                                 else (Float -1 0, Float 1 0))"
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1160
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1161
lemma floor_int:
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1162
  obtains k :: int where "real k = (floor_fl f)"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1163
proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1164
  assume *: "\<And> k :: int. real k = (floor_fl f) \<Longrightarrow> thesis"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1165
  obtain m e where fl: "Float m e = floor_fl f" by (cases "floor_fl f", auto)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1166
  from floor_pos_exp[OF this]
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1167
  have "real (m* 2^(nat e)) = (floor_fl f)"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1168
    by (auto simp add: fl[symmetric] real_of_float_def pow2_def)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1169
  from *[OF this] show thesis by blast
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1170
qed
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1171
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1172
lemma float_remove_real_numeral[simp]: "(number_of k :: float) = (number_of k :: real)"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1173
proof -
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1174
  have "(number_of k :: float) = real k"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1175
    unfolding number_of_float_def real_of_float_def pow2_def by auto
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1176
  also have "\<dots> = (number_of k :: int)"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1177
    by (simp add: number_of_is_id)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1178
  finally show ?thesis by auto
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1179
qed
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1180
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1181
lemma cos_periodic_nat[simp]: fixes n :: nat shows "cos (x + n * (2 * pi)) = cos x"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1182
proof (induct n arbitrary: x)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1183
  case (Suc n)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1184
  have split_pi_off: "x + (Suc n) * (2 * pi) = (x + n * (2 * pi)) + 2 * pi"
36778
739a9379e29b avoid using real-specific versions of generic lemmas
huffman
parents: 36531
diff changeset
  1185
    unfolding Suc_eq_plus1 real_of_nat_add real_of_one left_distrib by auto
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1186
  show ?case unfolding split_pi_off using Suc by auto
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1187
qed auto
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1188
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1189
lemma cos_periodic_int[simp]: fixes i :: int shows "cos (x + i * (2 * pi)) = cos x"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1190
proof (cases "0 \<le> i")
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1191
  case True hence i_nat: "real i = nat i" by auto
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1192
  show ?thesis unfolding i_nat by auto
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1193
next
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1194
  case False hence i_nat: "i = - real (nat (-i))" by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1195
  have "cos x = cos (x + i * (2 * pi) - i * (2 * pi))" by auto
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1196
  also have "\<dots> = cos (x + i * (2 * pi))"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1197
    unfolding i_nat mult_minus_left diff_minus_eq_add by (rule cos_periodic_nat)
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1198
  finally show ?thesis by auto
29805
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1199
qed
a5da150bd0ab Add approximation method
hoelzl
parents:
diff changeset
  1200
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1201
lemma bnds_cos: "\<forall> (x::real) lx ux. (l, u) = bnds_cos prec lx ux \<and> x \<in> {lx .. ux} \<longrightarrow> l \<le> cos x \<and> cos x \<le> u"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1202
proof ((rule allI | rule impI | erule conjE) +)
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1203
  fix x :: real fix lx ux
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset
  1204
  assume bnds: "(l, u) = bnds_cos prec lx ux" and x: "x \<in> {lx .. ux}"
31467
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1205
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1206
  let ?lpi = "round_down prec (lb_pi prec)"
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1207
  let ?upi = "round_up prec (ub_pi prec)"
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1208
  let ?k = "floor_fl (float_divr prec (lx + ?lpi) (2 * ?lpi))"
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1209
  let ?lx = "lx - ?k * 2 * (if ?k < 0 then ?lpi else ?upi)"
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1210
  let ?ux = "ux - ?k * 2 * (if ?k < 0 then ?upi else ?lpi)"
f7d2aa438bee Approximation: Implemented argument reduction for cosine. Sinus is now implemented in terms of cosine. Sqrt computes on the entire real numbers
hoelzl
parents: 31148
diff changeset
  1211
40881
e84f82418e09 Use coercions in Approximation (by Dmitriy Traytel).
hoelzl
parents: 39556
diff changeset