# Theory Sqrt_Bounds

```chapter ‹Square Root Upper and Lower Bounds›

theory Sqrt_Bounds
imports Bounds_Lemmas
"HOL-Library.Sum_of_Squares"

begin

text‹Covers all bounds used in sqrt-upper.ax, sqrt-lower.ax and sqrt-extended.ax.›

section‹Upper bounds›

primrec sqrtu :: "[real,nat] ⇒ real" where
"sqrtu x 0 = (x+1)/2"
| "sqrtu x (Suc n) = (sqrtu x n + x/sqrtu x n) / 2"

lemma sqrtu_upper: "x ≤ sqrtu x n ^ 2"
proof (induction n)
case 0 show ?case
apply (sos "(((A<0 * R<1) + (R<1 * (R<1 * [~1*x + 1]^2))))")
done
next
case (Suc n)
have xy: "⋀y. ⟦x ≤ y * y; y ≠ 0⟧ ⟹ x * (2 * (y * y)) ≤ x * x + y * (y * (y * y))"
by (sos "(((((A<0 * A<1) * R<1) + ((A<0 * R<1) * (R<1 * [~1*y^2 + x]^2)))) &
((((A<0 * A<1) * R<1) + ((A<0 * R<1) * (R<1 * [~1*y^2 + x]^2)))))")
show ?case using Suc
by (auto simp: power2_eq_square algebra_simps divide_simps xy)
qed

lemma sqrtu_numeral:
"sqrtu x (numeral n) = (sqrtu x (pred_numeral n) + x/sqrtu x (pred_numeral n)) / 2"

lemma sqrtu_gt_0: "x ≥ 0 ⟹ sqrtu x n > 0"
apply (induct n)
apply (auto simp: field_simps)

theorem gen_sqrt_upper: "0 ≤ x ⟹ sqrt x ≤ sqrtu x n"
using real_sqrt_le_mono [OF sqrtu_upper [of x n]]
by auto (metis abs_of_nonneg dual_order.strict_iff_order sqrtu_gt_0)

lemma sqrt_upper_bound_0:
assumes "x ≥ 0" shows "sqrt x ≤ (x+1)/2" (is "_ ≤ ?rhs")
proof -
have "sqrt x ≤ sqrtu x 0"
by (metis assms gen_sqrt_upper)
also have "... = ?rhs"
finally show ?thesis .
qed

lemma sqrt_upper_bound_1:
assumes "x ≥ 0" shows "sqrt x ≤ (1/4)*(x^2+6*x+1) / (x+1)" (is "_ ≤ ?rhs")
proof -
have "sqrt x ≤ sqrtu x 1"
by (metis assms gen_sqrt_upper)
also have "... = ?rhs"
finally show ?thesis .
qed

lemma sqrtu_2_eq:
"sqrtu x 2 = (1/8)*(x^4 + 28*x^3 + 70*x^2 + 28*x + 1) / ((x + 1)*(x^2 + 6*x + 1))"
by (simp add: sqrtu_numeral divide_simps) algebra

lemma sqrt_upper_bound_2:
assumes "x ≥ 0"
shows "sqrt x ≤ (1/8)*(x^4 + 28*x^3 + 70*x^2 + 28*x + 1) / ((x + 1)*(x^2 + 6*x + 1))"
by (metis assms gen_sqrt_upper sqrtu_2_eq)

lemma sqrtu_4_eq:
"x ≥ 0 ⟹
sqrtu x 4 = (1/32)*(225792840*x^6+64512240*x^5+601080390*x^8+471435600*x^7+496*x+1+35960*x^2+906192*x^3+10518300*x^4+x^16+906192*x^13+10518300*x^12+496*x^15+35960*x^14+471435600*x^9+225792840*x^10+64512240*x^11)
/ ((x+1)*(x^2+6*x+1)*(x^4+28*x^3+70*x^2+28*x+1)*(1820*x^6+8008*x^5+x^8+120*x^7+120*x+1+1820*x^2+8008*x^3+12870*x^4))"

lemma sqrt_upper_bound_4:
assumes "x ≥ 0"
shows "sqrt x ≤ (1/32)*(225792840*x^6+64512240*x^5+601080390*x^8+471435600*x^7+496*x+1+35960*x^2+906192*x^3+10518300*x^4+x^16+906192*x^13+10518300*x^12+496*x^15+35960*x^14+471435600*x^9+225792840*x^10+64512240*x^11)
/ ((x+1)*(x^2+6*x+1)*(x^4+28*x^3+70*x^2+28*x+1)*(1820*x^6+8008*x^5+x^8+120*x^7+120*x+1+1820*x^2+8008*x^3+12870*x^4))"
by (metis assms gen_sqrt_upper sqrtu_4_eq)

lemma gen_sqrt_upper_scaled:
assumes "0 ≤ x" "0 < u"
shows "sqrt x ≤ sqrtu (x*u^2) n / u"
proof -
have "sqrt x = sqrt x * sqrt (u^2) / u"
using assms
by simp
also have "... = sqrt (x*u^2) / u"
by (metis real_sqrt_mult)
also have "... ≤ sqrtu (x*u^2) n / u"
using assms
by (simp add: divide_simps) (metis gen_sqrt_upper zero_le_mult_iff zero_le_power2)
finally show ?thesis .
qed

lemma sqrt_upper_bound_2_small:
assumes "0 ≤ x"
shows "sqrt x ≤ (1/32)*(65536*x^4 + 114688*x^3 + 17920*x^2 + 448*x + 1) / ((16*x + 1)*(256*x^2 + 96*x + 1))"
apply (rule order_trans [OF gen_sqrt_upper_scaled [of x 4 2] eq_refl])
using assms
apply (auto simp: sqrtu_2_eq )
apply algebra
done

lemma sqrt_upper_bound_2_large:
assumes "0 ≤ x"
shows "sqrt x ≤ (1/32)*(65536 + 114688*x + 17920*x^2 + 448*x^3 + x^4) / ((x + 16)*(256 + 96*x + x^2))"
apply (rule order_trans [OF gen_sqrt_upper_scaled [of x "1/4" 2] eq_refl])
using assms
apply (auto simp: sqrtu_2_eq )
apply algebra
done

section‹Lower bounds›

lemma sqrt_lower_bound_id:
assumes "0 ≤ x" "x ≤ 1"
shows "x ≤ sqrt x"
proof -
have "x^2 ≤ x" using assms
by (metis one_le_numeral power_decreasing power_one_right)
then show ?thesis
by (metis real_le_rsqrt)
qed

definition sqrtl :: "[real,nat] ⇒ real" where
"sqrtl x n = x/sqrtu x n"

lemma sqrtl_lower: "0 ≤ x ⟹ sqrtl x n ^ 2 ≤ x"
unfolding sqrtl_def using sqrtu_upper [of x n]
by (auto simp: power2_eq_square divide_simps mult_left_mono)

theorem gen_sqrt_lower:  "0 ≤ x ⟹ sqrtl x n ≤ sqrt x"
using real_sqrt_le_mono [OF sqrtl_lower [of x n]]
by auto

lemma sqrt_lower_bound_0:
assumes "x ≥ 0" shows "2*x/(x+1) ≤ sqrt x" (is "?lhs ≤ _")
proof -
have "?lhs = sqrtl x 0"
also have "... ≤ sqrt x"
by (metis assms gen_sqrt_lower)
finally show ?thesis .
qed

lemma sqrt_lower_bound_1:
assumes "x ≥ 0" shows "4*x*(x+1) / (x^2+6*x+1) ≤ sqrt x" (is "?lhs ≤ _")
proof -
have "?lhs = sqrtl x 1" using assms
by (simp add: sqrtl_def power_eq_if algebra_simps divide_simps)
also have "... ≤ sqrt x"
by (metis assms gen_sqrt_lower)
finally show ?thesis .
qed

lemma sqrtl_2_eq: "sqrtl x 2 =
8*x*(x + 1)*(x^2 + 6*x + 1) / (x^4 + 28*x^3 + 70*x^2 + 28*x + 1)"
using sqrtu_gt_0 [of x 2]

lemma sqrt_lower_bound_2:
assumes "x ≥ 0"
shows "8*x*(x + 1)*(x^2 + 6*x + 1) / (x^4 + 28*x^3 + 70*x^2 + 28*x + 1) ≤ sqrt x"
by (metis assms sqrtl_2_eq gen_sqrt_lower)

lemma sqrtl_4_eq: "x ≥ 0 ⟹
sqrtl x 4
= (32*x*(x+1)*(x^2+6*x+1)*(x^4+28*x^3+70*x^2+28*x+1)*(1820*x^6+8008*x^5+x^8+120*x^7+120*x+1+1820*x^2+8008*x^3+12870*x^4))
/ (225792840*x^6+64512240*x^5+601080390*x^8+471435600*x^7+496*x+1+35960*x^2+906192*x^3+10518300*x^4+x^16+906192*x^13+10518300*x^12+496*x^15+35960*x^14+471435600*x^9+225792840*x^10+64512240*x^11)"
using sqrtu_gt_0 [of x 4]

lemma sqrt_lower_bound_4:
assumes "x ≥ 0"
shows "(32*x*(x+1)*(x^2+6*x+1)*(x^4+28*x^3+70*x^2+28*x+1)*(1820*x^6+8008*x^5+x^8+120*x^7+120*x+1+1820*x^2+8008*x^3+12870*x^4))
/ (225792840*x^6+64512240*x^5+601080390*x^8+471435600*x^7+496*x+1+35960*x^2+906192*x^3+10518300*x^4+x^16+906192*x^13+10518300*x^12+496*x^15+35960*x^14+471435600*x^9+225792840*x^10+64512240*x^11)
≤ sqrt x"
by (metis assms sqrtl_4_eq gen_sqrt_lower)

lemma gen_sqrt_lower_scaled:
assumes "0 ≤ x" "0 < u"
shows "sqrtl (x*u^2) n / u ≤ sqrt x"
proof -
have "sqrtl (x*u^2) n / u ≤ sqrt (x*u^2) / u"
using assms
by (simp add: divide_simps) (metis gen_sqrt_lower zero_le_mult_iff zero_le_power2)
also have "... = sqrt x * sqrt (u^2) / u"
by (metis real_sqrt_mult)
also have "... = sqrt x"
using assms
by simp
finally show ?thesis .
qed

lemma sqrt_lower_bound_2_small:
assumes "0 ≤ x"
shows "32*x*(16*x + 1)*(256*x^2 + 96*x + 1) / (65536*x^4 + 114688*x^3 + 17920*x^2 + 448*x + 1) ≤ sqrt x"
apply (rule order_trans [OF eq_refl gen_sqrt_lower_scaled [of x 4 2]])
using assms
apply (auto simp: sqrtl_2_eq )
apply algebra
done

lemma sqrt_lower_bound_2_large:
assumes "0 ≤ x"
shows "32*x*(x + 16)*(x^2 + 96*x + 256) / (x^4 + 448*x^3 + 17920*x^2 + 114688*x + 65536) ≤ sqrt x"
apply (rule order_trans [OF eq_refl gen_sqrt_lower_scaled [of x "1/4" 2]])
using assms
apply (auto simp: sqrtl_2_eq)