Skip to content

Ind2x mvec #32

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
mattuntergassmair wants to merge 4 commits into
base: master
Choose a base branch
from
Open

Ind2x mvec #32

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
b26d129
added benchmarking function for ind2x
Mar 5, 2020
4b0554e
use mvec for fixed size vector in ind2x
Mar 5, 2020
a7ac2aa
tests fix
Mar 5, 2020
16403a7
ind2x minor fixes
Mar 6, 2020
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion Project.toml
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
name = "GridInterpolations"
uuid = "bb4c363b-b914-514b-8517-4eb369bc008a"
version = "1.1.2"
version = "1.1.3"

[deps]
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
Expand Down
52 changes: 36 additions & 16 deletions bench/benchmark.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,20 +2,20 @@ using GridInterpolations
using Statistics
using BenchmarkTools

function benchmark_interpolate(grid, data, rand_points)
function _benchmark_interpolate(grid, data, rand_points)
for point in eachcol(rand_points)
interpolate(grid, data, point)
end
end

function benchmark(GridType::Type, n_dims::Int; n_points_per_dim::Int=15, n_rand_points::Int=100, verbosity::Int=1)
function benchmark_interpolate(GridType::Type, n_dims::Int; n_points_per_dim::Int=15, n_rand_points::Int=100, verbosity::Int=1)

# construct grid
grid = GridType((range(0, 1, length=n_points_per_dim) for _ in 1:n_dims)...)
rand_points = rand(n_dims, n_rand_points)
rand_data = rand(length(grid))

trial = @benchmark benchmark_interpolate($grid, $rand_data, $rand_points) samples=10000
trial = @benchmark _benchmark_interpolate($grid, $rand_data, $rand_points) samples=10000

if verbosity > 0
println("Benchmarking interpolation on $(label(grid))")
Expand All @@ -33,21 +33,41 @@ function benchmark(GridType::Type, n_dims::Int; n_points_per_dim::Int=15, n_rand
return trial
end

verbosity = 1
for ndims in 1:6
function _benchmark_ind2x(grid, indices)
x = zeros(ndims(grid))
for i in indices
x += ind2x(grid, i)
end
end

function benchmark_ind2x(n_dims=6, n_points_per_dim=10)
grid = RectangleGrid((range(0, 1, length=n_points_per_dim) for _ in 1:n_dims)...)
trial = @benchmark _benchmark_ind2x($grid, $(1:length(grid))) samples=1000
end

function benchmark(verbosity = 1)

for ndims in 1:6

println("\n##################################################")
println("### Benchmark for ndims = $(ndims)")
println("##################################################")
println("\n##################################################")
println("### Benchmark for ndims = $(ndims)")
println("##################################################")

trial_rectangle = benchmark(RectangleGrid, ndims; verbosity=verbosity)
trial_simplex = benchmark(SimplexGrid, ndims; verbosity=verbosity)
trial_rectangle = benchmark_interpolate(RectangleGrid, ndims; verbosity=verbosity)
trial_simplex = benchmark_interpolate(SimplexGrid, ndims; verbosity=verbosity)

# println(ratio(trial_estimate_rectangle, trial_estimate_simplex))
println(
"RectangleGrid vs SimplexGrid: ",
judge(median(trial_rectangle), median(trial_simplex))
)
println("##################################################\n")
# println(ratio(trial_estimate_rectangle, trial_estimate_simplex))
println(
"RectangleGrid vs SimplexGrid: ",
judge(median(trial_rectangle), median(trial_simplex))
)

println("ind2x:", benchmark_ind2x())
println("##################################################\n")


end

end

benchmark()
7 changes: 3 additions & 4 deletions src/GridInterpolations.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -108,11 +108,10 @@ function Base.show(io::IO, grid::AbstractGrid)
end
end

function ind2x(grid::AbstractGrid, ind::Int)
ndims = dimensions(grid)
x = Array{Float64}(undef, ndims)
function ind2x(grid::AbstractGrid{D}, ind::Int) where D
x::MVector{D, Float64} = @MVector zeros(D)

@MaximeBouton MaximeBouton Mar 6, 2020

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

if we want to reduce allocation to zero, adding a container for x in grid might do it. @zsunberg @mattuntergassmair

@mattuntergassmair mattuntergassmair Mar 6, 2020

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You mean adding a MVector{D, Float64} as a field to the RectangleGrid struct, so we can simply use that to store the result instead of allocating - did I understand the suggestion correctly?
It would help with the allocations, but I think what would happen is that if you call ind2x() twice, and you are still holding on to the result computed in the first call, then it gets overwritten in the second call without the user noticing, right?

@mattuntergassmair mattuntergassmair Mar 6, 2020

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I guess in that case users should just pre-allocate memory on their side and use ind2x!() rather than ind2x()

@zsunberg zsunberg Mar 6, 2020 •
edited
Loading

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

x is typically pretty low-dimensional, right? I think ind2x should return an SVector. I think that should get rid of all the allocations and save a large amount of time.

In fact, I think the calls between ind2x and ind2x! should be switched, i.e. ind2x should contain the logic, and ind2x! should be implemented as

function ind2x!(grid, ind, x)
x[:] = ind2x(grid, ind)
end

ind2x!(grid, ind, x)
x::Array{Float64}
return x
end

function ind2x!(grid::AbstractGrid, ind::Int, x::AbstractArray)
Expand Down
2 changes: 1 addition & 1 deletion test/runtests.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -99,7 +99,7 @@ function compareToGrid(testType::Symbol=:random, numDims::Int=3, pointsPerDim::I

end

function getFractionalIndexes(g::AbstractGrid, s::Array)
function getFractionalIndexes(g::AbstractGrid, s::AbstractArray)
# Returns the fractional index of sprime within the grid-defined discretization.

fracInd = Array{Int64}(undef, length(g.cutPoints))
Expand Down