feat: Thermodynamics (and some minor changes to StatisticalMechanics)

PhysLean.lean

@@ -383,6 +383,10 @@ import PhysLean.StringTheory.FTheory.SU5.Quanta.TenQuanta
 import PhysLean.Thermodynamics.Basic
 import PhysLean.Thermodynamics.IdealGas.Basic
 import PhysLean.Thermodynamics.Temperature.Basic
+import PhysLean.Thermodynamics.Temperature.Calculus
+import PhysLean.Thermodynamics.Temperature.Convergence
+import PhysLean.Thermodynamics.Temperature.Regularity
+import PhysLean.Thermodynamics.Temperature.TemperatureScales
 import PhysLean.Thermodynamics.Temperature.TemperatureUnits
 import PhysLean.Units.Basic
 import PhysLean.Units.Dimension

PhysLean/StatisticalMechanics/BoltzmannConstant.lean

@@ -1,5 +1,5 @@
 /-
-Copyright (c) 2025 Joseph Tooby-Smith. All rights reserved.
+Copyright (c) 2026 Joseph Tooby-Smith. All rights reserved.
 Released under Apache 2.0 license as described in the file LICENSE.
 Authors: Joseph Tooby-Smith
 -/
@@ -25,13 +25,13 @@ def kBAx : {p : ℝ | 0 < p} := ⟨1.380649e-23, by norm_num⟩
 
 /-- The Boltzmann constant in a given but arbitrary set of units.
   Boltzman's constant has dimension equivalent to `Energy/Temperature`. -/
-noncomputable def kB : ℝ := kBAx.1
+noncomputable def kB : ℝ := kBAx.val
 
 /-- The Boltzmann constant is positive. -/
-lemma kB_pos : 0 < kB := kBAx.2
+lemma kB_pos : 0 < kB := kBAx.property
 
 /-- The Boltzmann constant is non-negative. -/
-lemma kB_nonneg : 0 ≤ kB := le_of_lt kBAx.2
+lemma kB_nonneg : 0 ≤ kB := le_of_lt kBAx.property
 
 /-- The Boltzmann constant is not equal to zero. -/
 lemma kB_ne_zero : kB ≠ 0 := by

PhysLean/StatisticalMechanics/CanonicalEnsemble/Basic.lean

@@ -370,7 +370,8 @@ lemma μBolt_ne_zero_of_μ_ne_zero (T : Temperature) (h : 𝓒.μ ≠ 0) :
   simp [μBolt] at ⊢ h
   rw [Measure.ext_iff'] at ⊢ h
   simp only [Measure.coe_zero, Pi.zero_apply]
-  have hs : {x | ENNReal.ofReal (rexp (-(↑T.β * 𝓒.energy x))) ≠ 0} = Set.univ := by
+  have hs : {x | ENNReal.ofReal
+                (rexp (-(T.toReal⁻¹ * Constants.kB⁻¹ * 𝓒.energy x))) ≠ 0} = Set.univ := by
     ext i
     simp only [ne_eq, ENNReal.ofReal_eq_zero, not_le, Set.mem_setOf_eq, Set.mem_univ, iff_true]
     exact exp_pos _
@@ -454,7 +455,7 @@ lemma mathematicalPartitionFunction_eq_zero_iff (T : Temperature) [IsFiniteMeasu
   have h : s = Set.univ := by
     ext i
     simp [s]
-    exact exp_pos (-(T.β * 𝓒.energy i))
+    exact exp_pos (-(T.toReal⁻¹ * Constants.kB⁻¹ * 𝓒.energy i))
   change 𝓒.μ s = 0 ↔ 𝓒.μ = 0
   rw [h]
   simp only [Measure.measure_univ_eq_zero]

PhysLean/StatisticalMechanics/CanonicalEnsemble/Finite.lean

@@ -293,8 +293,18 @@ lemma sum_probability_eq_one
   have hZdef := mathematicalPartitionFunction_of_fintype (𝓒:=𝓒) T
   have hZpos := mathematicalPartitionFunction_pos_finite (𝓒:=𝓒) (T:=T)
   have hZne : 𝓒.mathematicalPartitionFunction T ≠ 0 := hZpos.ne'
-  simp [hZdef]
-  simp_all only [neg_mul, ne_eq, not_false_eq_true]
+  have hZdef' : 𝓒.mathematicalPartitionFunction T =
+      ∑ x, rexp (-(T.toReal⁻¹ * kB⁻¹ * 𝓒.energy x)) := by
+    simpa [Temperature.β, one_div, mul_comm, mul_left_comm, mul_assoc] using hZdef
+  have hZne' : ∑ x, rexp (-(T.toReal⁻¹ * kB⁻¹ * 𝓒.energy x)) ≠ 0 := by
+    rw [← hZdef']
+    exact hZne
+  have hnum :
+      ∑ x, rexp (-↑T.β * 𝓒.energy x) = ∑ x, rexp (-(T.toReal⁻¹ * kB⁻¹ * 𝓒.energy x)) := by
+    simp [Temperature.β, one_div, mul_comm, mul_assoc]
+  rw [hZdef']
+  rw [hnum]
+  field_simp [hZne']
 
 /-- The entropy of a finite canonical ensemble (Shannon entropy) is non-negative. -/
 lemma entropy_nonneg [MeasurableSingletonClass ι] [IsFinite 𝓒] [Nonempty ι] (T : Temperature) :
@@ -409,11 +419,29 @@ lemma meanEnergy_Beta_eq_finite [MeasurableSingletonClass ι] [IsFinite 𝓒] (b
     𝓒.meanEnergyBeta b = 𝓒.meanEnergyBetaReal b := by
   let T := Temperature.ofβ (Real.toNNReal b)
   have hT_beta : (T.β : ℝ) = b := by
-    simp [T, Real.toNNReal_of_nonneg hb.le]
+    change ((Temperature.ofβ (Real.toNNReal b)).β : ℝ) = b
+    simpa [Real.toNNReal_of_nonneg hb.le] using
+      congrArg (fun x : NNReal => (x : ℝ)) (Temperature.β_ofβ (Real.toNNReal b))
+  have hT_beta' : T.toReal⁻¹ * kB⁻¹ = b := by
+    simpa [Temperature.β, one_div, mul_comm, mul_left_comm, mul_assoc] using hT_beta
   rw [meanEnergyBeta, meanEnergy_of_fintype 𝓒 T, meanEnergyBetaReal]
   refine Finset.sum_congr rfl fun i _ => ?_
-  simp [CanonicalEnsemble.probability, probabilityBetaReal,
-        mathematicalPartitionFunction_of_fintype, mathematicalPartitionFunctionBetaReal, hT_beta]
+  have hden : 𝓒.mathematicalPartitionFunction T = 𝓒.mathematicalPartitionFunctionBetaReal b := by
+    rw [mathematicalPartitionFunction_of_fintype, mathematicalPartitionFunctionBetaReal]
+    refine Finset.sum_congr rfl ?_
+    intro j hj
+    simp [hT_beta']
+  have hprob :
+      rexp (-(T.toReal⁻¹ * kB⁻¹ * 𝓒.energy i)) / 𝓒.mathematicalPartitionFunction T
+        = rexp (-(b * 𝓒.energy i)) / 𝓒.mathematicalPartitionFunctionBetaReal b := by
+    calc
+      rexp (-(T.toReal⁻¹ * kB⁻¹ * 𝓒.energy i)) / 𝓒.mathematicalPartitionFunction T
+          = rexp (-(b * 𝓒.energy i)) / 𝓒.mathematicalPartitionFunction T := by
+              simp [hT_beta']
+      _ = rexp (-(b * 𝓒.energy i)) / 𝓒.mathematicalPartitionFunctionBetaReal b := by
+            rw [hden]
+  simpa [CanonicalEnsemble.probability, probabilityBetaReal] using
+    congrArg (fun p => 𝓒.energy i * p) hprob
 
 lemma differentiable_meanEnergyBetaReal
     [Nonempty ι] : Differentiable ℝ 𝓒.meanEnergyBetaReal := by
@@ -568,7 +596,7 @@ lemma derivWithin_meanEnergy_Beta_eq_neg_variance
     [MeasurableSingletonClass ι][𝓒.IsFinite] (T : Temperature) (hT_pos : 0 < T.val) :
     derivWithin 𝓒.meanEnergyBeta (Set.Ioi 0) (T.β : ℝ) = - 𝓒.energyVariance T := by
   let β₀ := (T.β : ℝ)
-  have hβ₀_pos : 0 < β₀ := beta_pos T hT_pos
+  have hβ₀_pos : 0 < β₀ := β_pos T hT_pos
   have h_eq_on : Set.EqOn 𝓒.meanEnergyBeta 𝓒.meanEnergyBetaReal (Set.Ioi 0) := by
     intro b hb; exact meanEnergy_Beta_eq_finite 𝓒 b hb
   rw [derivWithin_congr h_eq_on (h_eq_on hβ₀_pos)]
@@ -578,8 +606,11 @@ lemma derivWithin_meanEnergy_Beta_eq_neg_variance
   rw [deriv_meanEnergyBetaReal 𝓒 β₀]
   have h_U_eq : 𝓒.meanEnergyBetaReal β₀ = 𝓒.meanEnergy T := by
     rw [← meanEnergy_Beta_eq_finite 𝓒 β₀ hβ₀_pos]
-    simp [meanEnergyBeta]
-    simp_all only [NNReal.coe_pos, toNNReal_coe, ofβ_β, β₀]
+    change 𝓒.meanEnergy (Temperature.ofβ (Real.toNNReal β₀)) = 𝓒.meanEnergy T
+    have hβ₀_toNNReal : Real.toNNReal β₀ = T.β := by
+      change Real.toNNReal ((T.β : ℝ)) = T.β
+      simpa using (show Real.toNNReal ((T.β : ℝ)) = T.β from Real.toNNReal_coe)
+    rw [hβ₀_toNNReal, Temperature.ofβ_β]
   have h_prob_eq (i : ι) : 𝓒.probabilityBetaReal β₀ i = 𝓒.probability T i := by
     unfold probabilityBetaReal CanonicalEnsemble.probability
     congr 1
@@ -594,7 +625,7 @@ lemma derivWithin_meanEnergy_Beta_eq_neg_variance
 theorem fluctuation_dissipation_theorem_finite
     [MeasurableSingletonClass ι] [𝓒.IsFinite] (T : Temperature) (hT_pos : 0 < T.val) :
     𝓒.heatCapacity T = 𝓒.energyVariance T / (kB * (T.val : ℝ)^2) := by
-  have hβ₀_pos : 0 < (T.β : ℝ) := beta_pos T hT_pos
+  have hβ₀_pos : 0 < (T.β : ℝ) := β_pos T hT_pos
   let β₀ := (T.β : ℝ)
   have h_diff_U_beta : DifferentiableWithinAt ℝ 𝓒.meanEnergyBeta (Set.Ioi 0) β₀ := by
     have h_eq_on : Set.EqOn 𝓒.meanEnergyBeta 𝓒.meanEnergyBetaReal (Set.Ioi 0) := by

PhysLean/StatisticalMechanics/CanonicalEnsemble/Lemmas.lean

@@ -4,6 +4,7 @@
 Authors: Matteo Cipollina, Joseph Tooby-Smith
 -/
 import PhysLean.StatisticalMechanics.CanonicalEnsemble.Basic
+import PhysLean.Thermodynamics.Temperature.Calculus
 /-!
 # Canonical Ensemble: Thermodynamic Identities and Relations
 
@@ -120,7 +121,7 @@
   unfold probability
   simp [Real.log_div, hZpos.ne', Real.log_exp, sub_eq_add_neg]
 
 /-- Auxiliary identity: `kB · β = 1 / T`.
 `β` is defined as `1 / (kB · T)` (see `Temperature.β`). -/
 @[simp]
 lemma kB_mul_beta (T : Temperature) (hT : 0 < T.val) :
@@ -266,7 +267,10 @@
       _ = _ := by
               simp [add_comm, sub_eq_add_neg, mul_comm, mul_left_comm, mul_assoc]
   have hkβT : T.val * (kB * (T.β : ℝ)) = 1 := by
-    simp [hkβ, hTne]
+    rw [hkβ]
+    field_simp [hTne]
+  have hkβT' : T.val * (kB * (kB⁻¹ * T.toReal⁻¹)) = 1 := by
+    simpa [Temperature.β, one_div, mul_comm, mul_left_comm, mul_assoc] using hkβT
   have h_rhs :
       𝓒.meanEnergy T - T.val * 𝓒.thermodynamicEntropy T
         = -kB * T.val *
@@ -289,7 +293,7 @@
       _ = 𝓒.meanEnergy T - 1 * 𝓒.meanEnergy T
             - T.val * kB * Real.log (𝓒.mathematicalPartitionFunction T)
             + T.val * kB * 𝓒.dof * Real.log 𝓒.phaseSpaceunit := by
-              simp [hkβT, mul_comm, mul_assoc]
+              simp [hkβT', mul_comm, mul_assoc]
       _ = -kB * T.val *
             (Real.log (𝓒.mathematicalPartitionFunction T)
               - (𝓒.dof : ℝ) * Real.log 𝓒.phaseSpaceunit) := by
@@ -518,12 +522,11 @@
           (fun β : ℝ => Real.log (∫ i, Real.exp (-β * 𝓒.energy i) ∂𝓒.μ))
           (Set.Ioi 0) (T.β : ℝ)) := by
   set f : ℝ → ℝ := fun β => ∫ i, Real.exp (-β * 𝓒.energy i) ∂𝓒.μ
-  have hβ_pos : 0 < (T.β : ℝ) := beta_pos T hT_pos
+  have hβ_pos : 0 < (T.β : ℝ) := β_pos T hT_pos
   have hZpos : 0 < f (T.β : ℝ) := by
     have hZ := mathematicalPartitionFunction_pos (𝓒:=𝓒) (T:=T)
-    have hEq : f (T.β : ℝ) = 𝓒.mathematicalPartitionFunction T := by
-      simp [f, mathematicalPartitionFunction_eq_integral (𝓒:=𝓒) (T:=T)]
-    simpa [hEq] using hZ
+    simpa [f, mathematicalPartitionFunction_eq_integral (𝓒:=𝓒) (T:=T),
+      Temperature.β, one_div, mul_comm, mul_left_comm, mul_assoc] using hZ
   have h_log :
       HasDerivWithinAt
         (fun β : ℝ => Real.log (f β))
@@ -548,14 +551,12 @@
         = (1 / f (T.β : ℝ)) *
             (- ∫ i, 𝓒.energy i * Real.exp (-(T.β : ℝ) * 𝓒.energy i) ∂𝓒.μ) :=
     h_log.derivWithin hUD
-  have h_f_eval :
-      f (T.β : ℝ) = ∫ i, Real.exp (-(T.β : ℝ) * 𝓒.energy i) ∂𝓒.μ := rfl
   have h_ratio :
       (∫ i, 𝓒.energy i * Real.exp (-(T.β : ℝ) * 𝓒.energy i) ∂𝓒.μ) /
           (∫ i, Real.exp (-(T.β : ℝ) * 𝓒.energy i) ∂𝓒.μ)
         = (1 / f (T.β : ℝ)) *
             (∫ i, 𝓒.energy i * Real.exp (-(T.β : ℝ) * 𝓒.energy i) ∂𝓒.μ) := by
-    simp [h_f_eval, div_eq_mul_inv, mul_comm]
+    simp [f, div_eq_mul_inv, mul_comm, mul_left_comm, mul_assoc]
   calc
     𝓒.meanEnergy T = _ := h_mean_ratio
     _ = (1 / f (T.β : ℝ)) *
@@ -611,16 +612,29 @@
       Real.log (𝓒.partitionFunction (Temperature.ofβ (Real.toNNReal β))) =
         -((𝓒.dof : ℝ) * Real.log 𝓒.phaseSpaceunit)
           + Real.log (∫ i, Real.exp (-β * 𝓒.energy i) ∂ 𝓒.μ) := by
-    have h_integral_pos : 0 < ∫ i, Real.exp (-β * 𝓒.energy i) ∂ 𝓒.μ := by
-      have h_eq : ∫ i, Real.exp (-β * 𝓒.energy i) ∂ 𝓒.μ =
-        ∫ i, Real.exp (-(Real.toNNReal β).val * 𝓒.energy i) ∂ 𝓒.μ := by
-        simp [hβnn]
-      rw [h_eq]
-      simp [mathematicalPartitionFunction_eq_integral
-        (𝓒:=𝓒) (T:=Temperature.ofβ (Real.toNNReal β))] at hZpos
-      simp [hZpos]
     have h_beta_eq : (Temperature.ofβ (Real.toNNReal β)).β = Real.toNNReal β := by
-      simp_all only [gt_iff_lt, mem_Ioi, coe_toNNReal', sup_eq_left, log_pow, neg_mul, β_ofβ]
+      simp_all only [gt_iff_lt, mem_Ioi, coe_toNNReal', sup_eq_left, log_pow, β_ofβ]
+    have h_integral_pos : 0 < ∫ i, Real.exp (-β * 𝓒.energy i) ∂ 𝓒.μ := by
+      have h_int_eq0 :
+          ∫ i, Real.exp (-((Temperature.ofβ (Real.toNNReal β)).β : ℝ) * 𝓒.energy i) ∂ 𝓒.μ =
+            𝓒.mathematicalPartitionFunction (Temperature.ofβ (Real.toNNReal β)) := by
+        exact (mathematicalPartitionFunction_eq_integral
+          (𝓒:=𝓒) (T:=Temperature.ofβ (Real.toNNReal β))).symm
+      have h_int_eq : ∫ i, Real.exp (-β * 𝓒.energy i) ∂ 𝓒.μ =
+          𝓒.mathematicalPartitionFunction (Temperature.ofβ (Real.toNNReal β)) := by
+        have h_int_eq1 :
+            ∫ i, Real.exp (-(β * (kB * (kB⁻¹ * 𝓒.energy i)))) ∂ 𝓒.μ =
+              𝓒.mathematicalPartitionFunction (Temperature.ofβ (Real.toNNReal β)) := by
+          simpa [Temperature.β, hβnn, one_div, mul_comm, mul_left_comm, mul_assoc] using h_int_eq0
+        have h_int_eq2 :
+            ∫ i, Real.exp (-(β * (kB * (kB⁻¹ * 𝓒.energy i)))) ∂ 𝓒.μ =
+              ∫ i, Real.exp (-β * 𝓒.energy i) ∂ 𝓒.μ := by
+          refine integral_congr_ae ?_
+          filter_upwards with i
+          field_simp [kB_ne_zero]
+        exact h_int_eq2.symm.trans h_int_eq1
+      rw [h_int_eq]
+      exact hZpos
     rw [partitionFunction_def,
         mathematicalPartitionFunction_eq_integral (𝓒:=𝓒) (T:=Temperature.ofβ (Real.toNNReal β)),
         h_beta_eq,
@@ -663,7 +677,7 @@
   have h_eq' :
       Set.EqOn F_phys (fun β => F_math β - C) (Set.Ioi (0:ℝ)) := by
     simpa [F_phys, F_math] using h_eq
-  have h_mem : (T.β : ℝ) ∈ Set.Ioi (0:ℝ) := beta_pos T hT_pos
+  have h_mem : (T.β : ℝ) ∈ Set.Ioi (0:ℝ) := β_pos T hT_pos
   have h_congr :
       derivWithin F_phys (Set.Ioi 0) (T.β : ℝ)
         = derivWithin (fun β => F_math β - C) (Set.Ioi 0) (T.β : ℝ) := by
@@ -769,12 +783,11 @@
       = (derivWithin (𝓒.meanEnergyBeta) (Set.Ioi 0) (T.β : ℝ))
         * (-1 / (kB * (T.val : ℝ)^2)) := by
   unfold heatCapacity meanEnergy_T
-  have h_U_eq_comp : (𝓒.meanEnergy_T) = fun t : ℝ => (𝓒.meanEnergyBeta) (betaFromReal t) := by
+  have h_U_eq_comp : (𝓒.meanEnergy_T) = fun t : ℝ => (𝓒.meanEnergyBeta) (βFromReal t) := by
     funext t
-    dsimp [meanEnergy_T, meanEnergyBeta, betaFromReal]
-    simp
+    simp only [meanEnergy_T, meanEnergyBeta, βFromReal, Real.toNNReal_coe, Temperature.ofβ_β]
   let dUdβ := derivWithin (𝓒.meanEnergyBeta) (Set.Ioi 0) (T.β : ℝ)
-  have h_chain := chain_rule_T_beta (F:=𝓒.meanEnergyBeta) (F':=dUdβ) T hT_pos hU_deriv
+  have h_chain := chain_rule_T_β (F:=𝓒.meanEnergyBeta) (F':=dUdβ) T hT_pos hU_deriv
   have h_UD :
     UniqueDiffWithinAt ℝ (Set.Ioi (0 : ℝ)) (T.val : ℝ) :=
     (isOpen_Ioi : IsOpen (Set.Ioi (0 : ℝ))).uniqueDiffWithinAt hT_pos

PhysLean/Thermodynamics/Basic.lean

@@ -1,7 +1,7 @@
 /-
-Copyright (c) 2025 Joseph Tooby-Smith. All rights reserved.
+Copyright (c) 2026 Joseph Tooby-Smith. All rights reserved.
 Released under Apache 2.0 license as described in the file LICENSE.
-Authors: Joseph Tooby-Smith
+Authors: Trong-Nghia Be, Tan-Phuoc-Hung Le, Joseph Tooby-Smith
 -/
 /-!