JuliaMolSim / Molly.jl

Molecular simulation in Julia
Other
376 stars 51 forks source link

Force calculation suggestion #30

Open BradenDKelly opened 3 years ago

BradenDKelly commented 3 years ago

I was just looking through the Lennard-Jones force calculations and noticed that every time you calculate the force between two atoms, you are calculating the mixing rules

σ = sqrt(atom_i.σ * atom_j.σ) 
ϵ = sqrt(atom_i.ϵ * atom_j.ϵ)

Given the simplicity of the force calculation, these two square roots likely are adding a considerable amount of flops to the simulation. Typically simulation codes will precalculate the mixing rules and then simply look up the necessary sigma and epsilon, i.e.,

σ = vdw_table.σ(atom_i.type,  atom_j.type) 
ϵ = vdw_table.ϵ(atom_i.type,  atom_j.type)

I do not know if type is available, but, incorporating it into the atom struct so that you can use this parameter lookup table, will likely speed up the calculation a reasonable amount. The type would be an integer, so it should be isbits true and so shouldn't interfere with CUDA/GPU. This will also make the structures lighter since only type needs to be saved per atom, and not epsilon and sigma (also, other force-fields, such as EXP-6, have more parameters so you save even more with them).

Also, usually the size parameters use the Lorenz mixing rule i.e., are

σ = (atom_i.σ + atom_j.σ) / 2

It would be good to add the feature allowing the user to select the mixing rules to use for each parameter type.

jgreener64 commented 3 years ago

Thanks for the advice, I'll add it to the to-do list.

BradenDKelly commented 3 years ago

It is probably worth noting that the degrees of freedom are only 3*N, where N is number of particles, for the Andersen thermostat. DF are different for stochastic or deterministic thermostats, and this rippled to kinetic energy. The difference is small, but nonetheless, worth accounting for.

jgreener64 commented 3 years ago

Thanks, let me know as you see anything else, it's useful.

ehsanirani commented 2 years ago

I think @BradenDKelly suggestion would enhance the performance and simplicity of the api. The suggested table can also have a r_cutoff field for each pair of particle types. Then something like negative values can be set to it to exclude that pair of particle types from force calculation (or neighbor searching). That's the method that HOOMD-Blue uses.

jgreener64 commented 2 years ago

I think it's a good idea but will wait to see how the GPU kernel work (https://github.com/JuliaMolSim/Molly.jl/issues/60) progresses as this may require a refactor of other parts of the API too.