diff --git a/docs/src/DB.md b/docs/src/DB.md
index 7b63406..9a0071f 100644
--- a/docs/src/DB.md
+++ b/docs/src/DB.md
@@ -25,12 +25,13 @@ entry.energy
 
 Add a benchmark using `TwoBody.put!`:
 
-```@repl db-usage
-hamiltonian = Hamiltonian(
-    Kinetic(ℏ = 1.0, m = 1.0),
-    Coulomb(coefficient = -1.0),
-)
-TwoBody.put!(:example, hamiltonian, -0.5)
+```julia-repl
+julia> hamiltonian = Hamiltonian(
+           Kinetic(hbar = 1.0, m = 1.0),
+           Coulomb(coefficient = -1.0),
+       )
+
+julia> TwoBody.put!(:example, hamiltonian, -0.5)
 ```
 
 Duplicate keys are rejected, and Hamiltonians are copied on registration and
diff --git a/docs/src/FDM.md b/docs/src/FDM.md
index 0e24ec9..5248bac 100644
--- a/docs/src/FDM.md
+++ b/docs/src/FDM.md
@@ -38,7 +38,7 @@ Define the [Hamiltoninan](@ref Hamiltonian). This is an example for the non-rela
 ```
 ```@example example
 H = Hamiltonian(
-  Kinetic(ℏ = 1, m = 1),
+  Kinetic(hbar = 1, m = 1),
   Coulomb(coefficient = -1),
 )
 nothing # hide
diff --git a/docs/src/Rayleigh-Ritz.md b/docs/src/Rayleigh-Ritz.md
index e04f37b..8bb2536 100644
--- a/docs/src/Rayleigh-Ritz.md
+++ b/docs/src/Rayleigh-Ritz.md
@@ -35,7 +35,7 @@ Define the [Hamiltoninan](@ref Hamiltonian). This is an example for the non-rela
 
 ```@example example
 H = Hamiltonian(
-  Kinetic(ℏ = 1, m = 1),
+  Kinetic(hbar = 1, m = 1),
   Coulomb(coefficient = -1),
 )
 nothing # hide
diff --git a/docs/src/VMC.md b/docs/src/VMC.md
index 6f9db2a..2e637c5 100644
--- a/docs/src/VMC.md
+++ b/docs/src/VMC.md
@@ -51,7 +51,7 @@ using TwoBody
 
 # Hamiltonian
 H = Hamiltonian(
-  Kinetic(ℏ=1, m=1),
+  Kinetic(hbar=1, m=1),
   Coulomb(coefficient=-1),
 )
 
diff --git a/docs/src/index.md b/docs/src/index.md
index cdb459c..57f8f39 100644
--- a/docs/src/index.md
+++ b/docs/src/index.md
@@ -31,7 +31,7 @@ Define the [Hamiltoninan](@ref Hamiltonian). This is an example for the non-rela
 ```
 ```@example index
 H = Hamiltonian(
-  Kinetic(ℏ = 1, m = 1),
+  Kinetic(hbar = 1, m = 1),
   Coulomb(coefficient = -1),
 )
 nothing # hide
diff --git a/src/DB.jl b/src/DB.jl
index aede81a..c381b73 100644
--- a/src/DB.jl
+++ b/src/DB.jl
@@ -38,7 +38,7 @@ put!(key::AbstractString, hamiltonian::Hamiltonian, energy::Real) =
 put!(
     :hydrogen,
     Hamiltonian(
-        Kinetic(ℏ = 1.0, m = 1.0),
+        Kinetic(hbar = 1.0, m = 1.0),
         Coulomb(coefficient = -1.0),
     ),
     -0.5,
@@ -47,7 +47,7 @@ put!(
 put!(
     :positronium,
     Hamiltonian(
-        Kinetic(ℏ = 1.0, m = 0.5),
+        Kinetic(hbar = 1.0, m = 0.5),
         Coulomb(coefficient = -1.0),
     ),
     -0.25,
@@ -56,7 +56,7 @@ put!(
 put!(
     :harmonic_oscillator,
     Hamiltonian(
-        Kinetic(ℏ = 1.0, m = 1.0),
+        Kinetic(hbar = 1.0, m = 1.0),
         PowerLaw(coefficient = 0.5, exponent = 2.0),
     ),
     1.5,
diff --git a/src/FDM.jl b/src/FDM.jl
index 42ca116..e8b8cec 100644
--- a/src/FDM.jl
+++ b/src/FDM.jl
@@ -37,7 +37,7 @@ end
 function matrix(o::Kinetic, method::FiniteDifferenceMethod)
   D  = FiniteDifferenceMatrices.fdmatrix(Int64(length(method.R)), n=1, m=2, d=method.direction, h=method.Δr, t=typeof(method.Δr))
   D² = FiniteDifferenceMatrices.fdmatrix(Int64(length(method.R)), n=2, m=2, d=method.direction, h=method.Δr, t=typeof(method.Δr))
-  return -o.ℏ^2/2/o.m * (D² + SparseArrays.spdiagm(2 ./ method.R) * D - method.l*(method.l+1) * SparseArrays.spdiagm(1 ./ method.R .^ 2))
+  return -o.hbar^2/2/o.m * (D² + SparseArrays.spdiagm(2 ./ method.R) * D - method.l*(method.l+1) * SparseArrays.spdiagm(1 ./ method.R .^ 2))
 end
 
 function matrix(o::PotentialTerm, method::FiniteDifferenceMethod)
diff --git a/src/Hamiltonian.jl b/src/Hamiltonian.jl
index 0983cef..14133d5 100644
--- a/src/Hamiltonian.jl
+++ b/src/Hamiltonian.jl
@@ -18,7 +18,7 @@ Base.@kwdef struct Laplacian <: KineticTerm
 end
 
 Base.@kwdef struct Kinetic <: KineticTerm
-  ℏ = 1
+  hbar = 1
   m = 1
 end
 
@@ -135,7 +135,7 @@ The Hamiltonian is the input for each solver. This is an example for the non-rel
 
 ```@example
 H = Hamiltonian(
-  Kinetic(ℏ = 1, m = 1),
+  Kinetic(hbar = 1, m = 1),
   Coulomb(coefficient = -1),
 )
 ```
@@ -152,7 +152,7 @@ H = Hamiltonian(
 """ Laplacian
 
 @doc raw"""
-`Kinetic(ℏ=1, m=1)`
+`Kinetic(hbar=1, m=1)`
 ```math
 -\frac{\hbar^2}{2m} \nabla^2
 ```
diff --git a/src/Rayleigh-Ritz.jl b/src/Rayleigh-Ritz.jl
index acccf02..e735cfd 100644
--- a/src/Rayleigh-Ritz.jl
+++ b/src/Rayleigh-Ritz.jl
@@ -389,7 +389,7 @@ function element(o::Laplacian, SGB1::SimpleGaussianBasis, SGB2::SimpleGaussianBa
 end
 
 function element(o::Kinetic, SGB1::SimpleGaussianBasis, SGB2::SimpleGaussianBasis)
-  return o.ℏ^2/(2*o.m) * 6*π^(3/2)*SGB1.a*SGB2.a/(SGB1.a+SGB2.a)^(5/2)
+  return o.hbar^2/(2*o.m) * 6*π^(3/2)*SGB1.a*SGB2.a/(SGB1.a+SGB2.a)^(5/2)
 end
 
 function element(o::Constant, SGB1::SimpleGaussianBasis, SGB2::SimpleGaussianBasis)
diff --git a/src/VMC.jl b/src/VMC.jl
index 00faa19..45d49d0 100644
--- a/src/VMC.jl
+++ b/src/VMC.jl
@@ -56,7 +56,7 @@ function _local_energy(term::Laplacian, wavefunction::Function, position, ψ)
 end
 
 function _local_energy(term::Kinetic, wavefunction::Function, position, ψ)
-  return -term.ℏ^2 / (2 * term.m) * _laplacian(wavefunction, position) / ψ
+  return -term.hbar^2 / (2 * term.m) * _laplacian(wavefunction, position) / ψ
 end
 
 _local_energy(term::RestEnergy, wavefunction::Function, position, ψ) = term.m * term.c^2
diff --git a/test/DB.jl b/test/DB.jl
index 670e44a..d43ae5f 100644
--- a/test/DB.jl
+++ b/test/DB.jl
@@ -19,7 +19,7 @@
 
     key = :temporary_test_problem
     hamiltonian = Hamiltonian(
-        Kinetic(ℏ = 1.0, m = 1.0),
+        Kinetic(hbar = 1.0, m = 1.0),
         Coulomb(coefficient = -2.0),
     )
 
diff --git a/test/FDM.jl b/test/FDM.jl
index 95f2064..9a29164 100644
--- a/test/FDM.jl
+++ b/test/FDM.jl
@@ -3,7 +3,7 @@
   # Testing Results
 
   H = Hamiltonian(
-    Kinetic(ℏ = 1, m = 1),
+    Kinetic(hbar = 1, m = 1),
     Coulomb(coefficient = -1),
   )
   FDM = FiniteDifferenceMethod(
diff --git a/test/Rayleigh-Ritz.jl b/test/Rayleigh-Ritz.jl
index 3b69b66..18eb4d1 100644
--- a/test/Rayleigh-Ritz.jl
+++ b/test/Rayleigh-Ritz.jl
@@ -3,7 +3,7 @@
   # Testing Results
 
   H = Hamiltonian(
-    Kinetic(ℏ = 1, m = 1),
+    Kinetic(hbar = 1, m = 1),
     Coulomb(coefficient = -1),
   )
   BS = BasisSet(
@@ -75,7 +75,7 @@
   @testset "element()" begin
 
     # kinetic terms
-    o = Kinetic(ℏ=1, m=1)
+    o = Kinetic(hbar=1, m=1)
     println(o)
     println("  i  j\tnumerical  \tanalytical")
     l = 0
@@ -85,7 +85,7 @@
       φⱼ(r) = TwoBody.φ(BS.basis[j], r)
       dφⱼ(r) = ForwardDiff.derivative(r -> φⱼ(r), r)
       d²φⱼ(r) = ForwardDiff.derivative(r -> dφⱼ(r), r)
-      numerical  = 4*π*QuadGK.quadgk(r -> r^2 * φᵢ(r) * (-o.ℏ^2/2/o.m * (d²φⱼ(r) + 2/r*dφⱼ(r) - l*(l+1)/r^2*φⱼ(r))), 0, Inf, maxevals=10^3)[1]
+      numerical  = 4*π*QuadGK.quadgk(r -> r^2 * φᵢ(r) * (-o.hbar^2/2/o.m * (d²φⱼ(r) + 2/r*dφⱼ(r) - l*(l+1)/r^2*φⱼ(r))), 0, Inf, maxevals=10^3)[1]
       analytical = TwoBody.element(o, BS.basis[i], BS.basis[j])
       error = isinf(analytical) ? 0.0 : abs((numerical-analytical)/analytical)
       acceptance = error < 1e-5
diff --git a/test/VMC.jl b/test/VMC.jl
index 5398f62..faaa28e 100644
--- a/test/VMC.jl
+++ b/test/VMC.jl
@@ -13,7 +13,7 @@
 
   @testset "local energy" begin
     H = Hamiltonian(
-      Kinetic(ℏ=1, m=1),
+      Kinetic(hbar=1, m=1),
       PowerLaw(coefficient=1 / 2, exponent=2),
     )
     ψ(r) = exp(-sum(abs2, r) / 2)
@@ -28,7 +28,7 @@
 
   @testset "harmonic oscillator" begin
     H = Hamiltonian(
-      Kinetic(ℏ=1, m=1),
+      Kinetic(hbar=1, m=1),
       PowerLaw(coefficient=1 / 2, exponent=2),
     )
     ψ(r) = exp(-sum(abs2, r) / 2)
@@ -53,7 +53,7 @@
     @test mean_radius_squared ≈ 1.5 atol=0.5
 
     singular_H = Hamiltonian(
-      Kinetic(ℏ=1, m=1),
+      Kinetic(hbar=1, m=1),
       Custom(f=r -> r < 1 ? Inf : r^2 / 2),
     )
     singular_result = solve(
@@ -70,7 +70,7 @@
 
   @testset "multiple walkers" begin
     H = Hamiltonian(
-      Kinetic(ℏ=1, m=1),
+      Kinetic(hbar=1, m=1),
       PowerLaw(coefficient=1 / 2, exponent=2),
     )
     ψ(r) = exp(-sum(abs2, r) / 2)
@@ -103,7 +103,7 @@
 
   @testset "Coulomb singularity" begin
     H = Hamiltonian(
-      Kinetic(ℏ=1, m=1),
+      Kinetic(hbar=1, m=1),
       Coulomb(coefficient=-1),
     )
     α = 0.2829