Erdos124.lean

import Mathlib.Data.List.Sort
import Mathlib.Data.List.Basic
import Mathlib.Data.Nat.Digits.Defs
import Mathlib.Algebra.BigOperators.Group.Finset.Defs
import Mathlib.Analysis.SpecialFunctions.Log.Basic
import Mathlib.Algebra.Order.Archimedean.Basic
import Mathlib.Data.Set.Finite.Basic
import Mathlib.Data.Finset.Sort
import Mathlib.Algebra.Order.BigOperators.Group.Finset
import Mathlib.Algebra.Order.BigOperators.Group.List
import Mathlib.Data.List.Induction
import Mathlib.Algebra.Field.GeomSum
import Mathlib.Data.Rat.Defs
import Mathlib.Algebra.Field.Rat
import Mathlib.Data.Finset.Basic
import Mathlib.Data.Finset.Card
import Mathlib.Data.Nat.Basic
import Mathlib.Data.List.Permutation
import Mathlib.Order.Basic
import Mathlib.Data.Nat.Order.Lemmas
import Mathlib.Tactic.Linarith
import Mathlib.Tactic.Ring
import Mathlib.Data.List.MinMax
import Mathlib.Data.Multiset.Basic
import Mathlib.Algebra.BigOperators.Group.Multiset.Basic
import Mathlib.Algebra.Ring.GeomSum
import Mathlib.Algebra.BigOperators.Ring.Finset
import Mathlib.Algebra.Order.Ring.Defs
import Mathlib.Algebra.Order.Ring.Rat
import Mathlib.Data.Nat.ModEq
import Mathlib.Data.Nat.Digits.Lemmas

set_option linter.style.longLine false


open BigOperators
open scoped BigOperators

theorem Sorted_of_sorted_append {l : List ℕ} {v : ℕ} (h : (l ++ [v]).Sorted (· ≤ ·)) : l.Sorted (· ≤ ·) := by
  have sub : List.Sublist l (l ++ [v]) := List.sublist_append_left l [v]
  exact h.sublist sub

theorem erdos_val_diverges {k : ℕ} (d : Fin k → ℕ)
    (hd : ∀ i, 2 ≤ d i)
    (hk : 0 < k) :
    ∀ M, ∃ s : Finset (Fin k × ℕ), M ≤ ∑ p ∈ s, d p.1 ^ p.2 := by
  classical
  intro M
  -- choose a fixed index i₀ : Fin k, using hk : 0 < k
  let i₀ : Fin k := ⟨0, hk⟩
  -- take the finite set of pairs (i₀, j) for j < M
  let s : Finset (Fin k × ℕ) := (Finset.range M).image (fun j => (i₀, j))
  refine ⟨s, ?_⟩
  -- show that M ≤ the sum over s of d p.1 ^ p.2
  have hM : M = ∑ _ ∈ s, (1 : ℕ) := by
    -- compute this sum as the cardinality of s
    calc
      M = (Finset.range M).card := by simp
      _ = s.card := by
        -- s is the image of (fun j => (i₀, j)), which is injective
        have hinj : Function.Injective (fun j : ℕ => (i₀, j)) := by
          intro x y hxy
          -- equality of pairs implies equality of seconds
          exact congrArg Prod.snd hxy
        -- use the injectivity to identify the cardinalities
        have h :=
          Finset.card_image_of_injective (s := Finset.range M)
            (f := fun j : ℕ => (i₀, j)) hinj
        -- h : ((Finset.range M).image (fun j => (i₀, j))).card = (Finset.range M).card
        -- rewrite with the definition of s and symmetrize
        have h' : (Finset.range M).card = s.card := by
          simpa [s] using h.symm
        exact h'
      _ = ∑ _ ∈ s, (1 : ℕ) := by simp
  -- now compare the two sums pointwise: 1 ≤ d p.1 ^ p.2
  have hle : ∑ _ ∈ s, (1 : ℕ) ≤ ∑ p ∈ s, d p.1 ^ p.2 := by
    apply Finset.sum_le_sum
    intro p hp
    -- from hd, we know 2 ≤ d p.1, hence 1 ≤ d p.1
    have h1le : 1 ≤ d p.1 :=
      le_trans (by decide : (1 : ℕ) ≤ 2) (hd p.1)
    -- hence every power is at least 1
    have hpow : 1 ≤ (d p.1) ^ p.2 := by
      -- write d p.1 as (m + 1)
      have := Nat.one_le_pow' p.2 (d p.1 - 1)
      -- simplify the base (d p.1 - 1 + 1) to d p.1
      simpa [Nat.succ_eq_add_one, Nat.sub_add_cancel h1le] using this
    simpa using hpow
  -- combine the equalities and inequalities to get the desired bound
  -- rewrite using hM
  have : M ≤ ∑ p ∈ s, d p.1 ^ p.2 := by
    simpa [hM] using hle
  exact this

theorem exists_pow_ge (d : ℕ) (hd : d ≥ 2) (n : ℕ) : ∃ k, n ≤ d^k := by
  use n
  induction n with
  | zero => simp
  | succ n ih =>
    have h1 : 1 ≤ d^n := Nat.one_le_pow n d (by omega)
    calc n + 1 ≤ d^n + 1 := Nat.add_le_add_right ih 1
         _ ≤ d^n + d^n := Nat.add_le_add_left h1 (d^n)
         _ = 2 * d^n := by ring
         _ ≤ d * d^n := Nat.mul_le_mul_right (d^n) hd
         _ = d^(n + 1) := by ring

theorem filterMap_mask_sublist {α β : Type} (l : List (α × β)) (mask : List Bool) :
    List.Sublist ((l.zip mask).filterMap (fun (p, b) => if b then some p else none)) l := by
  induction l generalizing mask with
  | nil => exact List.Sublist.slnil
  | cons p ps ih =>
    cases mask with
    | nil => exact List.nil_sublist _
    | cons b bs =>
      cases b with
      | true =>
        change (p :: (ps.zip bs).filterMap (fun (p, b) => if b then some p else none)).Sublist (p :: ps)
        exact List.Sublist.cons_cons p (ih bs)
      | false =>
        change ((ps.zip bs).filterMap (fun (p, b) => if b then some p else none)).Sublist (p :: ps)
        exact List.Sublist.cons p (ih bs)

theorem gap_condition_set_to_list {ι : Type} [DecidableEq ι] (val : ι → ℕ)
  (h_fin : ∀ X, {i : ι | val i < X}.Finite)
  (h_gap : ∀ i : ι, val i ≤ 1 + ∑ j ∈ (h_fin (val i)).toFinset, val j)
  (X : ℕ)
  (vals : List ℕ)
  (hvals : vals.Perm (((h_fin X).toFinset.toList).map val))
  (hsorted : vals.Sorted (· ≤ ·)) :
  ∀ i : Fin (vals.length - 1),
    vals[i.val + 1] ≤ 1 + (vals.take (i.val + 1)).sum := by
  intro i
  set v := vals[i.val + 1] with hv_def
  have hi : i.val + 1 < vals.length := by omega
  have hv_mem : v ∈ vals := List.getElem_mem hi
  have hv_mapped : v ∈ (h_fin X).toFinset.toList.map val := by
    rw [← hvals.mem_iff]; exact hv_mem
  rw [List.mem_map] at hv_mapped
  obtain ⟨i₀, hi₀_mem, hi₀_eq⟩ := hv_mapped
  have h_gap_i₀ := h_gap i₀
  have hset_eq : (h_fin (val i₀)).toFinset = (h_fin v).toFinset := by rw [hi₀_eq]
  rw [Finset.mem_toList] at hi₀_mem

  have h_small : ∀ j ∈ (h_fin v).toFinset, val j < v := by
    intro j hj
    simp only [Set.Finite.mem_toFinset, Set.mem_setOf_eq] at hj
    exact hj

  have hv_lt_X : v < X := by
    rw [hi₀_eq.symm]
    simp only [Set.Finite.mem_toFinset, Set.mem_setOf_eq] at hi₀_mem
    exact hi₀_mem

  have h_subset : (h_fin v).toFinset ⊆ (h_fin X).toFinset := by
    intro j hj
    simp only [Set.Finite.mem_toFinset, Set.mem_setOf_eq] at hj ⊢
    exact lt_trans hj hv_lt_X

  have h_in_take : ∀ j ∈ (h_fin v).toFinset, val j ∈ vals.take (i.val + 1) := by
    intro j hj
    have hj_lt : val j < v := h_small j hj
    have hj_in_S : j ∈ (h_fin X).toFinset := h_subset hj
    have hvalj_mem : val j ∈ vals := by
      rw [hvals.mem_iff, List.mem_map]
      exact ⟨j, Finset.mem_toList.mpr hj_in_S, rfl⟩
    rw [List.mem_iff_getElem] at hvalj_mem
    obtain ⟨k, hk_lt, hk_eq⟩ := hvalj_mem
    have hk_le : k ≤ i.val := by
      by_contra hk_gt
      push_neg at hk_gt
      have hks : i.val + 1 ≤ k := hk_gt
      have hge : v ≤ vals[k] := by
        have h1 : (⟨i.val + 1, hi⟩ : Fin vals.length) ≤ ⟨k, hk_lt⟩ := hks
        exact hsorted.rel_get_of_le h1
      rw [hk_eq] at hge
      exact Nat.lt_irrefl _ (Nat.lt_of_lt_of_le hj_lt hge)
    rw [List.mem_take_iff_getElem]
    refine ⟨k, ?_, hk_eq⟩
    apply Nat.lt_min.mpr
    exact ⟨Nat.lt_succ_of_le hk_le, hk_lt⟩

  have h_sum_le : ∑ j ∈ (h_fin v).toFinset, val j ≤ (vals.take (i.val + 1)).sum := by
    have h_sub : (h_fin v).toFinset.val.map val ≤ (vals.take (i.val + 1) : Multiset ℕ) := by
      have h_multiset_sub : (h_fin v).toFinset.val ≤ (h_fin X).toFinset.val := by
        exact Finset.val_le_iff.mpr h_subset
      have h_map_sub : (h_fin v).toFinset.val.map val ≤ (h_fin X).toFinset.val.map val := by
        exact Multiset.map_le_map h_multiset_sub
      have h_val_eq : (h_fin X).toFinset.val.map val = ((h_fin X).toFinset.toList.map val : Multiset ℕ) := by
        rw [← Multiset.map_coe]
        congr 1
        exact (Finset.coe_toList _).symm
      have h_perm_eq : ((h_fin X).toFinset.toList.map val : Multiset ℕ) = (vals : Multiset ℕ) := by
        rw [Multiset.coe_eq_coe]
        exact hvals.symm
      have h_le_vals : (h_fin v).toFinset.val.map val ≤ (vals : Multiset ℕ) := by
        calc (h_fin v).toFinset.val.map val
            ≤ (h_fin X).toFinset.val.map val := h_map_sub
          _ = ((h_fin X).toFinset.toList.map val : Multiset ℕ) := h_val_eq
          _ = (vals : Multiset ℕ) := h_perm_eq
      rw [Multiset.le_iff_count]
      intro a
      by_cases ha : a < v
      · have h1 : Multiset.count a ((h_fin v).toFinset.val.map val) ≤ Multiset.count a (vals : Multiset ℕ) := by
          exact Multiset.count_le_of_le a h_le_vals
        have h2 : Multiset.count a (vals : Multiset ℕ) = Multiset.count a (vals.take (i.val + 1) : Multiset ℕ) := by
          rw [Multiset.coe_count, Multiset.coe_count]
          have h_idx : ∀ k : ℕ, (hk : k < vals.length) → vals[k]'hk = a → k < i.val + 1 := by
            intro k hk hka
            by_contra hk_ge
            push_neg at hk_ge
            have hge : v ≤ vals[k]'hk := by
              have h1 : (⟨i.val + 1, hi⟩ : Fin vals.length) ≤ ⟨k, hk⟩ := hk_ge
              exact hsorted.rel_get_of_le h1
            rw [hka] at hge
            exact Nat.not_lt.mpr hge ha
          rw [List.count_eq_countP, List.count_eq_countP]
          conv_lhs =>
            rw [← List.take_append_drop (i.val + 1) vals, List.countP_append]
          suffices h : (vals.drop (i.val + 1)).countP (· == a) = 0 by simp [h]
          rw [List.countP_eq_zero]
          intro x hx
          simp only [beq_iff_eq]
          intro hxa
          rw [List.mem_iff_getElem] at hx
          obtain ⟨k, hk_lt, hk_eq⟩ := hx
          rw [List.length_drop] at hk_lt
          have hk_vals_lt : (i.val + 1) + k < vals.length := by omega
          have hk_get : (vals.drop (i.val + 1))[k] = vals[(i.val + 1) + k] := by
            exact List.getElem_drop
          rw [hk_get] at hk_eq
          have h_contra := h_idx ((i.val + 1) + k) hk_vals_lt (hk_eq.trans hxa)
          omega
        exact le_trans h1 (le_of_eq h2)
      · have h_count_zero : Multiset.count a ((h_fin v).toFinset.val.map val) = 0 := by
          rw [Multiset.count_eq_zero]
          intro hcontra
          rw [Multiset.mem_map] at hcontra
          obtain ⟨j, hj_mem, hj_eq⟩ := hcontra
          have : val j < v := by
            simp only [Finset.mem_val] at hj_mem
            exact h_small j hj_mem
          rw [hj_eq] at this
          exact (Nat.not_lt.mpr (Nat.not_lt.mp ha)) this
        simp only [h_count_zero, Nat.zero_le]
    calc ∑ j ∈ (h_fin v).toFinset, val j
        = ((h_fin v).toFinset.val.map val).sum := rfl
      _ ≤ (vals.take (i.val + 1) : Multiset ℕ).sum := by
          obtain ⟨u, hu⟩ := Multiset.le_iff_exists_add.mp h_sub
          rw [hu, Multiset.sum_add]
          exact Nat.le_add_right _ _
      _ = (vals.take (i.val + 1)).sum := by rfl

  calc vals[i.val + 1] = v := rfl
    _ = val i₀ := hi₀_eq.symm
    _ ≤ 1 + ∑ j ∈ (h_fin (val i₀)).toFinset, val j := h_gap_i₀
    _ = 1 + ∑ j ∈ (h_fin v).toFinset, val j := by rw [hset_eq]
    _ ≤ 1 + (vals.take (i.val + 1)).sum := Nat.add_le_add_left h_sum_le _

theorem gap_last_append (init : List ℕ) (v : ℕ) (hne : init ≠ [])
    (hgap : ∀ i : Fin ((init ++ [v]).length - 1), (init ++ [v])[i.val + 1] ≤ 1 + ((init ++ [v]).take (i.val + 1)).sum) :
    v ≤ 1 + init.sum := by
  have hlen : init.length > 0 := List.length_pos_of_ne_nil hne
  let i : Fin ((init ++ [v]).length - 1) := ⟨init.length - 1, by simp; omega⟩
  have hidx : i.val + 1 = init.length := by simp [i]; omega

  specialize hgap i

  have h1 : (init ++ [v])[i.val + 1] = v := by
    rw [List.getElem_append_right]
    · simp only [hidx, Nat.sub_self]
      simp
    · omega

  have h2 : (init ++ [v]).take (i.val + 1) = init := by
    rw [hidx]
    simp

  rw [h1, h2] at hgap
  exact hgap

theorem gap_restrict_append (vals init : List ℕ) (v : ℕ) (h : vals = init ++ [v])
    (hgap : ∀ i : Fin (vals.length - 1), vals[i.val + 1] ≤ 1 + (vals.take (i.val + 1)).sum) :
    ∀ i : Fin (init.length - 1), init[i.val + 1] ≤ 1 + (init.take (i.val + 1)).sum := by
  subst h
  intro i
  let k := i.val + 1
  have hk : k < init.length := by
    have := i.isLt; dsimp [k]; omega
  have hlen : (init ++ [v]).length - 1 = init.length := by simp
  let j : Fin ((init ++ [v]).length - 1) := ⟨i.val, by rw [hlen]; exact (i.isLt.trans (by omega))⟩
  specialize hgap j
  dsimp [j] at hgap
  rw [List.getElem_append_left] at hgap
  swap
  · exact hk
  rw [List.take_append_of_le_length] at hgap
  swap
  · omega
  exact hgap

theorem complete_sequence_finite (vals : List ℕ) (hne : vals ≠ []) (hone : vals.head hne = 1) (hsorted : vals.Sorted (· ≤ ·))
    (hgap : ∀ i : Fin (vals.length - 1), vals[(i.val + 1)] ≤ 1 + (vals.take (i.val + 1)).sum) :
    ∀ n ≤ vals.sum, ∃ mask : List Bool, mask.length = vals.length ∧
      n = ((vals.zip mask).filterMap (fun (v, b) => if b then some v else none)).sum := by
  revert hne hone hsorted hgap
  induction vals using List.reverseRecOn with
  | nil =>
    intro hne
    contradiction
  | append_singleton init v ih =>
    intro hne hone hsorted hgap
    by_cases h_init : init = []
    · subst h_init
      simp at hone
      have hv : v = 1 := hone
      subst hv
      intro n hn
      simp at hn
      have : n = 0 ∨ n = 1 := by omega
      rcases this with rfl | rfl
      · use [false]
        simp
      · use [true]
        simp
    · have hne_init : init ≠ [] := h_init
      have hone_init : init.head hne_init = 1 := by
        have h1 : (init ++ [v]).head hne = init.head hne_init := by
          simp only [List.head_append]
          simp [List.isEmpty_iff, hne_init]
        rw [← h1]
        exact hone
      have hsorted_init : init.Sorted (· ≤ ·) := Sorted_of_sorted_append hsorted
      have hgap_init : ∀ i : Fin (init.length - 1), init[i.val + 1] ≤ 1 + (init.take (i.val + 1)).sum := by
        apply gap_restrict_append (init ++ [v]) init v rfl hgap
      specialize ih hne_init hone_init hsorted_init hgap_init
      have hv_gap : v ≤ 1 + init.sum := gap_last_append init v hne_init hgap
      intro n hn
      rw [List.sum_append, List.sum_cons, List.sum_nil, add_zero] at hn
      by_cases h_le : n ≤ init.sum
      · obtain ⟨mask, hlen, hsum⟩ := ih n h_le
        use mask ++ [false]
        constructor
        · simp [hlen]
        · rw [List.zip_append hlen.symm]
          rw [List.filterMap_append, List.sum_append]
          simp only [List.zip_cons_cons, List.filterMap_cons, Bool.false_eq_true, ↓reduceIte]
          simp [hsum]
      · have h_sub : n - v ≤ init.sum := by
          have : n ≤ init.sum + v := hn
          have : v ≤ 1 + init.sum := hv_gap
          have : init.sum < n := lt_of_not_ge h_le
          omega
        obtain ⟨mask, hlen, hsum⟩ := ih (n - v) h_sub
        use mask ++ [true]
        constructor
        · simp [hlen]
        · rw [List.zip_append hlen.symm]
          rw [List.filterMap_append, List.sum_append]
          simp only [List.zip_cons_cons, List.filterMap_cons]
          simp only [↓reduceIte]
          simp
          rw [←hsum]
          have : init.sum < n := lt_of_not_ge h_le
          have : v ≤ n := by
            have : v ≤ 1 + init.sum := hv_gap
            omega
          exact (Nat.sub_add_cancel this).symm

lemma helper_filterMap_zip_map {α β} (l : List (α × β)) (mask : List Bool) :
    ((l.map Prod.fst).zip mask).filterMap (fun (v, b) => if b then some v else none) =
    ((l.zip mask).filterMap (fun (p, b) => if b then some p else none)).map Prod.fst := by
  induction l generalizing mask with
  | nil => simp
  | cons p ps ih =>
    match mask with
    | [] => simp
    | b :: bs =>
      simp only [List.map_cons, List.zip_cons_cons, List.filterMap_cons]
      cases b <;> simp [ih]

lemma map_congr_local {α β} (f g : α → β) (l : List α) (h : ∀ x ∈ l, f x = g x) :
    l.map f = l.map g := by
  induction l with
  | nil => rfl
  | cons a l ih =>
    simp only [List.map_cons, List.cons.injEq]
    constructor
    · apply h a (@List.mem_cons_self _ a l)
    · apply ih
      intros x hx
      apply h x (List.mem_cons_of_mem a hx)

theorem complete_sequence_of_gap {ι : Type} [DecidableEq ι] (val : ι → ℕ)
    (h_pos : ∀ i, 0 < val i)
    (h_one : ∃ i, val i = 1)
    (h_fin : ∀ X, {i | val i < X}.Finite)
    (h_gap : ∀ i, val i ≤ 1 + ∑ j ∈ (h_fin (val i)).toFinset, val j)
    (h_div : ∀ M, ∃ s : Finset ι, M ≤ ∑ i ∈ s, val i) :
    ∀ n, ∃ s : Finset ι, n = ∑ i ∈ s, val i := by
  intro n
  obtain ⟨s0, hs0⟩ := h_div (n + 1)
  have hs0_nonempty : s0.Nonempty := by
    by_contra h_empty
    rw [Finset.not_nonempty_iff_eq_empty] at h_empty
    simp [h_empty] at hs0

  let B := s0.sup' hs0_nonempty val
  let S := (h_fin (B + 1)).toFinset
  have hs0_sub_S : s0 ⊆ S := by
    intro i hi
    rw [Set.Finite.mem_toFinset, Set.mem_setOf_eq]
    exact Nat.lt_succ_of_le (Finset.le_sup' val hi)

  let pairs := S.toList.map (fun i => (val i, i))

  let sort_le (p1 p2 : ℕ × ι) : Bool := decide (p1.1 ≤ p2.1)
  let sorted_pairs := pairs.mergeSort sort_le
  let vals := sorted_pairs.map Prod.fst

  have h_perm_pairs : sorted_pairs.Perm pairs := List.mergeSort_perm _ _

  have h_perm_vals : vals.Perm (S.toList.map val) := by
    have h1 : S.toList.map val = pairs.map Prod.fst := by simp [pairs, Function.comp]
    rw [h1]
    exact List.Perm.map Prod.fst h_perm_pairs

  have hsorted : vals.Sorted (· ≤ ·) := by
    have h_trans : ∀ a b c : ℕ × ι, sort_le a b → sort_le b c → sort_le a c := by
      intros a b c hab hbc
      simp only [sort_le, decide_eq_true_eq] at hab hbc ⊢
      exact Nat.le_trans hab hbc
    have h_total : ∀ a b : ℕ × ι, sort_le a b || sort_le b a := by
      intros a b
      simp only [sort_le, decide_eq_true_eq, Bool.or_eq_true]
      exact Nat.le_total a.1 b.1
    have h_sorted_pairs : List.Pairwise (fun a b => sort_le a b = true) sorted_pairs :=
      List.sorted_mergeSort h_trans h_total pairs
    rw [List.Sorted, List.pairwise_map]
    apply h_sorted_pairs.imp
    intros a b hab
    simp only [sort_le, decide_eq_true_eq] at hab
    exact hab

  obtain ⟨i1, hi1⟩ := h_one
  have hi1_S : i1 ∈ S := by
    rw [Set.Finite.mem_toFinset, Set.mem_setOf_eq, hi1]
    obtain ⟨i0, hi0⟩ := hs0_nonempty
    have : 1 ≤ B := Nat.le_trans (h_pos i0) (Finset.le_sup' val hi0)
    exact Nat.lt_succ_of_le this

  have hne : vals ≠ [] := by
    intro h
    have : 1 ∈ vals := by
      rw [h_perm_vals.mem_iff, List.mem_map]
      exact ⟨i1, Finset.mem_toList.mpr hi1_S, hi1⟩
    rw [h] at this
    have : 1 ∈ ([] : List ℕ) := this
    contradiction

  have hone : vals.head hne = 1 := by
    have h1_in : 1 ∈ vals := by
      rw [h_perm_vals.mem_iff, List.mem_map]
      exact ⟨i1, Finset.mem_toList.mpr hi1_S, hi1⟩
    apply Nat.le_antisymm
    · match h_vals : vals with
      | [] => contradiction
      | v :: vs =>
        simp at hsorted h1_in ⊢
        rcases h1_in with rfl | h_in_rest
        · exact Nat.le_refl 1
        · exact hsorted.1 1 h_in_rest
    · have h_head_in : vals.head hne ∈ vals := List.head_mem hne
      rw [h_perm_vals.mem_iff, List.mem_map] at h_head_in
      rcases h_head_in with ⟨i_head, _, h_eq⟩
      rw [← h_eq]
      exact h_pos i_head

  have hgap_vals : ∀ i : Fin (vals.length - 1), vals[i.val + 1] ≤ 1 + (vals.take (i.val + 1)).sum := by
    apply gap_condition_set_to_list val h_fin h_gap (B + 1) vals h_perm_vals hsorted

  have h_n_le : n ≤ vals.sum := by
    rw [h_perm_vals.sum_eq, Finset.sum_map_toList]
    apply le_of_lt
    calc n < n + 1 := Nat.lt_succ_self n
       _ ≤ ∑ i ∈ s0, val i := hs0
       _ ≤ ∑ i ∈ S, val i := Finset.sum_le_sum_of_subset_of_nonneg hs0_sub_S (fun _ _ _ => Nat.zero_le _)

  obtain ⟨mask, hmask_len, hmask_sum⟩ := complete_sequence_finite vals hne hone hsorted hgap_vals n h_n_le

  let chosen_pairs := (sorted_pairs.zip mask).filterMap (fun (p, b) => if b then some p else none)
  let chosen_indices := chosen_pairs.map Prod.snd

  use chosen_indices.toFinset

  have h_sub_pairs : List.Sublist chosen_pairs sorted_pairs := filterMap_mask_sublist sorted_pairs mask

  have h_nodup_pairs : sorted_pairs.Nodup := by
    apply h_perm_pairs.nodup_iff.mpr
    apply List.Nodup.map_on
    · intros x hx y hy h_eq
      injection h_eq
    · exact Finset.nodup_toList S

  have h_nodup_chosen : chosen_pairs.Nodup := List.Nodup.sublist h_sub_pairs h_nodup_pairs

  have h_pairs_val : ∀ p ∈ sorted_pairs, p.1 = val p.2 := by
    intros p hp
    have : p ∈ pairs := h_perm_pairs.mem_iff.mp hp
    rw [List.mem_map] at this
    rcases this with ⟨i, _, rfl⟩
    rfl

  have h_nodup_indices : chosen_indices.Nodup := by
    apply List.Nodup.map_on _ h_nodup_chosen
    intros x hx y hy h_eq
    have hx_pairs := List.Sublist.subset h_sub_pairs hx
    have hy_pairs := List.Sublist.subset h_sub_pairs hy
    have vx_eq : x.1 = val x.2 := h_pairs_val x hx_pairs
    have vy_eq : y.1 = val y.2 := h_pairs_val y hy_pairs
    ext
    · simp [vx_eq, vy_eq, h_eq]
    · exact h_eq

  rw [List.sum_toFinset _ h_nodup_indices]

  have h_sum_eq : (chosen_indices.map val).sum = (chosen_pairs.map Prod.fst).sum := by
    simp only [chosen_indices, List.map_map]
    apply congr_arg
    apply map_congr_local
    intro p hp
    symm
    apply h_pairs_val
    exact List.Sublist.subset h_sub_pairs hp

  rw [h_sum_eq, ← helper_filterMap_zip_map, ← hmask_sum]

lemma div_le_div_pos {a b c : ℚ} (hc : 0 < c) (hab : a ≤ b) : a / c ≤ b / c := by
  rw [div_eq_mul_inv, div_eq_mul_inv]
  exact mul_le_mul_of_nonneg_right hab (le_of_lt (inv_pos.mpr hc))

theorem sum_geometric_approx_ge {d_i : ℕ} (hd : 2 ≤ d_i) (X : ℕ) :
    ∀ (hS : {e : ℕ | d_i ^ e < X}.Finite),
    ((↑X : ℚ) - 1) / ((↑d_i : ℚ) - 1) ≤ ∑ e ∈ hS.toFinset, (↑d_i : ℚ) ^ e := by
  intro hS
  -- 1. Identify N as the smallest power d_i^N >= X
  have h_exists : ∃ k, X ≤ d_i ^ k := exists_pow_ge d_i hd X
  let N := Nat.find h_exists
  have hN_ge : X ≤ d_i ^ N := Nat.find_spec h_exists
  have hN_min : ∀ k < N, d_i ^ k < X := by
    intro k hk
    have := Nat.find_min h_exists hk
    exact Nat.lt_of_not_ge this

  -- 2. Show the set S corresponds exactly to range N
  have h_range : hS.toFinset = Finset.range N := by
    ext e
    rw [Set.Finite.mem_toFinset, Finset.mem_range, Set.mem_setOf_eq]
    constructor
    · intro he
      by_contra! h_ge
      have h1 : 1 ≤ d_i := by linarith
      have : d_i ^ N ≤ d_i ^ e := Nat.pow_le_pow_right h1 h_ge
      linarith
    · intro he
      exact hN_min e he

  rw [h_range]

  have h_one_lt_d : 1 < (d_i : ℚ) := by exact_mod_cast (Nat.lt_of_succ_le hd)
  have h_denom_pos : 0 < (d_i : ℚ) - 1 := sub_pos.mpr h_one_lt_d
  have h_denom_ne_zero : (d_i : ℚ) - 1 ≠ 0 := ne_of_gt h_denom_pos

  -- Geometric sum identity
  have h_geom : ∑ e ∈ Finset.range N, (↑d_i : ℚ) ^ e = ((↑d_i : ℚ) ^ N - 1) / ((↑d_i : ℚ) - 1) := by
    rw [eq_div_iff h_denom_ne_zero]
    exact geom_sum_mul (d_i : ℚ) N

  rw [h_geom]
  apply div_le_div_pos h_denom_pos
  rw [sub_le_sub_iff_right]
  exact_mod_cast hN_ge

theorem erdos_gap_condition_helper {k : ℕ} {d : Fin k → ℕ} (hd : ∀ i, 2 ≤ d i) (X : ℕ) (hX : 1 ≤ X)
    (hsum : 1 ≤ ∑ i : Fin k, (1 : ℚ) / (↑(d i) - 1))
    (h_fin : ∀ i, {e : ℕ | d i ^ e < X}.Finite) :
    (↑X : ℚ) - 1 ≤ ∑ i : Fin k, ∑ e ∈ (h_fin i).toFinset, (↑(d i) : ℚ) ^ e := by
  -- 1. One-dimensional estimate for each base
  have h_i (i : Fin k) :
      ((↑X : ℚ) - 1) / ((↑(d i) : ℚ) - 1) ≤ ∑ e ∈ (h_fin i).toFinset, (↑(d i) : ℚ) ^ e := by
    apply sum_geometric_approx_ge (hd := hd i) X (hS := h_fin i)

  -- 2. Sum the inequalities
  have h_sum_le :
      ∑ i : Fin k, ((↑X : ℚ) - 1) / ((↑(d i) : ℚ) - 1)
        ≤ ∑ i : Fin k, ∑ e ∈ (h_fin i).toFinset, (↑(d i) : ℚ) ^ e := by
    apply Finset.sum_le_sum
    intro i _
    exact h_i i

  -- 3. Factor out (X - 1)
  have h_factor :
      ∑ i : Fin k, ((↑X : ℚ) - 1) / ((↑(d i) : ℚ) - 1)
        = ((↑X : ℚ) - 1) * ∑ i : Fin k, (1 : ℚ) / ((↑(d i) : ℚ) - 1) := by
    rw [Finset.mul_sum]
    apply Finset.sum_congr rfl
    intro i _
    rw [div_eq_mul_one_div]

  -- 4. Show (X - 1) is non-negative
  have hX_nonneg : 0 ≤ (↑X : ℚ) - 1 := by
    rw [sub_nonneg]
    norm_cast

  -- 5. Global lower bound
  have h_lower : (↑X : ℚ) - 1
      ≤ ((↑X : ℚ) - 1) * ∑ i : Fin k, (1 : ℚ) / ((↑(d i) : ℚ) - 1) := by
    nth_rewrite 1 [← mul_one ((↑X : ℚ) - 1)]
    apply mul_le_mul_of_nonneg_left hsum hX_nonneg

  -- 6. Combine
  calc
    (↑X : ℚ) - 1 ≤ ((↑X : ℚ) - 1) * ∑ i : Fin k, (1 : ℚ) / ((↑(d i) : ℚ) - 1) := h_lower
    _ = ∑ i : Fin k, ((↑X : ℚ) - 1) / ((↑(d i) : ℚ) - 1) := h_factor.symm
    _ ≤ ∑ i : Fin k, ∑ e ∈ (h_fin i).toFinset, (↑(d i) : ℚ) ^ e := h_sum_le

theorem sum_over_S_by_base {k : ℕ} {d : Fin k → ℕ} (X : ℕ)
    (S : Set (Fin k × ℕ))
    (hS_def : ∀ q : Fin k × ℕ, q ∈ S ↔ d q.1 ^ q.2 < X)
    (hS_fin : S.Finite)
    (h_fin : ∀ i, {e : ℕ | d i ^ e < X}.Finite) :
    ∑ q ∈ hS_fin.toFinset, ((↑(d q.1) : ℚ) ^ q.2) =
    ∑ i : Fin k, ∑ e ∈ (h_fin i).toFinset, ((↑(d i) : ℚ) ^ e) := by
  classical
  -- We use the general lemma `Finset.sum_finset_product'` which rewrites a sum over pairs
  -- into an iterated sum over the first and second coordinates, given a characterization of
  -- membership in the finset of pairs.
  -- First, describe membership in `hS_fin.toFinset` in terms of the fibers `{e | d i ^ e < X}`.
  have h_mem :
      ∀ p : Fin k × ℕ,
        p ∈ hS_fin.toFinset ↔
          p.1 ∈ (Finset.univ : Finset (Fin k)) ∧ p.2 ∈ (h_fin p.1).toFinset := by
    intro p
    -- `p ∈ hS_fin.toFinset` iff `p ∈ S`, by definition of `toFinset`;
    -- `p ∈ S` iff `d p.1 ^ p.2 < X`, by `hS_def`;
    -- `p.1 ∈ univ` is always true; and `p.2 ∈ (h_fin p.1).toFinset`
    -- iff `d p.1 ^ p.2 < X`, by definition of `(h_fin p.1).toFinset`.
    -- So both sides are equivalent to `d p.1 ^ p.2 < X`.
    have := hS_def p
    -- Let `simp` turn all of these characterizations into the same predicate.
    -- `Set.Finite.mem_toFinset` rewrites membership in `toFinset` to membership in the set,
    -- and `Finset.mem_univ` simplifies the first coordinate condition.
    simp [Set.Finite.mem_toFinset, hS_def, Finset.mem_univ]  -- both sides reduce to the same
  -- Now apply `Finset.sum_finset_product'` with `γ = Fin k`, `α = ℕ`, and
  -- `r = hS_fin.toFinset`, `s = univ`, `t i = (h_fin i).toFinset`.
  have h_main :
      ∑ q ∈ hS_fin.toFinset, ((↑(d q.1) : ℚ) ^ q.2) =
        ∑ i ∈ (Finset.univ : Finset (Fin k)),
          ∑ e ∈ (h_fin i).toFinset, ((↑(d i) : ℚ) ^ e) := by
    -- Instantiate `sum_finset_product'` with the function
    -- `f i e = (↑(d i) : ℚ) ^ e`.
    simpa [h_mem] using
      (Finset.sum_finset_product'
        (r := hS_fin.toFinset)
        (s := (Finset.univ : Finset (Fin k)))
        (t := fun i : Fin k => (h_fin i).toFinset)
        (h := h_mem)
        (f := fun (i : Fin k) (e : ℕ) => ((↑(d i) : ℚ) ^ e)))
  -- Finally, rewrite the sum over `univ` as the unrestricted sum over `Fin k`.
  calc
    ∑ q ∈ hS_fin.toFinset, ((↑(d q.1) : ℚ) ^ q.2)
        = ∑ i ∈ (Finset.univ : Finset (Fin k)),
            ∑ e ∈ (h_fin i).toFinset, ((↑(d i) : ℚ) ^ e) := h_main
    _ = ∑ i : Fin k, ∑ e ∈ (h_fin i).toFinset, ((↑(d i) : ℚ) ^ e) := by
      -- The notation `∑ i : Fin k,` is by definition the sum over `Finset.univ`.
      -- So this is just definitional equality.
      -- We allow `simp` to unfold the big-operator notation.
      simp

theorem erdos_gap_condition {k : ℕ} {d : Fin k → ℕ}
    (hd : ∀ i, 2 ≤ d i)
    (hsum : 1 ≤ ∑ i : Fin k, (1 : ℚ) / (d i - 1))
    (p : Fin k × ℕ) :
    let S := {q : Fin k × ℕ | d q.1 ^ q.2 < d p.1 ^ p.2};
    ∀ (hS : S.Finite), d p.1 ^ p.2 ≤ 1 + ∑ q ∈ hS.toFinset, d q.1 ^ q.2 := by
  intro S hS_fin
  -- S is the set, hS_fin is the finiteness proof
  let X := d p.1 ^ p.2

  have hX : 1 ≤ X := by
    dsimp [X]
    apply Nat.one_le_pow
    apply le_trans _ (hd p.1)
    exact Nat.le_succ 1

  have h_fin (i : Fin k) : {e : ℕ | d i ^ e < X}.Finite := by
    let f : ℕ → Fin k × ℕ := fun e => (i, e)
    have hinj : Function.Injective f := fun x y h => (Prod.mk.inj h).2
    have h_pre : {e : ℕ | d i ^ e < X} = f ⁻¹' S := by
      ext e
      simp [f, X, S]
    rw [h_pre]
    exact Set.Finite.preimage (Set.injOn_of_injective hinj) hS_fin

  have hS_def : ∀ q : Fin k × ℕ, q ∈ S ↔ d q.1 ^ q.2 < X := by
    intro q
    simp [S, X]

  have h_reindex := sum_over_S_by_base X S hS_def hS_fin h_fin
  have h_gap := erdos_gap_condition_helper hd X hX hsum h_fin

  -- Convert goal to ℚ
  rw [← Nat.cast_le (α := ℚ)]
  simp only [Nat.cast_add, Nat.cast_one, Nat.cast_sum, Nat.cast_pow]

  -- Rewrite using h_reindex in h_gap
  rw [← h_reindex] at h_gap
  rw [sub_le_iff_le_add] at h_gap
  rw [add_comm] at h_gap
  -- Now h_gap : ↑X ≤ 1 + ∑ q ∈ hS_fin.toFinset, ↑(d q.1) ^ q.2
  -- Goal: ↑(d p.1) ^ p.2 ≤ 1 + ∑ x ∈ hS_fin.toFinset, ↑(d x.1) ^ x.2
  -- Since X = d p.1 ^ p.2, we have ↑X = (↑(d p.1))^p.2 by Nat.cast_pow
  have hX_eq : (↑X : ℚ) = (↑(d p.1) : ℚ) ^ p.2 := by
    simp only [X, Nat.cast_pow]
  rw [hX_eq] at h_gap
  exact h_gap

theorem sum_pow_lt_pow_bound {b : ℕ} (hb : 1 < b) (S : Finset ℕ) (k : ℕ) (hS : ∀ j ∈ S, j < k) :
    ∑ j ∈ S, b ^ j < b ^ k := by
  revert S
  induction k with
  | zero =>
    intro S hS
    have hEmpty : S = ∅ := by
      ext j
      simp only [Finset.notMem_empty, iff_false]
      intro hj
      exact Nat.not_lt_zero j (hS j hj)
    simp [hEmpty]
  | succ n ih =>
    intro S hS
    let S_lt := S.filter (· < n)
    let S_eq := S.filter (· = n)
    have h_decomp : S = S_lt ∪ S_eq := by
      ext x
      simp only [Finset.mem_union]
      constructor
      · intro hx
        have : x < n + 1 := hS x hx
        rw [Nat.lt_succ_iff] at this
        rcases Nat.lt_or_eq_of_le this with h | h
        · left; rw [Finset.mem_filter]; exact ⟨hx, h⟩
        · right; rw [Finset.mem_filter]; exact ⟨hx, h⟩
      · rintro (h1 | h2)
        · rw [Finset.mem_filter] at h1
          exact h1.1
        · rw [Finset.mem_filter] at h2
          exact h2.1
    have h_disj : Disjoint S_lt S_eq := by
      rw [Finset.disjoint_filter]
      intro x _ h1 h2
      rw [h2] at h1
      exact Nat.lt_irrefl n h1
    rw [h_decomp, Finset.sum_union h_disj]
    have h_bound_lt : ∑ j ∈ S_lt, b ^ j < b ^ n := by
      apply ih S_lt
      intro j hj
      rw [Finset.mem_filter] at hj
      exact hj.2

    cases Finset.eq_empty_or_nonempty S_eq with
    | inl h_eq_empty =>
      rw [h_eq_empty, Finset.sum_empty, add_zero]
      calc ∑ j ∈ S_lt, b ^ j < b ^ n := h_bound_lt
           _ < b ^ (n + 1) := Nat.pow_lt_pow_right hb (Nat.lt_succ_self n)
    | inr h_eq_nonempty =>
      have h_eq_singleton : S_eq = {n} := by
        ext x
        simp only [Finset.mem_singleton]
        constructor
        · intro h
          rw [Finset.mem_filter] at h
          exact h.2
        · intro h
          rw [h, Finset.mem_filter]
          constructor
          · obtain ⟨y, hy⟩ := h_eq_nonempty
            rw [Finset.mem_filter] at hy
            rw [hy.2] at hy
            exact hy.1
          · rfl
      rw [h_eq_singleton, Finset.sum_singleton]
      calc ∑ j ∈ S_lt, b ^ j + b ^ n < b ^ n + b ^ n := Nat.add_lt_add_right h_bound_lt _
           _ = 2 * b ^ n := by rw [two_mul]
           _ ≤ b * b ^ n := Nat.mul_le_mul_right _ (Nat.succ_le_of_lt hb)
           _ = b ^ (n + 1) := by rw [Nat.pow_succ', mul_comm]

theorem sum_pow_div_pow_eq {b : ℕ} (hb : 1 < b) (S : Finset ℕ) (k : ℕ) :
(∑ j ∈ S, b ^ j) / b ^ k = ∑ j ∈ S.filter (k ≤ ·), b ^ (j - k) := by
  classical
  let S_lt := S.filter (· < k)
  let S_ge := S.filter (k ≤ ·)

  have h_split : ∑ j ∈ S, b ^ j = (∑ j ∈ S_lt, b ^ j) + (∑ j ∈ S_ge, b ^ j) := by
    rw [←Finset.sum_filter_add_sum_filter_not S (· < k)]
    congr 1
    apply Finset.sum_congr
    · ext x
      simp only [Finset.mem_filter, not_lt]
      simp [S_ge]
    · intro _ _
      rfl

  have hA_lt : (∑ j ∈ S_lt, b ^ j) < b ^ k := by
    apply sum_pow_lt_pow_bound hb
    intro j hj
    simp [S_lt] at hj
    exact hj.2

  let X := ∑ j ∈ S_ge, b ^ (j - k)
  have hB_eq : (∑ j ∈ S_ge, b ^ j) = X * b ^ k := by
    rw [mul_comm]
    dsimp [X]
    rw [Finset.mul_sum]
    apply Finset.sum_congr rfl
    intro j hj
    simp [S_ge] at hj
    rw [←Nat.pow_add]
    congr 1
    rw [add_comm]
    exact (Nat.sub_add_cancel hj.2).symm

  rw [h_split, hB_eq]

  have b_pos : 0 < b := Nat.lt_trans Nat.zero_lt_one hb
  have bk_pos : 0 < b ^ k := Nat.pow_pos b_pos

  rw [Nat.add_mul_div_right _ _ bk_pos]
  rw [Nat.div_eq_of_lt hA_lt]
  rw [zero_add]

theorem sum_pow_mod_base_eq {b : ℕ} (hb : 1 < b) (S : Finset ℕ) (k : ℕ) :
(∑ j ∈ S.filter (k ≤ ·), b ^ (j - k)) % b = if k ∈ S then 1 else 0 := by
  classical
  have hb0 : 0 < b := lt_trans Nat.zero_lt_one hb
  -- term-wise behavior modulo b
  have hterm (j : ℕ) (hj : j ∈ S.filter (k ≤ ·)) :
      b ^ (j - k) % b = (if j = k then 1 else 0) := by
    have hjk : k ≤ j := (Finset.mem_filter.1 hj).2
    have hcases : k < j ∨ k = j := lt_or_eq_of_le hjk
    cases hcases with
    | inl hlt =>
        -- j > k: exponent j-k is positive, so b^(j-k) is a multiple of b
        have hpos : 0 < j - k := Nat.sub_pos_of_lt hlt
        obtain ⟨n, hn⟩ := Nat.exists_eq_succ_of_ne_zero (ne_of_gt hpos)
        have hdiv : b ∣ b ^ n.succ := by
          simp [Nat.pow_succ, Nat.mul_comm]
        have hmod0 : b ^ n.succ % b = 0 := Nat.mod_eq_zero_of_dvd hdiv
        have hneq : j ≠ k := by
          intro h
          have : j ≤ k := le_of_eq h
          exact (not_le_of_gt hlt) this
        calc
          b ^ (j - k) % b = b ^ n.succ % b := by simp [hn]
          _ = 0 := hmod0
          _ = (if j = k then 1 else 0) := by simp [hneq]
    | inr heq =>
        -- j = k: exponent is 0, so term is 1
        subst heq
        simp [Nat.mod_eq_of_lt hb]
  -- global relation between sum and termwise mod
  have hsum_mod :=
    Finset.sum_nat_mod (s := S.filter (k ≤ ·)) (n := b) (f := fun j => b ^ (j - k))
  have hinner :
      (∑ j ∈ S.filter (k ≤ ·), b ^ (j - k) % b)
        = ∑ j ∈ S.filter (k ≤ ·), (if j = k then 1 else 0) := by
    refine Finset.sum_congr rfl ?h
    intro j hj
    exact hterm j hj
  -- rewrite the modulo of the sum using hinner
  have hrewrite :
      (∑ j ∈ S.filter (k ≤ ·), b ^ (j - k)) % b =
        (∑ j ∈ S.filter (k ≤ ·), (if j = k then 1 else 0)) % b := by
    simpa [hinner] using hsum_mod
  -- evaluate the simple Kronecker-delta sum over the filtered set
  have hsum_ite_filter :
      (∑ j ∈ S.filter (k ≤ ·), (if j = k then 1 else 0) : ℕ)
        = (if k ∈ S.filter (k ≤ ·) then 1 else 0) := by
    -- prove by induction over the finset S.filter (k ≤ ·)
    refine Finset.induction_on (s := S.filter (k ≤ ·)) ?base ?step
    · simp
    · intro a s ha hrec
      by_cases hk : a = k
      · subst hk
        simp [ha]
      · simp
  -- translate membership in the filtered set back to membership in S
  have hsum_ite :
      (∑ j ∈ S.filter (k ≤ ·), (if j = k then 1 else 0) : ℕ)
        = if k ∈ S then 1 else 0 := by
    classical
    by_cases hkS : k ∈ S
    · simp [hkS]
    · simp [hkS]
  -- combine everything and analyze the two cases k∈S or k∉S
  have hrewrite2 :
      (∑ j ∈ S.filter (k ≤ ·), b ^ (j - k)) % b
        = (if k ∈ S then 1 else 0) % b := by
    simpa [hsum_ite] using hrewrite
  by_cases hkS : k ∈ S
  · have :
        (∑ j ∈ S.filter (k ≤ ·), b ^ (j - k)) % b = 1 % b := by
        simpa [hkS] using hrewrite2
    have :
        (∑ j ∈ S.filter (k ≤ ·), b ^ (j - k)) % b = 1 := by
      simpa [Nat.mod_eq_of_lt hb] using this
    simpa [hkS] using this
  · have :
        (∑ j ∈ S.filter (k ≤ ·), b ^ (j - k)) % b = 0 % b := by
        simpa [hkS] using hrewrite2
    have :
        (∑ j ∈ S.filter (k ≤ ·), b ^ (j - k)) % b = 0 := by
      simpa [Nat.mod_eq_of_lt hb0] using this
    simpa [hkS] using this

theorem digit_eq_div_mod_getD {b : ℕ} (hb : 1 < b) (n k : ℕ) :
    (n / b ^ k) % b = (Nat.digits b n).getD k 0 := by
  induction k generalizing n with
  | zero =>
    simp only [pow_zero, Nat.div_one]
    cases hn : n with
    | zero => simp [Nat.digits_zero]
    | succ m =>
      rw [Nat.digits_def' hb (Nat.succ_pos m)]
      simp only [List.getD_cons_zero]
  | succ k ih =>
    cases hn : n with
    | zero =>
      simp only [Nat.zero_div, Nat.zero_mod, Nat.digits_zero, List.getD_nil]
    | succ m =>
      rw [Nat.digits_def' hb (Nat.succ_pos m)]
      -- need to show: m.succ / b ^ (k+1) % b = ((m.succ % b) :: digits b (m.succ / b)).getD (k+1) 0
      -- RHS = (digits b (m.succ / b)).getD k 0
      simp only [List.getD_cons_succ]
      -- Now need: m.succ / b ^ (k+1) % b = (digits b (m.succ / b)).getD k 0
      -- By IH with n := m.succ / b:
      -- (m.succ / b) / b ^ k % b = (digits b (m.succ / b)).getD k 0
      rw [← ih (m.succ / b)]
      -- Now need: m.succ / b ^ (k+1) % b = (m.succ / b) / b ^ k % b
      -- i.e., m.succ / b ^ (k+1) = (m.succ / b) / b ^ k
      congr 1
      rw [pow_succ', Nat.div_div_eq_div_mul, mul_comm]

theorem digits_sum_distinct_pow {b : ℕ} (hb : 1 < b) (S : Finset ℕ) :
    (Nat.digits b (∑ j ∈ S, b ^ j)).toFinset ⊆ {0, 1} := by
  classical
  -- We start by working with the number represented by the sum of powers.
  set n : ℕ := ∑ j ∈ S, b ^ j with hn_def
  intro x hx
  -- Convert membership in the finset of digits to membership in the digits list.
  have hx_list : x ∈ Nat.digits b n := by
    simpa [hn_def] using (List.mem_toFinset.mp hx)
  -- Step 1: closed form for each digit of n in base b.
  have h_digit_form : ∀ k, (Nat.digits b n).getD k 0 = (if k ∈ S then 1 else 0) := by
    intro k
    have hdiv : ((∑ j ∈ S, b ^ j) / b ^ k) = ∑ j ∈ S.filter (k ≤ ·), b ^ (j - k) :=
      sum_pow_div_pow_eq (b := b) hb S k
    have hmod : (∑ j ∈ S.filter (k ≤ ·), b ^ (j - k)) % b = if k ∈ S then 1 else 0 :=
      sum_pow_mod_base_eq (b := b) hb S k
    have hmod' : (((∑ j ∈ S, b ^ j) / b ^ k) % b) = if k ∈ S then 1 else 0 := by
      simpa [hdiv] using hmod
    have hmod'' : (n / b ^ k) % b = if k ∈ S then 1 else 0 := by
      simpa [hn_def] using hmod'
    have hdig := digit_eq_div_mod_getD (b := b) hb n k
    simpa [hmod''] using hdig.symm
  -- Local helper: every element of a list occurs as some `getD` value.
  have h_getD_cover : ∀ (l : List ℕ) (d y : ℕ), y ∈ l → ∃ k, l.getD k d = y := by
    intro l
    induction l with
    | nil =>
      intro d y hy
      simp at hy
    | cons a l ih =>
      intro d y hy
      have hy' : y = a ∨ y ∈ l := by
        simpa [List.mem_cons] using hy
      rcases hy' with hya | hyl
      · subst hya
        refine ⟨0, ?_⟩
        rfl
      · obtain ⟨k, hk⟩ := ih d y hyl
        refine ⟨k.succ, ?_⟩
        -- We need to show (a :: l).getD k.succ d = y
        -- By definition, (a :: l).getD (k+1) d = l.getD k d when k < l.length
        -- Or by the getD definition: (a :: l)[k+1]?.getD d = l[k]?.getD d
        -- Since (a :: l)[k+1]? = l[k]? by List.getElem?_cons_succ
        -- And we know l.getD k d = y from hk
        show (a :: l).getD k.succ d = y
        have h1 : (a :: l).getD k.succ d = (a :: l)[k.succ]?.getD d := rfl
        have h2 : (a :: l)[k.succ]? = l[k]? := by simp [List.getElem?_cons_succ]
        have h3 : l.getD k d = l[k]?.getD d := rfl
        rw [h1, h2, ← h3, hk]
  -- Step 2: every digit of `n` is 0 or 1.
  have h_all_digits : ∀ y ∈ Nat.digits b n, y = 0 ∨ y = 1 := by
    intro y hy
    -- Find an index `k` such that `getD k 0` equals `y`.
    obtain ⟨k, hk⟩ := h_getD_cover (Nat.digits b n) 0 y hy
    -- Use the closed form for `getD` at this index.
    have hform := h_digit_form k
    -- Combine: y = getD k 0 = if k ∈ S then 1 else 0
    have heq : y = if k ∈ S then 1 else 0 := by rw [← hform, hk]
    -- Conclude that `y` is 0 or 1.
    by_cases hS : k ∈ S
    · right; simp [heq, hS]
    · left; simp [heq, hS]
  -- Apply the fact about all digits to `x`.
  have hx01 : x = 0 ∨ x = 1 := h_all_digits x hx_list
  -- Turn this into membership in the finset `{0, 1}`.
  have : x ∈ ({0, 1} : Finset ℕ) := by
    rcases hx01 with hx0 | hx1
    · subst hx0; simp
    · subst hx1; simp
  exact this

theorem erdos_124 (k : ℕ) (d : Fin k → ℕ)
    (hd : ∀ i, 2 ≤ d i) (hsum : 1 ≤ ∑ i : Fin k, (1 : ℚ) / (d i - 1)) :
    ∀ n, ∃ a : Fin k → ℕ,
    ∀ i, (Nat.digits (d i) (a i)).toFinset ⊆ {0, 1} ∧
    n = ∑ i, a i := by
  intro n
  -- Handle k = 0 case: hsum becomes 1 ≤ 0, contradiction
  by_cases hk : k = 0
  · subst hk
    simp only [Finset.univ_eq_empty, Finset.sum_empty] at hsum
    exact absurd hsum (by norm_num)
  · -- k > 0: apply the complete sequence theorem
    have hk_pos : 0 < k := Nat.pos_of_ne_zero hk
    -- Define val
    let val : Fin k × ℕ → ℕ := fun p => d p.1 ^ p.2
    -- Verify h_pos
    have h_pos : ∀ p, 0 < val p := fun p =>
      pow_pos (Nat.lt_of_lt_of_le zero_lt_two (hd p.1)) p.2
    -- Verify h_one: val (i, 0) = d i ^ 0 = 1
    have h_one : ∃ p, val p = 1 := ⟨(⟨0, hk_pos⟩, 0), pow_zero _⟩
    -- Verify h_fin: {p | d p.1 ^ p.2 < X} is finite
    have lt_two_pow : ∀ X : ℕ, X < 2 ^ X := fun X => by
      induction X with
      | zero => norm_num
      | succ n ih => calc n + 1 < 2 ^ n + 1 := by omega
        _ ≤ 2 ^ n + 2 ^ n := by omega
        _ = 2 ^ (n + 1) := by ring
    have h_fin : ∀ X, {p | val p < X}.Finite := by
      intro X
      have hbd : ∀ p : Fin k × ℕ, val p < X → p.2 < X := fun p hp => by
        by_contra h
        push_neg at h
        have h2 : 2 ≤ d p.1 := hd p.1
        have h3 : X < 2 ^ X := lt_two_pow X
        have h4 : 2 ^ X ≤ d p.1 ^ X := Nat.pow_le_pow_left h2 X
        have h5 : d p.1 ^ X ≤ d p.1 ^ p.2 := Nat.pow_le_pow_right (Nat.lt_of_lt_of_le zero_lt_two h2) h
        have : X < val p := calc X < 2 ^ X := h3
          _ ≤ d p.1 ^ X := h4
          _ ≤ d p.1 ^ p.2 := h5
        linarith
      have hsub : {p : Fin k × ℕ | val p < X} ⊆ (Set.univ : Set (Fin k)) ×ˢ Set.Iio X := by
        intro p hp
        simp only [Set.mem_prod, Set.mem_univ, Set.mem_Iio, true_and]
        exact hbd p hp
      exact Set.Finite.subset (Set.finite_univ.prod (Set.finite_Iio X)) hsub
    -- Verify h_gap (from erdos_gap_condition)
    have h_gap : ∀ p, val p ≤ 1 + ∑ q ∈ (h_fin (val p)).toFinset, val q := by
      intro p
      exact erdos_gap_condition hd hsum p (h_fin (val p))
    -- Verify h_div (from erdos_val_diverges)
    have h_div : ∀ M, ∃ s : Finset (Fin k × ℕ), M ≤ ∑ p ∈ s, val p :=
      erdos_val_diverges d hd hk_pos
    -- Apply complete_sequence_of_gap
    obtain ⟨s, hs⟩ := complete_sequence_of_gap val h_pos h_one h_fin h_gap h_div n
    -- Define a : Fin k → ℕ by grouping the sum by first coordinate
    let a : Fin k → ℕ := fun i => ∑ p ∈ s.filter (fun p => p.1 = i), d i ^ p.2
    use a
    intro i
    constructor
    · -- Show digits property: a i = sum of distinct powers of d i
      let exp_set := (s.filter (fun p => p.1 = i)).image Prod.snd
      have hinj : ∀ p1 ∈ s.filter (fun p => p.1 = i), ∀ p2 ∈ s.filter (fun p => p.1 = i),
          p1.2 = p2.2 → p1 = p2 := by
        intro p1 hp1 p2 hp2 heq
        simp only [Finset.mem_filter] at hp1 hp2
        exact Prod.ext (hp1.2.trans hp2.2.symm) heq
      have ha_eq : a i = ∑ j ∈ exp_set, d i ^ j := by
        simp only [a, exp_set]
        rw [Finset.sum_image hinj]
      rw [ha_eq]
      have hb : 1 < d i := Nat.lt_of_lt_of_le one_lt_two (hd i)
      exact digits_sum_distinct_pow hb exp_set
    · -- Show n = ∑ i, a i
      -- n = ∑ p ∈ s, val p = ∑ p ∈ s, d p.1 ^ p.2 = ∑ i, ∑ p ∈ s.filter (·.1 = i), d i ^ p.2 = ∑ i, a i
      have key : ∑ p ∈ s, d p.1 ^ p.2 = ∑ i : Fin k, ∑ p ∈ s.filter (fun p => p.1 = i), d p.1 ^ p.2 := by
        rw [← Finset.sum_fiberwise s Prod.fst (fun p => d p.1 ^ p.2)]
      calc n = ∑ p ∈ s, val p := hs
        _ = ∑ p ∈ s, d p.1 ^ p.2 := by rfl
        _ = ∑ i : Fin k, ∑ p ∈ s.filter (fun p => p.1 = i), d p.1 ^ p.2 := key
        _ = ∑ i : Fin k, ∑ p ∈ s.filter (fun p => p.1 = i), d i ^ p.2 := by
            congr 1 with i
            apply Finset.sum_congr rfl
            intro p hp
            simp only [Finset.mem_filter] at hp
            rw [hp.2]
        _ = ∑ i, a i := by rfl