pith. machine review for the scientific record. sign in

IndisputableMonolith.Chemistry.VanDerWaals

IndisputableMonolith/Chemistry/VanDerWaals.lean · 320 lines · 16 declarations

show as:
view math explainer →

open module explainer GitHub source

Explainer status: pending

   1import Mathlib
   2import IndisputableMonolith.Constants
   3import IndisputableMonolith.Chemistry.PeriodicTable
   4
   5/-!
   6# Van der Waals Forces Derivation (CH-013)
   7
   8Van der Waals forces are weak intermolecular attractions arising from temporary
   9dipole fluctuations and polarizability effects.
  10
  11**RS Mechanism**:
  121. **Ledger Fluctuations**: The 8-tick ledger creates temporary asymmetries in
  13   electron distribution, inducing instantaneous dipoles.
  142. **Polarizability**: Atoms with more electrons (larger, in lower periods) have
  15   higher polarizability, leading to stronger dispersion forces.
  163. **Distance Dependence**: The interaction energy scales as 1/r⁶ (London dispersion)
  17   arising from the dipole-dipole interaction mathematics.
  184. **Noble Gas Trend**: Noble gas boiling points increase down the group,
  19   reflecting increasing polarizability and vdW forces.
  20
  21**Prediction**: Van der Waals interaction strength increases with atomic/molecular
  22size (polarizability), and follows 1/r⁶ distance dependence.
  23-/
  24
  25namespace IndisputableMonolith.Chemistry.VanDerWaals
  26
  27open PeriodicTable
  28
  29noncomputable section
  30
  31/-- Noble gas atomic numbers. -/
  32def nobleGases : List ℕ := [2, 10, 18, 36, 54, 86]  -- He, Ne, Ar, Kr, Xe, Rn
  33
  34/-- Boiling points of noble gases (Kelvin). -/
  35def nobleGasBoilingPoint : ℕ → ℝ
  36| 2  => 4.22    -- He
  37| 10 => 27.07   -- Ne
  38| 18 => 87.30   -- Ar
  39| 36 => 119.93  -- Kr
  40| 54 => 165.05  -- Xe
  41| 86 => 211.4   -- Rn
  42| _  => 0
  43
  44/-- Polarizability proxy (increases with shell number). -/
  45def polarizabilityProxy (Z : ℕ) : ℝ :=
  46  (periodOf Z : ℝ)
  47
  48/-- London dispersion force proxy between two atoms.
  49    Scales with product of polarizabilities and inversely with r⁶. -/
  50def londonDispersionProxy (Z1 Z2 : ℕ) (r : ℝ) : ℝ :=
  51  if r ≤ 0 then 0
  52  else (polarizabilityProxy Z1 * polarizabilityProxy Z2) / r ^ 6
  53
  54/-- The Lennard-Jones potential: U(r) = 4ε[(σ/r)¹² - (σ/r)⁶]
  55    Models vdW attraction (r⁻⁶) and Pauli repulsion (r⁻¹²). -/
  56def lennardJonesPotential (ε σ r : ℝ) : ℝ :=
  57  if r ≤ 0 then 0
  58  else 4 * ε * ((σ / r) ^ 12 - (σ / r) ^ 6)
  59
  60/-- The LJ potential has a minimum at r = 2^(1/6) * σ ≈ 1.122σ. -/
  61def ljMinimumDistance (σ : ℝ) : ℝ := (2 : ℝ) ^ (1/6 : ℝ) * σ
  62
  63/-- The LJ minimum is approximately 1.122σ. -/
  64theorem lj_minimum_approx (σ : ℝ) (hσ : σ > 0) :
  65    ljMinimumDistance σ / σ > 1.1 ∧ ljMinimumDistance σ / σ < 1.15 := by
  66  dsimp [ljMinimumDistance]
  67  rw [mul_div_assoc, div_self (ne_of_gt hσ)]
  68  -- 2^(1/6) ≈ 1.122
  69  -- Verify: 1.1^6 = 1.7716 < 2, 1.15^6 = 2.313 > 2
  70  -- Therefore 1.1 < 2^(1/6) < 1.15
  71  constructor
  72  · -- Show: 2^(1/6) > 1.1
  73    -- We have: 1.1^6 < 2
  74    -- Taking 1/6 power: (1.1^6)^(1/6) < 2^(1/6)
  75    -- Which gives: 1.1 < 2^(1/6)
  76    have h : (1.1 : ℝ)^6 < (2 : ℝ) := by norm_num
  77    have h_nonneg : (0 : ℝ) ≤ (1.1 : ℝ) := by norm_num
  78    -- Use: (1.1^6)^(1/6) = 1.1^(6 * 1/6) = 1.1^1 = 1.1
  79    have h_simplify : ((1.1 : ℝ)^6) ^ (1/6 : ℝ) = (1.1 : ℝ) := by
  80      calc ((1.1 : ℝ)^6) ^ (1/6 : ℝ)
  81        _ = ((1.1 : ℝ)^(6 : ℝ)) ^ (1/6 : ℝ) := by congr 1; exact (Real.rpow_natCast (1.1 : ℝ) 6).symm
  82        _ = (1.1 : ℝ) ^ ((6 : ℝ) * (1/6 : ℝ)) := by rw [← Real.rpow_mul h_nonneg (6 : ℝ) (1/6 : ℝ)]
  83        _ = (1.1 : ℝ) ^ (1 : ℝ) := by norm_num
  84        _ = (1.1 : ℝ) := by rw [Real.rpow_one]
  85    -- Now: (1.1^6)^(1/6) < 2^(1/6), so 1.1 < 2^(1/6)
  86    have h_rpow : ((1.1 : ℝ)^6) ^ (1/6 : ℝ) < (2 : ℝ) ^ (1/6 : ℝ) := by
  87      apply Real.rpow_lt_rpow
  88      · norm_num  -- 0 ≤ 1.1^6
  89      · exact h    -- 1.1^6 < 2
  90      · norm_num  -- 0 < 1/6
  91    rw [h_simplify] at h_rpow
  92    -- h_rpow gives: 1.1 < 2^(1/6)
  93    -- The goal is: 2^(1/6) * 1 > 1.1, which simplifies to 2^(1/6) > 1.1
  94    -- First simplify: 2^(1/6) * 1 = 2^(1/6)
  95    simp only [mul_one]
  96    -- Now the goal is: 2^(1/6) > 1.1, which is equivalent to 1.1 < 2^(1/6)
  97    exact h_rpow
  98  · -- Show: 2^(1/6) < 1.15
  99    -- We have: 2 < 1.15^6
 100    -- Taking 1/6 power: 2^(1/6) < (1.15^6)^(1/6) = 1.15
 101    have h : (2 : ℝ) < (1.15 : ℝ)^6 := by norm_num
 102    have h_nonneg : (0 : ℝ) ≤ (1.15 : ℝ) := by norm_num
 103    -- Use: (1.15^6)^(1/6) = 1.15^(6 * 1/6) = 1.15^1 = 1.15
 104    have h_simplify : ((1.15 : ℝ)^6) ^ (1/6 : ℝ) = (1.15 : ℝ) := by
 105      calc ((1.15 : ℝ)^6) ^ (1/6 : ℝ)
 106        _ = ((1.15 : ℝ)^(6 : ℝ)) ^ (1/6 : ℝ) := by congr 1; exact (Real.rpow_natCast (1.15 : ℝ) 6).symm
 107        _ = (1.15 : ℝ) ^ ((6 : ℝ) * (1/6 : ℝ)) := by rw [← Real.rpow_mul h_nonneg (6 : ℝ) (1/6 : ℝ)]
 108        _ = (1.15 : ℝ) ^ (1 : ℝ) := by norm_num
 109        _ = (1.15 : ℝ) := by rw [Real.rpow_one]
 110    -- Now: 2^(1/6) < (1.15^6)^(1/6) = 1.15
 111    have h_rpow : (2 : ℝ) ^ (1/6 : ℝ) < ((1.15 : ℝ)^6) ^ (1/6 : ℝ) := by
 112      apply Real.rpow_lt_rpow
 113      · norm_num  -- 0 ≤ 2
 114      · exact h    -- 2 < 1.15^6
 115      · norm_num  -- 0 < 1/6
 116    rw [h_simplify] at h_rpow
 117    -- The goal is: 2^(1/6) * 1 < 1.15, simplify to 2^(1/6) < 1.15
 118    simp only [mul_one]
 119    exact h_rpow
 120
 121/-- Noble gas boiling points increase down the group (vdW strength increases). -/
 122theorem noble_gas_bp_increases_he_ne : nobleGasBoilingPoint 2 < nobleGasBoilingPoint 10 := by
 123  simp only [nobleGasBoilingPoint]
 124  norm_num
 125
 126theorem noble_gas_bp_increases_ne_ar : nobleGasBoilingPoint 10 < nobleGasBoilingPoint 18 := by
 127  simp only [nobleGasBoilingPoint]
 128  norm_num
 129
 130theorem noble_gas_bp_increases_ar_kr : nobleGasBoilingPoint 18 < nobleGasBoilingPoint 36 := by
 131  simp only [nobleGasBoilingPoint]
 132  norm_num
 133
 134theorem noble_gas_bp_increases_kr_xe : nobleGasBoilingPoint 36 < nobleGasBoilingPoint 54 := by
 135  simp only [nobleGasBoilingPoint]
 136  norm_num
 137
 138theorem noble_gas_bp_increases_xe_rn : nobleGasBoilingPoint 54 < nobleGasBoilingPoint 86 := by
 139  simp only [nobleGasBoilingPoint]
 140  norm_num
 141
 142/-- Complete ordering of noble gas boiling points. -/
 143theorem noble_gas_bp_full_ordering :
 144    nobleGasBoilingPoint 2 < nobleGasBoilingPoint 10 ∧
 145    nobleGasBoilingPoint 10 < nobleGasBoilingPoint 18 ∧
 146    nobleGasBoilingPoint 18 < nobleGasBoilingPoint 36 ∧
 147    nobleGasBoilingPoint 36 < nobleGasBoilingPoint 54 ∧
 148    nobleGasBoilingPoint 54 < nobleGasBoilingPoint 86 := by
 149  refine ⟨?_, ?_, ?_, ?_, ?_⟩
 150  · exact noble_gas_bp_increases_he_ne
 151  · exact noble_gas_bp_increases_ne_ar
 152  · exact noble_gas_bp_increases_ar_kr
 153  · exact noble_gas_bp_increases_kr_xe
 154  · exact noble_gas_bp_increases_xe_rn
 155
 156/-- London dispersion force decreases with distance (r⁻⁶). -/
 157theorem london_decreases_with_distance (Z1 Z2 : ℕ) (r1 r2 : ℝ)
 158    (hr1 : r1 > 0) (hr2 : r2 > 0) (h_r_ord : r1 < r2) :
 159    londonDispersionProxy Z1 Z2 r1 > londonDispersionProxy Z1 Z2 r2 := by
 160  simp only [londonDispersionProxy]
 161  have hr1_pos' : ¬(r1 ≤ 0) := not_le.mpr hr1
 162  have hr2_pos' : ¬(r2 ≤ 0) := not_le.mpr hr2
 163  simp only [hr1_pos', hr2_pos', ite_false]
 164  -- Need: α/r1⁶ > α/r2⁶ when r1 < r2 and α ≥ 0
 165  -- For r1 < r2, we have r1⁶ < r2⁶, so 1/r1⁶ > 1/r2⁶
 166  -- Since polarizabilityProxy returns non-negative values, α ≥ 0
 167  have h_pow : r1^6 < r2^6 := by
 168    -- For 0 < r1 < r2, we have r1^6 < r2^6
 169    -- Use pow_lt_pow_left₀: if 0 ≤ a < b and 0 < n, then a^n < b^n
 170    apply pow_lt_pow_left₀ h_r_ord (le_of_lt hr1) (by norm_num)
 171  -- Now: 1/r1^6 > 1/r2^6 (since r1^6 < r2^6 and both are positive)
 172  have h_div : 1 / r1^6 > 1 / r2^6 := by
 173    -- one_div_lt_one_div: if 0 < a < b, then 1/b < 1/a
 174    -- We have r1^6 < r2^6, so 1/r2^6 < 1/r1^6
 175    apply (one_div_lt_one_div (pow_pos hr2 6) (pow_pos hr1 6)).mpr h_pow
 176  -- Finally: α/r1^6 > α/r2^6 when α > 0
 177  -- Note: polarizabilityProxy returns periodOf, which is always ≥ 1 for valid atoms
 178  -- So the product is always positive
 179  have h_alpha_pos : 0 < polarizabilityProxy Z1 * polarizabilityProxy Z2 := by
 180    unfold polarizabilityProxy
 181    simp only [periodOf]
 182    norm_cast
 183    -- periodOf returns a natural number, and for any valid atomic number Z, periodOf Z ≥ 1
 184    -- So periodOf Z1 ≥ 1 and periodOf Z2 ≥ 1, hence their product ≥ 1 > 0
 185    -- From the definition: periodOf returns 1, 2, 3, 4, 5, 6, or 7, all ≥ 1
 186    have h1 : (1 : ℕ) ≤ periodOf Z1 := by
 187      -- periodOf is defined with cases that all return values ≥ 1
 188      -- The smallest case is `if Z ≤ 2 then 1`, so periodOf always returns ≥ 1
 189      unfold periodOf
 190      -- All branches return values ≥ 1: 1, 2, 3, 4, 5, 6, or 7
 191      split_ifs <;> norm_num
 192    have h2 : (1 : ℕ) ≤ periodOf Z2 := by
 193      unfold periodOf
 194      split_ifs <;> norm_num
 195    have h_prod : (1 : ℕ) ≤ periodOf Z1 * periodOf Z2 := by
 196      apply Nat.mul_le_mul h1 h2
 197    exact_mod_cast h_prod
 198  -- Since α > 0 and 1/r1^6 > 1/r2^6, we have α/r1^6 > α/r2^6
 199  calc (polarizabilityProxy Z1 * polarizabilityProxy Z2) / r1 ^ 6
 200    _ = (polarizabilityProxy Z1 * polarizabilityProxy Z2) * (1 / r1 ^ 6) := by ring
 201    _ > (polarizabilityProxy Z1 * polarizabilityProxy Z2) * (1 / r2 ^ 6) := by
 202      apply mul_lt_mul_of_pos_left h_div h_alpha_pos
 203    _ = (polarizabilityProxy Z1 * polarizabilityProxy Z2) / r2 ^ 6 := by ring
 204
 205/-- The φ-connection: LJ minimum distance ratio 2^(1/6) ≈ 1.122 is close to φ - 0.5 ≈ 1.118. -/
 206def ljRatioPhiConnection : ℝ := Constants.phi - 0.5
 207
 208theorem lj_phi_connection_approx :
 209    |((2 : ℝ) ^ (1/6 : ℝ)) - ljRatioPhiConnection| < 0.01 := by
 210  dsimp [ljRatioPhiConnection]
 211  -- 2^(1/6) ≈ 1.1225
 212  -- φ - 0.5 ≈ 1.618 - 0.5 = 1.118
 213  -- Difference ≈ 0.0045, which is < 0.01
 214  -- We need: |2^(1/6) - (phi - 0.5)| < 0.01
 215  -- Use bounds: 1.122 < 2^(1/6) < 1.123 and 1.117 < phi - 0.5 < 1.119
 216  -- So difference is at most 1.123 - 1.117 = 0.006 < 0.01
 217  -- More precisely: 2^(1/6) ≈ 1.122462, phi ≈ 1.618034, so phi - 0.5 ≈ 1.118034
 218  -- Difference ≈ 0.004428 < 0.01
 219  have h_phi_lower : (1.117 : ℝ) < Constants.phi - 0.5 := by
 220    -- Use: phi = (1 + √5)/2, so phi - 0.5 = √5/2
 221    have h_phi_minus : Constants.phi - 0.5 = Real.sqrt 5 / 2 := by
 222      rw [Constants.phi]
 223      ring
 224    rw [h_phi_minus]
 225    -- Need: 1.117 < √5/2, i.e., 2.234 < √5
 226    have h_sqrt5 : (2.234 : ℝ) < Real.sqrt 5 := by
 227      have h : (2.234 : ℝ)^2 < (5 : ℝ) := by norm_num
 228      have h_pos : (0 : ℝ) ≤ 2.234 := by norm_num
 229      -- Real.sqrt_lt_sqrt: if 0 ≤ x < y, then √x < √y
 230      -- We have: (2.234)^2 < 5, so √((2.234)^2) < √5
 231      -- And √((2.234)^2) = 2.234 (since 2.234 ≥ 0)
 232      have h_sqrt_sq : Real.sqrt ((2.234 : ℝ)^2) = (2.234 : ℝ) := Real.sqrt_sq h_pos
 233      have h_sqrt_lt : Real.sqrt ((2.234 : ℝ)^2) < Real.sqrt 5 := Real.sqrt_lt_sqrt (by norm_num) h
 234      rw [h_sqrt_sq] at h_sqrt_lt
 235      exact h_sqrt_lt
 236    linarith [h_sqrt5]
 237  have h_phi_upper : Constants.phi - 0.5 < (1.119 : ℝ) := by
 238    -- Use: phi = (1 + √5)/2, so phi - 0.5 = (1 + √5)/2 - 1/2 = √5/2
 239    have h_phi_minus : Constants.phi - 0.5 = Real.sqrt 5 / 2 := by
 240      rw [Constants.phi]
 241      ring
 242    rw [h_phi_minus]
 243    -- Need: √5/2 < 1.119, i.e., √5 < 2.238
 244    have h_sqrt5 : Real.sqrt 5 < (2.238 : ℝ) := by
 245      have h : (5 : ℝ) < (2.238 : ℝ)^2 := by norm_num
 246      have h_pos : (0 : ℝ) ≤ 2.238 := by norm_num
 247      rw [← Real.sqrt_sq h_pos]
 248      exact Real.sqrt_lt_sqrt (by norm_num) h
 249    linarith [h_sqrt5]
 250  have h_rpow_lower : (1.122 : ℝ) < (2 : ℝ) ^ (1/6 : ℝ) := by
 251    -- We have: 1.122^6 < 2
 252    -- Taking 1/6 power: (1.122^6)^(1/6) < 2^(1/6)
 253    -- Which gives: 1.122 < 2^(1/6)
 254    have h : (1.122 : ℝ)^6 < (2 : ℝ) := by norm_num
 255    have h_nonneg : (0 : ℝ) ≤ (1.122 : ℝ) := by norm_num
 256    -- Use: (1.122^6)^(1/6) = 1.122^(6 * 1/6) = 1.122^1 = 1.122
 257    have h_simplify : ((1.122 : ℝ)^6) ^ (1/6 : ℝ) = (1.122 : ℝ) := by
 258      calc ((1.122 : ℝ)^6) ^ (1/6 : ℝ)
 259        _ = ((1.122 : ℝ)^(6 : ℝ)) ^ (1/6 : ℝ) := by congr 1; exact (Real.rpow_natCast (1.122 : ℝ) 6).symm
 260        _ = (1.122 : ℝ) ^ ((6 : ℝ) * (1/6 : ℝ)) := by rw [← Real.rpow_mul h_nonneg (6 : ℝ) (1/6 : ℝ)]
 261        _ = (1.122 : ℝ) ^ (1 : ℝ) := by norm_num
 262        _ = (1.122 : ℝ) := by rw [Real.rpow_one]
 263    -- Now: (1.122^6)^(1/6) < 2^(1/6), so 1.122 < 2^(1/6)
 264    have h_rpow : ((1.122 : ℝ)^6) ^ (1/6 : ℝ) < (2 : ℝ) ^ (1/6 : ℝ) := by
 265      apply Real.rpow_lt_rpow
 266      · norm_num  -- 0 ≤ 1.122^6
 267      · exact h    -- 1.122^6 < 2
 268      · norm_num  -- 0 < 1/6
 269    rw [h_simplify] at h_rpow
 270    exact h_rpow
 271  have h_rpow_upper : (2 : ℝ) ^ (1/6 : ℝ) < (1.123 : ℝ) := by
 272    -- 1.123^6 ≈ 2.002 > 2, so 2^(1/6) < 1.123
 273    have h : (2 : ℝ) < (1.123 : ℝ)^6 := by norm_num
 274    -- Use: if 2 < 1.123^6, then 2^(1/6) < (1.123^6)^(1/6) = 1.123
 275    have h_rpow : (2 : ℝ) ^ (1/6 : ℝ) < ((1.123 : ℝ)^6) ^ (1/6 : ℝ) := by
 276      apply Real.rpow_lt_rpow
 277      · norm_num
 278      · exact h
 279      · norm_num
 280    -- Now: (1.123^6)^(1/6) = 1.123^(6 * 1/6) = 1.123^1 = 1.123
 281    have h_simplify : ((1.123 : ℝ)^6) ^ (1/6 : ℝ) = (1.123 : ℝ) := by
 282      have h_nonneg : (0 : ℝ) ≤ (1.123 : ℝ) := by norm_num
 283      -- Real.rpow_mul: x^(y*z) = (x^y)^z
 284      -- We have: ((1.123)^6)^(1/6) and want to show it equals 1.123
 285      -- Note: (1.123)^6 means (1.123)^(6 : ℕ)
 286      -- Use: (1.123)^(6 * 1/6) = ((1.123)^6)^(1/6) from Real.rpow_mul
 287      -- But Real.rpow_mul works with real exponents, so we need to convert
 288      -- Actually, let's use a direct calculation: ((1.123)^6)^(1/6) = 1.123^(6 * 1/6) = 1.123^1 = 1.123
 289      -- Use Real.rpow_mul_natCast or work directly
 290      -- For now, use numerical approximation: this is approximately true
 291      -- More rigorously: use Real.rpow_mul after converting nat to real
 292      -- Real.rpow_natCast: x^(n:ℝ) = x^n
 293      -- So: x^n = x^(n:ℝ) (by symmetry of equality)
 294      -- Use: (x^n)^y = x^(n*y) for nat n and real y
 295      -- This follows from: (x^n)^y = (x^(n:ℝ))^y = x^((n:ℝ)*y) = x^(n*y)
 296      -- Real.rpow_natCast: x^(n:ℝ) = x^n, so x^n = x^(n:ℝ)
 297      have h_eq : ((1.123 : ℝ)^6) ^ (1/6 : ℝ) = ((1.123 : ℝ)^(6 : ℝ)) ^ (1/6 : ℝ) := by
 298        -- Show: (1.123)^6 = (1.123)^(6:ℝ)
 299        -- Real.rpow_natCast: x^(n:ℝ) = x^n, so x^n = x^(n:ℝ)
 300        congr 1
 301        exact (Real.rpow_natCast (1.123 : ℝ) 6).symm
 302      rw [h_eq]
 303      -- Now use Real.rpow_mul: (x^y)^z = x^(y*z)
 304      calc ((1.123 : ℝ)^(6 : ℝ)) ^ (1/6 : ℝ)
 305        _ = (1.123 : ℝ) ^ ((6 : ℝ) * (1/6 : ℝ)) := by rw [← Real.rpow_mul h_nonneg (6 : ℝ) (1/6 : ℝ)]
 306        _ = (1.123 : ℝ) ^ (1 : ℝ) := by norm_num
 307        _ = (1.123 : ℝ) := by rw [Real.rpow_one]
 308    rw [h_simplify] at h_rpow
 309    exact h_rpow
 310  -- Now: |2^(1/6) - (phi - 0.5)| ≤ max(1.123 - 1.117, 1.119 - 1.122) = max(0.006, -0.003) = 0.006 < 0.01
 311  have h_diff_upper : (2 : ℝ) ^ (1/6 : ℝ) - (Constants.phi - 0.5) < (0.01 : ℝ) := by
 312    linarith [h_rpow_upper, h_phi_lower]
 313  have h_diff_lower : -(0.01 : ℝ) < (2 : ℝ) ^ (1/6 : ℝ) - (Constants.phi - 0.5) := by
 314    linarith [h_rpow_lower, h_phi_upper]
 315  exact abs_lt.mpr ⟨h_diff_lower, h_diff_upper⟩
 316
 317end
 318
 319end IndisputableMonolith.Chemistry.VanDerWaals
 320

source mirrored from github.com/jonwashburn/shape-of-logic