diff --git a/Project.toml b/Project.toml index 3f87d6c..bf1a57d 100644 --- a/Project.toml +++ b/Project.toml @@ -12,6 +12,7 @@ StaticArrays = "0.9, 0.10, 0.11, 0.12, 1" julia = "1" [extras] +BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf" DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab" Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" @@ -19,4 +20,4 @@ Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] -test = ["Test", "DelimitedFiles", "Random", "Interpolations", "Statistics"] +test = ["Test", "BenchmarkTools", "DelimitedFiles", "Random", "Interpolations", "Statistics"] diff --git a/bench/benchmark.jl b/bench/benchmark.jl index 9d8fb96..fd50551 100644 --- a/bench/benchmark.jl +++ b/bench/benchmark.jl @@ -33,7 +33,7 @@ function benchmark(GridType::Type, n_dims::Int; n_points_per_dim::Int=15, n_rand return trial end -verbosity = 1 +verbosity = 2 for ndims in 1:6 println("\n##################################################") diff --git a/src/GridInterpolations.jl b/src/GridInterpolations.jl index bf516aa..991ae5e 100644 --- a/src/GridInterpolations.jl +++ b/src/GridInterpolations.jl @@ -164,14 +164,21 @@ function interpolants(grid::RectangleGrid, x::AbstractVector) weight = MVector{num_points, eltype(x)}(undef) weight2 = MVector{num_points, eltype(x)}(undef) - # Note: these values are set explicitly because we have not verified that the logic below is independent of the initial values. See discussion in PR #47. These can be removed if it can be proved that the logic is independent of the initial values. - index .= 1 - index2 .= 1 - weight .= zero(eltype(weight)) - weight2 .= zero(eltype(weight2)) + # # Note: these values are set explicitly because we have not verified that the logic below is independent of the initial values. See discussion in PR #47. These can be removed if it can be proved that the logic is independent of the initial values. + for i in 1:num_points + @inbounds index[i] = 1 + @inbounds index2[i] = 1 + @inbounds weight[i] = zero(eltype(x)) + @inbounds weight2[i] = zero(eltype(x)) + end + + # index = @MVector ones(Int, num_points) + # index2 = @MVector ones(Int, num_points) + # weight = @MVector zeros(eltype(x), num_points) + # weight2 = @MVector zeros(eltype(x), num_points) - weight[1] = one(eltype(weight)) - weight2[1] = one(eltype(weight2)) + @inbounds weight[1] = one(eltype(weight)) + @inbounds weight2[1] = one(eltype(weight2)) l = 1 subblock_size = 1 @@ -288,16 +295,18 @@ function interpolants(grid::SimplexGrid, x::AbstractVector) # get weight for i = 1:(length(x_p)+1) if i == 1 - weight[i] = 1 - x_p[i] + @inbounds weight[i] = 1 - x_p[i] elseif i == length(x_p) + 1 - weight[i] = x_p[i-1] + @inbounds weight[i] = x_p[i-1] else - weight[i] = x_p[i-1] - x_p[i] + @inbounds weight[i] = x_p[i-1] - x_p[i] end end # get indices - fill!(index, 0) + for i in 1:length(index) + @inbounds index[i] = 0 + end i_index = 0 for i = 1:(length(x_p)+1) siz = 1 @@ -310,17 +319,17 @@ function interpolants(grid::SimplexGrid, x::AbstractVector) onHi = ((i_index & good_count) > 0) good_count <<= 1 if onHi - index[i] += (ihi[k] - 1 - ct) * siz + @inbounds index[i] += (ihi[k] - 1 - ct) * siz else - index[i] += (ilo[k] - 1 - ct) * siz + @inbounds index[i] += (ilo[k] - 1 - ct) * siz end siz = siz * cut_counts[k] ct += cut_counts[k] end - index[i] += 1 + @inbounds index[i] += 1 end - weight = weight ./ sum(weight) + @inbounds weight = weight ./ sum(weight) return SVector(index), SVector(weight) end diff --git a/test/REQUIRE b/test/REQUIRE deleted file mode 100644 index 97ff54f..0000000 --- a/test/REQUIRE +++ /dev/null @@ -1 +0,0 @@ -Interpolations diff --git a/test/runtests.jl b/test/runtests.jl index 5cfc075..e34da65 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -3,6 +3,8 @@ using Test using LinearAlgebra using Random using DelimitedFiles +using BenchmarkTools +using StaticArrays # use import otherwise Interpolations.interpolate conflicts with GridInterpolations.interpolate import Interpolations: BSpline, Linear, OnGrid, Flat @@ -261,14 +263,20 @@ function test_rect_implemented() @test x == [5, 5] @test interpolants(rgrid, [1, 1]) == ([1], [1.0]) + @show @ballocated(interpolants($rgrid, $([1,1]))) + @test @ballocated(interpolants($rgrid, $([1,1]))) == 0 @test interpolants(rgrid, [2, 5]) == ([3], [1.0]) + @test @ballocated(interpolants($rgrid, $([2, 5]))) == 0 indices, weights = interpolants(rgrid, [1.5, 3]) @test indices == [1, 3] @test isapprox(weights, [0.666667, 0.333333], rtol=1e-5) @test interpolate(rgrid, [1, 2, 3, 4], [1, 1]) == 1.0 + @test @ballocated(interpolate($rgrid, $([1, 2, 3, 4]), $([1, 1]))) == 0 @test maskedInterpolate(rgrid, [1, 2, 3, 4], [1, 1], BitArray([false, false, false, false])) == 1.0 + @test @ballocated(maskedInterpolate($rgrid, $([1, 2, 3, 4]), $([1, 1]), $(BitArray([false, false, false, false])))) == 0 + @test isnan(maskedInterpolate(rgrid, [1, 2, 3, 4], [1, 1], BitArray([true, false, false, false]))) @test isapprox(interpolate(rgrid, [1, 2, 3, 4], [1.5, 3]), 1.66666666, rtol=1e-5) @@ -296,6 +304,8 @@ function test_simplex_implemented() @test x == [5, 5] indices, weights = interpolants(grid, [1, 1]) + @test @ballocated(interpolants($grid, $([1, 1]))) == 0 + full_weight_indices = findall(x -> x == 1, weights) @test length(full_weight_indices) == 1 @test weights[full_weight_indices[1]] == 1.0 @@ -309,7 +319,9 @@ function test_simplex_implemented() @test isapprox(weights[full_weight_indices[1]], 0.333333, rtol=1e-5) == true @test interpolate(grid, [1, 2, 3, 4], [1, 1]) == 1.0 + @test @ballocated(interpolate($grid, $([1, 2, 3, 4]), $([1, 1]))) == 0 @test maskedInterpolate(grid, [1, 2, 3, 4], [1, 1], BitArray([false, false, false, false])) == 1.0 + @test @ballocated(maskedInterpolate($grid, $([1, 2, 3, 4]), $([1, 1]), $(BitArray([false, false, false, false])))) == 0 @test isnan(maskedInterpolate(grid, [1, 2, 3, 4], [1, 1], BitArray([true, false, false, false]))) @test isapprox(interpolate(grid, [1, 2, 3, 4], [1.5, 3]), 1.66666666, rtol=1e-5)