m3g / CellListMap.jl

Flexible implementation of cell lists to map the calculations of particle-pair dependent functions, such as forces, energies, neighbor lists, etc.
https://m3g.github.io/CellListMap.jl/
MIT License
87 stars 4 forks source link

Does the 'update_lists = false' option only apply to periodic systems? #80

Closed maphdze closed 1 year ago

maphdze commented 1 year ago

In the chapter of periodic systems of the document it is said that update cell list could be disabled. I tried it on the non-periodic systems, and there is no 'update_lists' in the variables for map_pairwise!(). Attached is my code

using CellListMap
using StaticArrays
using BenchmarkTools
@inline function my_norm(v::AbstractVector{T}) where {T}
    n = zero(T)
    @simd for x in v
        n += x ^ 2
    end
    n = sqrt(n)
    return n
end
function cal_damage!(x,y,i,j,d2,u_points,points_vol,γ,λᵦ,damage)
    d_new = my_norm(x + u_points[i] - y -u_points[j])
    λ = d_new - sqrt(d2)
    pair_dam = 0.0
    if λ > λᵦ
        pair_dam = 1.0 #- exp(-γ * (λ - λᵦ))
    end
    damage[i] += pair_dam * points_vol[j]
    damage[j] += pair_dam * points_vol[i]
    return damage
end
function map_damage!(damage_map,box,cl,u_points,points_vol,γ,λᵦ)
    damage_map .= 0.0
    map_pairwise!((x,y,i,j,d2,damage_map) -> cal_damage!(x,y,i,j,d2,u_points,points_vol,γ,λᵦ,damage_map),damage_map,box,cl; update_lists=false)
end
N = 300000
ℓ = 0.05
γ = 1.0e4
λᵦ = 1.0e-6
points = rand(SVector{2,Float64},N)
box = Box(limits(points),ℓ)
cl = CellList(points,box)
u_points = 1e-4 * rand(SVector{2,Float64},N)
points_vol = rand(N)
damage_map = zeros(Float64,N)
map_pairwise!((x,y,i,j,d2,damage_map) -> cal_damage!(x,y,i,j,d2,u_points,points_vol,γ,λᵦ,damage_map),damage_map,box,cl)
@btime map_damage!(damage_map,box,cl,u_points,points_vol,γ,λᵦ)

And there is the error information:

ERROR: MethodError: no method matching map_pairwise!(::var"#5#6"{Vector{SVector{2, Float64}}, Vector{Float64}, Float64, Float64}, ::Vector{Float64}, ::Box{NonPeriodicCell, 2, Float64, Float64, 4, Float64}, ::CellList{2, Float64}; update_lists=false)
Closest candidates are:
  map_pairwise!(::F, ::Any, ::Box, ::CellList; parallel, output_threaded, reduce, show_progress) where F at C:\Users\rydmaize\.julia\packages\CellListMap\BemzV\src\CellListMap.jl:90 got unsupported keyword argument "update_lists"
  map_pairwise!(::F1, ::Any, ::Box, ::CellListMap.CellListPair{V, N, T, Swap}; parallel, output_threaded, reduce, show_progress) where {F1, V, N, T, Swap} at C:\Users\rydmaize\.julia\packages\CellListMap\BemzV\src\CellListMap.jl:116 got unsupported keyword argument "update_lists"
Stacktrace:
  [1] kwerr(::NamedTuple{(:update_lists,), Tuple{Bool}}, ::Function, ::Function, ::Vector{Float64}, ::Box{NonPeriodicCell, 2, Float64, Float64, 4, Float64}, ::CellList{2, Float64})
    @ Base .\error.jl:165
  [2] map_damage!(damage_map::Vector{Float64}, box::Box{NonPeriodicCell, 2, Float64, Float64, 4, Float64}, cl::CellList{2, Float64}, u_points::Vector{SVector{2, Float64}}, points_vol::Vector{Float64}, γ::Float64, λᵦ::Float64)
    @ Main d:\Code\CellListMap\Update_lists.jl:26
  [3] var"##core#321"()
    @ Main C:\Users\rydmaize\.julia\packages\BenchmarkTools\0owsb\src\execution.jl:489
  [4] var"##sample#322"(::Tuple{}, __params::BenchmarkTools.Parameters)
    @ Main C:\Users\rydmaize\.julia\packages\BenchmarkTools\0owsb\src\execution.jl:495
  [5] _run(b::BenchmarkTools.Benchmark, p::BenchmarkTools.Parameters; verbose::Bool, pad::String, kwargs::Base.Pairs{Symbol, Integer, NTuple{4, Symbol}, NamedTuple{(:samples, :evals, :gctrial, :gcsample), Tuple{Int64, Int64, Bool, Bool}}})
    @ BenchmarkTools C:\Users\rydmaize\.julia\packages\BenchmarkTools\0owsb\src\execution.jl:99
  [6] #invokelatest#2
    @ .\essentials.jl:731 [inlined]
  [7] #run_result#45
    @ C:\Users\rydmaize\.julia\packages\BenchmarkTools\0owsb\src\execution.jl:34 [inlined]
  [8] run(b::BenchmarkTools.Benchmark, p::BenchmarkTools.Parameters; progressid::Nothing, nleaves::Float64, ndone::Float64, kwargs::Base.Pairs{Symbol, Integer, NTuple{5, Symbol}, NamedTuple{(:verbose, :samples, :evals, :gctrial, :gcsample), Tuple{Bool, Int64, Int64, Bool, Bool}}})
    @ BenchmarkTools C:\Users\rydmaize\.julia\packages\BenchmarkTools\0owsb\src\execution.jl:117
  [9] #warmup#54
    @ C:\Users\rydmaize\.julia\packages\BenchmarkTools\0owsb\src\execution.jl:169 [inlined]
 [10] warmup(item::BenchmarkTools.Benchmark)
    @ BenchmarkTools C:\Users\rydmaize\.julia\packages\BenchmarkTools\0owsb\src\execution.jl:168
 [11] top-level scope
    @ C:\Users\rydmaize\.julia\packages\BenchmarkTools\0owsb\src\execution.jl:575

Is it possible to make this option also possible for non-periodic systems? Thanks for your help!

lmiq commented 1 year ago

In that case the list is built when you call CellList or if you explicitly call the UpdateCellList! function before the mapping.

Thus, there is no cell list updating in that case, in the call to map_pairwise.

lmiq commented 1 year ago

Take a look at these examples:

https://m3g.github.io/CellListMap.jl/stable/LowLevel/#Preallocating-the-cell-lists-and-cell-list-auxiliary-arrays

maphdze commented 1 year ago

Thanks for your reply! I got it.

maphdze commented 1 year ago

Take a look at these examples:

https://m3g.github.io/CellListMap.jl/stable/LowLevel/#Preallocating-the-cell-lists-and-cell-list-auxiliary-arrays

I will try. I guess this will improve the performance.

lmiq commented 1 year ago

If you are running celllistmap from within an interative scheme (like the steps of a simulation), that is the pattern you should follow. Anyway, you need to update the cell lists (with UpdateCellList!), whenever the coordinates of the particles changes.

The PeriodicSystems interface takes care of all that for you, in fact, and does all that in background. With update_lists=false it will not update the lists, but that would be wrong if the positions change, and that option should only be used if one wants to compute a different property from the same set of positions, which is already mapped into the cell list.