Theory Closest_Pair_Alternative

section "Closest Pair Algorithm 2"

theory Closest_Pair_Alternative
  imports Common
begin

text‹
  Formalization of a divide-and-conquer algorithm solving the Closest Pair Problem
  based on the presentation of Cormen \emph{et al.} cite‹"Introduction-to-Algorithms:2009"›.
›

subsection "Functional Correctness Proof"

subsubsection "Core Argument"

lemma core_argument:
  assumes "distinct (p0 # ps)" "sorted_snd (p0 # ps)" "0 ≤ δ" "set (p0 # ps) = psL ∪ psR"
  assumes "∀p ∈ set (p0 # ps). l - δ ≤ fst p ∧ fst p ≤ l + δ"
  assumes "∀p ∈ psL. fst p ≤ l" "∀p ∈ psR. l ≤ fst p"
  assumes "sparse δ psL" "sparse δ psR"
  assumes "p1 ∈ set ps" "dist p0 p1 < δ"
  shows "p1 ∈ set (take 7 ps)"
proof -
  define PS where "PS = p0 # ps"
  define R where "R = cbox (l - δ, snd p0) (l + δ, snd p0 + δ)"
  define RPS where "RPS = { p ∈ set PS. p ∈ R }"
  define LSQ where "LSQ = cbox (l - δ, snd p0) (l, snd p0 + δ)"
  define LSQPS where "LSQPS = { p ∈ psL. p ∈ LSQ }"
  define RSQ where "RSQ = cbox (l, snd p0) (l + δ, snd p0 + δ)"
  define RSQPS where "RSQPS = { p ∈ psR. p ∈ RSQ }"
  note defs = PS_def R_def RPS_def LSQ_def LSQPS_def RSQ_def RSQPS_def

  have "R = LSQ ∪ RSQ"
    using defs cbox_right_un by auto
  moreover have "∀p ∈ psL. p ∈ RSQ ⟶ p ∈ LSQ"
    using RSQ_def LSQ_def assms(6) by auto
  moreover have "∀p ∈ psR. p ∈ LSQ ⟶ p ∈ RSQ"
    using RSQ_def LSQ_def assms(7) by auto
  ultimately have "RPS = LSQPS ∪ RSQPS"
    using LSQPS_def RSQPS_def PS_def RPS_def assms(4) by blast

  have "sparse δ LSQPS"
    using assms(8) LSQPS_def sparse_def by simp
  hence CLSQPS: "card LSQPS ≤ 4"
    using max_points_square[of LSQPS "l - δ" "snd p0" δ] assms(3) LSQ_def LSQPS_def by auto

  have "sparse δ RSQPS"
    using assms(9) RSQPS_def sparse_def by simp
  hence CRSQPS: "card RSQPS ≤ 4"
    using max_points_square[of RSQPS l "snd p0" δ] assms(3) RSQ_def RSQPS_def by auto

  have CRPS: "card RPS ≤ 8"
    using CLSQPS CRSQPS card_Un_le[of LSQPS RSQPS] ‹RPS = LSQPS ∪ RSQPS› by auto

  have "RPS ⊆ set (take 8 PS)"
  proof (rule ccontr)
    assume "¬ (RPS ⊆ set (take 8 PS))"
    then obtain p where *: "p ∈ set PS" "p ∈ RPS" "p ∉ set (take 8 PS)" "p ∈ R"
      using RPS_def by auto

    have "∀p0 ∈ set (take 8 PS). ∀p1 ∈ set (drop 8 PS). snd p0 ≤ snd p1"
      using sorted_wrt_take_drop[of "λp0 p1. snd p0 ≤ snd p1" PS 8] assms(2) sorted_snd_def PS_def by fastforce
    hence "∀p' ∈ set (take 8 PS). snd p' ≤ snd p"
      using append_take_drop_id set_append Un_iff *(1,3) by metis
    moreover have "snd p ≤ snd p0 + δ"
      using ‹p ∈ R› R_def mem_cbox_2D by (metis (mono_tags, lifting) prod.case_eq_if)
    ultimately have "∀p ∈ set (take 8 PS). snd p ≤ snd p0 + δ"
      by fastforce
    moreover have "∀p ∈ set (take 8 PS). snd p0 ≤ snd p"
      using sorted_wrt_hd_less_take[of "λp0 p1. snd p0 ≤ snd p1" p0 ps 8] assms(2) sorted_snd_def PS_def by fastforce
    moreover have "∀p ∈ set (take 8 PS). l - δ ≤ fst p ∧ fst p ≤ l + δ"
      using assms(5) PS_def by (meson in_set_takeD)
    ultimately have "∀p ∈ set (take 8 PS). p ∈ R"
      using R_def mem_cbox_2D by fastforce

    hence "set (take 8 PS) ⊆ RPS"
      using RPS_def set_take_subset by fastforce
    hence NINE: "{ p } ∪ set (take 8 PS) ⊆ RPS"
      using * by simp

    have "8 ≤ length PS"
      using *(1,3) nat_le_linear by fastforce
    hence "length (take 8 PS) = 8"
      by simp

    have "finite { p }" "finite (set (take 8 PS))"
      by simp_all
    hence "card ({ p } ∪ set (take 8 PS)) = card ({ p }) + card (set (take 8 PS))"
      using *(3) card_Un_disjoint by blast
    hence "card ({ p } ∪ set (take 8 PS)) = 9"
      using assms(1) ‹length (take 8 PS) = 8› distinct_card[of "take 8 PS"] distinct_take[of PS] PS_def by fastforce
    moreover have "finite RPS"
      using RPS_def by simp
    ultimately have "9 ≤ card RPS"
      using NINE card_mono by metis
    thus False
      using CRPS by simp
  qed

  have "dist (snd p0) (snd p1) < δ"
    using assms(11) dist_snd_le le_less_trans by (metis (no_types, lifting) prod.case_eq_if snd_conv)
  hence "snd p1 ≤ snd p0 + δ"
    by (simp add: dist_real_def)
  moreover have "l - δ ≤ fst p1" "fst p1 ≤ l + δ"
    using assms(5,10) by auto
  moreover have "snd p0 ≤ snd p1"
    using sorted_snd_def assms(2,10) by auto
  ultimately have "p1 ∈ R"
    using mem_cbox_2D[of "l - δ" "fst p1" "l + δ" "snd p0" "snd p1" "snd p0 + δ"] defs
    by (simp add: ‹l - δ ≤ fst p1› ‹snd p0 ≤ snd p1› prod.case_eq_if)
  moreover have "p1 ∈ set PS"
    using PS_def assms(10) by simp
  ultimately have "p1 ∈ set (take 8 PS)"
    using RPS_def ‹RPS ⊆ set (take 8 PS)› by auto
  thus ?thesis
    using assms(1,10) PS_def by auto
qed

subsubsection "Combine step"

lemma find_closest_bf_dist_take_7:
  assumes "∃p1 ∈ set ps. dist p0 p1 < δ"
  assumes "distinct (p0 # ps)" "sorted_snd (p0 # ps)" "0 < length ps" "0 ≤ δ" "set (p0 # ps) = psL ∪ psR"
  assumes "∀p ∈ set (p0 # ps). l - δ ≤ fst p ∧ fst p ≤ l + δ"
  assumes "∀p ∈ psL. fst p ≤ l" "∀p ∈ psR. l ≤ fst p"
  assumes "sparse δ psL" "sparse δ psR"
  shows "∀p1 ∈ set ps. dist p0 (find_closest_bf p0 (take 7 ps)) ≤ dist p0 p1"
proof -
  have "dist p0 (find_closest_bf p0 ps) < δ"
    using assms(1) dual_order.strict_trans2 find_closest_bf_dist by blast
  moreover have "find_closest_bf p0 ps ∈ set ps"
    using assms(4) find_closest_bf_set by blast
  ultimately have "find_closest_bf p0 ps ∈ set (take 7 ps)"
    using core_argument[of p0 ps δ psL psR l "find_closest_bf p0 ps"] assms by blast
  moreover have "∀p1 ∈ set (take 7 ps). dist p0 (find_closest_bf p0 (take 7 ps)) ≤ dist p0 p1"
    using find_closest_bf_dist by blast
  ultimately have "∀p1 ∈ set ps. dist p0 (find_closest_bf p0 (take 7 ps)) ≤ dist p0 p1"
    using find_closest_bf_dist order.trans by blast
  thus ?thesis .
qed

fun find_closest_pair_tm :: "(point * point) ⇒ point list ⇒ (point × point) tm" where
  "find_closest_pair_tm (c0, c1) [] =1 return (c0, c1)"
| "find_closest_pair_tm (c0, c1) [_] =1 return (c0, c1)"
| "find_closest_pair_tm (c0, c1) (p0 # ps) =1 (
    do {
      ps' <- take_tm 7 ps;
      p1 <- find_closest_bf_tm p0 ps';
      if dist c0 c1 ≤ dist p0 p1 then
        find_closest_pair_tm (c0, c1) ps
      else
        find_closest_pair_tm (p0, p1) ps
    }
  )"

fun find_closest_pair :: "(point * point) ⇒ point list ⇒ (point * point)" where
  "find_closest_pair (c0, c1) [] = (c0, c1)"
| "find_closest_pair (c0, c1) [_] = (c0, c1)"
| "find_closest_pair (c0, c1) (p0 # ps) = (
    let p1 = find_closest_bf p0 (take 7 ps) in
    if dist c0 c1 ≤ dist p0 p1 then
      find_closest_pair (c0, c1) ps
    else
      find_closest_pair (p0, p1) ps
  )"

lemma find_closest_pair_eq_val_find_closest_pair_tm:
  "val (find_closest_pair_tm (c0, c1) ps) = find_closest_pair (c0, c1) ps"
  by (induction "(c0, c1)" ps arbitrary: c0 c1 rule: find_closest_pair.induct)
     (auto simp: Let_def find_closest_bf_eq_val_find_closest_bf_tm take_eq_val_take_tm)

lemma find_closest_pair_set:
  assumes "(C0, C1) = find_closest_pair (c0, c1) ps"
  shows "(C0 ∈ set ps ∧ C1 ∈ set ps) ∨ (C0 = c0 ∧ C1 = c1)"
  using assms
proof (induction "(c0, c1)" ps arbitrary: c0 c1 C0 C1 rule: find_closest_pair.induct)
  case (3 c0 c1 p0 p2 ps)
  define p1 where p1_def: "p1 = find_closest_bf p0 (take 7 (p2 # ps))"
  hence A: "p1 ∈ set (p2 # ps)"
    using find_closest_bf_set[of "take 7 (p2 # ps)"] in_set_takeD by fastforce
  show ?case
  proof (cases "dist c0 c1 ≤ dist p0 p1")
    case True
    obtain C0' C1' where C'_def: "(C0', C1') = find_closest_pair (c0, c1) (p2 # ps)"
      using prod.collapse by blast
    note defs = p1_def C'_def
    hence "(C0' ∈ set (p2 # ps) ∧ C1' ∈ set (p2 # ps)) ∨ (C0' = c0 ∧ C1' = c1)"
      using "3.hyps"(1) True p1_def by blast
    moreover have "C0 = C0'" "C1 = C1'"
      using defs True "3.prems" apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      by auto
  next
    case False
    obtain C0' C1' where C'_def: "(C0', C1') = find_closest_pair (p0, p1) (p2 # ps)"
      using prod.collapse by blast
    note defs = p1_def C'_def
    hence "(C0' ∈ set (p2 # ps) ∧ C1' ∈ set (p2 # ps)) ∨ (C0' = p0 ∧ C1' = p1)"
      using "3.hyps"(2) p1_def False by blast
    moreover have "C0 = C0'" "C1 = C1'"
      using defs False "3.prems" apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      using A by auto
  qed
qed auto

lemma find_closest_pair_c0_ne_c1:
  "c0 ≠ c1 ⟹ distinct ps ⟹ (C0, C1) = find_closest_pair (c0, c1) ps ⟹ C0 ≠ C1"
proof (induction "(c0, c1)" ps arbitrary: c0 c1 C0 C1 rule: find_closest_pair.induct)
  case (3 c0 c1 p0 p2 ps)
  define p1 where p1_def: "p1 = find_closest_bf p0 (take 7 (p2 # ps))"
  hence "p1 ∈ set (p2 # ps)"
    using find_closest_bf_set[of "take 7 (p2 # ps)"] in_set_takeD by fastforce
  hence A: "p0 ≠ p1"
    using "3.prems"(1,2) by auto
  show ?case
  proof (cases "dist c0 c1 ≤ dist p0 p1")
    case True
    obtain C0' C1' where C'_def: "(C0', C1') = find_closest_pair (c0, c1) (p2 # ps)"
      using prod.collapse by blast
    note defs = p1_def C'_def
    hence "C0' ≠ C1'"
      using "3.hyps"(1) "3.prems"(1,2) True p1_def by simp
    moreover have "C0 = C0'" "C1 = C1'"
      using defs True "3.prems"(3) apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      by simp
  next
    case False
    obtain C0' C1' where C'_def: "(C0', C1') = find_closest_pair (p0, p1) (p2 # ps)"
      using prod.collapse by blast
    note defs = p1_def C'_def
    hence "C0' ≠ C1'"
      using "3.hyps"(2) "3.prems"(2) A False p1_def by simp
    moreover have "C0 = C0'" "C1 = C1'"
      using defs False "3.prems"(3) apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      by simp
  qed
qed auto

lemma find_closest_pair_dist_mono:
  assumes "(C0, C1) = find_closest_pair (c0, c1) ps"
  shows "dist C0 C1 ≤ dist c0 c1"
  using assms
proof (induction "(c0, c1)" ps arbitrary: c0 c1 C0 C1 rule: find_closest_pair.induct)
  case (3 c0 c1 p0 p2 ps)
  define p1 where p1_def: "p1 = find_closest_bf p0 (take 7 (p2 # ps))"
  show ?case
  proof (cases "dist c0 c1 ≤ dist p0 p1")
    case True
    obtain C0' C1' where C'_def: "(C0', C1') = find_closest_pair (c0, c1) (p2 # ps)"
      using prod.collapse by blast
    note defs = p1_def C'_def
    hence "dist C0' C1' ≤ dist c0 c1"
      using "3.hyps"(1) True p1_def by simp
    moreover have "C0 = C0'" "C1 = C1'"
      using defs True "3.prems" apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      by simp
  next
    case False
    obtain C0' C1' where C'_def: "(C0', C1') = find_closest_pair (p0, p1) (p2 # ps)"
      using prod.collapse by blast
    note defs = p1_def C'_def
    hence "dist C0' C1' ≤ dist p0 p1"
      using "3.hyps"(2) False p1_def by blast
    moreover have "C0 = C0'" "C1 = C1'"
      using defs False "3.prems"(1) apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      using False by simp
  qed
qed auto

lemma find_closest_pair_dist:
  assumes "sorted_snd ps" "distinct ps" "set ps = psL ∪ psR" "0 ≤ δ"
  assumes "∀p ∈ set ps. l - δ ≤ fst p ∧ fst p ≤ l + δ"
  assumes "∀p ∈ psL. fst p ≤ l" "∀p ∈ psR. l ≤ fst p"
  assumes "sparse δ psL" "sparse δ psR" "dist c0 c1 ≤ δ"
  assumes "(C0, C1) = find_closest_pair (c0, c1) ps"
  shows "sparse (dist C0 C1) (set ps)"
  using assms
proof (induction "(c0, c1)" ps arbitrary: c0 c1 C0 C1 psL psR rule: find_closest_pair.induct)
  case (1 c0 c1)
  thus ?case unfolding sparse_def
    by simp
next
  case (2 c0 c1 uu)
  thus ?case unfolding sparse_def
    by (metis length_greater_0_conv length_pos_if_in_set set_ConsD)
next
  case (3 c0 c1 p0 p2 ps)
  define p1 where p1_def: "p1 = find_closest_bf p0 (take 7 (p2 # ps))"
  define PSL where PSL_def: "PSL = psL - { p0 }"
  define PSR where PSR_def: "PSR = psR - { p0 }"

  have assms: "sorted_snd (p2 # ps)" "distinct (p2 # ps)" "set (p2 # ps) = PSL ∪ PSR"
              "∀p ∈ set (p2 # ps). l - δ ≤ (fst p) ∧ (fst p) ≤ l + δ"
              "∀p ∈ PSL. fst p ≤ l" "∀p ∈ PSR. l ≤ fst p"
              "sparse δ PSL" "sparse δ PSR"
    using "3.prems"(1-9) sparse_def sorted_snd_def PSL_def PSR_def by auto

  show ?case
  proof cases
    assume C1: "∃p ∈ set (p2 # ps). dist p0 p < δ"
    hence A: "∀p ∈ set (p2 # ps). dist p0 p1 ≤ dist p0 p"
      using p1_def find_closest_bf_dist_take_7 "3.prems" by blast
    hence B: "dist p0 p1 < δ"
      using C1 by auto
    show ?thesis
    proof cases
      assume C2: "dist c0 c1 ≤ dist p0 p1"
      obtain C0' C1' where def: "(C0', C1') = find_closest_pair (c0, c1) (p2 # ps)"
        using prod.collapse by blast
      hence "sparse (dist C0' C1') (set (p2 # ps))"
        using "3.hyps"(1)[of p1 PSL PSR] C2 p1_def "3.prems"(4,10) assms by blast
      moreover have "dist C0' C1' ≤ dist c0 c1"
        using def find_closest_pair_dist_mono by blast
      ultimately have "sparse (dist C0' C1') (set (p0 # p2 # ps))"
        using A C2 sparse_identity[of "dist C0' C1'" "p2 # ps" p0] by fastforce
      moreover have "C0' = C0" "C1' = C1"
        using def "3.prems"(11) C2 p1_def apply (auto) by (metis prod.inject)+
      ultimately show ?thesis
        by simp
    next
      assume C2: "¬ dist c0 c1 ≤ dist p0 p1"
      obtain C0' C1' where def: "(C0', C1') = find_closest_pair (p0, p1) (p2 # ps)"
        using prod.collapse by blast
      hence "sparse (dist C0' C1') (set (p2 # ps))"
        using "3.hyps"(2)[of p1 PSL PSR] C2 p1_def "3.prems"(4) assms B by auto
      moreover have "dist C0' C1' ≤ dist p0 p1"
        using def find_closest_pair_dist_mono by blast
      ultimately have "sparse (dist C0' C1') (set (p0 # p2 # ps))"
        using A sparse_identity order_trans by blast
      moreover have "C0' = C0" "C1' = C1"
        using def "3.prems"(11) C2 p1_def apply (auto) by (metis prod.inject)+
      ultimately show ?thesis
        by simp
    qed
  next
    assume C1: "¬ (∃p ∈ set (p2 # ps). dist p0 p < δ)"
    show ?thesis
    proof cases
      assume C2: "dist c0 c1 ≤ dist p0 p1"
      obtain C0' C1' where def: "(C0', C1') = find_closest_pair (c0, c1) (p2 # ps)"
        using prod.collapse by blast
      hence "sparse (dist C0' C1') (set (p2 # ps))"
        using "3.hyps"(1)[of p1 PSL PSR] C2 p1_def "3.prems"(4,10) assms by blast
      moreover have "dist C0' C1' ≤ dist c0 c1"
        using def find_closest_pair_dist_mono by blast
      ultimately have "sparse (dist C0' C1') (set (p0 # p2 # ps))"
        using "3.prems"(10) C1 sparse_identity[of "dist C0' C1'" "p2 # ps" p0] by force
      moreover have "C0' = C0" "C1' = C1"
        using def "3.prems"(11) C2 p1_def apply (auto) by (metis prod.inject)+
      ultimately show ?thesis
        by simp
    next
      assume C2: "¬ dist c0 c1 ≤ dist p0 p1"
      obtain C0' C1' where def: "(C0', C1') = find_closest_pair (p0, p1) (p2 # ps)"
        using prod.collapse by blast
      hence "sparse (dist C0' C1') (set (p2 # ps))"
        using "3.hyps"(2)[of p1 PSL PSR] C2 p1_def "3.prems"(4,10) assms by auto
      moreover have "dist C0' C1' ≤ dist p0 p1"
        using def find_closest_pair_dist_mono by blast
      ultimately have "sparse (dist C0' C1') (set (p0 # p2 # ps))"
        using "3.prems"(10) C1 C2 sparse_identity[of "dist C0' C1'" "p2 # ps" p0] by force
      moreover have "C0' = C0" "C1' = C1"
        using def "3.prems"(11) C2 p1_def apply (auto) by (metis prod.inject)+
      ultimately show ?thesis
        by simp
    qed
  qed
qed

declare find_closest_pair.simps [simp del]

fun combine_tm :: "(point × point) ⇒ (point × point) ⇒ int ⇒ point list ⇒ (point × point) tm" where
  "combine_tm (p0L, p1L) (p0R, p1R) l ps =1 (
    let (c0, c1) = if dist p0L p1L < dist p0R p1R then (p0L, p1L) else (p0R, p1R) in
    do {
      ps' <- filter_tm (λp. dist p (l, snd p) < dist c0 c1) ps;
      find_closest_pair_tm (c0, c1) ps'
    }
  )"

fun combine :: "(point * point) ⇒ (point * point) ⇒ int ⇒ point list ⇒ (point * point)" where
  "combine (p0L, p1L) (p0R, p1R) l ps = (
    let (c0, c1) = if dist p0L p1L < dist p0R p1R then (p0L, p1L) else (p0R, p1R) in
    let ps' = filter (λp. dist p (l, snd p) < dist c0 c1) ps in
    find_closest_pair (c0, c1) ps'
  )"

lemma combine_eq_val_combine_tm:
  "val (combine_tm (p0L, p1L) (p0R, p1R) l ps) = combine (p0L, p1L) (p0R, p1R) l ps"
  by (auto simp: filter_eq_val_filter_tm find_closest_pair_eq_val_find_closest_pair_tm)

lemma combine_set:
  assumes "(c0, c1) = combine (p0L, p1L) (p0R, p1R) l ps"
  shows "(c0 ∈ set ps ∧ c1 ∈ set ps) ∨ (c0 = p0L ∧ c1 = p1L) ∨ (c0 = p0R ∧ c1 = p1R)"
proof -
  obtain C0' C1' where C'_def: "(C0', C1') = (if dist p0L p1L < dist p0R p1R then (p0L, p1L) else (p0R, p1R))"
    by metis
  define ps' where ps'_def: "ps' = filter (λp. dist p (l, snd p) < dist C0' C1') ps"
  obtain C0 C1 where C_def: "(C0, C1) = find_closest_pair (C0', C1') ps'"
    using prod.collapse by blast
  note defs = C'_def ps'_def C_def
  have "(C0 ∈ set ps' ∧ C1 ∈ set ps') ∨ (C0 = C0' ∧ C1 = C1')"
    using C_def find_closest_pair_set by blast+
  hence "(C0 ∈ set ps ∧ C1 ∈ set ps)∨ (C0 = C0' ∧ C1 = C1')"
    using ps'_def by auto
  moreover have "C0 = c0" "C1 = c1"
    using assms defs apply (auto split: if_splits prod.splits) by (metis Pair_inject)+
  ultimately show ?thesis
    using C'_def by (auto split: if_splits)
qed

lemma combine_c0_ne_c1:
  assumes "p0L ≠ p1L" "p0R ≠ p1R" "distinct ps"
  assumes "(c0, c1) = combine (p0L, p1L) (p0R, p1R) l ps"
  shows "c0 ≠ c1"
proof -
  obtain C0' C1' where C'_def: "(C0', C1') = (if dist p0L p1L < dist p0R p1R then (p0L, p1L) else (p0R, p1R))"
    by metis
  define ps' where ps'_def: "ps' = filter (λp. dist p (l, snd p) < dist C0' C1') ps"
  obtain C0 C1 where C_def: "(C0, C1) = find_closest_pair (C0', C1') ps'"
    using prod.collapse by blast
  note defs = C'_def ps'_def C_def
  have "C0 ≠ C1"
    using defs find_closest_pair_c0_ne_c1[of C0' C1' ps'] assms by (auto split: if_splits)
  moreover have "C0 = c0" "C1 = c1"
    using assms defs apply (auto split: if_splits prod.splits) by (metis Pair_inject)+
  ultimately show ?thesis
    by blast
qed

lemma combine_dist:
  assumes "distinct ps" "sorted_snd ps" "set ps = psL ∪ psR"
  assumes "∀p ∈ psL. fst p ≤ l" "∀p ∈ psR. l ≤ fst p"
  assumes "sparse (dist p0L p1L) psL" "sparse (dist p0R p1R) psR"
  assumes "(c0, c1) = combine (p0L, p1L) (p0R, p1R) l ps"
  shows "sparse (dist c0 c1) (set ps)"
proof -
  obtain C0' C1' where C'_def: "(C0', C1') = (if dist p0L p1L < dist p0R p1R then (p0L, p1L) else (p0R, p1R))"
    by metis
  define ps' where ps'_def: "ps' = filter (λp. dist p (l, snd p) < dist C0' C1') ps"
  define PSL where PSL_def: "PSL = { p ∈ psL. dist p (l, snd p) < dist C0' C1' }"
  define PSR where PSR_def: "PSR = { p ∈ psR. dist p (l, snd p) < dist C0' C1' }"
  obtain C0 C1 where C_def: "(C0, C1) = find_closest_pair (C0', C1') ps'"
    using prod.collapse by blast
  note defs = C'_def ps'_def PSL_def PSR_def C_def
  have EQ: "C0 = c0" "C1 = c1"
    using defs assms(8) apply (auto split: if_splits prod.splits) by (metis Pair_inject)+
  have ps': "ps' = filter (λp. l - dist C0' C1' < fst p ∧ fst p < l + dist C0' C1') ps"
    using ps'_def dist_transform by simp
  have psL: "sparse (dist C0' C1') psL"
    using assms(6,8) C'_def sparse_def apply (auto split: if_splits) by force+
  hence PSL: "sparse (dist C0' C1') PSL"
    using PSL_def by (simp add: sparse_def)
  have psR: "sparse (dist C0' C1') psR"
    using assms(5,7) C'_def sparse_def apply (auto split: if_splits) by force+
  hence PSR: "sparse (dist C0' C1') PSR"
    using PSR_def by (simp add: sparse_def)
  have "sorted_snd ps'"
    using ps'_def assms(2) sorted_snd_def sorted_wrt_filter by blast
  moreover have "distinct ps'"
    using ps'_def assms(1) distinct_filter by blast
  moreover have "set ps' = PSL ∪ PSR "
    using ps'_def PSL_def PSR_def assms(3) filter_Un by auto
  moreover have "0 ≤ dist C0' C1'"
    by simp
  moreover have "∀p ∈ set ps'. l - dist C0' C1' ≤ fst p ∧ fst p ≤ l + dist C0' C1'"
    using ps' by simp
  ultimately have *: "sparse (dist C0 C1) (set ps')"
    using find_closest_pair_dist[of ps' PSL PSR "dist C0' C1'" l C0' C1'] C_def PSL PSR
    by (simp add: PSL_def PSR_def assms(4,5))
  have "∀p0 ∈ set ps. ∀p1 ∈ set ps. p0 ≠ p1 ∧ dist p0 p1 < dist C0' C1' ⟶ p0 ∈ set ps' ∧ p1 ∈ set ps'"
    using set_band_filter ps' psL psR assms(3,4,5) by blast
  moreover have "dist C0 C1 ≤ dist C0' C1'"
    using C_def find_closest_pair_dist_mono by blast
  ultimately have "∀p0 ∈ set ps. ∀p1 ∈ set ps. p0 ≠ p1 ∧ dist p0 p1 < dist C0 C1 ⟶ p0 ∈ set ps' ∧ p1 ∈ set ps'"
    by simp
  hence "sparse (dist C0 C1) (set ps)"
    using sparse_def * by (meson not_less)
  thus ?thesis
    using EQ by blast
qed

declare combine.simps [simp del]
declare combine_tm.simps [simp del]

subsubsection "Divide and Conquer Algorithm"

declare split_at_take_drop_conv [simp add]

function closest_pair_rec_tm :: "point list ⇒ (point list × point × point) tm" where
  "closest_pair_rec_tm xs =1 (
    do {
      n <- length_tm xs;
      if n ≤ 3 then
        do {
          ys <- mergesort_tm snd xs;
          p <- closest_pair_bf_tm xs;
          return (ys, p)
        }
      else
        do {
          (xsL, xsR) <- split_at_tm (n div 2) xs;
          (ysL, p0L, p1L) <- closest_pair_rec_tm xsL;
          (ysR, p0R, p1R) <- closest_pair_rec_tm xsR;
          ys <- merge_tm snd ysL ysR;
          (p0, p1) <- combine_tm (p0L, p1L) (p0R, p1R) (fst (hd xsR)) ys;
          return (ys, p0, p1)
       }
    }
  )"
  by pat_completeness auto
termination closest_pair_rec_tm
  by (relation "Wellfounded.measure (λxs. length xs)")
     (auto simp add: length_eq_val_length_tm split_at_eq_val_split_at_tm)

function closest_pair_rec :: "point list ⇒ (point list * point * point)" where
  "closest_pair_rec xs = (
    let n = length xs in
    if n ≤ 3 then
      (mergesort snd xs, closest_pair_bf xs)
    else
      let (xsL, xsR) = split_at (n div 2) xs in
      let (ysL, p0L, p1L) = closest_pair_rec xsL in
      let (ysR, p0R, p1R) = closest_pair_rec xsR in
      let ys = merge snd ysL ysR in
      (ys, combine (p0L, p1L) (p0R, p1R) (fst (hd xsR)) ys)
  )"
  by pat_completeness auto
termination closest_pair_rec
  by (relation "Wellfounded.measure (λxs. length xs)")
     (auto simp: Let_def)

declare split_at_take_drop_conv [simp del]

lemma closest_pair_rec_simps:
  assumes "n = length xs" "¬ (n ≤ 3)"
  shows "closest_pair_rec xs = (
    let (xsL, xsR) = split_at (n div 2) xs in
    let (ysL, p0L, p1L) = closest_pair_rec xsL in
    let (ysR, p0R, p1R) = closest_pair_rec xsR in
    let ys = merge snd ysL ysR in
    (ys, combine (p0L, p1L) (p0R, p1R) (fst (hd xsR)) ys)
  )"
  using assms by (auto simp: Let_def)

declare closest_pair_rec.simps [simp del]

lemma closest_pair_rec_eq_val_closest_pair_rec_tm:
  "val (closest_pair_rec_tm xs) = closest_pair_rec xs"
proof (induction rule: length_induct)
  case (1 xs)
  define n where "n = length xs"
  obtain xsL xsR where xs_def: "(xsL, xsR) = split_at (n div 2) xs"
    by (metis surj_pair)
  note defs = n_def xs_def
  show ?case
  proof cases
    assume "n ≤ 3"
    then show ?thesis
      using defs
      by (auto simp: length_eq_val_length_tm mergesort_eq_val_mergesort_tm
                     closest_pair_bf_eq_val_closest_pair_bf_tm closest_pair_rec.simps)
  next
    assume asm: "¬ n ≤ 3"
    have "length xsL < length xs" "length xsR < length xs"
      using asm defs by (auto simp: split_at_take_drop_conv)
    hence "val (closest_pair_rec_tm xsL) = closest_pair_rec xsL"
          "val (closest_pair_rec_tm xsR) = closest_pair_rec xsR"
      using "1.IH" by blast+
    thus ?thesis
      using asm defs
      apply (subst closest_pair_rec.simps, subst closest_pair_rec_tm.simps)
      by (auto simp del: closest_pair_rec_tm.simps
               simp add: Let_def length_eq_val_length_tm merge_eq_val_merge_tm
                         split_at_eq_val_split_at_tm combine_eq_val_combine_tm
               split: prod.split)
  qed
qed

lemma closest_pair_rec_set_length_sorted_snd:
  assumes "(ys, p) = closest_pair_rec xs"
  shows "set ys = set xs ∧ length ys = length xs ∧ sorted_snd ys"
  using assms
proof (induction xs arbitrary: ys p rule: length_induct)
  case (1 xs)
  let ?n = "length xs"
  show ?case
  proof (cases "?n ≤ 3")
    case True
    thus ?thesis using "1.prems" sorted_snd_def
      by (auto simp: mergesort closest_pair_rec.simps)
  next
    case False

    obtain XSL XSR where XSLR_def: "(XSL, XSR) = split_at (?n div 2) xs"
      using prod.collapse by blast
    define L where "L = fst (hd XSR)"
    obtain YSL PL where YSPL_def: "(YSL, PL) = closest_pair_rec XSL"
      using prod.collapse by blast
    obtain YSR PR where YSPR_def: "(YSR, PR) = closest_pair_rec XSR"
      using prod.collapse by blast
    define YS where "YS = merge (λp. snd p) YSL YSR"
    define P where "P = combine PL PR L YS"
    note defs = XSLR_def L_def YSPL_def YSPR_def YS_def P_def

    have "length XSL < length xs" "length XSR < length xs"
      using False defs by (auto simp: split_at_take_drop_conv)
    hence IH: "set XSL = set YSL" "set XSR = set YSR"
              "length XSL = length YSL" "length XSR = length YSR"
              "sorted_snd YSL" "sorted_snd YSR"
      using "1.IH" defs by metis+

    have "set xs = set XSL ∪ set XSR"
      using defs by (auto simp: set_take_drop split_at_take_drop_conv)
    hence SET: "set xs = set YS"
      using set_merge IH(1,2) defs by fast

    have "length xs = length XSL + length XSR"
      using defs by (auto simp: split_at_take_drop_conv)
    hence LENGTH: "length xs = length YS"
      using IH(3,4) length_merge defs by metis

    have SORTED: "sorted_snd YS"
      using IH(5,6) by (simp add: defs sorted_snd_def sorted_merge)

    have "(YS, P) = closest_pair_rec xs"
      using False closest_pair_rec_simps defs by (auto simp: Let_def split: prod.split)
    hence "(ys, p) = (YS, P)"
      using "1.prems" by argo
    thus ?thesis
      using SET LENGTH SORTED by simp
  qed
qed

lemma closest_pair_rec_distinct:
  assumes "distinct xs" "(ys, p) = closest_pair_rec xs"
  shows "distinct ys"
  using assms
proof (induction xs arbitrary: ys p rule: length_induct)
  case (1 xs)
  let ?n = "length xs"
  show ?case
  proof (cases "?n ≤ 3")
    case True
    thus ?thesis using "1.prems"
      by (auto simp: mergesort closest_pair_rec.simps)
  next
    case False

    obtain XSL XSR where XSLR_def: "(XSL, XSR) = split_at (?n div 2) xs"
      using prod.collapse by blast
    define L where "L = fst (hd XSR)"
    obtain YSL PL where YSPL_def: "(YSL, PL) = closest_pair_rec XSL"
      using prod.collapse by blast
    obtain YSR PR where YSPR_def: "(YSR, PR) = closest_pair_rec XSR"
      using prod.collapse by blast
    define YS where "YS = merge (λp. snd p) YSL YSR"
    define P where "P = combine PL PR L YS"
    note defs = XSLR_def L_def YSPL_def YSPR_def YS_def P_def

    have "length XSL < length xs" "length XSR < length xs"
      using False defs by (auto simp: split_at_take_drop_conv)
    moreover have "distinct XSL" "distinct XSR"
      using "1.prems"(1) defs by (auto simp: split_at_take_drop_conv)
    ultimately have IH: "distinct YSL" "distinct YSR"
      using "1.IH" defs by blast+

    have "set XSL ∩ set XSR = {}"
      using "1.prems"(1) defs by (auto simp: split_at_take_drop_conv set_take_disj_set_drop_if_distinct)
    moreover have "set XSL = set YSL" "set XSR = set YSR"
      using closest_pair_rec_set_length_sorted_snd defs by blast+
    ultimately have "set YSL ∩ set YSR = {}"
      by blast
    hence DISTINCT: "distinct YS"
      using distinct_merge IH defs by blast

    have "(YS, P) = closest_pair_rec xs"
      using False closest_pair_rec_simps defs by (auto simp: Let_def split: prod.split)
    hence "(ys, p) = (YS, P)"
      using "1.prems" by argo
    thus ?thesis
      using DISTINCT by blast
  qed
qed

lemma closest_pair_rec_c0_c1:
  assumes "1 < length xs" "distinct xs" "(ys, c0, c1) = closest_pair_rec xs"
  shows "c0 ∈ set xs ∧ c1 ∈ set xs ∧ c0 ≠ c1"
  using assms
proof (induction xs arbitrary: ys c0 c1 rule: length_induct)
  case (1 xs)
  let ?n = "length xs"
  show ?case
  proof (cases "?n ≤ 3")
    case True
    hence "(c0, c1) = closest_pair_bf xs"
      using "1.prems"(3) closest_pair_rec.simps by simp
    thus ?thesis
      using "1.prems"(1,2) closest_pair_bf_c0_c1 by simp
  next
    case False

    obtain XSL XSR where XSLR_def: "(XSL, XSR) = split_at (?n div 2) xs"
      using prod.collapse by blast
    define L where "L = fst (hd XSR)"

    obtain YSL C0L C1L where YSC01L_def: "(YSL, C0L, C1L) = closest_pair_rec XSL"
      using prod.collapse by metis
    obtain YSR C0R C1R where YSC01R_def: "(YSR, C0R, C1R) = closest_pair_rec XSR"
      using prod.collapse by metis

    define YS where "YS = merge (λp. snd p) YSL YSR"
    obtain C0 C1 where C01_def: "(C0, C1) = combine (C0L, C1L) (C0R, C1R) L YS"
      using prod.collapse by metis
    note defs = XSLR_def L_def YSC01L_def YSC01R_def YS_def C01_def

    have "1 < length XSL" "length XSL < length xs" "distinct XSL"
      using False "1.prems"(2) defs by (auto simp: split_at_take_drop_conv)
    hence "C0L ∈ set XSL" "C1L ∈ set XSL" and IHL1: "C0L ≠ C1L"
      using "1.IH" defs by metis+
    hence IHL2: "C0L ∈ set xs" "C1L ∈ set xs"
      using split_at_take_drop_conv in_set_takeD fst_conv defs by metis+

    have "1 < length XSR" "length XSR < length xs" "distinct XSR"
      using False "1.prems"(2) defs by (auto simp: split_at_take_drop_conv)
    hence "C0R ∈ set XSR" "C1R ∈ set XSR" and IHR1: "C0R ≠ C1R"
      using "1.IH" defs by metis+
    hence IHR2: "C0R ∈ set xs" "C1R ∈ set xs"
      using split_at_take_drop_conv in_set_dropD snd_conv defs by metis+

    have *: "(YS, C0, C1) = closest_pair_rec xs"
      using False closest_pair_rec_simps defs by (auto simp: Let_def split: prod.split)
    have YS: "set xs = set YS" "distinct YS"
      using "1.prems"(2) closest_pair_rec_set_length_sorted_snd closest_pair_rec_distinct * by blast+

    have "C0 ∈ set xs" "C1 ∈ set xs"
      using combine_set IHL2 IHR2 YS defs by blast+
    moreover have "C0 ≠ C1"
      using combine_c0_ne_c1 IHL1(1) IHR1(1) YS defs by blast
    ultimately show ?thesis
      using "1.prems"(3) * by (metis Pair_inject)
  qed
qed

lemma closest_pair_rec_dist:
  assumes "1 < length xs" "distinct xs" "sorted_fst xs" "(ys, c0, c1) = closest_pair_rec xs"
  shows "sparse (dist c0 c1) (set xs)"
  using assms
proof (induction xs arbitrary: ys c0 c1 rule: length_induct)
  case (1 xs)
  let ?n = "length xs"
  show ?case
  proof (cases "?n ≤ 3")
    case True
    hence "(c0, c1) = closest_pair_bf xs"
      using "1.prems"(4) closest_pair_rec.simps by simp
    thus ?thesis
      using "1.prems"(1,4) closest_pair_bf_dist by metis
  next
    case False

    obtain XSL XSR where XSLR_def: "(XSL, XSR) = split_at (?n div 2) xs"
      using prod.collapse by blast
    define L where "L = fst (hd XSR)"

    obtain YSL C0L C1L where YSC01L_def: "(YSL, C0L, C1L) = closest_pair_rec XSL"
      using prod.collapse by metis
    obtain YSR C0R C1R where YSC01R_def: "(YSR, C0R, C1R) = closest_pair_rec XSR"
      using prod.collapse by metis

    define YS where "YS = merge (λp. snd p) YSL YSR"
    obtain C0 C1 where C01_def: "(C0, C1) = combine (C0L, C1L) (C0R, C1R) L YS"
      using prod.collapse by metis
    note defs = XSLR_def L_def YSC01L_def YSC01R_def YS_def C01_def

    have XSLR: "XSL = take (?n div 2) xs" "XSR = drop (?n div 2) xs"
      using defs by (auto simp: split_at_take_drop_conv)

    have "1 < length XSL" "length XSL < length xs"
      using False XSLR by simp_all
    moreover have "distinct XSL" "sorted_fst XSL"
      using "1.prems"(2,3) XSLR by (auto simp: sorted_fst_def sorted_wrt_take)
    ultimately have L: "sparse (dist C0L C1L) (set XSL)"
                       "set XSL = set YSL"
      using 1 closest_pair_rec_set_length_sorted_snd closest_pair_rec_c0_c1
              YSC01L_def by blast+
    hence IHL: "sparse (dist C0L C1L) (set YSL)"
      by argo

    have "1 < length XSR" "length XSR < length xs"
      using False XSLR by simp_all
    moreover have "distinct XSR" "sorted_fst XSR"
      using "1.prems"(2,3) XSLR by (auto simp: sorted_fst_def sorted_wrt_drop)
    ultimately have R: "sparse (dist C0R C1R) (set XSR)"
                       "set XSR = set YSR"
      using 1 closest_pair_rec_set_length_sorted_snd closest_pair_rec_c0_c1
              YSC01R_def by blast+
    hence IHR: "sparse (dist C0R C1R) (set YSR)"
      by argo

    have *: "(YS, C0, C1) = closest_pair_rec xs"
      using False closest_pair_rec_simps defs by (auto simp: Let_def split: prod.split)

    have "set xs = set YS" "distinct YS" "sorted_snd YS"
      using "1.prems"(2) closest_pair_rec_set_length_sorted_snd closest_pair_rec_distinct * by blast+
    moreover have "∀p ∈ set YSL. fst p ≤ L"
      using False "1.prems"(3) XSLR L_def L(2) sorted_fst_take_less_hd_drop by simp
    moreover have "∀p ∈ set YSR. L ≤ fst p"
      using False "1.prems"(3) XSLR L_def R(2) sorted_fst_hd_drop_less_drop by simp
    moreover have "set YS = set YSL ∪ set YSR"
      using set_merge defs by fast
    moreover have "(C0, C1) = combine (C0L, C1L) (C0R, C1R) L YS"
      by (auto simp add: defs)
    ultimately have "sparse (dist C0 C1) (set xs)"
      using combine_dist IHL IHR by auto
    moreover have "(YS, C0, C1) = (ys, c0, c1)"
      using "1.prems"(4) * by simp
    ultimately show ?thesis
      by blast
  qed
qed

fun closest_pair_tm :: "point list ⇒ (point * point) tm" where
  "closest_pair_tm [] =1 return undefined"
| "closest_pair_tm [_] =1 return undefined"
| "closest_pair_tm ps =1 (
    do {
      xs <- mergesort_tm fst ps;
      (_, p) <- closest_pair_rec_tm xs;
      return p
    }
  )"

fun closest_pair :: "point list ⇒ (point * point)" where
  "closest_pair [] = undefined"
| "closest_pair [_] = undefined"
| "closest_pair ps = (let (_, c0, c1) = closest_pair_rec (mergesort fst ps) in (c0, c1))"

lemma closest_pair_eq_val_closest_pair_tm:
  "val (closest_pair_tm ps) = closest_pair ps"
  by (induction ps rule: induct_list012)
     (auto simp del: closest_pair_rec_tm.simps mergesort_tm.simps
           simp add: closest_pair_rec_eq_val_closest_pair_rec_tm mergesort_eq_val_mergesort_tm
           split: prod.split)

lemma closest_pair_simps:
  "1 < length ps ⟹ closest_pair ps = (let (_, c0, c1) = closest_pair_rec (mergesort fst ps) in (c0, c1))"
  by (induction ps rule: induct_list012) auto

declare closest_pair.simps [simp del]

theorem closest_pair_c0_c1:
  assumes "1 < length ps" "distinct ps" "(c0, c1) = closest_pair ps"
  shows "c0 ∈ set ps" "c1 ∈ set ps" "c0 ≠ c1"
  using assms closest_pair_rec_c0_c1[of "mergesort fst ps"]
  by (auto simp: mergesort closest_pair_simps split: prod.splits)

theorem closest_pair_dist:
  assumes "1 < length ps" "distinct ps" "(c0, c1) = closest_pair ps"
  shows "sparse (dist c0 c1) (set ps)"
  using assms closest_pair_rec_dist[of "mergesort fst ps"] closest_pair_rec_c0_c1[of "mergesort fst ps"]
  by (auto simp: sorted_fst_def mergesort closest_pair_simps split: prod.splits)


subsection "Time Complexity Proof"

subsubsection "Combine Step"

lemma time_find_closest_pair_tm:
  "time (find_closest_pair_tm (c0, c1) ps) ≤ 17 * length ps + 1"
proof (induction ps rule: find_closest_pair_tm.induct)
  case (3 c0 c1 p0 p2 ps)
  let ?ps = "p2 # ps"
  let ?p1 = "val (find_closest_bf_tm p0 (val (take_tm 7 ?ps)))"
  have *: "length (val (take_tm 7 ?ps)) ≤ 7"
    by (subst take_eq_val_take_tm, simp)
  show ?case
  proof cases
    assume C1: "dist c0 c1 ≤ dist p0 ?p1"
    hence "time (find_closest_pair_tm (c0, c1) (p0 # ?ps)) = 1 + time (take_tm 7 ?ps) +
           time (find_closest_bf_tm p0 (val (take_tm 7 ?ps))) + time (find_closest_pair_tm (c0, c1) ?ps)"
      by (auto simp: time_simps)
    also have "... ≤ 17 + time (find_closest_pair_tm (c0, c1) ?ps)"
      using time_take_tm[of 7 ?ps] time_find_closest_bf_tm[of p0 "val (take_tm 7 ?ps)"] * by auto
    also have "... ≤ 17 + 17 * (length ?ps) + 1"
      using "3.IH"(1) C1 by simp
    also have "... = 17 * length (p0 # ?ps) + 1"
      by simp
    finally show ?thesis .
  next
    assume C1: "¬ dist c0 c1 ≤ dist p0 ?p1"
    hence "time (find_closest_pair_tm (c0, c1) (p0 # ?ps)) = 1 + time (take_tm 7 ?ps) +
           time (find_closest_bf_tm p0 (val (take_tm 7 ?ps))) + time (find_closest_pair_tm (p0, ?p1) ?ps)"
      by (auto simp: time_simps)
    also have "... ≤ 17 + time (find_closest_pair_tm (p0, ?p1) ?ps)"
      using time_take_tm[of 7 ?ps] time_find_closest_bf_tm[of p0 "val (take_tm 7 ?ps)"] * by auto
    also have "... ≤ 17 + 17 * (length ?ps) + 1"
      using "3.IH"(2) C1 by simp
    also have "... = 17 * length (p0 # ?ps) + 1"
      by simp
    finally show ?thesis .
  qed
qed (auto simp: time_simps)

lemma time_combine_tm:
  fixes ps :: "point list"
  shows "time (combine_tm (p0L, p1L) (p0R, p1R) l ps) ≤ 3 + 18 * length ps"
proof -
  obtain c0 c1 where c_def:
    "(c0, c1) = (if dist p0L p1L < dist p0R p1R then (p0L, p1L) else (p0R, p1R))" by metis
  let ?P = "(λp. dist p (l, snd p) < dist c0 c1)"
  define ps' where "ps' = val (filter_tm ?P ps)"
  note defs = c_def ps'_def
  hence "time (combine_tm (p0L, p1L) (p0R, p1R) l ps) = 1 + time (filter_tm ?P ps) +
        time (find_closest_pair_tm (c0, c1) ps')"
    by (auto simp: combine_tm.simps Let_def time_simps split: prod.split)
  also have "... = 2 + length ps + time (find_closest_pair_tm (c0, c1) ps')"
    using time_filter_tm by auto
  also have "... ≤ 3 + length ps + 17 * length ps'"
    using defs time_find_closest_pair_tm by simp
  also have "... ≤ 3 + 18 * length ps"
     unfolding ps'_def by (subst filter_eq_val_filter_tm, simp)
  finally show ?thesis
    by blast
qed

subsubsection "Divide and Conquer Algorithm"

lemma time_closest_pair_rec_tm_simps_1:
  assumes "length xs ≤ 3"
  shows "time (closest_pair_rec_tm xs) = 1 + time (length_tm xs) + time (mergesort_tm snd xs) + time (closest_pair_bf_tm xs)"
  using assms by  (auto simp: time_simps length_eq_val_length_tm)

lemma time_closest_pair_rec_tm_simps_2:
  assumes "¬ (length xs ≤ 3)"
  shows "time (closest_pair_rec_tm xs) = 1 + (
    let (xsL, xsR) = val (split_at_tm (length xs div 2) xs) in
    let (ysL, pL) = val (closest_pair_rec_tm xsL) in
    let (ysR, pR) = val (closest_pair_rec_tm xsR) in
    let ys = val (merge_tm (λp. snd p) ysL ysR) in
    time (length_tm xs) + time (split_at_tm (length xs div 2) xs) + time (closest_pair_rec_tm xsL) +
    time (closest_pair_rec_tm xsR) + time (merge_tm (λp. snd p) ysL ysR) + time (combine_tm pL pR (fst (hd xsR)) ys)
  )"
  using assms
  apply (subst closest_pair_rec_tm.simps)
  by (auto simp del: closest_pair_rec_tm.simps
           simp add: time_simps length_eq_val_length_tm
              split: prod.split)

function closest_pair_recurrence :: "nat ⇒ real" where
  "n ≤ 3 ⟹ closest_pair_recurrence n = 3 + n + mergesort_recurrence n + n * n"
| "3 < n ⟹ closest_pair_recurrence n = 7 + 21 * n + closest_pair_recurrence (nat ⌊real n / 2⌋) +
    closest_pair_recurrence (nat ⌈real n / 2⌉)"
  by force simp_all
termination by akra_bazzi_termination simp_all

lemma closest_pair_recurrence_nonneg[simp]:
  "0 ≤ closest_pair_recurrence n"
  by (induction n rule: closest_pair_recurrence.induct) auto

lemma time_closest_pair_rec_conv_closest_pair_recurrence:
  "time (closest_pair_rec_tm ps) ≤ closest_pair_recurrence (length ps)"
proof (induction ps rule: length_induct)
  case (1 ps)
  let ?n = "length ps"
  show ?case
  proof (cases "?n ≤ 3")
    case True
    hence "time (closest_pair_rec_tm ps) = 1 + time (length_tm ps) + time (mergesort_tm snd ps) + time (closest_pair_bf_tm ps)"
      using time_closest_pair_rec_tm_simps_1 by simp
    moreover have "closest_pair_recurrence ?n = 3 + ?n + mergesort_recurrence ?n + ?n * ?n"
      using True by simp
    moreover have "time (length_tm ps) ≤ 1 + ?n" "time (mergesort_tm snd ps) ≤ mergesort_recurrence ?n"
                  "time (closest_pair_bf_tm ps) ≤ 1 + ?n * ?n"
      using time_length_tm[of ps] time_mergesort_conv_mergesort_recurrence[of snd ps] time_closest_pair_bf_tm[of ps] by auto
    ultimately show ?thesis
      by linarith
  next
    case False

    obtain XSL XSR where XS_def: "(XSL, XSR) = val (split_at_tm (?n div 2) ps)"
      using prod.collapse by blast
    obtain YSL C0L C1L where CPL_def: "(YSL, C0L, C1L) = val (closest_pair_rec_tm XSL)"
      using prod.collapse by metis
    obtain YSR C0R C1R where CPR_def: "(YSR, C0R, C1R) = val (closest_pair_rec_tm XSR)"
      using prod.collapse by metis
    define YS where "YS = val (merge_tm (λp. snd p) YSL YSR)"
    obtain C0 C1 where C01_def: "(C0, C1) = val (combine_tm (C0L, C1L) (C0R, C1R) (fst (hd XSR)) YS)"
      using prod.collapse by metis
    note defs = XS_def CPL_def CPR_def YS_def C01_def

    have XSLR: "XSL = take (?n div 2) ps" "XSR = drop (?n div 2) ps"
      using defs by (auto simp: split_at_take_drop_conv split_at_eq_val_split_at_tm)
    hence "length XSL = ?n div 2" "length XSR = ?n - ?n div 2"
      by simp_all
    hence *: "(nat ⌊real ?n / 2⌋) = length XSL" "(nat ⌈real ?n / 2⌉) = length XSR"
      by linarith+
    have "length XSL = length YSL" "length XSR = length YSR"
      using defs closest_pair_rec_set_length_sorted_snd closest_pair_rec_eq_val_closest_pair_rec_tm by metis+
    hence L: "?n = length YSL + length YSR"
      using defs XSLR by fastforce

    have "1 < length XSL" "length XSL < length ps"
      using False XSLR by simp_all
    hence "time (closest_pair_rec_tm XSL) ≤ closest_pair_recurrence (length XSL)"
      using "1.IH" by simp
    hence IHL: "time (closest_pair_rec_tm XSL) ≤ closest_pair_recurrence (nat ⌊real ?n / 2⌋)"
      using * by simp

    have "1 < length XSR" "length XSR < length ps"
      using False XSLR by simp_all
    hence "time (closest_pair_rec_tm XSR) ≤ closest_pair_recurrence (length XSR)"
      using "1.IH" by simp
    hence IHR: "time (closest_pair_rec_tm XSR) ≤ closest_pair_recurrence (nat ⌈real ?n / 2⌉)"
      using * by simp

    have "(YS, C0, C1) = val (closest_pair_rec_tm ps)"
      using False closest_pair_rec_simps defs by (auto simp: Let_def length_eq_val_length_tm split!: prod.split)
    hence "length ps = length YS"
      using closest_pair_rec_set_length_sorted_snd closest_pair_rec_eq_val_closest_pair_rec_tm by auto
    hence combine_bound: "time (combine_tm (C0L, C1L) (C0R, C1R) (fst (hd XSR)) YS) ≤ 3 + 18 * ?n"
      using time_combine_tm by simp
    have "time (closest_pair_rec_tm ps) = 1 + time (length_tm ps) + time (split_at_tm (?n div 2) ps) +
              time (closest_pair_rec_tm XSL) + time (closest_pair_rec_tm XSR) + time (merge_tm (λp. snd p) YSL YSR) +
              time (combine_tm (C0L, C1L) (C0R, C1R) (fst (hd XSR)) YS)"
      using time_closest_pair_rec_tm_simps_2[OF False] defs
      by (auto simp del: closest_pair_rec_tm.simps simp add: Let_def split: prod.split)
    also have "... ≤ 7 + 21 * ?n + time (closest_pair_rec_tm XSL) + time (closest_pair_rec_tm XSR)"
      using time_merge_tm[of "(λp. snd p)" YSL YSR] L combine_bound by (simp add: time_length_tm time_split_at_tm)
    also have "... ≤ 7 + 21 * ?n + closest_pair_recurrence (nat ⌊real ?n / 2⌋) +
              closest_pair_recurrence (nat ⌈real ?n / 2⌉)"
      using IHL IHR by simp
    also have "... = closest_pair_recurrence (length ps)"
      using False by simp
    finally show ?thesis
      by simp
  qed
qed

theorem closest_pair_recurrence:
  "closest_pair_recurrence ∈ Θ(λn. n * ln n)"
  by (master_theorem) auto

theorem time_closest_pair_rec_bigo:
  "(λxs. time (closest_pair_rec_tm xs)) ∈ O[length going_to at_top]((λn. n * ln n) o length)"
proof -
  have 0: "⋀ps. time (closest_pair_rec_tm ps) ≤ (closest_pair_recurrence o length) ps"
    unfolding comp_def using time_closest_pair_rec_conv_closest_pair_recurrence by auto
  show ?thesis
    using bigo_measure_trans[OF 0] bigthetaD1[OF closest_pair_recurrence] of_nat_0_le_iff by blast
qed

definition closest_pair_time :: "nat ⇒ real" where
  "closest_pair_time n = 1 + mergesort_recurrence n + closest_pair_recurrence n"

lemma time_closest_pair_conv_closest_pair_recurrence:
  "time (closest_pair_tm ps) ≤ closest_pair_time (length ps)"
  unfolding closest_pair_time_def
proof (induction rule: induct_list012)
  case (3 x y zs)
  let ?ps = "x # y # zs"
  define xs where "xs = val (mergesort_tm fst ?ps)"
  have *: "length xs = length ?ps"
    using xs_def mergesort(3)[of fst ?ps] mergesort_eq_val_mergesort_tm by metis
  have "time (closest_pair_tm ?ps) = 1 + time (mergesort_tm fst ?ps) + time (closest_pair_rec_tm xs)"
    using xs_def by (auto simp del: mergesort_tm.simps closest_pair_rec_tm.simps simp add: time_simps split: prod.split)
  also have "... ≤ 1 + mergesort_recurrence (length ?ps) + time (closest_pair_rec_tm xs)"
    using time_mergesort_conv_mergesort_recurrence[of fst ?ps] by simp
  also have "... ≤ 1 + mergesort_recurrence (length ?ps) + closest_pair_recurrence (length ?ps)"
    using time_closest_pair_rec_conv_closest_pair_recurrence[of xs] * by auto
  finally show ?case
    by blast
qed (auto simp: time_simps)

corollary closest_pair_time:
  "closest_pair_time ∈ O(λn. n * ln n)"
  unfolding closest_pair_time_def
  using mergesort_recurrence closest_pair_recurrence sum_in_bigo(1) const_1_bigo_n_ln_n by blast

corollary time_closest_pair_bigo:
  "(λps. time (closest_pair_tm ps)) ∈ O[length going_to at_top]((λn. n * ln n) o length)"
proof -
  have 0: "⋀ps. time (closest_pair_tm ps) ≤ (closest_pair_time o length) ps"
    unfolding comp_def using time_closest_pair_conv_closest_pair_recurrence by auto
  show ?thesis
    using bigo_measure_trans[OF 0] closest_pair_time by simp
qed


subsection "Code Export"

subsubsection "Combine Step"

fun find_closest_pair_code :: "(int * point * point) ⇒ point list ⇒ (int * point * point)" where
  "find_closest_pair_code (δ, c0, c1) [] = (δ, c0, c1)"
| "find_closest_pair_code (δ, c0, c1) [p] = (δ, c0, c1)"
| "find_closest_pair_code (δ, c0, c1) (p0 # ps) = (
    let (δ', p1) = find_closest_bf_code p0 (take 7 ps) in
    if δ ≤ δ' then
      find_closest_pair_code (δ, c0, c1) ps
    else
      find_closest_pair_code (δ', p0, p1) ps
  )"

lemma find_closest_pair_code_dist_eq:
  assumes "δ = dist_code c0 c1" "(Δ, C0, C1) = find_closest_pair_code (δ, c0, c1) ps"
  shows "Δ = dist_code C0 C1"
  using assms
proof (induction "(δ, c0, c1)" ps arbitrary: δ c0 c1 Δ C0 C1 rule: find_closest_pair_code.induct)
  case (3 δ c0 c1 p0 p2 ps)
  obtain δ' p1 where δ'_def: "(δ', p1) = find_closest_bf_code p0 (take 7 (p2 # ps))"
    by (metis surj_pair)
  hence A: "δ' = dist_code p0 p1"
    using find_closest_bf_code_dist_eq[of "take 7 (p2 # ps)"] by simp
  show ?case
  proof (cases "δ ≤ δ'")
    case True
    obtain Δ' C0' C1' where Δ'_def: "(Δ', C0', C1') = find_closest_pair_code (δ, c0, c1) (p2 # ps)"
      by (metis prod_cases4)
    note defs = δ'_def Δ'_def
    hence "Δ' = dist_code C0' C1'"
      using "3.hyps"(1)[of "(δ', p1)" δ' p1] "3.prems"(1) True δ'_def by blast
    moreover have "Δ = Δ'" "C0 = C0'" "C1 = C1'"
      using defs True "3.prems"(2) apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      by simp
  next
    case False
    obtain Δ' C0' C1' where Δ'_def: "(Δ', C0', C1') = find_closest_pair_code (δ', p0, p1) (p2 # ps)"
      by (metis prod_cases4)
    note defs = δ'_def Δ'_def
    hence "Δ' = dist_code C0' C1'"
      using "3.hyps"(2)[of "(δ', p1)" δ' p1] A False δ'_def by blast
    moreover have "Δ = Δ'" "C0 = C0'" "C1 = C1'"
      using defs False "3.prems"(2) apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      by simp
  qed
qed auto

declare find_closest_pair.simps [simp add]

lemma find_closest_pair_code_eq:
  assumes "δ = dist c0 c1" "δ' = dist_code c0 c1"
  assumes "(C0, C1) = find_closest_pair (c0, c1) ps"
  assumes "(Δ', C0', C1') = find_closest_pair_code (δ', c0, c1) ps"
  shows "C0 = C0' ∧ C1 = C1'"
  using assms
proof (induction "(c0, c1)" ps arbitrary: δ δ' c0 c1 C0 C1 Δ' C0' C1' rule: find_closest_pair.induct)
  case (3 c0 c1 p0 p2 ps)
  obtain p1 δp' p1' where δp_def: "p1 = find_closest_bf p0 (take 7 (p2 # ps))"
    "(δp', p1') = find_closest_bf_code p0 (take 7 (p2 # ps))"
    by (metis surj_pair)
  hence A: "δp' = dist_code p0 p1'"
    using find_closest_bf_code_dist_eq[of "take 7 (p2 # ps)"] by simp
  have B: "p1 = p1'"
    using "3.prems"(1,2,3) δp_def find_closest_bf_code_eq by auto
  show ?case
  proof (cases "δ ≤ dist p0 p1")
    case True
    hence C: "δ' ≤ δp'"
      by (simp add: "3.prems"(1,2) A B dist_eq_dist_code_le)
    obtain C0i C1i Δi' C0i' C1i' where Δi_def:
      "(C0i, C1i) = find_closest_pair (c0, c1) (p2 # ps)"
      "(Δi', C0i', C1i') = find_closest_pair_code (δ', c0, c1) (p2 # ps)"
      by (metis prod_cases3)
    note defs = δp_def Δi_def
    have "C0i = C0i' ∧ C1i = C1i'"
      using "3.hyps"(1)[of p1] "3.prems" True defs by blast
    moreover have "C0 = C0i" "C1 = C1i"
      using defs(1,3) True "3.prems"(1,3) apply (auto split: prod.splits) by (metis Pair_inject)+
    moreover have "Δ' = Δi'" "C0' = C0i'" "C1' = C1i'"
      using defs(2,4) C "3.prems"(4) apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      by simp
  next
    case False
    hence C: "¬ δ' ≤ δp'"
      by (simp add: "3.prems"(1,2) A B dist_eq_dist_code_le)
    obtain C0i C1i Δi' C0i' C1i' where Δi_def:
      "(C0i, C1i) = find_closest_pair (p0, p1) (p2 # ps)"
      "(Δi', C0i', C1i') = find_closest_pair_code (δp', p0, p1') (p2 # ps)"
      by (metis prod_cases3)
    note defs = δp_def Δi_def
    have "C0i = C0i' ∧ C1i = C1i'"
      using "3.prems" "3.hyps"(2)[of p1] A B False defs by blast
    moreover have "C0 = C0i" "C1 = C1i"
      using defs(1,3) False "3.prems"(1,3) apply (auto split: prod.splits) by (metis Pair_inject)+
    moreover have "Δ' = Δi'" "C0' = C0i'" "C1' = C1i'"
      using defs(2,4) C "3.prems"(4) apply (auto split: prod.splits) by (metis Pair_inject)+
    ultimately show ?thesis
      by simp
  qed
qed auto

fun combine_code :: "(int * point * point) ⇒ (int * point * point) ⇒ int ⇒ point list ⇒ (int * point * point)" where
  "combine_code (δL, p0L, p1L) (δR, p0R, p1R) l ps = (
    let (δ, c0, c1) = if δL < δR then (δL, p0L, p1L) else (δR, p0R, p1R) in
    let ps' = filter (λp. (fst p - l)2 < δ) ps in
    find_closest_pair_code (δ, c0, c1) ps'
  )"

lemma combine_code_dist_eq:
  assumes "δL = dist_code p0L p1L" "δR = dist_code p0R p1R"
  assumes "(δ, c0, c1) = combine_code (δL, p0L, p1L) (δR, p0R, p1R) l ps"
  shows "δ = dist_code c0 c1"
  using assms by (auto simp: find_closest_pair_code_dist_eq split: if_splits)

lemma combine_code_eq:
  assumes "δL' = dist_code p0L p1L" "δR' = dist_code p0R p1R"
  assumes "(c0, c1) = combine (p0L, p1L) (p0R, p1R) l ps"
  assumes "(δ', c0', c1') = combine_code (δL', p0L, p1L) (δR', p0R, p1R) l ps"
  shows "c0 = c0' ∧ c1 = c1'"
proof -
  obtain C0i C1i Δi' C0i' C1i' where Δi_def:
    "(C0i, C1i) = (if dist p0L p1L < dist p0R p1R then (p0L, p1L) else (p0R, p1R))"
    "(Δi', C0i', C1i') = (if δL' < δR' then (δL', p0L, p1L) else (δR', p0R, p1R))"
    by metis
  define ps' ps'' where ps'_def:
    "ps' = filter (λp. dist p (l, snd p) < dist C0i C1i) ps"
    "ps'' = filter (λp. (fst p - l)2 < Δi') ps"
  obtain C0 C1 Δ' C0' C1' where Δ_def:
    "(C0, C1) = find_closest_pair (C0i, C1i) ps'"
    "(Δ', C0', C1') = find_closest_pair_code (Δi', C0i', C1i') ps''"
    by (metis prod_cases3)
  note defs = Δi_def ps'_def Δ_def
  have *: "C0i = C0i'" "C1i = C1i'" "Δi' = dist_code C0i' C1i'"
    using Δi_def assms(1,2,3,4) dist_eq_dist_code_lt by (auto split: if_splits)
  hence "⋀p. ¦fst p - l¦ < dist C0i C1i ⟷ (fst p - l)2 < Δi'"
    using dist_eq_dist_code_abs_lt by (metis (mono_tags) of_int_abs)
  hence "ps' = ps''"
    using ps'_def dist_fst_abs by auto
  hence "C0 = C0'" "C1 = C1'"
    using * find_closest_pair_code_eq Δ_def by blast+
  moreover have "C0 = c0" "C1 = c1"
    using assms(3) defs(1,3,5) apply (auto simp: combine.simps split: prod.splits) by (metis Pair_inject)+
  moreover have "C0' = c0'" "C1' = c1'"
    using assms(4) defs(2,4,6) apply (auto split: prod.splits) by (metis prod.inject)+
  ultimately show ?thesis
    by blast
qed

subsubsection "Divide and Conquer Algorithm"

function closest_pair_rec_code :: "point list ⇒ (point list * int * point * point)" where
  "closest_pair_rec_code xs = (
    let n = length xs in
    if n ≤ 3 then
      (mergesort snd xs, closest_pair_bf_code xs)
    else
      let (xsL, xsR) = split_at (n div 2) xs in
      let l = fst (hd xsR) in

      let (ysL, pL) = closest_pair_rec_code xsL in
      let (ysR, pR) = closest_pair_rec_code xsR in

      let ys = merge snd ysL ysR in
      (ys, combine_code pL pR l ys)
  )"
  by pat_completeness auto
termination closest_pair_rec_code
  by (relation "Wellfounded.measure (λxs. length xs)")
     (auto simp: split_at_take_drop_conv Let_def)

lemma closest_pair_rec_code_simps:
  assumes "n = length xs" "¬ (n ≤ 3)"
  shows "closest_pair_rec_code xs = (
    let (xsL, xsR) = split_at (n div 2) xs in
    let l = fst (hd xsR) in
    let (ysL, pL) = closest_pair_rec_code xsL in
    let (ysR, pR) = closest_pair_rec_code xsR in
    let ys = merge snd ysL ysR in
    (ys, combine_code pL pR l ys)
  )"
  using assms by (auto simp: Let_def)

declare combine.simps combine_code.simps closest_pair_rec_code.simps [simp del]

lemma closest_pair_rec_code_dist_eq:
  assumes "1 < length xs" "(ys, δ, c0, c1) = closest_pair_rec_code xs"
  shows "δ = dist_code c0 c1"
  using assms
proof (induction xs arbitrary: ys δ c0 c1 rule: length_induct)
  case (1 xs)
  let ?n = "length xs"
  show ?case
  proof (cases "?n ≤ 3")
    case True
    hence "(δ, c0, c1) = closest_pair_bf_code xs"
      using "1.prems"(2) closest_pair_rec_code.simps by simp
    thus ?thesis
      using "1.prems"(1) closest_pair_bf_code_dist_eq by simp
  next
    case False

    obtain XSL XSR where XSLR_def: "(XSL, XSR) = split_at (?n div 2) xs"
      using prod.collapse by blast
    define L where "L = fst (hd XSR)"

    obtain YSL ΔL C0L C1L where YSC01L_def: "(YSL, ΔL, C0L, C1L) = closest_pair_rec_code XSL"
      using prod.collapse by metis
    obtain YSR ΔR C0R C1R where YSC01R_def: "(YSR, ΔR, C0R, C1R) = closest_pair_rec_code XSR"
      using prod.collapse by metis

    define YS where "YS = merge (λp. snd p) YSL YSR"
    obtain Δ C0 C1 where C01_def: "(Δ, C0, C1) = combine_code (ΔL, C0L, C1L) (ΔR, C0R, C1R) L YS"
      using prod.collapse by metis
    note defs = XSLR_def L_def YSC01L_def YSC01R_def YS_def C01_def

    have "1 < length XSL" "length XSL < length xs"
      using False "1.prems"(1) defs by (auto simp: split_at_take_drop_conv)
    hence IHL: "ΔL = dist_code C0L C1L"
      using "1.IH" defs by metis+

    have "1 < length XSR" "length XSR < length xs"
      using False "1.prems"(1) defs by (auto simp: split_at_take_drop_conv)
    hence IHR: "ΔR = dist_code C0R C1R"
      using "1.IH" defs by metis+

    have *: "(YS, Δ, C0, C1) = closest_pair_rec_code xs"
      using False closest_pair_rec_code_simps defs by (auto simp: Let_def split: prod.split)
    moreover have "Δ = dist_code C0 C1"
      using combine_code_dist_eq IHL IHR C01_def by blast
    ultimately show ?thesis
      using "1.prems"(2) * by (metis Pair_inject)
  qed
qed

lemma closest_pair_rec_ys_eq:
  assumes "1 < length xs"
  assumes "(ys, c0, c1) = closest_pair_rec xs"
  assumes "(ys', δ', c0', c1') = closest_pair_rec_code xs"
  shows "ys = ys'"
  using assms
proof (induction xs arbitrary: ys c0 c1 ys' δ' c0' c1' rule: length_induct)
  case (1 xs)
  let ?n = "length xs"
  show ?case
  proof (cases "?n ≤ 3")
    case True
    hence "ys = mergesort snd xs"
      using "1.prems"(2) closest_pair_rec.simps by simp
    moreover have "ys' = mergesort snd xs"
      using "1.prems"(3) closest_pair_rec_code.simps by (simp add: True)
    ultimately show ?thesis
      using "1.prems"(1) by simp
  next
    case False

    obtain XSL XSR where XSLR_def: "(XSL, XSR) = split_at (?n div 2) xs"
      using prod.collapse by blast
    define L where "L = fst (hd XSR)"

    obtain YSL C0L C1L YSL' ΔL' C0L' C1L' where YSC01L_def:
      "(YSL, C0L, C1L) = closest_pair_rec XSL"
      "(YSL', ΔL', C0L', C1L') = closest_pair_rec_code XSL"
      using prod.collapse by metis
    obtain YSR C0R C1R YSR' ΔR' C0R' C1R' where YSC01R_def:
      "(YSR, C0R, C1R) = closest_pair_rec XSR"
      "(YSR', ΔR', C0R', C1R') = closest_pair_rec_code XSR"
      using prod.collapse by metis

    define YS YS' where YS_def:
      "YS = merge (λp. snd p) YSL YSR"
      "YS' = merge (λp. snd p) YSL' YSR'"
    obtain C0 C1 Δ' C0' C1' where C01_def:
      "(C0, C1) = combine (C0L, C1L) (C0R, C1R) L YS"
      "(Δ', C0', C1') = combine_code (ΔL', C0L', C1L') (ΔR', C0R', C1R') L YS'"
      using prod.collapse by metis
    note defs = XSLR_def L_def YSC01L_def YSC01R_def YS_def C01_def

    have "1 < length XSL" "length XSL < length xs"
      using False "1.prems"(1) defs by (auto simp: split_at_take_drop_conv)
    hence IHL: "YSL = YSL'"
      using "1.IH" defs by metis

    have "1 < length XSR" "length XSR < length xs"
      using False "1.prems"(1) defs by (auto simp: split_at_take_drop_conv)
    hence IHR: "YSR = YSR'"
      using "1.IH" defs by metis

    have "(YS, C0, C1) = closest_pair_rec xs"
      using False closest_pair_rec_simps defs(1,2,3,5,7,9)
      by (auto simp: Let_def split: prod.split)
    moreover have "(YS', Δ', C0', C1') = closest_pair_rec_code xs"
      using False closest_pair_rec_code_simps defs(1,2,4,6,8,10)
      by (auto simp: Let_def split: prod.split)
    moreover have "YS = YS'"
      using IHL IHR YS_def by simp
    ultimately show ?thesis
      by (metis "1.prems"(2,3) Pair_inject)
  qed
qed

lemma closest_pair_rec_code_eq:
  assumes "1 < length xs"
  assumes "(ys, c0, c1) = closest_pair_rec xs"
  assumes "(ys', δ', c0', c1') = closest_pair_rec_code xs"
  shows "c0 = c0' ∧ c1 = c1'"
  using assms
proof (induction xs arbitrary: ys c0 c1 ys' δ' c0' c1' rule: length_induct)
  case (1 xs)
  let ?n = "length xs"
  show ?case
  proof (cases "?n ≤ 3")
    case True
    hence "(c0, c1) = closest_pair_bf xs"
      using "1.prems"(2) closest_pair_rec.simps by simp
    moreover have "(δ', c0', c1') = closest_pair_bf_code xs"
      using "1.prems"(3) closest_pair_rec_code.simps by (simp add: True)
    ultimately show ?thesis
      using "1.prems"(1) closest_pair_bf_code_eq by simp
  next
    case False

    obtain XSL XSR where XSLR_def: "(XSL, XSR) = split_at (?n div 2) xs"
      using prod.collapse by blast
    define L where "L = fst (hd XSR)"

    obtain YSL C0L C1L YSL' ΔL' C0L' C1L' where YSC01L_def:
      "(YSL, C0L, C1L) = closest_pair_rec XSL"
      "(YSL', ΔL', C0L', C1L') = closest_pair_rec_code XSL"
      using prod.collapse by metis
    obtain YSR C0R C1R YSR' ΔR' C0R' C1R' where YSC01R_def:
      "(YSR, C0R, C1R) = closest_pair_rec XSR"
      "(YSR', ΔR', C0R', C1R') = closest_pair_rec_code XSR"
      using prod.collapse by metis

    define YS YS' where YS_def:
      "YS = merge (λp. snd p) YSL YSR"
      "YS' = merge (λp. snd p) YSL' YSR'"
    obtain C0 C1 Δ' C0' C1' where C01_def:
      "(C0, C1) = combine (C0L, C1L) (C0R, C1R) L YS"
      "(Δ', C0', C1') = combine_code (ΔL', C0L', C1L') (ΔR', C0R', C1R') L YS'"
      using prod.collapse by metis
    note defs = XSLR_def L_def YSC01L_def YSC01R_def YS_def C01_def

    have "1 < length XSL" "length XSL < length xs"
      using False "1.prems"(1) defs by (auto simp: split_at_take_drop_conv)
    hence IHL: "C0L = C0L'" "C1L = C1L'"
      using "1.IH" defs by metis+

    have "1 < length XSR" "length XSR < length xs"
      using False "1.prems"(1) defs by (auto simp: split_at_take_drop_conv)
    hence IHR: "C0R = C0R'" "C1R = C1R'"
      using "1.IH" defs by metis+

    have "YS = YS'"
      using defs ‹1 < length XSL› ‹1 < length XSR› closest_pair_rec_ys_eq by blast
    moreover have "ΔL' = dist_code C0L' C1L'" "ΔR' = dist_code C0R' C1R'"
      using defs ‹1 < length XSL› ‹1 < length XSR› closest_pair_rec_code_dist_eq by blast+
    ultimately have "C0 = C0'" "C1 = C1'"
      using combine_code_eq IHL IHR C01_def by blast+
    moreover have "(YS, C0, C1) = closest_pair_rec xs"
      using False closest_pair_rec_simps defs(1,2,3,5,7,9)
      by (auto simp: Let_def split: prod.split)
    moreover have "(YS', Δ', C0', C1') = closest_pair_rec_code xs"
      using False closest_pair_rec_code_simps defs(1,2,4,6,8,10)
      by (auto simp: Let_def split: prod.split)
    ultimately show ?thesis
      using "1.prems"(2,3) by (metis Pair_inject)
  qed
qed

declare closest_pair.simps [simp add]

fun closest_pair_code :: "point list ⇒ (point * point)" where
  "closest_pair_code [] = undefined"
| "closest_pair_code [_] = undefined"
| "closest_pair_code ps = (let (_, _, c0, c1) = closest_pair_rec_code (mergesort fst ps) in (c0, c1))"

lemma closest_pair_code_eq:
  "closest_pair ps = closest_pair_code ps"
proof (induction ps rule: induct_list012)
  case (3 x y zs)
  obtain ys c0 c1 ys' δ' c0' c1' where *:
    "(ys, c0, c1) = closest_pair_rec (mergesort fst (x # y # zs))"
    "(ys', δ', c0', c1') = closest_pair_rec_code (mergesort fst (x # y # zs))"
    by (metis prod_cases3)
  moreover have "1 < length (mergesort fst (x # y # zs))"
    using length_mergesort[of fst "x # y # zs"] by simp
  ultimately have "c0 = c0'" "c1 = c1'"
    using closest_pair_rec_code_eq by blast+
  thus ?case
    using * by (auto split: prod.splits)
qed auto

export_code closest_pair_code in OCaml
  module_name Verified

end