Fix #483: Rename neg_energy → neg_kinetic_energy

src/hamiltonian.jl

@@ -67,7 +67,7 @@ function ∂H∂r(h::Hamiltonian{<:DenseEuclideanMetric,<:GaussianKinetic}, r::A
     return M⁻¹ * r
 end
 
-# TODO (kai) make the order of θ and r consistent with neg_energy
+# TODO (kai) make the order of θ and r consistent with neg_kinetic_energy
 # TODO (kai) add stricter types to block hamiltonian.jl#L37 from working on unknown metric/kinetic
 # The gradient of a position-dependent Hamiltonian system depends on both θ and r. 
 ∂H∂θ(h::Hamiltonian, θ::AbstractVecOrMat, r::AbstractVecOrMat) = ∂H∂θ(h, θ)
@@ -101,7 +101,7 @@ function Base.similar(z::PhasePoint{<:AbstractVecOrMat{T}}) where {T<:AbstractFl
 end
 
 function phasepoint(
-    h::Hamiltonian, θ::T, r::T; ℓπ=∂H∂θ(h, θ), ℓκ=DualValue(neg_energy(h, r, θ), ∂H∂r(h, r))
+    h::Hamiltonian, θ::T, r::T; ℓπ=∂H∂θ(h, θ), ℓκ=DualValue(neg_kinetic_energy(h, r, θ), ∂H∂r(h, r))
 ) where {T<:AbstractVecOrMat}
     return PhasePoint(θ, r, ℓπ, ℓκ)
 end
@@ -115,7 +115,7 @@ function phasepoint(
     _r::T2;
     r=safe_rsimilar(θ, _r),
     ℓπ=∂H∂θ(h, θ),
-    ℓκ=DualValue(neg_energy(h, r, θ), ∂H∂r(h, r)),
+    ℓκ=DualValue(neg_kinetic_energy(h, r, θ), ∂H∂r(h, r)),
 ) where {T1<:AbstractVecOrMat,T2<:AbstractVecOrMat}
     return PhasePoint(θ, r, ℓπ, ℓκ)
 end
@@ -140,38 +140,39 @@ neg_energy(h::Hamiltonian, θ::AbstractVecOrMat) = h.ℓπ(θ)
 
 # GaussianKinetic
 
-function neg_energy(
+function neg_kinetic_energy(
     h::Hamiltonian{<:UnitEuclideanMetric,<:GaussianKinetic}, r::T, θ::T
 ) where {T<:AbstractVector}
     return -sum(abs2, r) / 2
 end
 
-function neg_energy(
+function neg_kinetic_energy(
     h::Hamiltonian{<:UnitEuclideanMetric,<:GaussianKinetic}, r::T, θ::T
 ) where {T<:AbstractMatrix}
     return -vec(sum(abs2, r; dims=1)) / 2
 end
 
-function neg_energy(
+function neg_kinetic_energy(
     h::Hamiltonian{<:DiagEuclideanMetric,<:GaussianKinetic}, r::T, θ::T
 ) where {T<:AbstractVector}
     return -sum(abs2.(r) .* h.metric.M⁻¹) / 2
 end
 
-function neg_energy(
+function neg_kinetic_energy(
     h::Hamiltonian{<:DiagEuclideanMetric,<:GaussianKinetic}, r::T, θ::T
 ) where {T<:AbstractMatrix}
     return -vec(sum(abs2.(r) .* h.metric.M⁻¹; dims=1)) / 2
 end
 
-function neg_energy(
+function neg_kinetic_energy(
     h::Hamiltonian{<:DenseEuclideanMetric,<:GaussianKinetic}, r::T, θ::T
 ) where {T<:AbstractVecOrMat}
     mul!(h.metric._temp, h.metric.M⁻¹, r)
     return -dot(r, h.metric._temp) / 2
 end
 
 energy(args...) = -neg_energy(args...)
+energy(h::Hamiltonian, r::AbstractVecOrMat, θ::AbstractVecOrMat) = -neg_kinetic_energy(h, r, θ)
 
 ####
 #### Momentum refreshment

src/riemannian/hamiltonian.jl

@@ -103,7 +103,7 @@ function step(
     return res
 end
 
-# TODO Make the order of θ and r consistent with neg_energy
+# TODO Make the order of θ and r consistent with neg_kinetic_energy
 ∂H∂θ(h::Hamiltonian, θ::AbstractVecOrMat, r::AbstractVecOrMat) = ∂H∂θ(h, θ)
 ∂H∂r(h::Hamiltonian, θ::AbstractVecOrMat, r::AbstractVecOrMat) = ∂H∂r(h, r)
 
@@ -221,7 +221,7 @@ end
 
 ### hamiltonian.jl
 
-import AdvancedHMC: phasepoint, neg_energy, ∂H∂θ, ∂H∂r
+import AdvancedHMC: phasepoint, neg_kinetic_energy, ∂H∂θ, ∂H∂r
 using LinearAlgebra: logabsdet, tr
 
 # QUES Do we want to change everything to position dependent by default?
@@ -231,14 +231,14 @@ function phasepoint(
     θ::T,
     r::T;
     ℓπ=∂H∂θ(h, θ),
-    ℓκ=DualValue(neg_energy(h, r, θ), ∂H∂r(h, θ, r)),
+    ℓκ=DualValue(neg_kinetic_energy(h, r, θ), ∂H∂r(h, θ, r)),
 ) where {T<:AbstractVecOrMat}
     return PhasePoint(θ, r, ℓπ, ℓκ)
 end
 
 # Negative kinetic energy
 #! Eq (13) of Girolami & Calderhead (2011)
-function neg_energy(
+function neg_kinetic_energy(
     h::Hamiltonian{<:DenseRiemannianMetric}, r::T, θ::T
 ) where {T<:AbstractVecOrMat}
     G = h.metric.map(h.metric.G(θ))

test/hamiltonian.jl

@@ -60,19 +60,19 @@ end
             r_init = randn(T, D)
 
             h = Hamiltonian(UnitEuclideanMetric(T, D), ℓπ, ∂ℓπ∂θ)
-            @test -AdvancedHMC.neg_energy(h, r_init, θ_init) == sum(abs2, r_init) / 2
+            @test -AdvancedHMC.neg_kinetic_energy(h, r_init, θ_init) == sum(abs2, r_init) / 2
             @test AdvancedHMC.∂H∂r(h, r_init) == r_init
 
             M⁻¹ = ones(T, D) + abs.(randn(T, D))
             h = Hamiltonian(DiagEuclideanMetric(M⁻¹), ℓπ, ∂ℓπ∂θ)
-            @test -AdvancedHMC.neg_energy(h, r_init, θ_init) ≈
+            @test -AdvancedHMC.neg_kinetic_energy(h, r_init, θ_init) ≈
                 r_init' * diagm(0 => M⁻¹) * r_init / 2
             @test AdvancedHMC.∂H∂r(h, r_init) == M⁻¹ .* r_init
 
             m = randn(T, D, D)
             M⁻¹ = m' * m
             h = Hamiltonian(DenseEuclideanMetric(M⁻¹), ℓπ, ∂ℓπ∂θ)
-            @test -AdvancedHMC.neg_energy(h, r_init, θ_init) ≈ r_init' * M⁻¹ * r_init / 2
+            @test -AdvancedHMC.neg_kinetic_energy(h, r_init, θ_init) ≈ r_init' * M⁻¹ * r_init / 2
             @test AdvancedHMC.∂H∂r(h, r_init) == M⁻¹ * r_init
         end
     end
@@ -86,15 +86,15 @@ end
             r_init = ComponentArray(; a=randn(T, D), b=randn(T, D))
 
             h = Hamiltonian(UnitEuclideanMetric(T, 2 * D), ℓπ, ∂ℓπ∂θ)
-            @test -AdvancedHMC.neg_energy(h, r_init, θ_init) == sum(abs2, r_init) / 2
+            @test -AdvancedHMC.neg_kinetic_energy(h, r_init, θ_init) == sum(abs2, r_init) / 2
             @test AdvancedHMC.∂H∂r(h, r_init) == r_init
             @test typeof(AdvancedHMC.∂H∂r(h, r_init)) == typeof(r_init)
 
             M⁻¹ = ComponentArray(;
                 a=ones(T, D) + abs.(randn(T, D)), b=ones(T, D) + abs.(randn(T, D))
             )
             h = Hamiltonian(DiagEuclideanMetric(M⁻¹), ℓπ, ∂ℓπ∂θ)
-            @test -AdvancedHMC.neg_energy(h, r_init, θ_init) ≈
+            @test -AdvancedHMC.neg_kinetic_energy(h, r_init, θ_init) ≈
                 r_init' * diagm(0 => M⁻¹) * r_init / 2
             @test AdvancedHMC.∂H∂r(h, r_init) == M⁻¹ .* r_init
             @test typeof(AdvancedHMC.∂H∂r(h, r_init)) == typeof(r_init)
@@ -103,7 +103,7 @@ end
             ax = getaxes(r_init)[1]
             M⁻¹ = ComponentArray(m' * m, ax, ax)
             h = Hamiltonian(DenseEuclideanMetric(M⁻¹), ℓπ, ∂ℓπ∂θ)
-            @test -AdvancedHMC.neg_energy(h, r_init, θ_init) ≈ r_init' * M⁻¹ * r_init / 2
+            @test -AdvancedHMC.neg_kinetic_energy(h, r_init, θ_init) ≈ r_init' * M⁻¹ * r_init / 2
             @test all(AdvancedHMC.∂H∂r(h, r_init) .== M⁻¹ * r_init)
             @test typeof(AdvancedHMC.∂H∂r(h, r_init)) == typeof(r_init)
         end

test/riemannian.jl

@@ -6,7 +6,7 @@ include("../src/riemannian_hmc_utility.jl")
 using FiniteDiff:
     finite_difference_gradient, finite_difference_hessian, finite_difference_jacobian
 using Distributions: MvNormal
-using AdvancedHMC: neg_energy, energy
+using AdvancedHMC: neg_kinetic_energy, energy
 
 # Taken from https://github.com/JuliaDiff/FiniteDiff.jl/blob/master/test/finitedifftests.jl
 δ(a, b) = maximum(abs.(a - b))
@@ -45,7 +45,7 @@ using AdvancedHMC: neg_energy, energy
                 all(iszero, x) # or for x==0 that I know it's PD
                 @testset "Kinetic energy" begin
                     Σ = hamiltonian.metric.map(hamiltonian.metric.G(x))
-                    @test neg_energy(hamiltonian, r, x) ≈ logpdf(MvNormal(zeros(D), Σ), r)
+                    @test neg_kinetic_energy(hamiltonian, r, x) ≈ logpdf(MvNormal(zeros(D), Σ), r)
                 end
             end