Implement external potentials on a grid

[Luthaf]

In Monte-Carlo simulations, when working with a fixed framework a nice way to compute the energy is to pre-compute the potential (electrostatic or non-bonded) on a grid, and then interpolate the energy for each particle given its position.

For electrostatic potential, the idea to compute the electrostatic potential U on the grid, and then get the energy as q U. For Lennard-Jones potential, we need one grid for each atom type of the moving molecule, and we can precompute the energy by placing a particle on each point of the grid and using the standard energy computing method.

Then, we can save the pre-computed grids to the disk, and reload them in each simulation using this grid. This is very useful for adsorption in a fixed framework, where you can reuse the same pre-computed grid in multiple simulations, and you don't need to simulate the framework particles.

Implementation strategy

We can have two steps here: use a new kind of Simulation which compute the grid and output it to the disk; and add a new kind of GlobalPotential that read the cube file and can be used in the [potential] section of the input files.

I think the Gaussian cube file format will be nice to use here, as it represent volumetric data, and can be visualized with standard tools. The only downside is that we need more that one file for a simulation.

[g-bauer]

Interesting. For adsorption, you would still need the grid positions to compute the interactions with the adsorbent, right?

A simple alternative for adsorption (if we'd have particle insertion) would be to add the grid/solid structure, but don't add any moves and then add specific moves for the adsorbed species. Would the "grid method" be more useful in that case?

[Luthaf]

For adsorption, you would still need the grid positions to compute the interactions with the adsorbent, right?

I am not sure what you mean here. You don't need the positions of the atoms in the framework, because they are already accounted for in the pre-computed grid.

The point here is that once the grid is computed, computing the interaction of a particle with the grid (so with the rigid structure) is a O(1) operation, while computing the interaction with all the atoms in the structure is still an O(n) operation.

[g-bauer]

Ah, I see. You would take the adsorbing molecule and move it through the whole grid beforehand. Wouldn't it be necessary to also sample multiple intra molecular configurations as well as orientations? I guess the initial computation could be very time consuming then.

[Luthaf]

Yes, the initial step is pretty long, that why you run it as a separated simulation, and you save the result to the disk.

You don't need to sample orientation, because this work at the particle level, not at the molecule level. Which means you need one grid for each particle type.

[g-bauer]

Note: There is an explanation of the method and references to other papers for this algorithm in this paper.