fix 0.7 deprecations and improve performance for "normal" arrays and views (#91)

* split up evaluate to enable better codegen for Arrays

* fix for 0.7

Author: KristofferC

REQUIRE

@@ -1 +1,2 @@
 julia 0.6
+Compat 0.54.0

benchmark/benchmarks.jl

@@ -53,7 +53,7 @@ SUITE["colwise"] = BenchmarkGroup()
 function evaluate_colwise(dist, x, y)
     n = size(x, 2)
     T = typeof(evaluate(dist, x[:, 1], y[:, 1]))
-    r = Vector{T}(n)
+    r = Vector{T}(uninitialized, n)
     for j = 1:n
         r[j] = evaluate(dist, x[:, j], y[:, j])
     end

src/Distances.jl

@@ -2,6 +2,9 @@ __precompile__()
 
 module Distances
 
+using Compat
+using Compat.LinearAlgebra
+
 export
     # generic types/functions
     PreMetric,

src/common.jl

@@ -101,7 +101,7 @@ end
 function sumsq_percol(a::AbstractMatrix{T}) where {T}
     m = size(a, 1)
     n = size(a, 2)
-    r = Vector{T}(n)
+    r = Vector{T}(uninitialized, n)
     for j = 1:n
         aj = view(a, :, j)
         r[j] = dot(aj, aj)
@@ -113,7 +113,7 @@ function wsumsq_percol(w::AbstractArray{T1}, a::AbstractMatrix{T2}) where {T1, T
     m = size(a, 1)
     n = size(a, 2)
     T = typeof(one(T1) * one(T2))
-    r = Vector{T}(n)
+    r = Vector{T}(uninitialized, n)
     for j = 1:n
         aj = view(a, :, j)
         s = zero(T)

src/generic.jl

@@ -32,26 +32,26 @@ result_type(::PreMetric, ::AbstractArray, ::AbstractArray) = Float64
 function colwise!(r::AbstractArray, metric::PreMetric, a::AbstractVector, b::AbstractMatrix)
     n = size(b, 2)
     length(r) == n || throw(DimensionMismatch("Incorrect size of r."))
-    for j = 1:n
-        @inbounds r[j] = evaluate(metric, a, view(b, :, j))
+    @inbounds for j = 1:n
+        r[j] = evaluate(metric, a, view(b, :, j))
     end
     r
 end
 
 function colwise!(r::AbstractArray, metric::PreMetric, a::AbstractMatrix, b::AbstractVector)
     n = size(a, 2)
     length(r) == n || throw(DimensionMismatch("Incorrect size of r."))
-    for j = 1:n
-        @inbounds r[j] = evaluate(metric, view(a, :, j), b)
+    @inbounds for j = 1:n
+        r[j] = evaluate(metric, view(a, :, j), b)
     end
     r
 end
 
 function colwise!(r::AbstractArray, metric::PreMetric, a::AbstractMatrix, b::AbstractMatrix)
     n = get_common_ncols(a, b)
     length(r) == n || throw(DimensionMismatch("Incorrect size of r."))
-    for j = 1:n
-        @inbounds r[j] = evaluate(metric, view(a, :, j), view(b, :, j))
+    @inbounds for j = 1:n
+        r[j] = evaluate(metric, view(a, :, j), view(b, :, j))
     end
     r
 end
@@ -62,19 +62,19 @@ end
 
 function colwise(metric::PreMetric, a::AbstractMatrix, b::AbstractMatrix)
     n = get_common_ncols(a, b)
-    r = Vector{result_type(metric, a, b)}(n)
+    r = Vector{result_type(metric, a, b)}(uninitialized, n)
     colwise!(r, metric, a, b)
 end
 
 function colwise(metric::PreMetric, a::AbstractVector, b::AbstractMatrix)
     n = size(b, 2)
-    r = Vector{result_type(metric, a, b)}(n)
+    r = Vector{result_type(metric, a, b)}(uninitialized, n)
     colwise!(r, metric, a, b)
 end
 
 function colwise(metric::PreMetric, a::AbstractMatrix, b::AbstractVector)
     n = size(a, 2)
-    r = Vector{result_type(metric, a, b)}(n)
+    r = Vector{result_type(metric, a, b)}(uninitialized, n)
     colwise!(r, metric, a, b)
 end
 
@@ -85,10 +85,10 @@ function pairwise!(r::AbstractMatrix, metric::PreMetric, a::AbstractMatrix, b::A
     na = size(a, 2)
     nb = size(b, 2)
     size(r) == (na, nb) || throw(DimensionMismatch("Incorrect size of r."))
-    for j = 1:size(b, 2)
+    @inbounds for j = 1:size(b, 2)
         bj = view(b, :, j)
         for i = 1:size(a, 2)
-            @inbounds r[i, j] = evaluate(metric, view(a, :, i), bj)
+            r[i, j] = evaluate(metric, view(a, :, i), bj)
         end
     end
     r
@@ -101,14 +101,14 @@ end
 function pairwise!(r::AbstractMatrix, metric::SemiMetric, a::AbstractMatrix)
     n = size(a, 2)
     size(r) == (n, n) || throw(DimensionMismatch("Incorrect size of r."))
-    for j = 1:n
+    @inbounds for j = 1:n
         aj = view(a, :, j)
         for i = (j + 1):n
-            @inbounds r[i, j] = evaluate(metric, view(a, :, i), aj)
+            r[i, j] = evaluate(metric, view(a, :, i), aj)
         end
-        @inbounds r[j, j] = 0
+        r[j, j] = 0
         for i = 1:(j - 1)
-            @inbounds r[i, j] = r[j, i]   # leveraging the symmetry of SemiMetric
+            r[i, j] = r[j, i]   # leveraging the symmetry of SemiMetric
         end
     end
     r
@@ -117,12 +117,12 @@ end
 function pairwise(metric::PreMetric, a::AbstractMatrix, b::AbstractMatrix)
     m = size(a, 2)
     n = size(b, 2)
-    r = Matrix{result_type(metric, a, b)}(m, n)
+    r = Matrix{result_type(metric, a, b)}(uninitialized, m, n)
     pairwise!(r, metric, a, b)
 end
 
 function pairwise(metric::PreMetric, a::AbstractMatrix)
     n = size(a, 2)
-    r = Matrix{result_type(metric, a, a)}(n, n)
+    r = Matrix{result_type(metric, a, a)}(uninitialized, n, n)
     pairwise!(r, metric, a)
 end

src/metrics.jl

@@ -145,25 +145,48 @@ SqEuclidean() = SqEuclidean(0)
 #
 ###########################################################
 
-function evaluate(d::UnionMetrics, a::AbstractArray, b::AbstractArray)
-    if length(a) != length(b)
+const ArraySlice{T} = SubArray{T,1,Array{T,2},Tuple{Base.Slice{Base.OneTo{Int}},Int},true}
+
+# Specialized for Arrays and avoids a branch on the size
+@inline Base.@propagate_inbounds function evaluate(d::UnionMetrics, a::Union{Array, ArraySlice}, b::Union{Array, ArraySlice})
+    @boundscheck if length(a) != length(b)
         throw(DimensionMismatch("first array has length $(length(a)) which does not match the length of the second, $(length(b))."))
     end
     if length(a) == 0
         return zero(result_type(d, a, b))
     end
-    s = eval_start(d, a, b)
-    if size(a) == size(b)
+    @inbounds begin
+        s = eval_start(d, a, b)
         @simd for I in eachindex(a, b)
-            @inbounds ai = a[I]
-            @inbounds bi = b[I]
+            ai = a[I]
+            bi = b[I]
             s = eval_reduce(d, s, eval_op(d, ai, bi))
         end
-    else
-        for (Ia, Ib) in zip(eachindex(a), eachindex(b))
-            @inbounds ai = a[Ia]
-            @inbounds bi = b[Ib]
-            s = eval_reduce(d, s, eval_op(d, ai, bi))
+        return eval_end(d, s)
+    end
+end
+
+@inline function evaluate(d::UnionMetrics, a::AbstractArray, b::AbstractArray)
+    @boundscheck if length(a) != length(b)
+        throw(DimensionMismatch("first array has length $(length(a)) which does not match the length of the second, $(length(b))."))
+    end
+    if length(a) == 0
+        return zero(result_type(d, a, b))
+    end
+    @inbounds begin
+        s = eval_start(d, a, b)
+        if size(a) == size(b)
+            @simd for I in eachindex(a, b)
+                ai = a[I]
+                bi = b[I]
+                s = eval_reduce(d, s, eval_op(d, ai, bi))
+            end
+        else
+            for (Ia, Ib) in zip(eachindex(a), eachindex(b))
+                ai = a[Ia]
+                bi = b[Ib]
+                s = eval_reduce(d, s, eval_op(d, ai, bi))
+            end
         end
     end
     return eval_end(d, s)
@@ -200,7 +223,7 @@ cityblock(a::T, b::T) where {T <: Number} = evaluate(Cityblock(), a, b)
 @inline eval_op(::Chebyshev, ai, bi) = abs(ai - bi)
 @inline eval_reduce(::Chebyshev, s1, s2) = max(s1, s2)
 # if only NaN, will output NaN
-@inline eval_start(::Chebyshev, a::AbstractArray, b::AbstractArray) = abs(a[1] - b[1])
+@inline Base.@propagate_inbounds eval_start(::Chebyshev, a::AbstractArray, b::AbstractArray) = abs(a[1] - b[1])
 chebyshev(a::AbstractArray, b::AbstractArray) = evaluate(Chebyshev(), a, b)
 chebyshev(a::T, b::T) where {T <: Number} = evaluate(Chebyshev(), a, b)
 
@@ -218,7 +241,7 @@ hamming(a::AbstractArray, b::AbstractArray) = evaluate(Hamming(), a, b)
 hamming(a::T, b::T) where {T <: Number} = evaluate(Hamming(), a, b)
 
 # Cosine dist
-function eval_start(::CosineDist, a::AbstractArray{T}, b::AbstractArray{T}) where {T <: Real}
+@inline function eval_start(::CosineDist, a::AbstractArray{T}, b::AbstractArray{T}) where {T <: Real}
     zero(T), zero(T), zero(T)
 end
 @inline eval_op(::CosineDist, ai, bi) = ai * bi, ai * ai, bi * bi
@@ -236,6 +259,8 @@ cosine_dist(a::AbstractArray, b::AbstractArray) = evaluate(CosineDist(), a, b)
 # Correlation Dist
 _centralize(x::AbstractArray) = x .- mean(x)
 evaluate(::CorrDist, a::AbstractArray, b::AbstractArray) = cosine_dist(_centralize(a), _centralize(b))
+# Ambiguity resolution
+evaluate(::CorrDist, a::Array, b::Array) = cosine_dist(_centralize(a), _centralize(b))
 corr_dist(a::AbstractArray, b::AbstractArray) = evaluate(CorrDist(), a, b)
 result_type(::CorrDist, a::AbstractArray, b::AbstractArray) = result_type(CosineDist(), a, b)
 
@@ -255,7 +280,7 @@ kl_divergence(a::AbstractArray, b::AbstractArray) = evaluate(KLDivergence(), a,
 gkl_divergence(a::AbstractArray, b::AbstractArray) = evaluate(GenKLDivergence(), a, b)
 
 # RenyiDivergence
-function eval_start(::RenyiDivergence, a::AbstractArray{T}, b::AbstractArray{T}) where {T <: Real}
+@inline Base.@propagate_inbounds function eval_start(::RenyiDivergence, a::AbstractArray{T}, b::AbstractArray{T}) where {T <: Real}
     zero(T), zero(T), T(sum(a)), T(sum(b))
 end
 
@@ -316,7 +341,7 @@ end
 js_divergence(a::AbstractArray, b::AbstractArray) = evaluate(JSDivergence(), a, b)
 
 # SpanNormDist
-function eval_start(::SpanNormDist, a::AbstractArray, b::AbstractArray)
+@inline Base.@propagate_inbounds function eval_start(::SpanNormDist, a::AbstractArray, b::AbstractArray)
     a[1] - b[1], a[1] - b[1]
 end
 @inline eval_op(::SpanNormDist, ai, bi)  = ai - bi

src/wmetrics.jl

@@ -45,37 +45,39 @@ end
 function result_type(dist::UnionWeightedMetrics, ::AbstractArray{T1}, ::AbstractArray{T2}) where {T1, T2}
     typeof(evaluate(dist, one(T1), one(T2)))
 end
-function eval_start(d::UnionWeightedMetrics, a::AbstractArray, b::AbstractArray)
+@inline function eval_start(d::UnionWeightedMetrics, a::AbstractArray, b::AbstractArray)
     zero(result_type(d, a, b))
 end
 eval_end(d::UnionWeightedMetrics, s) = s
 
 
 
-function evaluate(d::UnionWeightedMetrics, a::AbstractArray, b::AbstractArray)
-    if length(a) != length(b)
+@inline function evaluate(d::UnionWeightedMetrics, a::AbstractArray, b::AbstractArray)
+    @boundscheck if length(a) != length(b)
         throw(DimensionMismatch("first array has length $(length(a)) which does not match the length of the second, $(length(b))."))
     end
-    if length(a) != length(d.weights)
+    @boundscheck if length(a) != length(d.weights)
         throw(DimensionMismatch("arrays have length $(length(a)) but weights have length $(length(d.weights))."))
     end
     if length(a) == 0
         return zero(result_type(d, a, b))
     end
-    s = eval_start(d, a, b)
-    if size(a) == size(b)
-        @simd for I in eachindex(a, b, d.weights)
-            @inbounds ai = a[I]
-            @inbounds bi = b[I]
-            @inbounds wi = d.weights[I]
-            s = eval_reduce(d, s, eval_op(d, ai, bi, wi))
-        end
-    else
-        for (Ia, Ib, Iw) in zip(eachindex(a), eachindex(b), eachindex(d.weights))
-            @inbounds ai = a[Ia]
-            @inbounds bi = b[Ib]
-            @inbounds wi = d.weights[Iw]
-            s = eval_reduce(d, s, eval_op(d, ai, bi, wi))
+    @inbounds begin
+        s = eval_start(d, a, b)
+        if size(a) == size(b)
+            @simd for I in eachindex(a, b, d.weights)
+                ai = a[I]
+                bi = b[I]
+                wi = d.weights[I]
+                s = eval_reduce(d, s, eval_op(d, ai, bi, wi))
+            end
+        else
+            for (Ia, Ib, Iw) in zip(eachindex(a), eachindex(b), eachindex(d.weights))
+                ai = a[Ia]
+                bi = b[Ib]
+                wi = d.weights[Iw]
+                s = eval_reduce(d, s, eval_op(d, ai, bi, wi))
+            end
         end
     end
     return eval_end(d, s)

test/F64.jl

@@ -1,19 +1,25 @@
 # dummy type wrapping a Float64 used in tests
-struct F64 <: Real
+struct F64 <: AbstractFloat
     x::Float64
 end
+F64(x::F64) = x
 
 # operations
-for op in (:+, :-)
+for op in (:+, :-, :sin, :cos, :asin, :acos)
     @eval Base.$op(a::F64) = F64($op(a.x))
 end
-for op in (:+, :-, :*, :/)
+for op in (:+, :-, :*, :/, :atan2)
     @eval Base.$op(a::F64, b::F64) = F64($op(a.x, b.x))
 end
-for op in (:zero, :one)
+for op in (:zero, :one,)
     @eval Base.$op(::Type{F64}) = F64($op(Float64))
 end
-Base.rand(rng::AbstractRNG, ::Type{F64}) = F64(rand())
+
+if VERSION.minor >= 7
+    Random.rand(rng::AbstractRNG, ::Random.SamplerTrivial{Random.CloseOpen01{F64}}) = F64(rand(rng))
+else
+    Base.rand(rng::AbstractRNG, ::Type{F64}) = F64(rand())
+end
 Base.sqrt(a::F64) = F64(sqrt(a.x))
 Base.:^(a::F64, b::Number) = F64(a.x^b)
 Base.:^(a::F64, b::Int) = F64(a.x^b)
@@ -32,8 +38,9 @@ Base.eps(::Type{F64}) = eps(Float64)
 # promotion
 Base.promote_type(::Type{Float32}, ::Type{F64}) = Float64 # for eig
 Base.promote_type(::Type{Float64}, ::Type{F64}) = Float64 # for vecnorm
-Base.promote(a::F64, b::T) where {T <: Number} = a, F64(b)
-Base.promote(a::T, b::F64) where {T <: Number} = F64(a), b
+Base.promote(a::F64, b::T) where {T <: Number} = a, F64(float(b))
+Base.promote(a::T, b::F64) where {T <: Number} = F64(float(a)), b
+
 Base.convert(::Type{F64}, a::F64) = a
 Base.convert(::Type{Float64}, a::F64) = a.x
 Base.convert(::Type{F64}, a::T) where {T <: Number} = F64(a)

test/runtests.jl

@@ -1,5 +1,7 @@
 using Distances
 using Compat.Test
+using Compat.LinearAlgebra
+using Compat.Random
 
 include("F64.jl")
 include("test_dists.jl")

test/test_dists.jl

@@ -158,7 +158,7 @@ end
             w = ones(4)
             @test sqeuclidean(x, y) ≈ wsqeuclidean(x, y, w)
 
-            w = rand(size(x))
+            w = rand(Float64, size(x))
             @test wsqeuclidean(x, y, w) ≈ dot((x - vec(y)).^2, w)
             @test weuclidean(x, y, w) == sqrt(wsqeuclidean(x, y, w))
             @test wcityblock(x, y, w) ≈ dot(abs.(x - vec(y)), w)