Browse chapters

Imports

import CLRSLean.FourthEdition.Chapter_03.Section_03_1_Asymptotic_Notation
import CLRSLean.FourthEdition.Chapter_07.Section_07_1_Description_Of_Quicksort
import CLRSLean.Probability.FiniteExpectation
import Mathlib
import Mathlib.NumberTheory.Harmonic.Bounds

7.3. Randomized Quicksort

This section defines the expected-comparison recurrence for randomized quicksort (CLRS equation (7.4)) and proves its closed-form solution, giving the O(n log n) average-case bound for the first time in CLRS-Lean.

The expected number of comparisons expectedComparisons n = E[T(n)] satisfies:

T(0) = 0, T(1) = 0
For n >= 1: T(n) = n-1 + (2/n) * sum_{k=0}^{n-1} T(k)

The closed form is T(n) = 2(n+1)H_n - 4n where H_n is the n-th harmonic number. This yields T(n) <= 2n H_n and T(n) <= n^2 (quadratic fallback).

Main results:

Lemma harmonic_succ: recurrence for harmonic numbers
Lemma harmonic_le_n: H_n <= n
Lemma sum_mul_harmonic_eq: sum_{k=1}^{n} k H_k = n(n+1)/2 H_n - n(n-1)/4
Lemma sum_expectedComparisons_eq: closed form of sum_{k=0}^{n-1} T(k)
Theorem expectedComparisons_closed_form: named CLRS closed-form formula
Theorem expectedComparisons_recurrence: closed form satisfies CLRS (7.4)
Theorem expectedComparisons_telescope: (n+1)T(n+1) = (n+2)T(n) + 2n
Theorem expectedComparisons_clrs_harmonic_bound: T(n) <= 2(n+1)H_n
Theorem expectedComparisons_harmonic_bound: T(n) <= 2n H_n
Theorem expectedComparisons_quadratic: T(n) <= n^2
Theorem expectedComparisons_monotone: T(n) <= T(n+1)

Implementation details

The detailed probability proof remains available outside the main sidebar:

7.3 Randomized Quicksort: Comparison Probability

Notation conventions:

harmonic n : H_n, the n-th harmonic number in Q
expectedComparisons n : T(n), expected number of comparisons for randomized quicksort on n distinct elements

              namespace CLRSnamespace Chapter07open Chapter07

Harmonic numbers

The n-th harmonic number as a rational. H_0 = 0, H_{n+1} = H_n + 1/(n+1).

def harmonic : Nat → Rat
  | 0 => 0
  | n+1 => harmonic n + 1 / ((n+1 : Nat) : Rat)

@[simp]
theorem harmonic_zero : harmonic 0 = 0 := rfl

@[simp]
theorem harmonic_one : harmonic 1 = 1 := by
  simp [harmonic]

Recurrence for harmonic numbers: H_{n+1} = H_n + 1/(n+1).

theorem harmonic_succ (n : Nat) : harmonic (n+1) = harmonic n + (1 : Rat) / ((n+1 : Nat) : Rat) :=
  rfl

Harmonic numbers are nonnegative.


theorem harmonic_nonneg (n : Nat) : 0 ≤ harmonic n := by
  induction n with
  | zero => simp
  | succ n ih =>
      rw [harmonic_succ]
      have hpos : 0 ≤ (1 : Rat) / ((n+1 : Nat) : Rat) := by
        positivity
      nlinarith

The harmonic number is bounded by its index: H_n <= n for all n.

This trivial bound is enough for many estimates.


theorem harmonic_le_n (n : Nat) : harmonic n ≤ (n : Rat) := by
  induction n with
  | zero => simp
  | succ n ih =>
      rw [harmonic_succ]
      push_cast
      have hdiv : (1 : Rat) / ((n : Rat) + 1) ≤ 1 :=
        (div_le_one (by positivity)).mpr (by nlinarith)
      nlinarith

Expected comparisons: closed form

Expected number of comparisons in randomized quicksort on n distinct elements, given by the closed-form solution of CLRS recurrence (7.4):

T(n) = 2(n+1)H_n - 4n

where H_n is the n-th harmonic number. This is a computable deterministic rational function; the expectation is folded into the recurrence coefficients in this core file. The companion ExplicitRandomness modules construct the finite permutation sample space and prove its executable pair-trace-count expectation equals this function.

def expectedComparisons (n : Nat) : Rat :=
  2 * ((n : Rat) + 1) * harmonic n - 4 * (n : Rat)

Named CLRS closed form for randomized-quicksort expected comparisons.

theorem expectedComparisons_closed_form (n : Nat) :
    expectedComparisons n = 2 * ((n : Rat) + 1) * harmonic n - 4 * (n : Rat) :=
  rfl

@[simp]
theorem expectedComparisons_zero : expectedComparisons 0 = 0 := by
  simp [expectedComparisons, harmonic]

@[simp]
theorem expectedComparisons_one : expectedComparisons 1 = 0 := by
  simp [expectedComparisons, harmonic]
  ring

Explicit formula for expectedComparisons (n+1) in terms of harmonic (n+1).

theorem expectedComparisons_succ (n : Nat) :
    expectedComparisons (n+1) = 2 * ((n+1 : Rat) + 1) * harmonic (n+1) - 4 * ((n+1 : Rat)) := by
  simp [expectedComparisons]

Key combinatorial identity - sum of k times harmonic k

Central combinatorial identity for the expected-quicksort closed form:

sum_{k=1}^{n} k * H_k = (n(n+1)/2) * H_n - n(n-1)/4

This is proved by induction on n using the harmonic recurrence to express H_n in terms of H_{n+1} in the inductive step.


theorem sum_mul_harmonic_eq (n : Nat) :
    (∑ k ∈ Finset.Icc 1 n, ((k : Rat) * harmonic k)) =
    (((n : Rat) * ((n : Rat) + 1)) / 2) * harmonic n - ((n : Rat) * ((n : Rat) - 1) / 4) := by
  induction n with
  | zero =>
      simp [harmonic]
  | succ n ih =>
      rw [Finset.sum_Icc_succ_top (by omega) (fun k => (k : Rat) * harmonic k)]
      rw [ih]
      -- Now: (n(n+1)/2)*H_n - n(n-1)/4 + (n+1)*H_{n+1} = ((n+1)(n+2)/2)*H_{n+1} - (n+1)n/4
      -- Use H_n = H_{n+1} - 1/(n+1)
      have hH_n : harmonic n = harmonic (n+1) - (1 : Rat) / ((n+1 : Nat) : Rat) := by
        rw [harmonic_succ]
        ring
      rw [hH_n]
      push_cast
      ring_nf
      have hpos : ((n : Nat) : Rat) + 1 ≠ 0 := by
        intro hzero
        have hsum : ((n+1 : Nat) : Rat) = 0 := by push_cast; simpa using hzero
        exact Nat.succ_ne_zero n (by exact_mod_cast hsum)
      field_simp [hpos]
      ring

Sum of expected comparisons

Closed form for the sum of expected comparisons up to n-1:

sum_{k=0}^{n-1} T(k) = n(n+1)*H_n - (5 n^2 - n)/2


theorem sum_expectedComparisons_eq (n : Nat) :
    (∑ k ∈ Finset.range n, expectedComparisons k) =
    ((n : Rat) * ((n : Rat) + 1)) * harmonic n - ((5 : Rat) * (n : Rat) * (n : Rat) - (n : Rat)) / 2 := by
  induction n with
  | zero => simp
  | succ n ih =>
      rw [Finset.sum_range_succ, expectedComparisons, ih]
      have hH_succ : harmonic (n+1) = harmonic n + (1 : Rat) / ((n+1 : Nat) : Rat) := harmonic_succ n
      rw [hH_succ]
      push_cast
      ring_nf
      have hpos : ((n : Nat) : Rat) + 1 ≠ 0 := by
        intro hzero
        have hsum : ((n+1 : Nat) : Rat) = 0 := by push_cast; simpa using hzero
        exact Nat.succ_ne_zero n (by exact_mod_cast hsum)
      field_simp [hpos]
      ring

Recurrence verification

The closed-form expectedComparisons satisfies the CLRS expected-comparison recurrence (7.4): for n >= 1,

T(n) = n-1 + (2/n) * sum_{k=0}^{n-1} T(k).

The proof multiplies through by n and uses the closed form of the sum.


theorem expectedComparisons_recurrence (n : Nat) (hn : n ≥ 1) :
    expectedComparisons n = ((n : Rat) - 1) + (2 / (n : Rat)) *
      (∑ k ∈ Finset.range n, expectedComparisons k) := by
  have hnpos : (n : Rat) ≠ 0 := by
    intro hzero
    have : n = 0 := by exact_mod_cast hzero
    omega
  -- Clear denominator by multiplying both sides by n
  field_simp [hnpos]
  -- Goal: n * T(n) = n * (n-1) + 2 * S(n)
  rw [sum_expectedComparisons_eq n]
  rw [expectedComparisons]
  ring

Alternative form of the recurrence, clearing denominators:

(n+1) * T(n+1) = (n+2) * T(n) + 2n for all n >= 0.

This telescoping identity is the key to the closed form and is used in the inductive proofs below.


theorem expectedComparisons_telescope (n : Nat) :
    ((n+1 : Nat) : Rat) * expectedComparisons (n+1) =
    (((n : Rat) + 2)) * expectedComparisons n + 2 * (n : Rat) := by
  rw [expectedComparisons, expectedComparisons]
  have hH_succ : harmonic (n+1) = harmonic n + (1 : Rat) / ((n+1 : Nat) : Rat) := harmonic_succ n
  rw [hH_succ]
  push_cast
  ring_nf
  have hpos : ((n : Nat) : Rat) + 1 ≠ 0 := by
    intro hzero
    have hsum : ((n+1 : Nat) : Rat) = 0 := by push_cast; simpa using hzero
    exact Nat.succ_ne_zero n (by exact_mod_cast hsum)
  field_simp [hpos]
  ring

Expected comparisons: nonnegativity

Expected comparisons are nonnegative.


theorem expectedComparisons_nonneg (n : Nat) : 0 ≤ expectedComparisons n := by
  induction n with
  | zero => simp
  | succ n ih =>
      have ht := expectedComparisons_telescope n
      -- ht: (n+1)*T(n+1) = (n+2)*T(n) + 2n
      -- RHS >= 0 since T(n) >= 0 and n >= 0, and (n+1) > 0 so T(n+1) >= 0
      have hpos_denom : ((n+1 : Nat) : Rat) ≠ 0 :=
        Nat.cast_ne_zero.mpr (Nat.succ_ne_zero n)
      have hnum_nonneg : 0 ≤ (((n : Rat) + 2)) * expectedComparisons n + 2 * (n : Rat) := by
        nlinarith
      -- From ht: T(n+1) = numerator / (n+1)
      have hT_expr : expectedComparisons (n+1) =
          ((((n : Rat) + 2)) * expectedComparisons n + 2 * (n : Rat)) / ((n+1 : Nat) : Rat) :=
        (eq_div_iff_mul_eq hpos_denom).mpr (by
          -- Need: T(n+1) * (n+1) = numerator
          -- ht gives: (n+1) * T(n+1) = numerator
          simpa [mul_comm] using ht)
      rw [hT_expr]
      refine div_nonneg hnum_nonneg (by positivity)

Bounds

Harmonic upper bound. The expected number of comparisons in randomized quicksort is at most 2 n * H_n.

Since H_n = Theta(log n), this gives T(n) = O(n log n).


theorem expectedComparisons_harmonic_bound (n : Nat) :
    expectedComparisons n ≤ 2 * (n : Rat) * harmonic n := by
  have hle : harmonic n ≤ (n : Rat) := harmonic_le_n n
  rw [expectedComparisons]
  nlinarith

CLRS-facing harmonic upper bound using the closed-form scale 2(n+1)H_n.


theorem expectedComparisons_clrs_harmonic_bound (n : Nat) :
    expectedComparisons n ≤ 2 * ((n : Rat) + 1) * harmonic n := by
  rw [expectedComparisons_closed_form]
  have hn : 0 ≤ (4 : Rat) * (n : Rat) := by positivity
  nlinarith

Quadratic upper bound. On any input of length n, the expected number of comparisons is at most n^2.

The proof uses induction with the telescope identity: T(n+1) = ((n+2)T(n) + 2n)/(n+1). The inductive hypothesis T(n) <= n^2 and a simple polynomial inequality n^2 + n + 1 >= 0 close the step.


theorem expectedComparisons_quadratic (n : Nat) :
    expectedComparisons n ≤ (n : Rat) * (n : Rat) := by
  induction n with
  | zero => simp
  | succ n ih =>
      have ht := expectedComparisons_telescope n
      -- ht: (n+1)*T(n+1) = (n+2)*T(n) + 2n
      have hpos : ((n+1 : Nat) : Rat) ≠ 0 :=
        Nat.cast_ne_zero.mpr (Nat.succ_ne_zero n)
      -- From ht: T(n+1) = ((n+2)*T(n) + 2n) / (n+1)
      have hT_succ : expectedComparisons (n+1) =
          ((((n : Rat) + 2)) * expectedComparisons n + 2 * (n : Rat)) / ((n+1 : Nat) : Rat) :=
        (eq_div_iff_mul_eq hpos).mpr (by
          simpa [mul_comm] using ht)
      rw [hT_succ]
      -- Need: ((n+2)*T(n) + 2n) / (n+1) <= (n+1)^2
      -- First, bound the numerator using ih: T(n) <= n^2
      have hnum_bound : (((n : Rat) + 2)) * expectedComparisons n + 2 * (n : Rat) ≤
          ((n : Rat) + 1) * ((n : Rat) + 1) * ((n : Rat) + 1) := by
        -- (n+2)*T(n) + 2n <= (n+2)*n^2 + 2n = n^3 + 2n^2 + 2n
        -- <= n^3 + 3n^2 + 3n + 1 = (n+1)^3  (since n^2 + n + 1 >= 0)
        nlinarith
      -- Apply the division lemma: if a <= b and c > 0, then a/c <= b/c
      refine le_trans (div_le_div_of_nonneg_right hnum_bound (by positivity)) ?_
      -- Now need: (n+1)^3 / (n+1) <= (n+1)^2
      -- Since (n+1)^3 / (n+1) = (n+1)^2 exactly, this is equality
      push_cast
      have h_eq : ((n : Rat) + 1) * ((n : Rat) + 1) * ((n : Rat) + 1) / ((n : Rat) + 1) =
          ((n : Rat) + 1) * ((n : Rat) + 1) := by
        field_simp [show ((n : Rat) + 1) ≠ 0 from by positivity]
      exact h_eq.le

Monotonicity. The expected comparison count is non-decreasing: T(n) <= T(n+1).

From the telescope identity, T(n+1) - T(n) = (T(n) + 2n)/(n+1) >= 0.


theorem expectedComparisons_monotone (n : Nat) : expectedComparisons n ≤ expectedComparisons (n+1) := by
  have ht := expectedComparisons_telescope n
  -- ht: (n+1)*T(n+1) = (n+2)*T(n) + 2n
  -- Rearranged: (n+1)*(T(n+1) - T(n)) = T(n) + 2n
  -- Since T(n) >= 0, RHS >= 0, so T(n+1) - T(n) >= 0
  have hpos : ((n+1 : Nat) : Rat) ≠ 0 :=
    Nat.cast_ne_zero.mpr (Nat.succ_ne_zero n)
  have hnonneg : 0 ≤ expectedComparisons n := expectedComparisons_nonneg n
  have hdiff : expectedComparisons (n+1) - expectedComparisons n =
      (expectedComparisons n + 2 * (n : Rat)) / ((n+1 : Nat) : Rat) :=
    (eq_div_iff_mul_eq hpos).mpr (by
      -- Need: (T(n+1) - T(n)) * (n+1) = T(n) + 2n
      -- Start from ht: (n+1)*T(n+1) = (n+2)*T(n) + 2n
      calc
        (expectedComparisons (n+1) - expectedComparisons n) * ((n+1 : Nat) : Rat)
            = ((n+1 : Nat) : Rat) * expectedComparisons (n+1) -
              ((n+1 : Nat) : Rat) * expectedComparisons n := by ring
        _ = (((n : Rat) + 2) * expectedComparisons n + 2 * (n : Rat)) -
              ((n+1 : Nat) : Rat) * expectedComparisons n := by rw [ht]
        _ = expectedComparisons n + 2 * (n : Rat) := by push_cast; ring
      )
  have hdiff_nonneg : 0 ≤ expectedComparisons (n+1) - expectedComparisons n := by
    rw [hdiff]
    refine div_nonneg ?_ (by positivity)
    nlinarith
  linarith

Asymptotic Θ(n log n) bound

We now lift the harmonic upper bound to the textbook asymptotic statement T(n) = Θ(n log n) using the standard harmonic bounds log(n+1) ≤ H_n ≤ 1 + log n from Mathlib.

open Chapter03

The rational harmonic number defined in this section equals Mathlib's global harmonic number after casting to ℝ.


theorem harmonic_eq_mathlib_harmonic (n : ℕ) : (harmonic n : ℝ) = (_root_.harmonic n : ℝ) := by
  induction n with
  | zero => simp [harmonic, _root_.harmonic]
  | succ n ih =>
    rw [harmonic_succ, _root_.harmonic_succ]
    push_cast
    rw [ih]
    simp

Expected comparisons cast to ℝ, for use with the Chapter 3 asymptotic wrappers.

noncomputable def expectedComparisonsReal (n : ℕ) : ℝ := (expectedComparisons n : ℝ)

Cast of the harmonic upper bound to ℝ.

theorem expectedComparisons_harmonic_bound_real (n : ℕ) :
    expectedComparisonsReal n ≤ 2 * (n : ℝ) * (harmonic n : ℝ) := by
  have h := expectedComparisons_harmonic_bound n
  dsimp [expectedComparisonsReal]
  exact_mod_cast h

Lower bound. For n ≥ 1, the expected number of comparisons is at least n * H_n - 4n.


theorem expectedComparisons_lower_bound_real (n : ℕ) (_hn : 1 ≤ n) :
    (n : ℝ) * (harmonic n : ℝ) - 4 * (n : ℝ) ≤ expectedComparisonsReal n := by
  dsimp [expectedComparisonsReal, expectedComparisons]
  push_cast
  have h_nonneg : 0 ≤ (harmonic n : ℝ) := by exact mod_cast harmonic_nonneg n
  nlinarith

The local harmonic after casting is bounded by 1 + log n.


theorem harmonic_le_one_add_log' (n : ℕ) : (harmonic n : ℝ) ≤ 1 + Real.log (n : ℝ) := by
  rw [harmonic_eq_mathlib_harmonic n]
  exact harmonic_le_one_add_log n

The local harmonic after casting is bounded below by log (n+1).


theorem log_add_one_le_harmonic' (n : ℕ) : Real.log ((n : ℝ) + 1) ≤ (harmonic n : ℝ) := by
  rw [harmonic_eq_mathlib_harmonic n]
  simpa [Nat.cast_add, Nat.cast_one] using log_add_one_le_harmonic n

Randomized quicksort is O(n log n). The expected number of comparisons satisfies T(n) = O(n log n).


theorem expectedComparisons_isBigO_nlogn :
    isBigO expectedComparisonsReal (fun n : ℕ => (n : ℝ) * Real.log (n : ℝ)) := by
  rw [isBigO_iff]
  have h_harm_le : ∀ n : ℕ, (harmonic n : ℝ) ≤ 1 + Real.log (n : ℝ) := harmonic_le_one_add_log'
  -- Real.log → ∞, so eventually log n ≥ 1
  have h_log_eventually : ∀ᶠ (n : ℕ) in Filter.atTop, (1 : ℝ) ≤ Real.log (n : ℝ) :=
    (Real.tendsto_log_atTop.comp tendsto_natCast_atTop_atTop) (Filter.eventually_ge_atTop (1 : ℝ))
  rcases Filter.eventually_atTop.mp h_log_eventually with ⟨n₁, hn₁⟩
  refine ⟨4, by norm_num, max 2 n₁, fun n hn => ?_⟩
  have hn2 : 2 ≤ n := le_trans (le_max_left _ _) hn
  have hn_log_ge_one : (1 : ℝ) ≤ Real.log (n : ℝ) := hn₁ n (le_trans (le_max_right _ _) hn)
  have hn_pos : 1 ≤ n := by omega
  have hn_real_pos : 1 ≤ (n : ℝ) := by exact_mod_cast hn_pos
  have hlog_nonneg : 0 ≤ Real.log (n : ℝ) := Real.log_nonneg hn_real_pos
  have hT_nonneg : 0 ≤ expectedComparisonsReal n := by
    dsimp [expectedComparisonsReal]; exact mod_cast expectedComparisons_nonneg n
  have hmul_nonneg : 0 ≤ (n : ℝ) * Real.log (n : ℝ) := by positivity
  rw [abs_of_nonneg hT_nonneg, abs_of_nonneg hmul_nonneg]
  calc
    expectedComparisonsReal n ≤ 2 * (n : ℝ) * (harmonic n : ℝ) := expectedComparisons_harmonic_bound_real n
    _ ≤ 2 * (n : ℝ) * (1 + Real.log (n : ℝ)) := by
      have h_nonneg : 0 ≤ 2 * (n : ℝ) := by positivity
      gcongr
      exact h_harm_le n
    _ = 2 * (n : ℝ) + 2 * ((n : ℝ) * Real.log (n : ℝ)) := by ring
    _ ≤ 4 * ((n : ℝ) * Real.log (n : ℝ)) := by
      -- 2n ≤ 2n*log n  when log n ≥ 1, so 2n + 2n*log n ≤ 4n*log n
      have h : 2 * (n : ℝ) ≤ 2 * ((n : ℝ) * Real.log (n : ℝ)) := by
        have hn_nonneg : 0 ≤ (n : ℝ) := Nat.cast_nonneg _
        calc
          2 * (n : ℝ) = 2 * (n : ℝ) * (1 : ℝ) := by ring
          _ ≤ 2 * (n : ℝ) * Real.log (n : ℝ) := by gcongr
          _ = 2 * ((n : ℝ) * Real.log (n : ℝ)) := by ring
      nlinarith

Randomized quicksort is Ω(n log n). The expected number of comparisons satisfies T(n) = Ω(n log n).


theorem expectedComparisons_isBigOmega_nlogn :
    isBigOmega expectedComparisonsReal (fun n : ℕ => (n : ℝ) * Real.log (n : ℝ)) := by
  rw [isBigOmega_iff]
  -- Use log(n+1) ≤ H_n and T(n) ≥ n*H_n - 4n
  have h_harm_lower : ∀ n : ℕ, Real.log ((n : ℝ) + 1) ≤ (harmonic n : ℝ) :=
    log_add_one_le_harmonic'
  -- Real.log → ∞, so eventually log n ≥ 8
  have h_log_eventually : ∀ᶠ (n : ℕ) in Filter.atTop, (8 : ℝ) ≤ Real.log (n : ℝ) :=
    (Real.tendsto_log_atTop.comp tendsto_natCast_atTop_atTop) (Filter.eventually_ge_atTop (8 : ℝ))
  rcases Filter.eventually_atTop.mp h_log_eventually with ⟨n₀₁, hn₀₁⟩
  -- Also need log(n+1) ≥ (1/2)*log n for large n
  -- Since log(n+1)/log n → 1, for log n ≥ 8 we have log(n+1) ≥ (3/4)*log n
  -- Actually log(n+1) ≥ log n ≥ (1/2)*log n trivially
  let n₀ := max n₀₁ 8
  refine ⟨1/8, by norm_num, n₀, fun n hn => ?_⟩
  have hn₁ : n₀₁ ≤ n := le_trans (le_max_left _ _) hn
  have hn_pos : 8 ≤ n := le_trans (le_max_right _ _) hn
  have hn_real_pos : 0 < (n : ℝ) := by
    have : 0 < n := by omega
    exact_mod_cast this
  have hT_nonneg : 0 ≤ expectedComparisonsReal n := by
    dsimp [expectedComparisonsReal]; exact mod_cast expectedComparisons_nonneg n
  have hn1pos : 1 ≤ n := by omega
  have hn1real : 1 ≤ (n : ℝ) := by exact_mod_cast hn1pos
  have hlog_nonneg : 0 ≤ Real.log (n : ℝ) := Real.log_nonneg hn1real
  have hmul_nonneg : 0 ≤ (n : ℝ) * Real.log (n : ℝ) := by positivity
  rw [abs_of_nonneg hT_nonneg, abs_of_nonneg hmul_nonneg]
  have h_log_ge_eight : (8 : ℝ) ≤ Real.log (n : ℝ) := hn₀₁ n hn₁
  -- log(n+1) ≥ log n ≥ 8 for n ≥ n₀
  have h_log_succ_ge : Real.log (n : ℝ) ≤ Real.log ((n : ℝ) + 1) :=
    Real.log_le_log (by positivity) (by nlinarith)
  -- T(n) ≥ n*H_n - 4n ≥ n*log(n+1) - 4n ≥ n*log n - 4n
  -- Since log n ≥ 8, we have log n/8 ≥ 1, so 4n ≤ (log n/2)*n = n*log n/2
  -- Thus n*log n - 4n ≥ n*log n/2 ≥ n*log n/8
  calc
    (1/8 : ℝ) * ((n : ℝ) * Real.log (n : ℝ)) = ((n : ℝ) * Real.log (n : ℝ)) / 8 := by ring
    _ ≤ ((n : ℝ) * Real.log (n : ℝ)) - 4 * (n : ℝ) := by
      -- Need: (n*log n)/8 ≤ n*log n - 4n  ⇔  4n ≤ (7/8)*n*log n  ⇔  32/7 ≤ log n ≈ 4.57
      -- Since log n ≥ 8, this holds.
      have h : 4 * (n : ℝ) ≤ (7/8 : ℝ) * ((n : ℝ) * Real.log (n : ℝ)) := by
        calc
          4 * (n : ℝ) = (n : ℝ) * 4 := by ring
          _ ≤ (n : ℝ) * ((7/8 : ℝ) * Real.log (n : ℝ)) := by
            nlinarith [h_log_ge_eight]
          _ = (7/8 : ℝ) * ((n : ℝ) * Real.log (n : ℝ)) := by ring
      nlinarith
    _ ≤ (n : ℝ) * Real.log ((n : ℝ) + 1) - 4 * (n : ℝ) := by nlinarith
    _ ≤ (n : ℝ) * (harmonic n : ℝ) - 4 * (n : ℝ) := by nlinarith [h_harm_lower n]
    _ ≤ expectedComparisonsReal n := expectedComparisons_lower_bound_real n hn1pos

Randomized quicksort is Θ(n log n). The expected number of comparisons satisfies T(n) = Θ(n log n).

theorem expectedComparisons_isBigTheta_nlogn :
    isBigTheta expectedComparisonsReal (fun n : ℕ => (n : ℝ) * Real.log (n : ℝ)) :=
  ⟨expectedComparisons_isBigO_nlogn, expectedComparisons_isBigOmega_nlogn⟩

Bridge: probability model to closed form

We connect the random-permutation pairwise comparison probability (compared_prob, CLRS Theorem 7.3) to the deterministic closed form expectedComparisons n and the Θ(n log n) asymptotic.

open CLRS.Probability

Additive recurrence: T(n+1) = T(n) + 2*(H_{n+1} - 1).


theorem expectedComparisons_succ_add_two (n : ℕ) :
    expectedComparisons (n+1) = expectedComparisons n + 2 * (harmonic (n+1) - 1) := by
  have ht := expectedComparisons_telescope n
  have hpos : ((n+1 : ℕ) : ℚ) ≠ 0 := Nat.cast_ne_zero.mpr (Nat.succ_ne_zero n)
  have hT_succ : expectedComparisons (n+1) =
      (((n : ℚ) + 2) * expectedComparisons n + 2 * (n : ℚ)) / ((n+1 : ℕ) : ℚ) :=
    (eq_div_iff_mul_eq hpos).mpr (by simpa [mul_comm] using ht)
  rw [hT_succ]
  rw [show expectedComparisons n = 2 * ((n : ℚ) + 1) * harmonic n - 4 * (n : ℚ) from rfl]
  have hH_succ : harmonic (n+1) = harmonic n + (1 : ℚ) / ((n+1 : ℕ) : ℚ) := harmonic_succ n
  rw [hH_succ]
  push_cast
  field_simp [show ((n : ℚ) + 1) ≠ 0 from by positivity]
  ring

The double sum of pairwise comparison probabilities 2/(j-i+1) over all 0 ≤ i < j < n equals the expected-comparison closed form.


theorem sum_compared_prob_eq_expectedComparisons (n : ℕ) :
    (∑ i ∈ Finset.range n, ∑ j ∈ Finset.range n,
      if i < j then (2 : ℚ) / ((j - i + 1 : ℕ) : ℚ) else 0) =
    (expectedComparisons n : ℚ) := by
  induction n with
  | zero => simp [expectedComparisons]
  | succ n ih =>
    -- S(n+1) = S(n) + A(n), where A(n) = Σ_{i<n} 2/(n-i+1)
    -- Split the outer sum: i=n contributes nothing (n < j never holds in range (n+1))
    rw [Finset.sum_range_succ]
    have h_last_row_zero : (∑ j ∈ Finset.range (n+1),
        if (n : ℕ) < j then (2 : ℚ) / ((j - n + 1 : ℕ) : ℚ) else 0) = 0 := by
      apply Finset.sum_eq_zero; intro j hj
      rw [Finset.mem_range] at hj
      simp [show ¬ (n : ℕ) < j from by omega]
    rw [h_last_row_zero, add_zero]
    -- For i < n, split inner sum at j = n (the new column)
    have h_inner_split : (∑ i ∈ Finset.range n, ∑ j ∈ Finset.range (n+1),
        if i < j then (2 : ℚ) / ((j - i + 1 : ℕ) : ℚ) else 0) =
        (∑ i ∈ Finset.range n, ∑ j ∈ Finset.range n,
          if i < j then (2 : ℚ) / ((j - i + 1 : ℕ) : ℚ) else 0) +
        (∑ i ∈ Finset.range n, (2 : ℚ) / (((n : ℕ) - i + 1 : ℕ) : ℚ)) := by
      calc
        (∑ i ∈ Finset.range n, ∑ j ∈ Finset.range (n+1),
          if i < j then (2 : ℚ) / ((j - i + 1 : ℕ) : ℚ) else 0)
            = (∑ i ∈ Finset.range n,
                ((∑ j ∈ Finset.range n, if i < j then (2 : ℚ) / ((j - i + 1 : ℕ) : ℚ) else 0) +
                  (if i < n then (2 : ℚ) / (((n : ℕ) - i + 1 : ℕ) : ℚ) else 0))) := by
          refine Finset.sum_congr rfl (fun i hi => ?_)
          rw [Finset.sum_range_succ]
        _ = (∑ i ∈ Finset.range n, ∑ j ∈ Finset.range n,
              if i < j then (2 : ℚ) / ((j - i + 1 : ℕ) : ℚ) else 0) +
            (∑ i ∈ Finset.range n,
              (if i < n then (2 : ℚ) / (((n : ℕ) - i + 1 : ℕ) : ℚ) else 0)) := by
          rw [Finset.sum_add_distrib]
        _ = (∑ i ∈ Finset.range n, ∑ j ∈ Finset.range n,
              if i < j then (2 : ℚ) / ((j - i + 1 : ℕ) : ℚ) else 0) +
            (∑ i ∈ Finset.range n, (2 : ℚ) / (((n : ℕ) - i + 1 : ℕ) : ℚ)) := by
          congr 1
          apply Finset.sum_congr rfl; intro i hi
          have hi_lt_n : i < n := Finset.mem_range.1 hi
          simp [hi_lt_n]
    rw [h_inner_split, ih]
    rw [expectedComparisons_succ_add_two n]
    congr 1
    -- Prove A(n) = 2*(H_{n+1} - 1) using the same recurrence
    -- A(0) = 0, A(n+1) = A(n) + 2/(n+2)
    -- Both sides satisfy this recurrence
    have hA_recurrence : ∀ m, (∑ i ∈ Finset.range m, (2 : ℚ) / (((m : ℕ) - i + 1 : ℕ) : ℚ)) =
        2 * (harmonic (m+1) - 1) := by
      intro m
      induction m with
      | zero => simp [harmonic]
      | succ m ih =>
        -- A(m+1) = Σ_{i∈range(m+1)} 2/((m+1)-i+1)
        -- Decompose: i=0 term = 2/(m+2), remaining shifted by i↦i+1
        have h_decomp : (Finset.range (m+1) : Finset ℕ) =
            ({0} : Finset ℕ) ∪ ((Finset.range m).map ⟨(· + 1), Nat.succ_injective⟩) := by
          ext i; constructor
          · intro hi
            have hi_val : i < m+1 := Finset.mem_range.1 hi
            rcases Nat.eq_zero_or_pos i with (rfl | hpos)
            · apply Finset.mem_union_left; simp
            · apply Finset.mem_union_right
              apply Finset.mem_map.mpr
              have h_bound : i - 1 < m := by omega
              refine ⟨i-1, Finset.mem_range.2 h_bound, ?_⟩
              have h_one_le : 1 ≤ i := Nat.one_le_of_lt hpos
              dsimp
              rw [Nat.sub_add_cancel h_one_le]
          · intro hi
            rcases Finset.mem_union.1 hi with (h | h)
            · rcases Finset.mem_singleton.1 h with rfl
              exact Finset.mem_range.2 (by
                have : 0 < m+1 := Nat.zero_lt_succ m
                exact this)
            · rcases Finset.mem_map.1 h with ⟨j, hj, rfl⟩
              have hj_val : j < m := Finset.mem_range.1 hj
              have : j+1 < m+1 := Nat.add_lt_add_right hj_val 1
              exact Finset.mem_range.2 this
        have h_disjoint : Disjoint ({0} : Finset ℕ) ((Finset.range m).map ⟨(· + 1), Nat.succ_injective⟩) := by
          refine Finset.disjoint_singleton_left.mpr (fun h => ?_)
          rcases Finset.mem_map.1 h with ⟨j, hj, h⟩
          have : j + 1 = 0 := h
          omega
        rw [h_decomp, Finset.sum_union h_disjoint,
          Finset.sum_singleton, Finset.sum_map]
        -- Now: 2/((m+1)-0+1) + Σ_{j∈range m} 2/((m+1)-(j+1)+1)
        -- = 2/(m+2) + Σ_{j∈range m} 2/(m-j+1)
        -- = 2/(m+2) + A(m)
        simp only [Function.Embedding.coeFn_mk]
        have h0 : ((m+1 : ℕ) - 0 + 1 : ℕ) = (m+2 : ℕ) := by omega
        have h_shift : ∀ j, ((m+1 : ℕ) - (j+1) + 1 : ℕ) = ((m : ℕ) - j + 1 : ℕ) := by
          intro j; omega
        simp_rw [h0, h_shift]
        rw [ih]
        rw [harmonic_succ (m+1)]
        push_cast; ring
    exact hA_recurrence n

end Chapter07end CLRS

Definitions and proofs

CLRSLean.FourthEdition.Chapter_07.Section_07_3_Randomized_Quicksort.Comparison_Probability

open Finsetopen Classical

Random Permutation Symmetry Lemma

Core: for a uniform random permutation of Fin n, each element of a subset has equal probability 1/|S| of appearing first. Proved via the transposition bijection π ↦ swap(a,b) * π.

              namespace CLRSnamespace Chapter07variable {n : ℕ}

"Position" of element x in permutation π = π.symm x (the index mapping to x).

def pos (π : Equiv.Perm (Fin n)) (x : Fin n) : Fin n := π.symm x

IsFirstIn S x π means x has minimal position among S in π.

def IsFirstIn (S : Finset (Fin n)) (x : Fin n) (π : Equiv.Perm (Fin n)) : Prop :=
  x ∈ S ∧ ∀ y ∈ S, pos π x ≤ pos π y

The swap bijection



lemma pos_swap_comp {a b : Fin n} (π : Equiv.Perm (Fin n)) (x : Fin n) :
    pos ((Equiv.swap a b) * π) x = pos π ((Equiv.swap a b) x) := by
  dsimp [pos]
  change ((Equiv.swap a b) * π)⁻¹ x = π⁻¹ ((Equiv.swap a b) x)
  rw [mul_inv_rev]
  rw [Equiv.swap_inv]
  rfl

Composing with swap a b on the LEFT is a bijection on Perm(Fin n).

lemma swapComp_bijective {a b : Fin n} :
    Function.Bijective (fun (π : Equiv.Perm (Fin n)) => (Equiv.swap a b) * π) := by
  constructor
  · intro π₁ π₂ h
    apply_fun (fun φ => (Equiv.swap a b).symm * φ) at h
    simpa [mul_assoc] using h
  · intro π
    refine ⟨(Equiv.swap a b).symm * π, ?_⟩
    simp

If a is first in S under π, then b is first under swap a b ∘ π.


lemma IsFirstIn_swap {S : Finset (Fin n)} {a b : Fin n} (Variable name `hne` is not explicitly referenced.

The binding can be removed (if unused) or named `_` (if used implicitly).

Note: This linter can be disabled with `set_option linter.unusedVariables false`hne : a ≠ b)
    (π : Equiv.Perm (Fin n)) (ha : IsFirstIn S a π) (hbS : b ∈ S) :
    IsFirstIn S b ((Equiv.swap a b) * π) := by
  rcases ha with ⟨haS, ha_min⟩
  let σ := Equiv.swap a b
  have hσa : σ a = b := by simp [σ]
  have hσb : σ b = a := by simp [σ]
  refine ⟨hbS, ?_⟩
  intro y hyS
  rw [pos_swap_comp π b, pos_swap_comp π y, hσb]
  -- Need: pos π a ≤ pos π (σ y)
  -- Since σ y ∈ S (σ permutes S), ha_min applies
  have hσy_S : σ y ∈ S := by
    -- σ swaps a and b, both in S, fixes others
    by_cases hya : y = a
    · subst y; simpa [σ]
    · by_cases hyb : y = b
      · subst y; simpa [σ]
      · have : σ y = y := Equiv.swap_apply_of_ne_of_ne hya hyb
        rw [this]; exact hyS
  exact ha_min (σ y) hσy_S

Equal cardinality via bijection

The sets of permutations where a is first vs where b is first have equal cardinality.


lemma card_firstSet_eq {S : Finset (Fin n)} {a b : Fin n}
    (haS : a ∈ S) (hbS : b ∈ S) (hne : a ≠ b) :
    ((Finset.univ : Finset (Equiv.Perm (Fin n))).filter (IsFirstIn S a)).card =
    ((Finset.univ : Finset (Equiv.Perm (Fin n))).filter (IsFirstIn S b)).card := by
  let σ := Equiv.swap a b
  -- The map f(π) = σ * π is a bijection that maps firstSet S a to firstSet S b
  apply Finset.card_bij (fun π _ => σ * π) ?_ ?_ ?_
  · -- f maps firstSet a into firstSet b
    intro π hπ
    rw [Finset.mem_filter] at hπ
    rcases hπ with ⟨hπu, ha⟩
    refine Finset.mem_filter.mpr ⟨Finset.mem_univ _, ?_⟩
    exact IsFirstIn_swap hne π ha hbS
  · -- f is injective (follows from bijectivity)
    intro π₁ _ π₂ _ h
    exact (swapComp_bijective (a:=a) (b:=b)).1 h
  · -- f is surjective onto firstSet b
    intro π hπ
    rw [Finset.mem_filter] at hπ
    rcases hπ with ⟨hπu, hb⟩
    -- preimage: σ * π (since σ⁻¹ = σ, left-multiplying by σ maps back)
    have h_pre : IsFirstIn S a (σ * π) := by
      -- Apply IsFirstIn_swap with a↔b swapped
      have h := IsFirstIn_swap (Ne.symm hne) π hb haS
      simpa [σ, Equiv.swap_comm] using h
    refine ⟨σ * π, Finset.mem_filter.mpr ⟨Finset.mem_univ _, h_pre⟩, ?_⟩
    simp [σ]

Partition: each permutation has exactly one first element

For a nonempty S, each permutation has a unique first element in S.


lemma existsUnique_firstInSet (S : Finset (Fin n)) (hSne : S.Nonempty) (π : Equiv.Perm (Fin n)) :
    ∃! s, s ∈ S ∧ IsFirstIn S s π := by
  -- The set {pos π s | s ∈ S} is a nonempty finite set of Fin n, so it has a minimum
  -- The s that maps to that minimum is unique (since π.symm is injective)
  let positions : Finset (Fin n) := S.image (pos π)
  have hpos_ne : positions.Nonempty := by
    rcases hSne with ⟨s, hs⟩
    exact ⟨pos π s, Finset.mem_image.mpr ⟨s, hs, rfl⟩⟩
  let p := positions.min' hpos_ne
  have hp_mem : p ∈ positions := Finset.min'_mem _ hpos_ne
  rcases Finset.mem_image.mp hp_mem with ⟨s, hsS, hsp⟩
  have h_first : IsFirstIn S s π := by
    refine ⟨hsS, ?_⟩
    intro y hyS
    have hy_pos : pos π y ∈ positions := Finset.mem_image.mpr ⟨y, hyS, rfl⟩
    have hp_le : p ≤ pos π y := Finset.min'_le _ _ hy_pos
    rw [← hsp] at hp_le
    exact hp_le
  refine ⟨s, ⟨hsS, h_first⟩, ?_⟩
  intro s' ⟨hs'S, hs'_first⟩
  have hpos_le : pos π s' ≤ pos π s := hs'_first.2 s hsS
  have hpos_ge : pos π s ≤ pos π s' := h_first.2 s' hs'S
  have hpos_eq : pos π s' = pos π s := le_antisymm hpos_le hpos_ge
  apply π.symm.injective
  exact hpos_eq

Main symmetry theorem

Symmetry Lemma. For a nonempty set S ⊆ Fin n and s ∈ S, under a uniform random permutation, P(s is first in S) = 1 / |S|.

Proof: All elements of S have equal probability by card_firstSet_eq. Since the sum over t in S of P(t first) = 1 (each π has exactly one first element), we get |S| × P(s first) = 1, hence P = 1 / |S|.


theorem isFirst_prob (S : Finset (Fin n)) (hSne : S.Nonempty) (s : Fin n) (hsS : s ∈ S) :
    ((Finset.filter (IsFirstIn S s) Finset.univ).card : ℝ) / (Nat.factorial n : ℝ) =
    1 / (S.card : ℝ) := by
  -- Let A_t = {π | IsFirstIn S t π}
  -- Step 1: All A_t have equal cardinality
  have h_eq_card : ∀ t ∈ S, ((Finset.filter (IsFirstIn S t) Finset.univ).card : ℝ) =
      ((Finset.filter (IsFirstIn S s) Finset.univ).card : ℝ) := by
    intro t htS
    by_cases hts : t = s
    · subst t; rfl
    · have hcard := card_firstSet_eq htS hsS hts
      -- This gives Nat equality; cast to ℝ
      exact_mod_cast hcard
  -- Step 2: The sets A_t for t ∈ S are pairwise disjoint
  have h_disjoint : ∀ t₁ t₂, t₁ ∈ S → t₂ ∈ S → t₁ ≠ t₂ →
      Disjoint (Finset.filter (IsFirstIn S t₁) Finset.univ)
               (Finset.filter (IsFirstIn S t₂) Finset.univ) := by
    intro t₁ t₂ ht₁ ht₂ hne
    apply Finset.disjoint_filter.2
    intro π _ h₁ h₂
    -- If π is in both filters, then both t₁ and t₂ are first in S for π
    -- This contradicts uniqueness from existsUnique_firstInSet
    rcases existsUnique_firstInSet S hSne π with ⟨t, ⟨htS, ht_first⟩, hunique⟩
    have heq : t₁ = t₂ := by
      have h1 := hunique t₁ ⟨ht₁, h₁⟩
      have h2 := hunique t₂ ⟨ht₂, h₂⟩
      exact h1.trans h2.symm
    exact hne heq
  -- Step 3: The sets cover Finset.univ
  have h_cover : (Finset.biUnion S (fun t => Finset.filter (IsFirstIn S t) Finset.univ)) =
      Finset.univ := by
    apply Finset.Subset.antisymm
    · exact Finset.subset_univ _
    · intro π hπ
      have hπu : π ∈ Finset.univ := Finset.mem_univ _
      rcases existsUnique_firstInSet S hSne π with ⟨t, ⟨htS, ht_first⟩, _⟩
      apply Finset.mem_biUnion.mpr
      exact ⟨t, htS, Finset.mem_filter.mpr ⟨hπu, ht_first⟩⟩
  -- Step 4: Sum of cardinalities = |univ| = n!
  have h_total_nat : (Finset.univ : Finset (Equiv.Perm (Fin n))).card = Nat.factorial n := by
    simp [Fintype.card_perm]
  have h_pairwise : (S : Set (Fin n)).PairwiseDisjoint
      (fun t => Finset.filter (IsFirstIn S t) Finset.univ) := by
    intro t₁ ht₁ t₂ ht₂ hne
    exact h_disjoint t₁ t₂ ht₁ ht₂ hne
  have h_sum_card : (∑ t ∈ S, (Finset.filter (IsFirstIn S t) Finset.univ).card) =
      Nat.factorial n := by
    calc
      (∑ t ∈ S, (Finset.filter (IsFirstIn S t) Finset.univ).card)
          = (Finset.biUnion S (fun t => Finset.filter (IsFirstIn S t) Finset.univ)).card := by
        rw [Finset.card_biUnion (h := h_pairwise)]
      _ = (Finset.univ : Finset (Equiv.Perm (Fin n))).card := by rw [h_cover]
      _ = Nat.factorial n := h_total_nat
  -- Step 5: Since all |A_t| are equal, S.card * |A_s| = n!
  -- Hence |A_s| / n! = 1 / S.card
  have h_total_real : (∑ t ∈ S, ((Finset.filter (IsFirstIn S t) Finset.univ).card : ℝ)) =
      (Nat.factorial n : ℝ) := by exact_mod_cast h_sum_card
  have h_all_eq : (∑ t ∈ S, ((Finset.filter (IsFirstIn S t) Finset.univ).card : ℝ)) =
      (S.card : ℝ) * ((Finset.filter (IsFirstIn S s) Finset.univ).card : ℝ) := by
    calc
      (∑ t ∈ S, ((Finset.filter (IsFirstIn S t) Finset.univ).card : ℝ))
          = (∑ t ∈ S, ((Finset.filter (IsFirstIn S s) Finset.univ).card : ℝ)) :=
        Finset.sum_congr rfl (fun t ht => by rw [h_eq_card t ht])
      _ = (S.card : ℝ) * ((Finset.filter (IsFirstIn S s) Finset.univ).card : ℝ) := by
        simp [Finset.sum_const, nsmul_eq_mul]
  rw [h_all_eq] at h_total_real
  have hS_card_ne_zero : (S.card : ℝ) ≠ 0 := by
    have hpos : 0 < S.card := Finset.card_pos.mpr hSne
    positivity
  have h_nfac_ne_zero : (Nat.factorial n : ℝ) ≠ 0 := by positivity
  field_simp [hS_card_ne_zero, h_nfac_ne_zero]
  linarith

Application to quicksort comparison probability

Convert i..j (all < n) to Finset (Fin n) via an embedding.

def rangeFin (n i j : ℕ) (_hij : i ≤ j) (hjn : j < n) : Finset (Fin n) :=
  ((Finset.Icc i j).attach).map
    ⟨fun ⟨k, hk⟩ =>
        ⟨k, lt_of_le_of_lt (Finset.mem_Icc.mp hk).2 hjn⟩,
     fun ⟨a, ha⟩ ⟨b, hb⟩ h => by
        apply Subtype.ext
        simpa using congrArg Fin.val h⟩

ranks i,j are compared iff first of i,...,j is i or j.

def comparedInQuicksort (n i j : ℕ) (hij : i < j) (hjn : j < n) (π : Equiv.Perm (Fin n)) : Prop :=
  have hi : i < n := lt_trans hij hjn
  let si : Fin n := ⟨i, hi⟩
  let sj : Fin n := ⟨j, hjn⟩
  let S := rangeFin n i j (Nat.le_of_lt hij) hjn
  IsFirstIn S si π ∨ IsFirstIn S sj π

|rangeFin n i j| = j-i+1.


lemma card_rangeFin (n i j : ℕ) (_hij : i ≤ j) (hjn : j < n) :
    (rangeFin n i j _hij hjn).card = j - i + 1 := by
  unfold rangeFin
  rw [Finset.card_map]
  simp
  omega

CLRS Theorem 7.3. P(compared) = 2/(j-i+1).


theorem compared_prob (n i j : ℕ) (hij : i < j) (hjn : j < n) :
    ((Finset.filter (comparedInQuicksort n i j hij hjn) Finset.univ).card : ℝ) /
      (Nat.factorial n : ℝ) = (2 : ℝ) / ((j - i + 1 : ℕ) : ℝ) := by
  have hi : i < n := lt_trans hij hjn
  have h_le : i ≤ j := Nat.le_of_lt hij
  let si : Fin n := ⟨i, hi⟩
  let sj : Fin n := ⟨j, hjn⟩
  let S := rangeFin n i j h_le hjn
  have hcard_S : S.card = j - i + 1 := card_rangeFin n i j h_le hjn
  have hS_si : si ∈ S := by
    dsimp [S, rangeFin]
    have hi_mem : i ∈ Finset.Icc i j := Finset.mem_Icc.mpr ⟨le_refl i, h_le⟩
    refine Finset.mem_map.mpr ⟨⟨i, hi_mem⟩, Finset.mem_attach _ _, ?_⟩
    ext; rfl
  have hS_sj : sj ∈ S := by
    dsimp [S, rangeFin]
    have hj_mem : j ∈ Finset.Icc i j := Finset.mem_Icc.mpr ⟨h_le, le_refl j⟩
    refine Finset.mem_map.mpr ⟨⟨j, hj_mem⟩, Finset.mem_attach _ _, ?_⟩
    ext; rfl
  have hSne : S.Nonempty := ⟨si, hS_si⟩
  have hne : si ≠ sj := by
    intro h
    have hval : i = j := Fin.ext_iff.mp h
    omega
  have h_disjoint : Disjoint
      (Finset.filter (IsFirstIn S si) Finset.univ)
      (Finset.filter (IsFirstIn S sj) Finset.univ) := by
    apply Finset.disjoint_filter.2
    intro π _ hsi hsj
    rcases hsi with ⟨_, hsi_min⟩; rcases hsj with ⟨_, hsj_min⟩
    have hle1 : pos π si ≤ pos π sj := hsi_min sj hS_sj
    have hle2 : pos π sj ≤ pos π si := hsj_min si hS_si
    have heq : pos π si = pos π sj := le_antisymm hle1 hle2
    apply hne; apply π.symm.injective; exact heq
  have h_si_prob : ((Finset.filter (IsFirstIn S si) Finset.univ).card : ℝ) /
      (Nat.factorial n : ℝ) = 1 / (S.card : ℝ) :=
    isFirst_prob S hSne si hS_si
  have h_sj_prob : ((Finset.filter (IsFirstIn S sj) Finset.univ).card : ℝ) /
      (Nat.factorial n : ℝ) = 1 / (S.card : ℝ) :=
    isFirst_prob S hSne sj hS_sj
  have h_union : (Finset.filter (comparedInQuicksort n i j hij hjn) Finset.univ) =
      (Finset.filter (IsFirstIn S si) Finset.univ) ∪
      (Finset.filter (IsFirstIn S sj) Finset.univ) := by
    ext π; constructor
    · intro h; rcases Finset.mem_filter.mp h with ⟨hu, hc⟩
      unfold comparedInQuicksort at hc
      rcases hc with (h | h)
      · exact Finset.mem_union_left _ (Finset.mem_filter.mpr ⟨hu, h⟩)
      · exact Finset.mem_union_right _ (Finset.mem_filter.mpr ⟨hu, h⟩)
    · intro h; rcases Finset.mem_union.mp h with (h' | h')
      · rcases Finset.mem_filter.mp h' with ⟨hu, h⟩
        refine Finset.mem_filter.mpr ⟨hu, ?_⟩
        unfold comparedInQuicksort
        exact Or.inl h
      · rcases Finset.mem_filter.mp h' with ⟨hu, h⟩
        refine Finset.mem_filter.mpr ⟨hu, ?_⟩
        unfold comparedInQuicksort
        exact Or.inr h
  rw [h_union, Finset.card_union_of_disjoint h_disjoint, Nat.cast_add]
  have h_nfac_ne_zero : (Nat.factorial n : ℝ) ≠ 0 := by positivity
  have h_card_ne_zero : (S.card : ℝ) ≠ 0 := by
    have hpos : 0 < S.card := Finset.card_pos.mpr hSne
    positivity
  -- Goal: (|A_si| + |A_sj|) / n! = 2 / (j-i+1)
  -- From h_si_prob: |A_si|/n! = 1/S.card
  -- From h_sj_prob: |A_sj|/n! = 1/S.card
  -- Add them: (|A_si| + |A_sj|)/n! = 2/S.card = 2/(j-i+1)
  have h_sum : ((Finset.filter (IsFirstIn S si) Finset.univ).card : ℝ) +
      ((Finset.filter (IsFirstIn S sj) Finset.univ).card : ℝ) =
      (2 : ℝ) * (Nat.factorial n : ℝ) / (S.card : ℝ) := by
    field_simp [h_nfac_ne_zero, h_card_ne_zero] at h_si_prob h_sj_prob ⊢
    linarith [h_si_prob, h_sj_prob]
  calc
    (((Finset.filter (IsFirstIn S si) Finset.univ).card : ℝ) +
      ((Finset.filter (IsFirstIn S sj) Finset.univ).card : ℝ)) / (Nat.factorial n : ℝ)
        = ((2 : ℝ) * (Nat.factorial n : ℝ) / (S.card : ℝ)) / (Nat.factorial n : ℝ) := by rw [h_sum]
    _ = (2 : ℝ) / (S.card : ℝ) := by
      field_simp [h_nfac_ne_zero]
    _ = (2 : ℝ) / ((j - i + 1 : ℕ) : ℝ) := by rw [hcard_S]

Layer boundary

This file stops at random-permutation symmetry and the pairwise comparison probability. The imported §7.3 and §7.4 layers close the two former follow-up items: sum_compared_prob_eq_expectedComparisons identifies the pairwise sum with the expected comparison count, and expectedRunningTime_isBigTheta_nlogn proves the textbook Θ(n log n) result. The ExplicitRandomness companion now closes the former probability-space gap: a uniform priority permutation feeds the executable first-pivot quicksort, its natural pair-trace counter is a sum of these same events, and its finite expectation is proved equal to the closed form. Pointwise refinement to the separate operational counter remains explicitly tracked by the Chapter 7 guide.

end Chapter07end CLRS

CLRSLean.FourthEdition.Chapter_07.Section_07_3_Randomized_Quicksort.ExplicitRandomness.Bridge

Expectation bridge for executable random-priority quicksort

The natural-number comparison counter from Execution is averaged over the explicit finite sample space Equiv.Perm (Fin n). Its pointwise indicator decomposition, finite linearity of expectation, and the proved permutation symmetry for each pair yield the existing CLRS pairwise-probability sum. The already-proved algebraic sum theorem then supplies the closed form and the Theta(n log n) result.

              namespace CLRSnamespace Chapter07open CLRS.Probabilityopen Chapter03

Expected comparison count of the executable random-priority model.

noncomputable def explicitRandomizedQuicksortExpectedComparisons (n : Nat) : Real :=
  fintypeExpect (fun priority : Equiv.Perm (Fin n) =>
    (randomizedQuicksortComparisonCount priority : Real))

The expectation of the executable comparison counter is the CLRS sum of pairwise comparison probabilities. No independence assumption is used: finite linearity of expectation suffices.


theorem explicitRandomizedQuicksortExpectedComparisons_eq_pairSum (n : Nat) :
    explicitRandomizedQuicksortExpectedComparisons n =
      ∑ i ∈ Finset.range n, ∑ j ∈ Finset.range n,
        if i < j then (2 : Real) / ((j - i + 1 : Nat) : Real) else 0 := by
  unfold explicitRandomizedQuicksortExpectedComparisons
  rw [show (fun priority : Equiv.Perm (Fin n) =>
      (randomizedQuicksortComparisonCount priority : Real)) =
      (fun priority : Equiv.Perm (Fin n) =>
        ∑ i : Fin n, ∑ j : Fin n,
          if hij : i.val < j.val then
            @indicator (comparedInQuicksort n i.val j.val hij j.isLt priority)
              (Classical.propDecidable _)
          else 0) by
    funext priority
    exact randomizedQuicksortComparisonCount_cast priority]
  calc
    fintypeExpect (fun priority : Equiv.Perm (Fin n) =>
        ∑ i : Fin n, ∑ j : Fin n,
          if hij : i.val < j.val then
            @indicator (comparedInQuicksort n i.val j.val hij j.isLt priority)
              (Classical.propDecidable _)
          else 0) =
        ∑ i : Fin n, fintypeExpect (fun priority : Equiv.Perm (Fin n) =>
          ∑ j : Fin n,
            if hij : i.val < j.val then
              @indicator (comparedInQuicksort n i.val j.val hij j.isLt priority)
                (Classical.propDecidable _)
            else 0) := by
      simpa using
        (fintypeExpect_sum (Ω := Equiv.Perm (Fin n))
          (Finset.univ : Finset (Fin n))
          (fun i priority =>
            ∑ j : Fin n,
              if hij : i.val < j.val then
                @indicator (comparedInQuicksort n i.val j.val hij j.isLt priority)
                  (Classical.propDecidable _)
              else 0))
    _ = ∑ i : Fin n, ∑ j : Fin n,
          fintypeExpect (fun priority : Equiv.Perm (Fin n) =>
            if hij : i.val < j.val then
              @indicator (comparedInQuicksort n i.val j.val hij j.isLt priority)
                (Classical.propDecidable _)
            else 0) := by
      apply Finset.sum_congr rfl
      intro i _
      simpa using
        (fintypeExpect_sum (Ω := Equiv.Perm (Fin n))
          (Finset.univ : Finset (Fin n))
          (fun j priority =>
            if hij : i.val < j.val then
              @indicator (comparedInQuicksort n i.val j.val hij j.isLt priority)
                (Classical.propDecidable _)
            else 0))
    _ = ∑ i : Fin n, ∑ j : Fin n,
          if i.val < j.val then
            (2 : Real) / ((j.val - i.val + 1 : Nat) : Real)
          else 0 := by
      apply Finset.sum_congr rfl
      intro i _
      apply Finset.sum_congr rfl
      intro j _
      by_cases hij : i.val < j.val
      · simp only [dif_pos hij, if_pos hij]
        exact comparedIndicator_expectation n i j hij
      · simp [hij, fintypeExpect]
    _ = ∑ i ∈ Finset.range n, ∑ j ∈ Finset.range n,
          if i < j then (2 : Real) / ((j - i + 1 : Nat) : Real) else 0 := by
      have hinner (i : Fin n) :
          (∑ j : Fin n,
            if i.val < j.val then
              (2 : Real) / ((j.val - i.val + 1 : Nat) : Real)
            else 0) =
          ∑ j ∈ Finset.range n,
            if i.val < j then
              (2 : Real) / ((j - i.val + 1 : Nat) : Real)
            else 0 := by
        rw [Fin.sum_univ_eq_sum_range
          (fun j : Nat =>
            if i.val < j then
              (2 : Real) / ((j - i.val + 1 : Nat) : Real)
            else 0) n]
      simp_rw [hinner]
      rw [Fin.sum_univ_eq_sum_range
        (fun i : Nat =>
          ∑ j ∈ Finset.range n,
            if i < j then
              (2 : Real) / ((j - i + 1 : Nat) : Real)
            else 0) n]

The explicit finite execution has exactly the previously proved CLRS expected-comparison closed form.


theorem explicitRandomizedQuicksortExpectedComparisons_eq (n : Nat) :
    explicitRandomizedQuicksortExpectedComparisons n = expectedComparisonsReal n := by
  rw [explicitRandomizedQuicksortExpectedComparisons_eq_pairSum]
  have h := sum_compared_prob_eq_expectedComparisons n
  dsimp [expectedComparisonsReal]
  rw [← h]
  push_cast
  apply Finset.sum_congr rfl
  intro i _
  apply Finset.sum_congr rfl
  intro j _
  by_cases hij : i < j
  · simp [hij]
  · simp [hij]

The executable finite model inherits the existing textbook Theta(n log n) theorem through the exact expectation bridge.


theorem explicitRandomizedQuicksortExpectedComparisons_isBigTheta_nlogn :
    isBigTheta explicitRandomizedQuicksortExpectedComparisons
      (fun n : Nat => (n : Real) * Real.log (n : Real)) := by
  have hfun : explicitRandomizedQuicksortExpectedComparisons = expectedComparisonsReal := by
    funext n
    exact explicitRandomizedQuicksortExpectedComparisons_eq n
  rw [hfun]
  exact expectedComparisons_isBigTheta_nlogn

end Chapter07end CLRS

CLRSLean.FourthEdition.Chapter_07.Section_07_3_Randomized_Quicksort.ExplicitRandomness.Execution

Executable random-priority quicksort model

A sample priority : Equiv.Perm (Fin n) is read as the input order of the distinct ranks 0, ..., n-1. Running the existing first-pivot functional quicksort on that list is therefore the standard random-permutation implementation of randomized quicksort.

The comparison trace uses CLRS's exact pair characterization: ranks i < j are compared precisely when i or j has minimum priority among the ranks in [i,j]. This predicate is executable after unfolding its finite quantifier, so the total trace cardinality is a natural-number random variable on the same permutation sample.

              namespace CLRSnamespace Chapter07open CLRS.Probability

Turn a priority permutation into the concrete list of distinct ranks fed to the existing first-pivot quicksort.

def randomizedQuicksortInput {n : Nat} (priority : Equiv.Perm (Fin n)) : List Nat :=
  (List.finRange n).map (fun position => (priority position).val)

@[simp] theorem randomizedQuicksortInput_length {n : Nat}
    (priority : Equiv.Perm (Fin n)) :
    (randomizedQuicksortInput priority).length = n := by
  simp [randomizedQuicksortInput]

The concrete randomized input contains every rank exactly once.

theorem randomizedQuicksortInput_perm_range {n : Nat}
    (priority : Equiv.Perm (Fin n)) :
    (randomizedQuicksortInput priority).Perm (List.range n) := by
  have hperm := (Equiv.Perm.map_finRange_perm priority).map (fun rank => rank.val)
  simpa [randomizedQuicksortInput, List.map_map, Function.comp_def] using hperm

The output of the existing executable first-pivot quicksort on the sampled priority order.

def randomizedQuicksortOutput {n : Nat} (priority : Equiv.Perm (Fin n)) : List Nat :=
  quickSort (randomizedQuicksortInput priority)

Definitional refinement to the executable Chapter 7 quicksort. This theorem deliberately concerns the returned list; identifying the abstract CLRS pair trace below with quickSortComparisons is a separate operational counter-refinement obligation.

theorem randomizedQuicksortOutput_eq_quickSort {n : Nat}
    (priority : Equiv.Perm (Fin n)) :
    randomizedQuicksortOutput priority =
      quickSort (randomizedQuicksortInput priority) := rfl

Exact relation to the executable Chapter 7 quicksort: the sampled run is ordered and preserves precisely the ranks 0, ..., n-1.

theorem randomizedQuicksortOutput_correct {n : Nat}
    (priority : Equiv.Perm (Fin n)) :
    Ordered (randomizedQuicksortOutput priority) ∧
      (randomizedQuicksortOutput priority).Perm (List.range n) := by
  refine ⟨quickSort_ordered _, ?_⟩
  exact (quickSort_perm _).trans (randomizedQuicksortInput_perm_range priority)

A constructive decision procedure for the CLRS comparison event.

def comparedInQuicksortDecidable (n : Nat) (i j : Fin n)
    (hij : i.val < j.val) (priority : Equiv.Perm (Fin n)) :
    Decidable (comparedInQuicksort n i.val j.val hij j.isLt priority) := by
  unfold comparedInQuicksort IsFirstIn pos
  infer_instance

The executable 0/1 contribution of one ordered rank pair.

def randomizedQuicksortComparisonBit {n : Nat} (i j : Fin n)
    (hij : i.val < j.val) (priority : Equiv.Perm (Fin n)) : Nat :=
  @ite Nat (comparedInQuicksort n i.val j.val hij j.isLt priority)
    (comparedInQuicksortDecidable n i j hij priority) 1 0

Total cardinality of the CLRS rank-pair comparison trace induced by the same priority sample that supplies randomizedQuicksortInput.

def randomizedQuicksortComparisonCount {n : Nat}
    (priority : Equiv.Perm (Fin n)) : Nat :=
  ∑ i : Fin n, ∑ j : Fin n,
    if hij : i.val < j.val then
      randomizedQuicksortComparisonBit i j hij priority
    else 0

Casting an executable comparison bit to Real gives the corresponding indicator random variable.

theorem randomizedQuicksortComparisonBit_cast {n : Nat} (i j : Fin n)
    (hij : i.val < j.val) (priority : Equiv.Perm (Fin n)) :
    (randomizedQuicksortComparisonBit i j hij priority : Real) =
      @indicator (comparedInQuicksort n i.val j.val hij j.isLt priority)
        (Classical.propDecidable _) := by
  unfold randomizedQuicksortComparisonBit indicator
  split <;> simp_all

The natural comparison counter is pointwise equal, after casting, to the finite sum of its pairwise indicators.

theorem randomizedQuicksortComparisonCount_cast {n : Nat}
    (priority : Equiv.Perm (Fin n)) :
    (randomizedQuicksortComparisonCount priority : Real) =
      ∑ i : Fin n, ∑ j : Fin n,
        if hij : i.val < j.val then
          @indicator (comparedInQuicksort n i.val j.val hij j.isLt priority)
            (Classical.propDecidable _)
        else 0 := by
  unfold randomizedQuicksortComparisonCount
  push_cast
  apply Finset.sum_congr rfl
  intro i _
  apply Finset.sum_congr rfl
  intro j _
  split
  · apply randomizedQuicksortComparisonBit_cast
  · norm_num

end Chapter07end CLRS

CLRSLean.FourthEdition.Chapter_07.Section_07_3_Randomized_Quicksort.ExplicitRandomness.OperationalBridge.ExecutionToBST

From executable quicksort to binary-search-tree depth

This module proves the operational half of the randomized-quicksort bridge. For a duplicate-free input, the recursive comparison counter is exactly the sum of the search depths in the binary search tree obtained by inserting the same keys in input order.

namespace CLRSnamespace Chapter07

Sum of the search depths of the keys in xs within t.

def bstDepthSumOn (xs : List Nat) (t : Chapter12.BSTree) : Nat :=
  (xs.map (fun x => Chapter12.BSTree.depth x t)).sum

On a duplicate-free pivot call, quicksort's ≤ pivot left partition is the BST construction's strict < pivot partition.


theorem partitionAround_left_eq_filter_lt_of_nodup
    (pivot : Nat) (tail : List Nat) (h : (pivot :: tail).Nodup) :
    (partitionAround pivot tail).1 =
      tail.filter (fun x => decide (x < pivot)) := by
  rw [partitionAround_left_eq_filter]
  apply List.filter_congr
  intro x hx
  have hne : x ≠ pivot := by
    intro hxp
    subst x
    exact (List.nodup_cons.mp h).1 hx
  simp only [decide_eq_decide]
  omega



private theorem bstDepthSumOn_node_left
    (xs : List Nat) (left right : Chapter12.BSTree) (pivot : Nat)
    (hall : ∀ x ∈ xs, x < pivot) :
    bstDepthSumOn xs (.node left pivot right) =
      xs.length + bstDepthSumOn xs left := by
  induction xs with
  | nil => simp [bstDepthSumOn]
  | cons x xs ih =>
      have hx : x < pivot := hall x (by simp)
      have hne : x ≠ pivot := Nat.ne_of_lt hx
      have hrest : ∀ y ∈ xs, y < pivot := by
        intro y hy
        exact hall y (by simp [hy])
      change Chapter12.BSTree.depth x (.node left pivot right) +
          bstDepthSumOn xs (.node left pivot right) =
        (x :: xs).length + bstDepthSumOn (x :: xs) left
      simp only [Chapter12.BSTree.depth, hne, hx, if_false, if_true]
      rw [ih hrest]
      simp only [List.length_cons, bstDepthSumOn, List.map_cons, List.sum_cons]
      omega



private theorem bstDepthSumOn_node_right
    (xs : List Nat) (left right : Chapter12.BSTree) (pivot : Nat)
    (hall : ∀ x ∈ xs, pivot < x) :
    bstDepthSumOn xs (.node left pivot right) =
      xs.length + bstDepthSumOn xs right := by
  induction xs with
  | nil => simp [bstDepthSumOn]
  | cons x xs ih =>
      have hx : pivot < x := hall x (by simp)
      have hne : x ≠ pivot := Nat.ne_of_gt hx
      have hnlt : ¬x < pivot := Nat.not_lt_of_ge (Nat.le_of_lt hx)
      have hrest : ∀ y ∈ xs, pivot < y := by
        intro y hy
        exact hall y (by simp [hy])
      change Chapter12.BSTree.depth x (.node left pivot right) +
          bstDepthSumOn xs (.node left pivot right) =
        (x :: xs).length + bstDepthSumOn (x :: xs) right
      simp only [Chapter12.BSTree.depth, hne, hnlt, if_false]
      rw [ih hrest]
      simp only [List.length_cons, bstDepthSumOn, List.map_cons, List.sum_cons]
      omega

The recursive comparison counter and the BST depth sum satisfy the same recurrence on every duplicate-free input whose length is covered by fuel.


theorem quickSortComparisonsFuel_eq_bstDepthSumOn
    (fuel : Nat) (xs : List Nat) (hlen : xs.length ≤ fuel)
    (hnodup : xs.Nodup) :
    quickSortComparisonsFuel fuel xs =
      bstDepthSumOn xs (Chapter12.BSTree.buildFromList xs) := by
  induction fuel generalizing xs with
  | zero =>
      have hnil : xs = [] :=
        List.eq_nil_of_length_eq_zero (Nat.eq_zero_of_le_zero hlen)
      subst xs
      simp [quickSortComparisonsFuel, bstDepthSumOn, Chapter12.BSTree.buildFromList]
  | succ fuel ih =>
      cases xs with
      | nil =>
          simp [quickSortComparisonsFuel, bstDepthSumOn, Chapter12.BSTree.buildFromList]
      | cons pivot tail =>
          let left := tail.filter (fun x => decide (x < pivot))
          let right := tail.filter (fun x => decide (pivot < x))
          have htailLength : tail.length ≤ fuel := by
            exact Nat.succ_le_succ_iff.mp (by simpa using hlen)
          have htailNodup : tail.Nodup := hnodup.tail
          have hleft : (partitionAround pivot tail).1 = left := by
            simpa [left] using
              partitionAround_left_eq_filter_lt_of_nodup pivot tail hnodup
          have hright : (partitionAround pivot tail).2 = right := by
            simpa [right] using partitionAround_right_eq_filter pivot tail
          have hleftLength : left.length ≤ fuel :=
            Nat.le_trans (List.length_filter_le _ _) htailLength
          have hrightLength : right.length ≤ fuel :=
            Nat.le_trans (List.length_filter_le _ _) htailLength
          have hleftNodup : left.Nodup := htailNodup.filter _
          have hrightNodup : right.Nodup := htailNodup.filter _
          have ihLeft := ih left hleftLength hleftNodup
          have ihRight := ih right hrightLength hrightNodup
          have hperm : (left ++ right).Perm tail := by
            simpa [hleft, hright] using partitionAround_perm pivot tail
          let leftTree := Chapter12.BSTree.buildFromList left
          let rightTree := Chapter12.BSTree.buildFromList right
          have hsplit :
              bstDepthSumOn tail (.node leftTree pivot rightTree) =
                bstDepthSumOn left (.node leftTree pivot rightTree) +
                  bstDepthSumOn right (.node leftTree pivot rightTree) := by
            have hsum :=
              (hperm.map (fun x => Chapter12.BSTree.depth x
                (.node leftTree pivot rightTree))).sum_eq
            simpa [bstDepthSumOn, List.map_append] using hsum.symm
          have hallLeft : ∀ x ∈ left, x < pivot := by
            intro x hx
            exact of_decide_eq_true (List.mem_filter.mp hx).2
          have hallRight : ∀ x ∈ right, pivot < x := by
            intro x hx
            exact of_decide_eq_true (List.mem_filter.mp hx).2
          have hleftDepth :=
            bstDepthSumOn_node_left left leftTree rightTree pivot hallLeft
          have hrightDepth :=
            bstDepthSumOn_node_right right leftTree rightTree pivot hallRight
          have hpartsLength : left.length + right.length = tail.length := by
            simpa [hleft, hright] using partitionAround_length_add pivot tail
          have hdepthRec :
              bstDepthSumOn (pivot :: tail)
                  (Chapter12.BSTree.buildFromList (pivot :: tail)) =
                tail.length + bstDepthSumOn left leftTree +
                  bstDepthSumOn right rightTree := by
            rw [Chapter12.BSTree.buildFromList_cons]
            change bstDepthSumOn (pivot :: tail)
                (.node leftTree pivot rightTree) = _
            rw [show bstDepthSumOn (pivot :: tail)
                (.node leftTree pivot rightTree) =
                  bstDepthSumOn tail (.node leftTree pivot rightTree) by
                simp [bstDepthSumOn, Chapter12.BSTree.depth]]
            rw [hsplit, hleftDepth, hrightDepth]
            omega
          calc
            quickSortComparisonsFuel (fuel + 1) (pivot :: tail) =
                tail.length + quickSortComparisonsFuel fuel left +
                  quickSortComparisonsFuel fuel right := by
              simp [quickSortComparisonsFuel, hleft, hright]
            _ = tail.length + bstDepthSumOn left leftTree +
                  bstDepthSumOn right rightTree := by
              simp only [leftTree, rightTree]
              rw [ihLeft, ihRight]
            _ = bstDepthSumOn (pivot :: tail)
                  (Chapter12.BSTree.buildFromList (pivot :: tail)) := hdepthRec.symm

Chapter 7 and Chapter 12 use definitionally identical encodings of a permutation as an insertion/input list.

theorem randomizedQuicksortInput_eq_permKeys {n : Nat}
    (priority : Equiv.Perm (Fin n)) :
    randomizedQuicksortInput priority = Chapter12.BSTree.permKeys priority := rfl

On a sampled permutation, executable quicksort comparisons equal total BST depth over all ranks.


theorem quickSortComparisons_randomizedInput_eq_totalDepth {n : Nat}
    (priority : Equiv.Perm (Fin n)) :
    quickSortComparisons (randomizedQuicksortInput priority) =
      ∑ j : Fin n,
        Chapter12.BSTree.depth j.val (Chapter12.BSTree.buildFromPerm priority) := by
  have hnodup : (randomizedQuicksortInput priority).Nodup := by
    exact (randomizedQuicksortInput_perm_range priority).nodup_iff.mpr
      List.nodup_range
  have hrun := quickSortComparisonsFuel_eq_bstDepthSumOn
    (randomizedQuicksortInput priority).length
    (randomizedQuicksortInput priority) (by omega) hnodup
  rw [show quickSortComparisons (randomizedQuicksortInput priority) =
      quickSortComparisonsFuel (randomizedQuicksortInput priority).length
        (randomizedQuicksortInput priority) from rfl]
  rw [hrun]
  have hperm := randomizedQuicksortInput_perm_range priority
  have hsum := (hperm.map (fun x =>
    Chapter12.BSTree.depth x (Chapter12.BSTree.buildFromPerm priority))).sum_eq
  rw [show Chapter12.BSTree.buildFromList (randomizedQuicksortInput priority) =
      Chapter12.BSTree.buildFromPerm priority by
        rw [randomizedQuicksortInput_eq_permKeys]
        rfl]
  unfold bstDepthSumOn
  rw [hsum]
  rw [Fin.sum_univ_eq_sum_range
    (fun j : Nat =>
      Chapter12.BSTree.depth j (Chapter12.BSTree.buildFromPerm priority)) n]
  have hrange :
      (Finset.range n).sum (fun j : Nat =>
          Chapter12.BSTree.depth j (Chapter12.BSTree.buildFromPerm priority)) =
        ((List.range n).map (fun j : Nat =>
          Chapter12.BSTree.depth j
            (Chapter12.BSTree.buildFromPerm priority))).sum := by
    rw [← List.toFinset_range n]
    exact List.sum_toFinset
      (fun j : Nat =>
        Chapter12.BSTree.depth j (Chapter12.BSTree.buildFromPerm priority))
      (l := List.range n)
      List.nodup_range
  exact hrange.symm

end Chapter07end CLRS

CLRSLean.FourthEdition.Chapter_07.Section_07_3_Randomized_Quicksort.ExplicitRandomness.OperationalBridge.PairTraceToBST

From the CLRS comparison trace to binary-search-tree depth

The pair characterization of randomized quicksort and the ancestor characterization of a randomly built BST both say that one endpoint is the first key in the same closed rank interval. This module makes that shared semantics explicit and converts the pair count into total BST depth.

              namespace CLRSnamespace Chapter07open CLRS.Probabilityprivate noncomputable instance ancestorDecidable
    (x y : Nat) (tree : Chapter12.BSTree) :
    Decidable (Chapter12.BSTree.isAncestorOf x y tree) :=
  Classical.propDecidable _

The embedded natural-number interval used by Chapter 7 is the ordinary closed interval of bounded ranks used by Chapter 12.

theorem rangeFin_eq_intervalIcc {n : Nat} (i j : Fin n) (hij : i.val < j.val) :
    rangeFin n i.val j.val (Nat.le_of_lt hij) j.isLt = Finset.Icc i j := by
  ext k
  simp [rangeFin, Finset.mem_Icc, Fin.ext_iff]



private theorem isFirstIn_Icc_left_iff_firstInInterval {n : Nat}
    (i j : Fin n) (hij : i.val < j.val) (priority : Equiv.Perm (Fin n)) :
    IsFirstIn (Finset.Icc i j) i priority ↔
      Chapter12.BSTree.firstInInterval priority i j := by
  have hle : i ≤ j := by exact Fin.mk_le_mk.mpr (Nat.le_of_lt hij)
  simp [IsFirstIn, pos, Chapter12.BSTree.firstInInterval,
    Finset.mem_Icc, hle]



private theorem isFirstIn_Icc_right_iff_firstInInterval {n : Nat}
    (i j : Fin n) (hij : i.val < j.val) (priority : Equiv.Perm (Fin n)) :
    IsFirstIn (Finset.Icc i j) j priority ↔
      Chapter12.BSTree.firstInInterval priority j i := by
  have hle : i ≤ j := by exact Fin.mk_le_mk.mpr (Nat.le_of_lt hij)
  simp [IsFirstIn, pos, Chapter12.BSTree.firstInInterval,
    Finset.mem_Icc, hle]

A rank pair is compared by quicksort exactly when one endpoint is an ancestor of the other in the BST built from the same priority permutation.


theorem comparedInQuicksort_iff_ancestor {n : Nat} (i j : Fin n)
    (hij : i.val < j.val) (priority : Equiv.Perm (Fin n)) :
    comparedInQuicksort n i.val j.val hij j.isLt priority ↔
      Chapter12.BSTree.isAncestorOf i.val j.val
          (Chapter12.BSTree.buildFromPerm priority) ∨
        Chapter12.BSTree.isAncestorOf j.val i.val
          (Chapter12.BSTree.buildFromPerm priority) := by
  unfold comparedInQuicksort
  dsimp only
  rw [rangeFin_eq_intervalIcc i j hij]
  have hi : (⟨i.val, Nat.lt_trans hij j.isLt⟩ : Fin n) = i := Fin.ext rfl
  have hj : (⟨j.val, j.isLt⟩ : Fin n) = j := Fin.ext rfl
  rw [hi, hj]
  rw [isFirstIn_Icc_left_iff_firstInInterval i j hij priority]
  rw [isFirstIn_Icc_right_iff_firstInInterval i j hij priority]
  rw [← Chapter12.BSTree.isAncestorOf_buildFromPerm_iff_firstInInterval]
  rw [← Chapter12.BSTree.isAncestorOf_buildFromPerm_iff_firstInInterval]

              private noncomputable def ancestorBit {n : Nat} (priority : Equiv.Perm (Fin n))
    (i j : Fin n) : Nat :=
  if Chapter12.BSTree.isAncestorOf i.val j.val
      (Chapter12.BSTree.buildFromPerm priority) then 1 else 0private theorem not_ancestor_both {n : Nat} (priority : Equiv.Perm (Fin n))
    {i j : Fin n} (hne : i ≠ j) :
    ¬(Chapter12.BSTree.isAncestorOf i.val j.val
          (Chapter12.BSTree.buildFromPerm priority) ∧
        Chapter12.BSTree.isAncestorOf j.val i.val
          (Chapter12.BSTree.buildFromPerm priority)) := by
  intro hboth
  have hfirstI :=
    (Chapter12.BSTree.isAncestorOf_buildFromPerm_iff_firstInInterval
      priority i j).mp hboth.1
  have hfirstJ :=
    (Chapter12.BSTree.isAncestorOf_buildFromPerm_iff_firstInInterval
      priority j i).mp hboth.2
  have hposIJ : (priority.symm i : Nat) ≤ (priority.symm j : Nat) :=
    hfirstI j (min_le_right _ _) (le_max_right _ _)
  have hposJI : (priority.symm j : Nat) ≤ (priority.symm i : Nat) :=
    hfirstJ i (min_le_right _ _) (le_max_right _ _)
  have hpos : priority.symm i = priority.symm j :=
    Fin.ext (Nat.le_antisymm hposIJ hposJI)
  exact hne (priority.symm.injective hpos)

private theorem randomizedQuicksortComparisonBit_eq_ancestorBits {n : Nat}
    (i j : Fin n) (hij : i.val < j.val)
    (priority : Equiv.Perm (Fin n)) :
    randomizedQuicksortComparisonBit i j hij priority =
      ancestorBit priority i j + ancestorBit priority j i := by
  have hne : i ≠ j := by
    intro h
    exact (Nat.ne_of_lt hij) (congrArg Fin.val h)
  have hiff := comparedInQuicksort_iff_ancestor i j hij priority
  unfold randomizedQuicksortComparisonBit ancestorBit
  simp only [hiff]
  by_cases hleft : Chapter12.BSTree.isAncestorOf i.val j.val
      (Chapter12.BSTree.buildFromPerm priority)
  · have hright : ¬Chapter12.BSTree.isAncestorOf j.val i.val
        (Chapter12.BSTree.buildFromPerm priority) := by
      intro hright
      exact not_ancestor_both priority hne ⟨hleft, hright⟩
    simp [hleft, hright]
  · by_cases hright : Chapter12.BSTree.isAncestorOf j.val i.val
        (Chapter12.BSTree.buildFromPerm priority)
    · simp [hleft, hright]
    · simp [hleft, hright]

The depth of one key is the number of its strict ancestors.


theorem depth_eq_sum_strictAncestorBits {n : Nat}
    (priority : Equiv.Perm (Fin n)) (j : Fin n) :
    Chapter12.BSTree.depth j.val
        (Chapter12.BSTree.buildFromPerm priority) =
      ∑ i : Fin n,
        if i ≠ j ∧ Chapter12.BSTree.isAncestorOf i.val j.val
            (Chapter12.BSTree.buildFromPerm priority) then 1 else 0 := by
  classical
  let tree := Chapter12.BSTree.buildFromPerm priority
  let allAncestors : Finset (Fin n) :=
    Finset.univ.filter (fun i : Fin n =>
      Chapter12.BSTree.isAncestorOf i.val j.val tree)
  let strictAncestors : Finset (Fin n) :=
    Finset.univ.filter (fun i : Fin n =>
      i ≠ j ∧ Chapter12.BSTree.isAncestorOf i.val j.val tree)
  have hself : Chapter12.BSTree.isAncestorOf j.val j.val tree := by
    simpa [tree] using
      Chapter12.BSTree.isAncestorOf_self_buildFromPerm priority j
  have hselfMem : j ∈ allAncestors := by
    simp [allAncestors, hself]
  have herase : allAncestors.erase j = strictAncestors := by
    ext i
    simp [allAncestors, strictAncestors]
  have hcard : strictAncestors.card + 1 = allAncestors.card := by
    simpa [herase] using Finset.card_erase_add_one hselfMem
  have hordered : Chapter12.BSTree.Ordered tree := by
    simpa [tree, Chapter12.BSTree.buildFromPerm] using
      Chapter12.BSTree.buildFromList_ordered
      (Chapter12.BSTree.permKeys priority)
  have hbounded : ∀ z, Chapter12.BSTree.InTree z tree → z < n := by
    intro z hz
    exact Chapter12.BSTree.InTree_buildFromPerm_lt priority (by simpa [tree] using hz)
  have hcount : Chapter12.BSTree.ancestorCount j.val tree = allAncestors.card := by
    change Chapter12.BSTree.ancestorCount j.val tree =
      (Finset.univ.filter (fun i : Fin n =>
        Chapter12.BSTree.isAncestorOf i.val j.val tree)).card
    exact Chapter12.BSTree.ancestorCount_eq_sum j.val tree hordered hbounded
  have hdepth : Chapter12.BSTree.ancestorCount j.val tree =
      Chapter12.BSTree.depth j.val tree + 1 :=
    Chapter12.BSTree.ancestorCount_eq_depth_add_one j.val tree hself
  have hsum :
      (∑ i : Fin n,
        if i ≠ j ∧ Chapter12.BSTree.isAncestorOf i.val j.val tree
          then 1 else 0) = strictAncestors.card := by
    simp [strictAncestors]
  change Chapter12.BSTree.depth j.val tree = _
  rw [hsum]
  omega

The executable CLRS pair-trace counter is total BST depth for the tree built from the same permutation.


theorem randomizedQuicksortComparisonCount_eq_totalDepth {n : Nat}
    (priority : Equiv.Perm (Fin n)) :
    randomizedQuicksortComparisonCount priority =
      ∑ j : Fin n,
        Chapter12.BSTree.depth j.val
          (Chapter12.BSTree.buildFromPerm priority) := by
  classical
  calc
    randomizedQuicksortComparisonCount priority =
        ∑ i : Fin n, ∑ j : Fin n,
          if hij : i.val < j.val then
            ancestorBit priority i j + ancestorBit priority j i else 0 := by
      unfold randomizedQuicksortComparisonCount
      apply Finset.sum_congr rfl
      intro i _
      apply Finset.sum_congr rfl
      intro j _
      split
      · exact randomizedQuicksortComparisonBit_eq_ancestorBits i j _ priority
      · rfl
    _ = ∑ i : Fin n, ∑ j : Fin n,
          if i ≠ j then ancestorBit priority i j else 0 := by
      exact (sum_offDiagonal_eq_sum_strictPairs (ancestorBit priority)).symm
    _ = ∑ j : Fin n, ∑ i : Fin n,
          if i ≠ j then ancestorBit priority i j else 0 := by
      rw [Finset.sum_comm]
    _ = ∑ j : Fin n,
          Chapter12.BSTree.depth j.val
            (Chapter12.BSTree.buildFromPerm priority) := by
      apply Finset.sum_congr rfl
      intro j _
      rw [depth_eq_sum_strictAncestorBits priority j]
      apply Finset.sum_congr rfl
      intro i _
      simp only [ancestorBit]
      by_cases hne : i ≠ j <;>
        by_cases hancestor : Chapter12.BSTree.isAncestorOf i.val j.val
          (Chapter12.BSTree.buildFromPerm priority) <;>
        simp [hne, hancestor]

end Chapter07end CLRS

CLRSLean.FourthEdition.Chapter_07.Section_07_3_Randomized_Quicksort.ExplicitRandomness.OperationalBridge.UnorderedPairs

Finite unordered-pair sums

This module isolates the finite reindexing used by the operational quicksort bridge. An off-diagonal sum over ordered pairs is rewritten as one sum over strict pairs, with the two orientations added at each strict pair.

Main result:

Theorem CLRS.Chapter07.sum_offDiagonal_eq_sum_strictPairs: off-diagonal and strict-pair sums agree.

namespace CLRSnamespace Chapter07

Summing both orientations of every strict pair is the same as summing an off-diagonal matrix.


theorem sum_offDiagonal_eq_sum_strictPairs {n : Nat}
    (f : Fin n → Fin n → Nat) :
    (∑ i : Fin n, ∑ j : Fin n, if i ≠ j then f i j else 0) =
      ∑ i : Fin n, ∑ j : Fin n,
        if i.val < j.val then f i j + f j i else 0 := by
  induction n with
  | zero => simp
  | succ n ih =>
      simp only [Fin.sum_univ_succ]
      have hsub :
          (∑ i : Fin n, ∑ j : Fin n,
              if i.succ ≠ j.succ then f i.succ j.succ else 0) =
            ∑ i : Fin n, ∑ j : Fin n,
              if i.val < j.val then
                f i.succ j.succ + f j.succ i.succ else 0 := by
        simpa using ih (fun i j => f i.succ j.succ)
      simp_rw [Finset.sum_add_distrib]
      rw [hsub]
      have hzero (i : Fin n) : (0 : Fin (n + 1)) ≠ i.succ :=
        (Fin.succ_ne_zero i).symm
      simp [Fin.succ_ne_zero, Finset.sum_add_distrib, add_assoc,
        add_left_comm, add_comm, hzero]

end Chapter07end CLRS

CLRSLean.FourthEdition.Chapter_07.Section_07_3_Randomized_Quicksort.ExplicitRandomness.OperationalBridge

Pointwise operational bridge for randomized quicksort

The abstract CLRS rank-pair trace and the recursive executable quicksort counter have independently been normalized to the total depth of the same randomly built BST. Their composition gives the pointwise operational refinement and transfers the exact expectation and asymptotic theorem to the actual recursive comparison counter.

              namespace CLRSnamespace Chapter07open CLRS.Probabilityopen Chapter03

The CLRS pair-trace counter is pointwise identical to the comparisons made by the recursive executable quicksort on the sampled permutation input.


theorem randomizedQuicksortComparisonCount_eq_quickSortComparisons {n : Nat}
    (priority : Equiv.Perm (Fin n)) :
    randomizedQuicksortComparisonCount priority =
      quickSortComparisons (randomizedQuicksortInput priority) := by
  rw [randomizedQuicksortComparisonCount_eq_totalDepth priority]
  exact (quickSortComparisons_randomizedInput_eq_totalDepth priority).symm

Expected comparison count of the actual recursive quicksort counter over the uniform finite permutation sample space.

noncomputable def operationalRandomizedQuicksortExpectedComparisons
    (n : Nat) : Real :=
  fintypeExpect (fun priority : Equiv.Perm (Fin n) =>
    (quickSortComparisons (randomizedQuicksortInput priority) : Real))

The operational and pair-trace expectations are equal because their counters agree for every individual permutation.

theorem operationalRandomizedQuicksortExpectedComparisons_eq_explicit
    (n : Nat) :
    operationalRandomizedQuicksortExpectedComparisons n =
      explicitRandomizedQuicksortExpectedComparisons n := by
  unfold operationalRandomizedQuicksortExpectedComparisons
  unfold explicitRandomizedQuicksortExpectedComparisons
  congr 1
  funext priority
  exact_mod_cast
    (randomizedQuicksortComparisonCount_eq_quickSortComparisons priority).symm

The actual recursive comparison counter has exactly the textbook closed form in expectation.


theorem operationalRandomizedQuicksortExpectedComparisons_eq (n : Nat) :
    operationalRandomizedQuicksortExpectedComparisons n =
      expectedComparisonsReal n := by
  rw [operationalRandomizedQuicksortExpectedComparisons_eq_explicit]
  exact explicitRandomizedQuicksortExpectedComparisons_eq n

The expected number of comparisons performed by the actual recursive quicksort implementation is Theta(n log n).


theorem operationalRandomizedQuicksortExpectedComparisons_isBigTheta_nlogn :
    isBigTheta operationalRandomizedQuicksortExpectedComparisons
      (fun n : Nat => (n : Real) * Real.log (n : Real)) := by
  have hfun : operationalRandomizedQuicksortExpectedComparisons =
      explicitRandomizedQuicksortExpectedComparisons := by
    funext n
    exact operationalRandomizedQuicksortExpectedComparisons_eq_explicit n
  rw [hfun]
  exact explicitRandomizedQuicksortExpectedComparisons_isBigTheta_nlogn

end Chapter07end CLRS

CLRSLean.FourthEdition.Chapter_07.Section_07_3_Randomized_Quicksort.ExplicitRandomness.Probability

Uniform random priorities for quicksort

A permutation of Fin n assigns a distinct random priority to every rank: pos pi x is the priority of rank x. The pivot of any recursive subproblem is the rank of minimum priority in that subproblem. This single finite sample space couples all recursive pivot choices.

The expectation proof below uses linearity, so it does not require pairwise independence of comparison indicators. The symmetry actually needed by the algorithm -- every member of a nonempty subproblem is equally likely to have minimum priority -- is proved by the transposition bijection in Comparison_Probability and exposed here in fintypeExpect form.

              namespace CLRSnamespace Chapter07open CLRS.Probability

The expectation of a uniform finite event indicator is its filtered-card ratio.

theorem uniformEvent_expectation_eq_filterRatio {Omega : Type}
    [Fintype Omega] [DecidableEq Omega] (event : Omega -> Prop)
    [DecidablePred event] :
    fintypeExpect (fun sample : Omega => indicator (event sample)) =
      (((Finset.univ : Finset Omega).filter event).card : Real) /
        (Fintype.card Omega : Real) := by
  unfold fintypeExpect indicator
  congr 1
  simp

In the uniform-priority model, every rank in a nonempty recursive subproblem is equally likely to be selected as its pivot.


theorem priorityPivot_uniform {n : Nat} (subproblem : Finset (Fin n))
    (hne : subproblem.Nonempty) (pivot : Fin n) (hpivot : pivot ∈ subproblem) :
    fintypeExpect (fun priority : Equiv.Perm (Fin n) =>
      @indicator (IsFirstIn subproblem pivot priority)
        (Classical.propDecidable _)) =
      1 / (subproblem.card : Real) := by
  rw [uniformEvent_expectation_eq_filterRatio]
  simpa [Fintype.card_perm, Fintype.card_fin] using
    isFirst_prob subproblem hne pivot hpivot

The indicator expectation for one rank pair is the comparison probability proved by the permutation-symmetry argument.


theorem comparedIndicator_expectation (n : Nat) (i j : Fin n)
    (hij : i.val < j.val) :
    fintypeExpect (fun priority : Equiv.Perm (Fin n) =>
      @indicator (comparedInQuicksort n i.val j.val hij j.isLt priority)
        (Classical.propDecidable _)) =
      (2 : Real) / ((j.val - i.val + 1 : Nat) : Real) := by
  rw [uniformEvent_expectation_eq_filterRatio]
  simpa [Fintype.card_perm, Fintype.card_fin] using
    compared_prob n i.val j.val hij j.isLt

end Chapter07end CLRS