diff --git a/src/LinearAlgebra.jl b/src/LinearAlgebra.jl index 2c76dcdc..e6f4f1de 100644 --- a/src/LinearAlgebra.jl +++ b/src/LinearAlgebra.jl @@ -130,6 +130,10 @@ export lyap, mul!, norm, + lpdist, + l1dist, + euclidean, + l2dist, normalize!, normalize, nullspace, diff --git a/src/generic.jl b/src/generic.jl index 2b03b249..520c1507 100644 --- a/src/generic.jl +++ b/src/generic.jl @@ -892,6 +892,89 @@ opnorm(v::TransposeAbsVec) = norm(v.parent) norm(v::AdjOrTrans, p::Real) = norm(v.parent, p) +""" + lpdist(x::T, y::T, p::Int) where {T <: AbstractVector} + lpdist(x::T, y::T, p::Real) where {T <: AbstractVector} + + Computes the Lₚ distance between two vectors `x` and `y`, defined as: + for 0 < p < Inf: Lₚ dist = ᵖ√(Σ|xᵢ - yᵢ|ᵖ) + for p = 0: count of mismatches for p = 0 + for p + Inf: L_∞ = maxᵢ(abs(xᵢ - yᵢ)) +""" + +function _lpunnorm(x::T, y::T, fn) where {T} + r = 0.0 + for i in eachindex(x,y) + @inbounds r += fn(x[i] - y[i]) + end + + return r +end + + +function lpdist(x::T, y::T, p::Int) where {T <: AbstractVector} + p < 0 && throw(DomainError("p must be non-negative")) + length(x) == length(y) || throw(DimensionMismatch("x and y have different lenghts")) + if p == 0 + @warn("Technically not a distance metric for p = 0") + fn = !iszero # simply count number of nonzeros + elseif p == 1 + fn = abs + elseif p == 2 + fn = abs2 + else + fn = u -> abs(u)^p + end + r = _lpunnorm(x, y, fn) + + p <= 1 && return r + p == 2 && return √r + + return r^(1/p) +end + +function lpdist(x::T, y::T, p::Real) where {T <: AbstractVector} + length(x) == length(y) || throw(DimensionMismatch("x and y have different lenghts")) + p < 0 && throw(DomainError("p must be non-negative")) + p < 1 && @warn("Technically not a distance metric for 0 < p < 1") + + # handle inf norm separatey + if p == Inf + r = 0.0 + for i in eachindex(x,y) + @inbounds r = max(abs(x[i] - y[i]), r) + end + return r + elseif iszero(p) + return lpdist(x, y, 0) + else + fn = u -> abs(u)^p + end + r = _lpunnorm(x, y, fn) + + return r^(1/p) +end + +""" + euclidean(x::T, y::T) + l2dist(x::T, y::T) + + Compute the L₂ distance between vectors x and y. + Basically just aliases for lpdist(x, y, 2) +""" +euclidean(x::T, y::T) where {T <: AbstractVector} = lpdist(x, y, 2) +l2dist(x::T, y::T) where {T <: AbstractVector} = lpdist(x, y, 2) + +""" + l1dist(x::T, y::T) + + Compute the L₁ distance between vectors x and y. + Basically just an alias for lpdist(x, y, 1) +""" +l1dist(x::T, y::T) where {T <: AbstractVector} = lpdist(x, y, 1) + + +## dot products """ dot(x, y) x ⋅ y diff --git a/test/generic.jl b/test/generic.jl index 6d11ec82..232d7be9 100644 --- a/test/generic.jl +++ b/test/generic.jl @@ -837,4 +837,20 @@ end end end +@testset "distance functions" begin + x = [1., 2., 3.]; y = [1., 2.1, 3.4] + @test lpdist(x, y, 1) ≈ norm(x .- y, 1) + @test lpdist(x, y, 2) ≈ norm(x .- y, 2) + @test lpdist(x, y, 4) ≈ norm(x .- y, 4) + @test lpdist(x, y, 3.5) ≈ norm(x .- y, 3.5) + @test euclidean(x, y) ≈ norm(x .- y, 2) + @test l2dist(x, y) ≈ norm(x .- y, 2) + @test l1dist(x, y) ≈ norm(x .- y, 1) + @test lpdist(x, y, 0) ≈ norm(x .- y, 0) + @test lpdist(x, y, Inf) ≈ norm(x .- y, Inf) + @test_throws DomainError lpdist(x, y, -1) + @test_throws DomainError lpdist(x, y, -0.5) + @test_throws DimensionMismatch lpdist(x, [1.5, 2.0, 2.0, 3.1], 1) +end + end # module TestGeneric