diff --git a/PhysLean.lean b/PhysLean.lean
index 67cc5431..283c464e 100644
--- a/PhysLean.lean
+++ b/PhysLean.lean
@@ -56,6 +56,7 @@ import PhysLean.Mathematics.Geometry.Metric.Riemannian.Defs
 import PhysLean.Mathematics.InnerProductSpace.Adjoint
 import PhysLean.Mathematics.InnerProductSpace.Basic
 import PhysLean.Mathematics.InnerProductSpace.Calculus
+import PhysLean.Mathematics.KroneckerDelta
 import PhysLean.Mathematics.LinearMaps
 import PhysLean.Mathematics.List
 import PhysLean.Mathematics.List.InsertIdx
diff --git a/PhysLean/Mathematics/KroneckerDelta.lean b/PhysLean/Mathematics/KroneckerDelta.lean
new file mode 100644
index 00000000..f3288474
--- /dev/null
+++ b/PhysLean/Mathematics/KroneckerDelta.lean
@@ -0,0 +1,26 @@
+/-
+Copyright (c) 2026 Gregory J. Loges. All rights reserved.
+Released under Apache 2.0 license as described in the file LICENSE.
+Authors: Gregory J. Loges
+-/
+import Mathlib.Algebra.Ring.Basic
+import PhysLean.Meta.TODO.Basic
+/-!
+
+# Kronecker delta
+
+This file defines compact notation for the Kronecker delta, `δ[i,j] ≔ if (i = j) then 1 else 0`.
+
+-/
+TODO "YVABB" "Build functionality for working with sums involving Kronecker deltas."
+
+
+namespace KroneckerDelta
+
+/-- The Kronecker delta function, `ite (i = j) 1 0`. -/
+def kronecker_delta [Ring R] (i j : Fin d) : R := if (i = j) then 1 else 0
+
+@[inherit_doc kronecker_delta]
+macro "δ[" i:term "," j:term "]" : term => `(kronecker_delta $i $j)
+
+end KroneckerDelta
diff --git a/PhysLean/QuantumMechanics/DDimensions/Operators/AngularMomentum.lean b/PhysLean/QuantumMechanics/DDimensions/Operators/AngularMomentum.lean
index 34dfa0f1..e86ea855 100644
--- a/PhysLean/QuantumMechanics/DDimensions/Operators/AngularMomentum.lean
+++ b/PhysLean/QuantumMechanics/DDimensions/Operators/AngularMomentum.lean
@@ -29,37 +29,31 @@ open ContDiff SchwartzMap
 
 /-- Component `i j` of the angular momentum operator is the continuous linear map
 from `𝓢(Space d, ℂ)` to itself defined by `𝐋ᵢⱼ ≔ 𝐱ᵢ∘𝐩ⱼ - 𝐱ⱼ∘𝐩ᵢ`. -/
-@[sorryful]
 def angularMomentumOperator {d : ℕ} (i j : Fin d) : 𝓢(Space d, ℂ) →L[ℂ] 𝓢(Space d, ℂ) :=
   𝐱[i] ∘L 𝐩[j] - 𝐱[j] ∘L 𝐩[i]
 
 @[inherit_doc angularMomentumOperator]
 macro "𝐋[" i:term "," j:term "]" : term => `(angularMomentumOperator $i $j)
 
-@[sorryful]
 lemma angularMomentumOperator_apply_fun {d : ℕ} (i j : Fin d) (ψ : 𝓢(Space d, ℂ)) :
     𝐋[i,j] ψ = 𝐱[i] (𝐩[j] ψ) - 𝐱[j] (𝐩[i] ψ) := rfl
 
-@[sorryful]
 lemma angularMomentumOperator_apply {d : ℕ} (i j : Fin d) (ψ : 𝓢(Space d, ℂ)) (x : Space d) :
     𝐋[i,j] ψ x = 𝐱[i] (𝐩[j] ψ) x - 𝐱[j] (𝐩[i] ψ) x := rfl
 
 /-- The square of the angular momentum operator, `𝐋² ≔ ½ ∑ᵢⱼ 𝐋ᵢⱼ∘𝐋ᵢⱼ`. -/
-@[sorryful]
 def angularMomentumOperatorSqr {d : ℕ} : 𝓢(Space d, ℂ) →L[ℂ] 𝓢(Space d, ℂ) :=
   ∑ i, ∑ j, (2 : ℂ)⁻¹ • 𝐋[i,j] ∘L 𝐋[i,j]
 
 @[inherit_doc angularMomentumOperatorSqr]
 notation "𝐋²" => angularMomentumOperatorSqr
 
-@[sorryful]
 lemma angularMomentumOperatorSqr_apply_fun {d : ℕ} (ψ : 𝓢(Space d, ℂ)) :
     𝐋² ψ = ∑ i, ∑ j, (2 : ℂ)⁻¹ • 𝐋[i,j] (𝐋[i,j] ψ) := by
   dsimp only [angularMomentumOperatorSqr]
   simp only [ContinuousLinearMap.coe_sum', ContinuousLinearMap.coe_smul',
     ContinuousLinearMap.coe_comp', Finset.sum_apply, Pi.smul_apply, Function.comp_apply]
 
-@[sorryful]
 lemma angularMomentumOperatorSqr_apply {d : ℕ} (ψ : 𝓢(Space d, ℂ)) (x : Space d) :
     𝐋² ψ x = ∑ i, ∑ j, (2 : ℂ)⁻¹ * 𝐋[i,j] (𝐋[i,j] ψ) x := by
   rw [angularMomentumOperatorSqr_apply_fun]
@@ -74,12 +68,10 @@ lemma angularMomentumOperatorSqr_apply {d : ℕ} (ψ : 𝓢(Space d, ℂ)) (x :
 -/
 
 /-- The angular momentum operator is antisymmetric, `𝐋ᵢⱼ = -𝐋ⱼᵢ` -/
-@[sorryful]
 lemma angularMomentumOperator_antisymm {d : ℕ} (i j : Fin d) : 𝐋[i,j] = - 𝐋[j,i] :=
   Eq.symm (neg_sub _ _)
 
 /-- Angular momentum operator components with repeated index vanish, `𝐋ᵢᵢ = 0`. -/
-@[sorryful]
 lemma angularMomentumOperator_eq_zero {d : ℕ} (i : Fin d) : 𝐋[i,i] = 0 := sub_self _
 
 /-
@@ -98,18 +90,15 @@ lemma angularMomentumOperator_eq_zero {d : ℕ} (i : Fin d) : 𝐋[i,i] = 0 := s
 -/
 
 /-- In one dimension the angular momentum operator is trivial. -/
-@[sorryful]
 lemma angularMomentumOperator1D_trivial : ∀ (i j : Fin 1), 𝐋[i,j] = 0 := by
   intro i j
   fin_cases i, j
   exact angularMomentumOperator_eq_zero 0
 
 /-- The angular momentum (pseudo)scalar operator in two dimensions, `𝐋 ≔ 𝐋₀₁`. -/
-@[sorryful]
 def angularMomentumOperator2D : 𝓢(Space 2, ℂ) →L[ℂ] 𝓢(Space 2, ℂ) := 𝐋[0,1]
 
 /-- The angular momentum (pseudo)vector operator in three dimension, `𝐋ᵢ ≔ ½ ∑ⱼₖ εᵢⱼₖ 𝐋ⱼₖ`. -/
-@[sorryful]
 def angularMomentumOperator3D (i : Fin 3) : 𝓢(Space 3, ℂ) →L[ℂ] 𝓢(Space 3, ℂ) :=
   match i with
     | 0 => 𝐋[1,2]
diff --git a/PhysLean/QuantumMechanics/DDimensions/Operators/Commutation.lean b/PhysLean/QuantumMechanics/DDimensions/Operators/Commutation.lean
index 29506242..6b9645d3 100644
--- a/PhysLean/QuantumMechanics/DDimensions/Operators/Commutation.lean
+++ b/PhysLean/QuantumMechanics/DDimensions/Operators/Commutation.lean
@@ -3,6 +3,7 @@ Copyright (c) 2026 Gregory J. Loges. All rights reserved.
 Released under Apache 2.0 license as described in the file LICENSE.
 Authors: Gregory J. Loges
 -/
+import PhysLean.Mathematics.KroneckerDelta
 import PhysLean.QuantumMechanics.DDimensions.Operators.AngularMomentum
 /-!
 
@@ -13,98 +14,226 @@ import PhysLean.QuantumMechanics.DDimensions.Operators.AngularMomentum
 namespace QuantumMechanics
 noncomputable section
 open Constants
-open SchwartzMap
+open KroneckerDelta
+open SchwartzMap ContinuousLinearMap
+
+private lemma ite_cond_symm (i j : Fin d) :
+    (if i = j then A else B) = (if j = i then A else B) :=
+  ite_cond_congr (Eq.propIntro Eq.symm Eq.symm)
+
+lemma lie_leibniz_left {d : ℕ} (A B C : 𝓢(Space d, ℂ) →L[ℂ] 𝓢(Space d, ℂ)) :
+    ⁅A ∘L B, C⁆ = A ∘L ⁅B, C⁆ + ⁅A, C⁆ ∘L B := by
+  dsimp only [Bracket.bracket]
+  simp only [ContinuousLinearMap.mul_def, comp_assoc, comp_sub, sub_comp, sub_add_sub_cancel]
+
+lemma lie_leibniz_right {d : ℕ} (A B C : 𝓢(Space d, ℂ) →L[ℂ] 𝓢(Space d, ℂ)) :
+    ⁅A, B ∘L C⁆ = B ∘L ⁅A, C⁆ + ⁅A, B⁆ ∘L C := by
+  dsimp only [Bracket.bracket]
+  simp only [ContinuousLinearMap.mul_def, comp_assoc, comp_sub, sub_comp, sub_add_sub_cancel']
 
 /-
 ## Position / position commutators
 -/
 
-/-- `[xᵢ, xⱼ] = 0` -/
-@[sorryful]
+/-- Position operators commute: `[xᵢ, xⱼ] = 0`. -/
 lemma position_commutation_position {d : ℕ} (i j : Fin d) : ⁅𝐱[i], 𝐱[j]⁆ = 0 := by
   dsimp only [Bracket.bracket]
   ext ψ x
-  simp only [ContinuousLinearMap.coe_sub', ContinuousLinearMap.coe_mul, Pi.sub_apply, sub_apply,
-    Function.comp_apply, ContinuousLinearMap.zero_apply, zero_apply, positionOperator_apply]
+  simp only [coe_sub', coe_mul, Pi.sub_apply, Function.comp_apply, SchwartzMap.sub_apply,
+    ContinuousLinearMap.zero_apply, SchwartzMap.zero_apply, positionOperator_apply]
   ring
 
 /-
 ## Momentum / momentum commutators
 -/
 
-/-- `[pᵢ, pⱼ] = 0` -/
-@[sorryful]
+/-- Momentum operators commute: `[pᵢ, pⱼ] = 0`. -/
 lemma momentum_commutation_momentum {d : ℕ} (i j : Fin d) : ⁅𝐩[i], 𝐩[j]⁆ = 0 := by
-  sorry
+  dsimp only [Bracket.bracket]
+  ext ψ x
+  simp only [coe_sub', coe_mul, Pi.sub_apply, Function.comp_apply, SchwartzMap.sub_apply,
+    ContinuousLinearMap.zero_apply, SchwartzMap.zero_apply, momentumOperator_apply_fun]
+  rw [Space.deriv_const_smul _ ?_, Space.deriv_const_smul _ ?_]
+  · rw [Space.deriv_commute _ (ψ.smooth _), sub_self]
+  · exact Space.deriv_differentiable (ψ.smooth _) i
+  · exact Space.deriv_differentiable (ψ.smooth _) j
 
-@[sorryful]
 lemma momentum_momentum_eq {d : ℕ} (i j : Fin d) : 𝐩[i] ∘L 𝐩[j] = 𝐩[j] ∘L 𝐩[i] := by
   rw [← sub_eq_zero]
   exact momentum_commutation_momentum i j
 
-/-- `[𝐩², 𝐩ᵢ] = 0` -/
-@[sorryful]
-lemma momentumSqr_commutation_momentum {d : ℕ} (i : Fin d) : 𝐩² ∘L 𝐩[i] - 𝐩[i] ∘L 𝐩² = 0 := by
-  dsimp only [momentumOperatorSqr]
-  simp only [ContinuousLinearMap.finset_sum_comp, ContinuousLinearMap.comp_finset_sum]
-  rw [sub_eq_zero]
-  congr
-  ext j ψ x
-  rw [ContinuousLinearMap.comp_assoc]
-  rw [momentum_momentum_eq]
-  rw [← ContinuousLinearMap.comp_assoc]
-  rw [momentum_momentum_eq]
-  rfl
+lemma momentumSqr_commutation_momentum {d : ℕ} (i : Fin d) :
+    ⁅momentumOperatorSqr (d := d), 𝐩[i]⁆ = 0 := by
+  dsimp only [Bracket.bracket, momentumOperatorSqr]
+  rw [Finset.mul_sum, Finset.sum_mul, ← Finset.sum_sub_distrib]
+  conv_lhs =>
+    enter [2, j]
+    simp only [ContinuousLinearMap.mul_def]
+    rw [comp_assoc]
+    rw [momentum_momentum_eq j i, ← comp_assoc]
+    rw [momentum_momentum_eq j i, comp_assoc]
+    rw [sub_self]
+  rw [Finset.sum_const_zero]
 
 /-
 ## Position / momentum commutators
 -/
 
-/-- `[xᵢ, pⱼ] = iℏ δᵢⱼ𝟙` -/
-@[sorryful]
+/-- The canonical commutation relations: `[xᵢ, pⱼ] = iℏ δᵢⱼ𝟙`. -/
 lemma position_commutation_momentum {d : ℕ} (i j : Fin d) : ⁅𝐱[i], 𝐩[j]⁆ =
-    (Complex.I * ℏ) • (if i = j then 1 else 0) • ContinuousLinearMap.id ℂ 𝓢(Space d, ℂ) := by
-  sorry
+    (Complex.I * ℏ * δ[i,j]) • ContinuousLinearMap.id ℂ 𝓢(Space d, ℂ) := by
+  dsimp only [Bracket.bracket, kronecker_delta]
+  ext ψ x
+  simp only [ContinuousLinearMap.smul_apply, SchwartzMap.smul_apply, coe_id', id_eq, smul_eq_mul,
+    coe_sub', coe_mul, Pi.sub_apply, Function.comp_apply, SchwartzMap.sub_apply]
+  rw [positionOperator_apply, momentumOperator_apply_fun]
+  rw [momentumOperator_apply, positionOperator_apply_fun]
+  simp only [neg_mul, Pi.smul_apply, smul_eq_mul, mul_neg, sub_neg_eq_add]
+
+  have h : (fun x ↦ ↑(x i) * ψ x) = (fun (x : Space d) ↦ x i) • ψ := rfl
+  rw [h]
+  rw [Space.deriv_smul (by fun_prop) (SchwartzMap.differentiableAt ψ)]
+  rw [Space.deriv_component, ite_cond_symm j i]
+  simp only [mul_add, Complex.real_smul, ite_smul, one_smul, zero_smul, mul_ite, mul_one, mul_zero,
+    ite_mul, zero_mul]
+  ring
+
+lemma position_position_commutation_momentum {d : ℕ} (i j k : Fin d) : ⁅𝐱[i] ∘L 𝐱[j], 𝐩[k]⁆ =
+    (Complex.I * ℏ * δ[i,k]) • 𝐱[j] + (Complex.I * ℏ * δ[j,k]) • 𝐱[i] := by
+  rw [lie_leibniz_left]
+  rw [position_commutation_momentum, position_commutation_momentum]
+  rw [ContinuousLinearMap.comp_smul, ContinuousLinearMap.smul_comp]
+  rw [id_comp, comp_id]
+  rw [add_comm]
+
+lemma position_commutation_momentum_momentum {d : ℕ} (i j k : Fin d) : ⁅𝐱[i], 𝐩[j] ∘L 𝐩[k]⁆ =
+    (Complex.I * ℏ * δ[i,k]) • 𝐩[j] + (Complex.I * ℏ * δ[i,j]) • 𝐩[k] := by
+  rw [lie_leibniz_right]
+  rw [position_commutation_momentum, position_commutation_momentum]
+  rw [ContinuousLinearMap.comp_smul, ContinuousLinearMap.smul_comp]
+  rw [id_comp, comp_id]
+
+lemma position_commutation_momentumSqr {d : ℕ} (i : Fin d) : ⁅𝐱[i], 𝐩²⁆ =
+    (2 * Complex.I * ℏ) • 𝐩[i] := by
+  unfold momentumOperatorSqr
+  rw [lie_sum]
+  simp only [position_commutation_momentum_momentum]
+  dsimp only [kronecker_delta]
+  simp only [mul_ite, mul_one, mul_zero, ite_smul, zero_smul, Finset.sum_add_distrib,
+    Finset.sum_ite_eq, Finset.mem_univ, ↓reduceIte]
+  ext ψ x
+  simp only [ContinuousLinearMap.add_apply, coe_smul', Pi.smul_apply, SchwartzMap.add_apply,
+    SchwartzMap.smul_apply, smul_eq_mul]
+  ring
 
 /-
 ## Angular momentum / position commutators
 -/
 
-@[sorryful]
 lemma angularMomentum_commutation_position {d : ℕ} (i j k : Fin d) : ⁅𝐋[i,j], 𝐱[k]⁆ =
-    (Complex.I * ℏ) • ((if i = k then 1 else 0) * 𝐱[j] - (if j = k then 1 else 0) * 𝐱[i]) := by
-  sorry
-
-@[sorryful]
-lemma angularMomentumSqr_commutation_position {d : ℕ} (i : Fin d) :
-    𝐋² ∘L 𝐱[i] - 𝐱[i] ∘L 𝐋² = 0 := by
-  sorry
+    (Complex.I * ℏ * δ[i,k]) • 𝐱[j] - (Complex.I * ℏ * δ[j,k]) • 𝐱[i] := by
+  unfold angularMomentumOperator
+  rw [sub_lie]
+  rw [lie_leibniz_left, lie_leibniz_left]
+  rw [position_commutation_position, position_commutation_position]
+  rw [← lie_skew, position_commutation_momentum]
+  rw [← lie_skew, position_commutation_momentum]
+  dsimp only [kronecker_delta]
+  rw [ite_cond_symm k i, ite_cond_symm k j]
+  simp only [ContinuousLinearMap.comp_neg, ContinuousLinearMap.comp_smul, comp_id, zero_comp,
+    add_zero, add_comm, sub_neg_eq_add, ← sub_eq_add_neg]
 
 /-
 ## Angular momentum / momentum commutators
 -/
 
-@[sorryful]
 lemma angularMomentum_commutation_momentum {d : ℕ} (i j k : Fin d) : ⁅𝐋[i,j], 𝐩[k]⁆ =
-    (Complex.I * ℏ) • ((if i = k then 1 else 0) * 𝐩[j] - (if j = k then 1 else 0) * 𝐩[i]) := by
-  sorry
+    (Complex.I * ℏ * δ[i,k]) • 𝐩[j] - (Complex.I * ℏ * δ[j,k]) • 𝐩[i] := by
+  unfold angularMomentumOperator
+  rw [sub_lie]
+  rw [lie_leibniz_left, lie_leibniz_left]
+  rw [momentum_commutation_momentum, momentum_commutation_momentum]
+  rw [position_commutation_momentum, position_commutation_momentum]
+  simp only [ContinuousLinearMap.smul_comp, id_comp, comp_zero, zero_add]
+
+lemma angularMomentum_commutation_momentumSqr {d : ℕ} (i j : Fin d) :
+    ⁅𝐋[i,j], momentumOperatorSqr (d := d)⁆ = 0 := by
+  unfold momentumOperatorSqr
+  conv_lhs =>
+    rw [lie_sum]
+    enter [2, k]
+    rw [lie_leibniz_right]
+    rw [angularMomentum_commutation_momentum]
+    simp only [comp_sub, comp_smulₛₗ, RingHom.id_apply, sub_comp, smul_comp]
+    rw [momentum_momentum_eq _ i, momentum_momentum_eq j _]
+  dsimp only [kronecker_delta]
+  simp only [Finset.sum_add_distrib, Finset.sum_sub_distrib, mul_ite, mul_zero, ite_smul,
+    zero_smul, Finset.sum_ite_eq, Finset.mem_univ, ↓reduceIte, sub_self, add_zero]
+
+lemma angularMomentumSqr_commutation_momentumSqr {d : ℕ} :
+    ⁅angularMomentumOperatorSqr (d := d), momentumOperatorSqr (d := d)⁆ = 0 := by
+  unfold angularMomentumOperatorSqr
+  conv_lhs =>
+    rw [sum_lie]
+    enter [2, i]
+    rw [sum_lie]
+    enter [2, j]
+    rw [smul_lie, lie_leibniz_left]
+    rw [angularMomentum_commutation_momentumSqr]
+    rw [comp_zero, zero_comp, add_zero, smul_zero]
+  simp only [Finset.sum_const_zero]
 
 /-
 ## Angular momentum / angular momentum commutators
 -/
 
-@[sorryful]
 lemma angularMomentum_commutation_angularMomentum {d : ℕ} (i j k l : Fin d) : ⁅𝐋[i,j], 𝐋[k,l]⁆ =
-    (Complex.I * ℏ) • ((if i = k then 1 else 0) * 𝐋[j,l]
-                      - (if i = l then 1 else 0) * 𝐋[j,k]
-                      - (if j = k then 1 else 0) * 𝐋[i,l]
-                      + (if j = l then 1 else 0) * 𝐋[i,k]) := by
-  sorry
+    (Complex.I * ℏ * δ[i,k]) • 𝐋[j,l] - (Complex.I * ℏ * δ[i,l]) • 𝐋[j,k]
+    - (Complex.I * ℏ * δ[j,k]) • 𝐋[i,l] + (Complex.I * ℏ * δ[j,l]) • 𝐋[i,k] := by
+  nth_rw 2 [angularMomentumOperator]
+  rw [lie_sub]
+  rw [lie_leibniz_right, lie_leibniz_right]
+  rw [angularMomentum_commutation_momentum, angularMomentum_commutation_position]
+  rw [angularMomentum_commutation_momentum, angularMomentum_commutation_position]
+  dsimp only [angularMomentumOperator, kronecker_delta]
+  simp only [ContinuousLinearMap.comp_sub, ContinuousLinearMap.sub_comp,
+    ContinuousLinearMap.comp_smul, ContinuousLinearMap.smul_comp]
+  ext ψ x
+  simp only [mul_ite, mul_one, mul_zero, ite_smul, zero_smul, coe_sub', Pi.sub_apply,
+    ContinuousLinearMap.add_apply, SchwartzMap.sub_apply, SchwartzMap.add_apply, smul_sub]
+  ring
 
-@[sorryful]
 lemma angularMomentumSqr_commutation_angularMomentum {d : ℕ} (i j : Fin d) :
-    𝐋² ∘L 𝐋[i,j] - 𝐋[i,j] ∘L 𝐋² = 0 := by
-  sorry
+    ⁅angularMomentumOperatorSqr (d := d), 𝐋[i,j]⁆ = 0 := by
+  unfold angularMomentumOperatorSqr
+  conv_lhs =>
+    rw [sum_lie]
+    enter [2, k]
+    rw [sum_lie]
+    enter [2, l]
+    simp only [smul_lie]
+    rw [lie_leibniz_left]
+    rw [angularMomentum_commutation_angularMomentum]
+  dsimp only [kronecker_delta]
+  simp only [comp_add, comp_sub, add_comp, sub_comp, comp_smul, smul_comp, smul_add, smul_sub,
+    smul_smul, mul_ite, mul_zero, mul_one, ← mul_assoc]
+  simp only [ite_smul, zero_smul]
+
+  -- Split into individual terms to do one of the sums, then recombine
+  simp only [Finset.sum_add_distrib, Finset.sum_sub_distrib, Finset.sum_ite_irrel,
+    Finset.sum_const_zero, Finset.sum_ite_eq', Finset.mem_univ, ↓reduceIte]
+  simp only [← Finset.sum_add_distrib, ← Finset.sum_sub_distrib]
+
+  ext ψ x
+  simp only [angularMomentumOperator_antisymm _ i, angularMomentumOperator_antisymm j _,
+    neg_comp, comp_neg, neg_neg, smul_neg, sub_neg_eq_add]
+  simp only [ContinuousLinearMap.sum_apply, ContinuousLinearMap.add_apply,
+    ContinuousLinearMap.sub_apply, ContinuousLinearMap.smul_apply, ContinuousLinearMap.comp_apply,
+    ContinuousLinearMap.neg_apply, ContinuousLinearMap.zero_apply, SchwartzMap.add_apply,
+    SchwartzMap.sum_apply, SchwartzMap.sub_apply, SchwartzMap.smul_apply, SchwartzMap.neg_apply,
+    SchwartzMap.zero_apply]
+  ring_nf
+  exact Fintype.sum_eq_zero _ (congrFun rfl)
 
 end
 end QuantumMechanics
diff --git a/PhysLean/QuantumMechanics/DDimensions/Operators/Position.lean b/PhysLean/QuantumMechanics/DDimensions/Operators/Position.lean
index 1bb54147..c7a4db2e 100644
--- a/PhysLean/QuantumMechanics/DDimensions/Operators/Position.lean
+++ b/PhysLean/QuantumMechanics/DDimensions/Operators/Position.lean
@@ -18,9 +18,80 @@ noncomputable section
 open Space
 open ContDiff SchwartzMap
 
+/-
+Some helper lemmas for computing the iterated derivatives of `fun (x : Space d) ↦ (x i : ℂ)`,
+generalized from `ℂ` to any `RCLike`.
+-/
+
+private lemma norm_fderiv_ofReal_coord {M} [RCLike M] {d : ℕ} {i : Fin d} {x : Space d}:
+    ‖fderiv ℝ (RCLike.ofRealCLM (K := M) ∘L coordCLM i) x‖ = 1 := by
+  rw [ContinuousLinearMap.fderiv]
+  rw [ContinuousLinearMap.norm_def]
+  simp only [ContinuousLinearMap.coe_comp', RCLike.ofRealCLM_apply, Function.comp_apply,
+    norm_algebraMap', Real.norm_eq_abs]
+  conv_lhs =>
+    enter [1, 1, c, 2, 2]
+    rw [Space.coordCLM_apply, Space.coord_apply]
+
+  refine csInf_eq_of_forall_ge_of_forall_gt_exists_lt ?_ ?_ ?_
+  · use 1
+    simp only [Set.mem_setOf_eq, zero_le_one, one_mul, true_and]
+    exact fun x => abs_eval_le_norm x i
+  · intro c
+    rw [Set.mem_setOf_eq, and_imp]
+    intro hc h
+    have h' : ∃ (x : Space d), |x i| = ‖x‖ ∧ 0 < |x i| := by
+      use (basis i)
+      rw [basis_self, abs_one, OrthonormalBasis.norm_eq_one]
+      simp only [zero_lt_one, and_self]
+    rcases h' with ⟨x2, hx2, hx2'⟩
+    apply one_le_of_le_mul_right₀ hx2'
+    apply le_of_le_of_eq (h x2) (congrArg (HMul.hMul c) (Eq.symm hx2))
+  · intro w hw
+    use 1
+    simp only [Set.mem_setOf_eq, zero_le_one, one_mul, true_and]
+    exact ⟨fun x => abs_eval_le_norm x i, hw⟩
+
+private lemma iteratedFDeriv_succ_succ_ofReal_coord_eq_zero {M} [RCLike M] {d : ℕ} {i : Fin d} :
+    ∀ n : ℕ, iteratedFDeriv ℝ n.succ.succ (RCLike.ofRealCLM (K := M) ∘L coordCLM i) = 0 := by
+  intro n
+  induction n with
+  | zero =>
+    ext x _
+    rw [iteratedFDeriv_succ_apply_right]
+    conv_lhs =>
+      enter [1, 1, 3, y]
+      change fderiv ℝ (RCLike.ofRealCLM ∘L coordCLM i) y
+      rw [ContinuousLinearMap.fderiv]
+    simp only [Nat.reduceAdd, iteratedFDeriv_one_apply, fderiv_fun_const, Pi.zero_apply,
+      Fin.isValue, ContinuousLinearMap.zero_apply, Fin.reduceLast, Nat.succ_eq_add_one,
+      ContinuousMultilinearMap.zero_apply]
+  | succ n hn =>
+    ext x _
+    rw [iteratedFDeriv_succ_apply_left]
+    rw [hn]
+    simp only [fderiv_zero, Pi.zero_apply, ContinuousLinearMap.zero_apply,
+      ContinuousMultilinearMap.zero_apply, Nat.succ_eq_add_one]
+
+private lemma norm_iteratedFDeriv_ofRealCLM_coord {M} [RCLike M] (n : ℕ) {d : ℕ} (i : Fin d)
+    (x : Space d) : ‖iteratedFDeriv ℝ n (RCLike.ofRealCLM (K := M) ∘L coordCLM i) x‖ =
+    if n = 0 then |x i| else if n = 1 then 1 else 0 := by
+  match n with
+  | 0 =>
+    simp only [ContinuousLinearMap.coe_comp', RCLike.ofRealCLM_apply, norm_iteratedFDeriv_zero,
+      Function.comp_apply, norm_algebraMap', Real.norm_eq_abs, ↓reduceIte]
+    rw [coordCLM_apply, coord_apply]
+  | 1 =>
+    simp only [ContinuousLinearMap.coe_comp', one_ne_zero, ↓reduceIte]
+    rw [← norm_iteratedFDeriv_fderiv, norm_iteratedFDeriv_zero]
+    exact norm_fderiv_ofReal_coord
+  | .succ (.succ n) =>
+    rw [iteratedFDeriv_succ_succ_ofReal_coord_eq_zero]
+    simp only [Nat.succ_eq_add_one, Pi.zero_apply, norm_zero, Nat.add_eq_zero_iff, one_ne_zero,
+      and_false, and_self, ↓reduceIte, Nat.add_eq_right]
+
 /-- Component `i` of the position operator is the continuous linear map
 from `𝓢(Space d, ℂ)` to itself which maps `ψ` to `xᵢψ`. -/
-@[sorryful]
 def positionOperator {d : ℕ} (i : Fin d) : 𝓢(Space d, ℂ) →L[ℂ] 𝓢(Space d, ℂ) := by
   refine SchwartzMap.mkCLM (fun ψ x ↦ x i * ψ x) ?hadd ?hsmul ?hsmooth ?hbound
   -- hadd
@@ -50,7 +121,6 @@ def positionOperator {d : ℕ} (i : Fin d) : 𝓢(Space d, ℂ) →L[ℂ] 𝓢(S
       * ‖iteratedFDeriv ℝ j (fun x ↦ (x i : ℂ)) x‖
       * ‖iteratedFDeriv ℝ (n - j) ψ x‖
     · apply (mul_le_mul_of_nonneg_left ?_ (pow_nonneg (norm_nonneg x) k))
-
       have hcd : ContDiff ℝ ∞ (fun (x : Space d) ↦ (x i : ℂ)) := by
         apply ContDiff.fun_comp
         · change ContDiff ℝ ∞ RCLike.ofRealCLM
@@ -58,35 +128,19 @@ def positionOperator {d : ℕ} (i : Fin d) : 𝓢(Space d, ℂ) →L[ℂ] 𝓢(S
         · fun_prop
       apply norm_iteratedFDeriv_mul_le (N := ∞) hcd (SchwartzMap.smooth ψ ⊤) x ENat.LEInfty.out
 
-    -- h0, h1 and hj are the analogues of `norm_iteratedFDeriv_ofRealCLM ℂ j`
-    -- but including a projection to the i-th component of x
-    have h0 : ‖iteratedFDeriv ℝ 0 (fun x ↦ (x i : ℂ)) x‖ = ‖x i‖ := by
-      simp only [norm_iteratedFDeriv_zero, Complex.norm_real, Real.norm_eq_abs]
-
-    have h1 : ‖iteratedFDeriv ℝ 1 (fun x ↦ (x i : ℂ)) x‖ = 1 := by
-      rw [← norm_iteratedFDeriv_fderiv, norm_iteratedFDeriv_zero]
-      sorry
-
-    have hj : ∀ (j : ℕ), ‖iteratedFDeriv ℝ (j + 2) (fun x ↦ (x i : ℂ)) x‖ = 0 := by
-      intro j
-      rw [← norm_iteratedFDeriv_fderiv, ← norm_iteratedFDeriv_fderiv]
-      sorry
-
-    have hproj : ∀ (j : ℕ), ‖iteratedFDeriv ℝ j (fun x ↦ (x i : ℂ)) x‖ =
-        if j = 0 then ‖x i‖ else if j = 1 then 1 else 0 := by
-      intro j
-      match j with
-        | 0 => rw [h0]; simp
-        | 1 => rw [h1]; simp
-        | k + 2 => rw [hj]; simp
-
+    have h' : (fun x ↦ ↑(x i)) = RCLike.ofRealCLM (K := ℂ) ∘L coordCLM i := by
+      ext x
+      simp only [ContinuousLinearMap.coe_comp', RCLike.ofRealCLM_apply, Complex.coe_algebraMap,
+        Function.comp_apply, Complex.ofReal_inj]
+      rw [Space.coordCLM_apply, Space.coord_apply]
+    rw [h']
     conv_lhs =>
-      enter [2, 2, j, 1, 2]
-      rw [hproj]
+      enter [2, 2, j]
+      rw [norm_iteratedFDeriv_ofRealCLM_coord]
 
     match n with
       | 0 =>
-        simp only [zero_add, Finset.range_one, Real.norm_eq_abs, mul_ite, mul_one, mul_zero,
+        simp only [zero_add, Finset.range_one, mul_ite, mul_one, mul_zero,
           ite_mul, zero_mul, Finset.sum_singleton, ↓reduceIte, Nat.choose_self, Nat.cast_one,
           one_mul, Nat.sub_zero, norm_iteratedFDeriv_zero, CharP.cast_eq_zero]
         trans (SchwartzMap.seminorm ℝ (k + 1) 0) ψ
@@ -102,7 +156,7 @@ def positionOperator {d : ℕ} (i : Fin d) : 𝓢(Space d, ℂ) →L[ℂ] 𝓢(S
         simp only [Nat.succ_eq_add_one, Nat.add_eq_zero_iff, one_ne_zero, and_false, and_self,
           ↓reduceIte, Nat.add_eq_right, mul_zero, zero_mul, Finset.sum_const_zero, zero_add,
           Nat.choose_one_right, Nat.cast_add, Nat.cast_one, mul_one, Nat.reduceAdd,
-          Nat.add_one_sub_one, Nat.choose_zero_right, Real.norm_eq_abs, one_mul, Nat.sub_zero,
+          Nat.add_one_sub_one, Nat.choose_zero_right, one_mul, Nat.sub_zero,
           add_tsub_cancel_right, ge_iff_le]
 
         trans (↑n + 1) * (‖x‖ ^ k * ‖iteratedFDeriv ℝ n ψ x‖)
@@ -145,11 +199,9 @@ def positionOperator {d : ℕ} (i : Fin d) : 𝓢(Space d, ℂ) →L[ℂ] 𝓢(S
 @[inherit_doc positionOperator]
 macro "𝐱[" i:term "]" : term => `(positionOperator $i)
 
-@[sorryful]
 lemma positionOperator_apply_fun {d : ℕ} (i : Fin d) (ψ : 𝓢(Space d, ℂ)) :
     𝐱[i] ψ = (fun x ↦ x i * ψ x) := rfl
 
-@[sorryful]
 lemma positionOperator_apply {d : ℕ} (i : Fin d) (ψ : 𝓢(Space d, ℂ)) (x : Space d) :
     𝐱[i] ψ x = x i * ψ x := rfl
 
diff --git a/PhysLean/SpaceAndTime/Space/Derivatives/Basic.lean b/PhysLean/SpaceAndTime/Space/Derivatives/Basic.lean
index 9bd94e7f..54f1c697 100644
--- a/PhysLean/SpaceAndTime/Space/Derivatives/Basic.lean
+++ b/PhysLean/SpaceAndTime/Space/Derivatives/Basic.lean
@@ -132,20 +132,29 @@ lemma deriv_coord_add (f1 f2 : Space d → EuclideanSpace ℝ (Fin d))
 
 -/
 
-/-- Scalar multiplication on space derivatives. -/
-lemma deriv_smul [NormedAddCommGroup M] [NormedSpace ℝ M]
-    (f : Space d → M) (k : ℝ) (hf : Differentiable ℝ f) :
-    ∂[u] (k • f) = (k • ∂[u] f) := by
+/-- Space derivatives on scalar product of functions. -/
+lemma deriv_smul [NormedAddCommGroup M] [NormedSpace ℝ M] [NontriviallyNormedField 𝕜]
+    [NormedAlgebra ℝ 𝕜] [NormedSpace 𝕜 M] {c : Space d → 𝕜} {f : Space d → M}
+    (hc : DifferentiableAt ℝ c x) (hf : DifferentiableAt ℝ f x) :
+    ∂[u] (c • f) x = c x • ∂[u] f x + ∂[u] c x • f x := by
+  unfold deriv
+  rw [fderiv_smul hc hf]
+  rfl
+
+/-- Space derivatives on scalar times function. -/
+lemma deriv_const_smul [NormedAddCommGroup M] [NormedSpace ℝ M] [Semiring R]
+    [Module R M] [SMulCommClass ℝ R M] [ContinuousConstSMul R M] {f : Space d → M} (c : R)
+    (h : Differentiable ℝ f) : ∂[u] (c • f) = c • ∂[u] f := by
   unfold deriv
   ext x
   rw [fderiv_const_smul]
-  rfl
+  rw [ContinuousLinearMap.coe_smul', Pi.smul_apply, Pi.smul_apply]
   fun_prop
 
 /-- Coordinate-wise scalar multiplication on space derivatives. -/
 lemma deriv_coord_smul (f : Space d → EuclideanSpace ℝ (Fin d)) (k : ℝ)
     (hf : Differentiable ℝ f) :
-    ∂[u] (fun x => k * f x i) x= (k • ∂[u] (fun x => f x i)) x:= by
+    ∂[u] (fun x => k * f x i) x = k * ∂[u] (fun x => f x i) x := by
   unfold deriv
   rw [fderiv_const_mul]
   simp only [ContinuousLinearMap.coe_smul', Pi.smul_apply, smul_eq_mul]
diff --git a/PhysLean/SpaceAndTime/Space/Derivatives/Curl.lean b/PhysLean/SpaceAndTime/Space/Derivatives/Curl.lean
index 41c9027e..3b29816a 100644
--- a/PhysLean/SpaceAndTime/Space/Derivatives/Curl.lean
+++ b/PhysLean/SpaceAndTime/Space/Derivatives/Curl.lean
@@ -129,7 +129,6 @@ lemma curl_smul (f : Space → EuclideanSpace ℝ (Fin 3)) (k : ℝ)
   fin_cases i <;>
   · simp only [Fin.isValue, Pi.smul_apply, PiLp.smul_apply, smul_eq_mul, Fin.zero_eta]
     rw [deriv_coord_smul, deriv_coord_smul, mul_sub]
-    simp only [Fin.isValue, Pi.smul_apply, smul_eq_mul]
     repeat fun_prop
 
 /-!
diff --git a/PhysLean/SpaceAndTime/Space/Derivatives/Grad.lean b/PhysLean/SpaceAndTime/Space/Derivatives/Grad.lean
index 3612186d..6eae8e3a 100644
--- a/PhysLean/SpaceAndTime/Space/Derivatives/Grad.lean
+++ b/PhysLean/SpaceAndTime/Space/Derivatives/Grad.lean
@@ -116,7 +116,7 @@ lemma grad_smul (f : Space d → ℝ) (k : ℝ)
   unfold grad
   ext x i
   simp only [Pi.smul_apply]
-  rw [deriv_smul]
+  rw [deriv_const_smul]
   rfl
   exact hf