Cell list algorithm is slower than double for loop

Just a note, it seems

    NeiCell = Array{Tuple{Int64,Int64}}(undef,3,1)

is unused. You could also hoist the creation of

    NeiCell = Array{Tuple{Int64,Int64}}(undef,4,1)

to outside the function and modify it in place instead. Will save some allocations.

I would recommend you to read the style guide (Style Guide · The Julia Language), especially the naming section because the way the code is written right now looks very different from most Julia code and it adds some “friction” to read it.

A small rewrite (with some other small tweaks) could look something like:

this
# fixed-radius nearest neighbor
#cell list algorithm

const Intx2 = Tuple{Int, Int}

function FRNN_main(num_nod::Int, dim::Int, A::Matrix{Float64}, s::Float64)
    poi_pair = Intx2[]
    cell = poi_cell(A,num_nod, s)
    nei_cell = Vector{Intx2}(undef, 4)
    for key in keys(cell)
        value = cell[key]
        inner_cell_pair!(value, poi_pair, A, s)
        loop_neighbor!(key, nei_cell, cell, poi_pair, A, s)
    end
    return poi_pair
end

function poi_cell(A::Matrix{Float64}, num_nod::Int, s::Float64)
    cell = Dict{Intx2,Vector{Int}}()
    for i = 1:1:num_nod
        x = fld(A[1,i],s)
        y = fld(A[2,i],s)
        z = (x,y)
        if z ∉ keys(cell)
            cell[z] = [i]
        else
            push!(cell[z],i)
        end
    end
    return cell
end

function get_forward_cell!(nei_cell::Vector{Intx2}, key::Intx2)
    nei_cell[1] = (key[1]+1,key[2])
    nei_cell[2] = (key[1],  key[2]-1)
    nei_cell[3] = (key[1]+1,key[2]-1)
    nei_cell[4] = (key[1]+1,key[2]+1)
    return nei_cell
end

function inner_cell_pair!(value::Array{Int}, poi_pair::Vector{Intx2}, A::Matrix{Float64}, s::Float64)
    len = length(value)
    if len >= 2
        for i = 1:len-1
            for j = i+1:len
                add_or_not!(value[i], value[j], A, s, poi_pair)
            end
        end
    end
    return poi_pair
end

function loop_neighbor!(key::Intx2, nei_cell::Vector{Intx2}, cell::Dict{Intx2,Vector{Int}}, poi_pair::Vector{Intx2}, A::Matrix{Float64},s ::Float64)
    get_forward_cell!(nei_cell, key)
    for i = 1:4
        if nei_cell[i] ∈ keys(cell)
            find_pair!(cell[key], cell[nei_cell[i]], poi_pair, A, s)
        end
    end
    return poi_pair
end

function find_pair!(value::Array{Int}, value1::Array{Int}, poi_pair::Vector{Intx2}, A::Matrix{Float64}, s::Float64)
    for P1 in value
        for P2 in value1
            add_or_not!(P1, P2, A, s, poi_pair)
        end
    end
    return poi_pair
end

function add_or_not!(P1::Int, P2::Int, A::Matrix{Float64}, s::Float64, poi_pair::Vector{Intx2})
    dist = (A[1,P1]-A[1,P2])^2 + (A[2,P1]-A[2,P2])^2
    if dist <= s^2
        push!(poi_pair, (P1,P2))
    end
    return poi_pair
end

num_nod = 10000
dim = 2
s = 0.1*100
println("For Loop")

println("cell List")
for i = collect(1:10)*num_nod
    A = rand(dim,i)*100
    @time FRNN_main(num_nod,dim,A,s)
end

Just a question, what does it mean when NumNod is not the same as the number of points in A?

GitHub - KristofferC/NearestNeighbors.jl: High performance nearest neighbor data structures and algorithms for Julia. might be interesting to look at. You could do your query with

tree = KDTree(a)
inrange(tree, A, S)
1 Like