Theory Ptolemys_Theorem

(* Author: Lukas Bulwahn <lukas.bulwahn-at-gmail.com> *)

section ‹Ptolemy's Theorem›

theory Ptolemys_Theorem
imports
  "HOL-Analysis.Multivariate_Analysis"
begin

subsection ‹Preliminaries›

subsubsection ‹Additions to Rat theory›

hide_const (open) normalize

subsubsection ‹Additions to Transcendental theory›

text ‹
Lemmas about @{const arcsin} and @{const arccos} commonly involve to show that their argument is
in the domain of those partial functions, i.e., the argument @{term y} is between @{term "-1::real"}
and @{term "1::real"}.
As the argumentation for @{term "(-1::real) ≤ y"} and @{term "y ≤ (1::real)"} is often very similar,
we prefer to prove @{term "¦y¦ ≤ (1::real)"} to the two goals above.

The lemma for rewriting the term @{term "cos (arccos y)"} is already provided in the Isabelle
distribution with name @{thm [source] cos_arccos_abs}. Here, we further provide the analogue on
@{term "arcsin"} for rewriting @{term "sin (arcsin y)"}.
›

lemma sin_arcsin_abs: "¦y¦ ≤ 1 ⟹ sin (arcsin y) = y"
  by (simp add: abs_le_iff)

text ‹
The further lemmas are the required variants from existing lemmas @{thm [source] arccos_lbound}
and @{thm [source] arccos_ubound}.
›

lemma arccos_lbound_abs [simp]:
  "¦y¦ ≤ 1 ⟹ 0 ≤ arccos y"
by (simp add: arccos_lbound)

lemma arccos_ubound_abs [simp]:
  "¦y¦ ≤ 1 ⟹ arccos y ≤ pi"
by (simp add: arccos_ubound)

text ‹
As we choose angles to be between @{term "0::real"} between @{term "2 * pi"},
we need some lemmas to reason about the sign of @{term "sin x"}
for angles @{term "x"}.
›

lemma sin_ge_zero_iff:
  assumes "0 ≤ x" "x < 2 * pi"
  shows "0 ≤ sin x ⟷ x ≤ pi"
proof
  assume "0 ≤ sin x"
  show "x ≤ pi"
  proof (rule ccontr)
    assume "¬ x ≤ pi"
    from this ‹x < 2 * pi› have "sin x < 0"
      using sin_lt_zero by auto
    from this ‹0 ≤ sin x› show False by auto
  qed
next
  assume "x ≤ pi"
  from this ‹0 ≤ x› show "0 ≤ sin x" by (simp add: sin_ge_zero)
qed

lemma sin_less_zero_iff:
  assumes "0 ≤ x" "x < 2 * pi"
  shows "sin x < 0 ⟷ pi < x"
using assms sin_ge_zero_iff by fastforce

subsubsection ‹Addition to Finite-Cartesian-Product theory›

text ‹
Here follow generally useful additions and specialised equations
for two-dimensional real-valued vectors.
›

lemma axis_nth_eq_0 [simp]:
  assumes "i ≠ j"
  shows "axis i x $ j = 0"
using assms unfolding axis_def by simp

lemma norm_axis:
  fixes x :: real
  shows "norm (axis i x) = abs x"
by (simp add: norm_eq_sqrt_inner inner_axis_axis)

lemma norm_eq_on_real_2_vec:
  fixes x :: "real ^ 2"
  shows "norm x = sqrt ((x $ 1) ^ 2 + (x $ 2) ^ 2)"
by (simp add: norm_eq_sqrt_inner inner_vec_def UNIV_2 power2_eq_square)

lemma dist_eq_on_real_2_vec:
  fixes a b :: "real ^ 2"
  shows "dist a b = sqrt ((a $ 1 - b $ 1) ^ 2 + (a $ 2 - b $ 2) ^ 2)"
unfolding dist_norm norm_eq_on_real_2_vec by simp

subsection ‹Polar Form of Two-Dimensional Real-Valued Vectors›

subsubsection ‹Definitions to Transfer to Polar Form and Back›

definition of_radiant :: "real ⇒ real ^ 2"
where
  "of_radiant ω = axis 1 (cos ω) + axis 2 (sin ω)"

definition normalize :: "real ^ 2 ⇒ real ^ 2"
where
  "normalize p = (if p = 0 then axis 1 1 else (1 / norm p) *⇩_R p)"

definition radiant_of :: "real ^ 2 ⇒ real"
where
  "radiant_of p = (THE ω. 0 ≤ ω ∧ ω < 2 * pi ∧ of_radiant ω = normalize p)"

text ‹
The vector @{term "of_radiant ω"} is the vector with length @{term "1::real"} and angle @{term "ω"}
to the first axis.
We normalize vectors to length @{term "1::real"} keeping their orientation with the normalize function.
Conversely, @{term "radiant_of p"} is the angle of vector @{term p} to the first axis, where we
choose @{term "radiant_of"} to return angles between @{term "0::real"} and @{term "2 * pi"},
following the usual high-school convention.
With these definitions, we can express the main result
@{term "norm p *⇩_R of_radiant (radiant_of p) = p"}.
Note that the main result holds for any definition of @{term "radiant_of 0"}.
So, we choose to define @{term "normalize 0"} and @{term "radiant_of 0"}, such that
@{term "radiant_of 0 = 0"}.
›

subsubsection ‹Lemmas on @{const of_radiant}›

lemma nth_of_radiant_1 [simp]:
  "of_radiant ω $ 1 = cos ω"
unfolding of_radiant_def by simp

lemma nth_of_radiant_2 [simp]:
  "of_radiant ω $ 2 = sin ω"
unfolding of_radiant_def by simp

lemma norm_of_radiant:
  "norm (of_radiant ω) = 1"
unfolding of_radiant_def norm_eq_on_real_2_vec by simp

lemma of_radiant_plus_2pi:
  "of_radiant (ω + 2 * pi) = of_radiant ω"
unfolding of_radiant_def by simp

lemma of_radiant_minus_2pi:
  "of_radiant (ω - 2 * pi) = of_radiant ω"
proof -
  have "of_radiant (ω - 2 * pi) = of_radiant (ω - 2 * pi + 2 * pi)"
    by (simp only: of_radiant_plus_2pi)
  also have "… = of_radiant ω" by simp
  finally show ?thesis .
qed

subsubsection ‹Lemmas on @{const normalize}›

lemma normalize_eq:
  "norm p *⇩_R normalize p = p"
unfolding normalize_def by simp

lemma norm_normalize:
  "norm (normalize p) = 1"
unfolding normalize_def by (auto simp add: norm_axis)

lemma nth_normalize [simp]:
  "¦normalize p $ i¦ ≤ 1"
using norm_normalize component_le_norm_cart by metis

lemma normalize_square:
  "(normalize p $ 1)⇧² + (normalize p $ 2)⇧² = 1"
using dot_square_norm[of "normalize p"]
by (simp add: inner_vec_def UNIV_2 power2_eq_square norm_normalize)

lemma nth_normalize_ge_zero_iff:
  "0 ≤ normalize p $ i ⟷ 0 ≤ p $ i"
proof
  assume "0 ≤ normalize p $ i"
  from this show "0 ≤ p $ i"
    unfolding normalize_def by (auto split: if_split_asm simp add: zero_le_divide_iff)
next
  assume "0 ≤ p $ i"
  have "0 ≤ axis 1 (1 :: real) $ i"
    using exhaust_2[of i] by auto
  from this ‹0 ≤ p $ i› show "0 ≤ normalize p $ i"
    unfolding normalize_def by auto
qed

lemma nth_normalize_less_zero_iff:
  "normalize p $ i < 0 ⟷ p $ i < 0"
using nth_normalize_ge_zero_iff leD leI by metis

lemma normalize_boundary_iff:
  "¦normalize p $ 1¦ = 1 ⟷ p $ 2 = 0"
proof
  assume "¦normalize p $ 1¦ = 1"
  from this have 1: "(p $ 1) ^ 2 = norm p ^ 2"
    unfolding normalize_def by (auto split: if_split_asm simp add: power2_eq_iff)
  moreover have "(p $ 1) ^ 2 + (p $ 2) ^ 2 = norm p ^ 2"
    using norm_eq_on_real_2_vec by auto
  ultimately show "p $ 2 = 0" by simp
next
  assume "p $ 2 = 0"
  from this have "¦p $ 1¦ = norm p"
    by (auto simp add: norm_eq_on_real_2_vec)
  from this show "¦normalize p $ 1¦ = 1"
    unfolding normalize_def by simp
qed

lemma between_normalize_if_distant_from_0:
  assumes "norm p ≥ 1"
  shows "between (0, p) (normalize p)"
using assms by (auto simp add: between_mem_segment closed_segment_def normalize_def)

lemma between_normalize_if_near_0:
  assumes "norm p ≤ 1"
  shows "between (0, normalize p) p"
proof -
  have "0 ≤ norm p" by simp
  from assms have "p = (norm p / norm p) *⇩_R p ∧ 0 ≤ norm p ∧ norm p ≤ 1" by auto
  from this have "∃u. p = (u / norm p) *⇩_R p ∧ 0 ≤ u ∧ u ≤ 1" by blast
  from this show ?thesis
    by (auto simp add: between_mem_segment closed_segment_def normalize_def)
qed

subsubsection ‹Lemmas on @{const radiant_of}›

lemma radiant_of:
  "0 ≤ radiant_of p ∧ radiant_of p < 2 * pi ∧ of_radiant (radiant_of p) = normalize p"
proof -
  let ?a = "if 0 ≤ p $ 2 then arccos (normalize p $ 1) else pi + arccos (- (normalize p $ 1))"
  have "0 ≤ ?a ∧ ?a < 2 * pi ∧ of_radiant ?a = normalize p"
  proof -
    have "0 ≤ ?a" by auto
    moreover have "?a < 2 * pi"
    proof cases
      assume "0 ≤ p $ 2"
      from this have "?a ≤ pi" by simp
      from this show ?thesis
        using pi_gt_zero by linarith
    next
      assume "¬ 0 ≤ p $ 2"
      have "arccos (- normalize p $ 1) < pi"
      proof -
        have "¦normalize p $ 1¦ ≠ 1"
          using ‹¬ 0 ≤ p $ 2› by (simp only: normalize_boundary_iff)
        from this have "arccos (- normalize p $ 1) ≠ pi"
          unfolding arccos_minus_1[symmetric] by (subst arccos_eq_iff) auto
        moreover have "arccos (- normalize p $ 1) ≤ pi" by simp
        ultimately show "arccos (- normalize p $ 1) < pi" by linarith
      qed
      from this ‹¬ 0 ≤ p $ 2› show ?thesis by simp
    qed
    moreover have "of_radiant ?a = normalize p"
    proof -
      have "of_radiant ?a $ i = normalize p $ i" for i
      proof -
        have "of_radiant ?a $ 1 = normalize p $ 1"
          unfolding of_radiant_def by (simp add: cos_arccos_abs)
        moreover have "of_radiant ?a $ 2 = normalize p $ 2"
        proof cases
          assume "0 ≤ p $ 2"
          have "sin (arccos (normalize p $ 1)) = sqrt (1 - (normalize p $ 1) ^ 2)"
            by (simp add: sin_arccos_abs)
          also have "… = normalize p $ 2"
          proof -
            have "1 - (normalize p $ 1)⇧² = (normalize p $ 2)⇧²"
              using normalize_square[of p] by auto
            from this ‹0 ≤ p $ 2› show ?thesis by (simp add: nth_normalize_ge_zero_iff)
          qed
          finally show ?thesis
            using ‹0 ≤ p $ 2› unfolding of_radiant_def by auto
        next
          assume "¬ 0 ≤ p $ 2"
          have "- sin (arccos (- normalize p $ 1)) = - sqrt (1 - (normalize p $ 1)⇧²)"
            by (simp add: sin_arccos_abs)
          also have "… = normalize p $ 2"
          proof -
            have "1 - (normalize p $ 1)⇧² = (normalize p $ 2)⇧²"
              using normalize_square[of p] by auto
            from this ‹¬ 0 ≤ p $ 2› show ?thesis
              using nth_normalize_ge_zero_iff by fastforce
          qed
          finally show ?thesis
            using ‹¬ 0 ≤ p $ 2› unfolding of_radiant_def by auto
        qed
        ultimately show ?thesis by (metis exhaust_2[of i])
      qed
      from this show ?thesis by (simp add: vec_eq_iff)
    qed
    ultimately show ?thesis by blast
  qed
  moreover {
    fix ω
    assume "0 ≤ ω ∧ ω < 2 * pi ∧ of_radiant ω = normalize p"
    from this have "0 ≤ ω" "ω < 2 * pi" "normalize p = of_radiant ω" by auto
    from this have "cos ω = normalize p $ 1" "sin ω = normalize p $ 2" by auto
    have "ω = ?a"
    proof cases
      assume "0 ≤ p $ 2"
      from this have "ω ≤ pi"
        using ‹0 ≤ ω› ‹ω < 2 * pi› ‹sin ω = normalize p $ 2›
        by (simp add: sin_ge_zero_iff[symmetric] nth_normalize_ge_zero_iff)
      from ‹0 ≤ ω› this have "ω = arccos (cos ω)" by (simp add: arccos_cos)
      from ‹cos ω = normalize p $ 1› this have "ω = arccos (normalize p $ 1)"
        by (simp add: arccos_eq_iff)
      from this show "ω = ?a" using ‹0 ≤ p $ 2› by auto
    next
      assume "¬ 0 ≤ p $ 2"
      from this have "ω > pi"
        using ‹0 ≤ ω› ‹ω < 2 * pi› ‹sin ω = normalize p $ 2›
        by (simp add: sin_less_zero_iff[symmetric] nth_normalize_less_zero_iff)
      from this ‹ω < 2 * pi› have "ω - pi = arccos (cos (ω - pi))"
        by (auto simp only: arccos_cos)
      from this ‹cos ω = normalize p $ 1› have "ω - pi = arccos (- normalize p $ 1)" by simp
      from this have "ω = pi + arccos (- normalize p $ 1)" by simp
      from this show "ω = ?a" using ‹¬ 0 ≤ p $ 2› by auto
    qed
  }
  ultimately show ?thesis
    unfolding radiant_of_def by (rule theI)
qed

lemma radiant_of_bounds [simp]:
  "0 ≤ radiant_of p" "radiant_of p < 2 * pi"
using radiant_of by auto

lemma radiant_of_weak_ubound [simp]:
  "radiant_of p ≤ 2 * pi"
using radiant_of_bounds(2)[of p] by linarith

subsubsection ‹Main Equations for Transforming to Polar Form›

lemma polar_form_eq:
  "norm p *⇩_R of_radiant (radiant_of p) = p"
using radiant_of normalize_eq by simp

lemma relative_polar_form_eq:
  "Q + dist P Q *⇩_R of_radiant (radiant_of (P - Q)) = P"
proof -
  have "norm (P - Q) *⇩_R of_radiant (radiant_of (P - Q)) = P - Q"
    unfolding polar_form_eq ..
  moreover have "dist P Q = norm (P - Q)" by (simp add: dist_norm)
  ultimately show ?thesis by (metis add.commute diff_add_cancel)
qed

subsection ‹Ptolemy's Theorem›

lemma dist_circle_segment:
  assumes "0 ≤ radius" "0 ≤ α" "α ≤ β" "β ≤ 2 * pi"
  shows "dist (center + radius *⇩_R of_radiant α) (center + radius *⇩_R of_radiant β) = 2 * radius * sin ((β - α) / 2)"
    (is "?lhs = ?rhs")
proof -
  have trigonometry: "(cos α - cos β)⇧² + (sin α - sin β)⇧² = (2 *  sin ((β - α) / 2))⇧²"
  proof -
    have sin_diff_minus: "sin ((α - β) / 2) = - sin ((β - α) / 2)"
      by (simp only: sin_minus[symmetric] minus_divide_left minus_diff_eq)
    have "(cos α - cos β)⇧² + (sin α - sin β)⇧² =
      (2 * sin ((α + β) / 2) * sin ((β - α) / 2))⇧² + (2 * sin ((α - β) / 2) * cos ((α + β) / 2))⇧²"
      by (simp only: cos_diff_cos sin_diff_sin)
    also have "… = (2 * sin ((β - α) / 2))⇧² * ((sin ((α + β) / 2))⇧² + (cos ((α + β) / 2))⇧²)"
      unfolding sin_diff_minus by algebra
    also have "… = (2 *  sin ((β - α) / 2))⇧²" by simp
    finally show ?thesis .
  qed
  from assms have "0 ≤ sin ((β - α) / 2)" by (simp add: sin_ge_zero)
  have "?lhs = sqrt (radius⇧² * ((cos α - cos β)⇧² + (sin α - sin β)⇧²))"
    unfolding dist_eq_on_real_2_vec by simp algebra
  also have "… = sqrt (radius⇧² *  (2 * sin ((β - α) / 2))⇧²)" by (simp add: trigonometry)
  also have "… = ?rhs"
    using ‹0 ≤ radius› ‹0 ≤ sin ((β - α) / 2)› by (simp add: real_sqrt_mult)
  finally show ?thesis .
qed

theorem ptolemy_trigonometric:
  fixes ω⇩₁ ω⇩₂ ω⇩₃ :: real
  shows "sin (ω⇩₁ + ω⇩₂) * sin (ω⇩₂ + ω⇩₃) = sin ω⇩₁ * sin ω⇩₃ + sin ω⇩₂ * sin (ω⇩₁ + ω⇩₂ + ω⇩₃)"
proof -
  have "sin (ω⇩₁ + ω⇩₂) * sin (ω⇩₂ + ω⇩₃) = ((sin ω⇩₂)⇧² + (cos ω⇩₂)⇧²) * sin ω⇩₁ * sin ω⇩₃ + sin ω⇩₂ * sin (ω⇩₁ + ω⇩₂ + ω⇩₃)"
    by (simp only: sin_add cos_add) algebra
  also have "… = sin ω⇩₁ * sin ω⇩₃ + sin ω⇩₂ * sin (ω⇩₁ + ω⇩₂ + ω⇩₃)" by simp
  finally show ?thesis .
qed

theorem ptolemy:
  fixes A B C D center :: "real ^ 2"
  assumes "dist center A = radius" and "dist center B = radius"
  assumes "dist center C = radius" and "dist center D = radius"
  assumes ordering_of_points:
    "radiant_of (A - center) ≤ radiant_of (B - center)"
    "radiant_of (B - center) ≤ radiant_of (C - center)"
    "radiant_of (C - center) ≤ radiant_of (D - center)"
  shows "dist A C * dist B D = dist A B * dist C D + dist A D * dist B C"
proof -
  from ‹dist center A = radius› have "0 ≤ radius" by auto
  define α β γ δ
    where "α = radiant_of (A - center)" and "β = radiant_of (B - center)"
    and "γ = radiant_of (C - center)" and "δ = radiant_of (D - center)"
  from ordering_of_points have angle_basics:
    "α ≤ β" "β ≤ γ" "γ ≤ δ"
    "0 ≤ α" "α ≤ 2 * pi" "0 ≤ β" "β ≤ 2 * pi"
    "0 ≤ γ" "γ ≤ 2 * pi" "0 ≤ δ" "δ ≤ 2 * pi"
    unfolding α_def β_def γ_def δ_def by auto
  from assms(1-4) have
    "A = center + radius *⇩_R of_radiant α" "B = center + radius *⇩_R of_radiant β"
    "C = center + radius *⇩_R of_radiant γ" "D = center + radius *⇩_R of_radiant δ"
    unfolding α_def β_def γ_def δ_def
    using relative_polar_form_eq dist_commute by metis+

  from this have dist_eqs:
    "dist A C = 2 * radius * sin ((γ - α) / 2)"
    "dist B D = 2 * radius * sin ((δ - β) / 2)"
    "dist A B = 2 * radius * sin ((β - α) / 2)"
    "dist C D = 2 * radius * sin ((δ - γ) / 2)"
    "dist A D = 2 * radius * sin ((δ - α) / 2)"
    "dist B C = 2 * radius * sin ((γ - β) / 2)"
    using angle_basics ‹radius ≥ 0› dist_circle_segment by (auto)

  have "dist A C * dist B D = 4 * radius ^ 2 * sin ((γ - α) / 2) * sin ((δ - β) / 2)"
    unfolding dist_eqs by (simp add: power2_eq_square)
  also have "… = 4 * radius ^ 2 * (sin ((β - α) / 2) * sin ((δ - γ) / 2) + sin ((γ - β) / 2) * sin ((δ - α) / 2))"
  proof -
    define ω⇩₁ ω⇩₂ ω⇩₃ where "ω⇩₁ = (β - α) / 2" and "ω⇩₂ = (γ - β) / 2" and "ω⇩₃ = (δ - γ) / 2"
    have "(γ - α) / 2 = ω⇩₁ + ω⇩₂" and "(δ - β) / 2 = ω⇩₂ + ω⇩₃" and "(δ - α) / 2 = ω⇩₁ + ω⇩₂ + ω⇩₃"
      unfolding ω⇩₁_def ω⇩₂_def ω⇩₃_def by (auto simp add: field_simps)
    have "sin ((γ - α) / 2) * sin ((δ - β) / 2) = sin (ω⇩₁ + ω⇩₂) * sin (ω⇩₂ + ω⇩₃)"
      using ‹(γ - α) / 2 = ω⇩₁ + ω⇩₂› ‹(δ - β) / 2 = ω⇩₂ + ω⇩₃› by (simp only:)
    also have "… = sin ω⇩₁ * sin ω⇩₃ + sin ω⇩₂ * sin (ω⇩₁ + ω⇩₂ + ω⇩₃)"
      by (rule ptolemy_trigonometric)
    also have "… = (sin ((β - α) / 2) * sin ((δ - γ) / 2) + sin ((γ - β) / 2) * sin ((δ - α) / 2))"
      using ω⇩₁_def ω⇩₂_def ω⇩₃_def ‹(δ - α) / 2 = ω⇩₁ + ω⇩₂ + ω⇩₃› by (simp only:)
    finally show ?thesis by simp
  qed
  also have "… = dist A B * dist C D + dist A D * dist B C"
    unfolding dist_eqs by (simp add: distrib_left power2_eq_square)
  finally show ?thesis .
qed

end