Electrostatics in Julia

Specific Domains Modelling & Simulations
21

Yes, that seems its natural home, let’s transfer it there.

I think making me a member of the org will let me make the transfer, could someone do that? Thanks.

1 Like
22

I’ve just given you “owner” privileges.

Note - the purpose of this organisation is two-fold. (1) just collect packages relevant to molecular simulation; (2) manage our own registry MolSim where we register our packages essentially without and constraints and quality-checks other than our own, so that we can have version control for our research projects. (This doesn’t prevent us from registering polished version with General as well.)

2 Likes
23

Great, I’ve made the switch.

4 Likes
24

I wrote an Ewald sums code here:
https://github.com/SimonEnsemble/PorousMaterials.jl/blob/master/src/electrostatics.jl
despite my best efforts to speed it up (resulting in very ugly code), it is too slow for production. also the default choices for the number of replications in k-space are not optimal. I tested it with some NIST data (see https://github.com/SimonEnsemble/PorousMaterials.jl/blob/master/test/electrostatics.jl). feel free to copy any of my code.
I would love to both use an EWalds.jl package and contribute to one.
we have a package for crystal structures here:
https://github.com/SimonEnsemble/Xtals.jl
and we simulate gas adsorption in MOFs (GCMC sims) with PorousMaterials.jl for our research:
https://github.com/SimonEnsemble/PorousMaterials.jl
(only handles LJ potentials and point charges).

Cory

5 Likes
25

A bit off topic, but regarding Xtals.jl: I’ve been thinking about something similar as a light-weight ASE to provide more convenience for starting calculations in DFTK. Is that the intention you have in mind with the package?

1 Like
26

we built the package over the last three years, out of pure selfishness-- whatever we need to get our research done in Julia!
that would be awesome if Xtals.jl turns into a more general purpose package for crystals, PRs welcome!
we have the capability to infer bonding graphs and work with those (another package on top of Xtals.jl coming soon). the only qualification is that there are speed considerations for how to represent atoms and charges for molecular simulations, in PorousMaterials.jl. so I wouldn’t want those to change, unless they were faster/as fast.
plz try out the package and let us know what you think!

2 Likes
27

Cool. Thanks. Will do :wink:

28

I’ve got some functionality for that in JuLIP and wouldn’t mind at all getting rid of it and using another oackage

29

Bumping this thread. Any updates?
I might implement a barebones particle mesh since it’s just couple lines of fourier transform. Probably similar to Michael’s DFTK.jl ewald.

30

The situation hasn’t changed from my end. If you or somebody else want to register a package and maintain it we’d be happy to have DFTK depend on it, and I would imagine a number of people are in the same situation.

31

Thanks I’ll take a look.

32

getting serious about this. we do molecular simulations of gases adsorbing in porous crystals such as MOFs with PorousMaterials.jl. need to compute electrostatic potential inside the MOF very quickly to enable simulation of molecular models with point charges on them.

to start a detailed discussion:

goal

a simple function to compute electrostatic potential at a point given in fractional coords xf, arising from charges in a periodic box.

phi = electrostatic_potential(xf::Array{Float64, 1}, charges::Charges, box::Box, ewald_settings::Settings)

is this what you are all looking for too?

but we need to call the above function hundreds of thousands of times, involving many loops over ~100 atoms—same charges different xf’s. for this function to be of reasonable speed for us, we need to pre-compute (and re-use) the k-vectors and pre-allocate memory for arrays that are used over and over in the for loop over the charges to compute e^{i n k ⋅ r}.
my code for this is here– very ugly. also has several tricks that take advantage of (a) symmetry to halve the computations and (b) recurrence relations to speed up e^{i n k ⋅ r}.

anyway, I propose the following for building blocks for this.

how to store the point charges

this is based on our Xtals.jl package here. can we use Xtals.jl or is there a more flexible way to represent charges for the community?

  • a Box to tell us about the periodicity of the system. see here.

  • a Charges structure to tell us where the point charges are in the box and the magnitude of the charges. we use this because we distinguish between fractional and cartesian coordinates in our code to avoid bugs (treat them as separate types). maybe not needed for ya’ll? if not, we can do more general:

struct Charges
    q::Array{Float64, 1} # values of the point charges
    xf::Array{Float64, 2} # fractional coordinates of charges in the box
end

and maybe put info about the box in there too.

in practice I found the above much faster than storing charges as charges::Array{Charge, 1} with

struct Charge
   q::Float64
   xf::Array{Float64, 1}
end

since the former allows for quick array operations through broadcasting (memory locality thing?).

for an Ewald.jl would we start over in defining these basic building blocks or use Xtals.jl?

Offsetarrays

I found OffsetArrays.jl to be very helpful for making the code match the math since it allows for k = -2:2 and to index arrays as x[k]. would this give a performance hit, though?

resources

Allen & Tildesley give a sample Ewald code that I found helpful.

these notes for understanding Ewald.

1 Like