Ewald

class Ewald : public feasst::VisitModel

The Ewald summation accounts for the long-range nature of the electrostatic interaction by applying a Gaussian screening charge, computing a Fourier-space long-range component, and then correcting for the various spurious terms that are included in the Fourier summation, such as self and intra-particle.

See “Computer simulation of liquids” by M. P. Allen and D. J. Tildesley. The LAMMPS and DL_POLY manuals also include thorough descriptions of the Ewald summation.

Ewald is not supported for use as a reference state for dual-cut.

Following the description in the classic DL-POLY user manual (version 1.9), if the real space basis vectors are given by \(\vec{a}, \vec{b}, \vec{c}\), then reciprocal space basis vectors are

\(V = \vec{a} \cdot \vec{b} \times \vec{c}\)

\(\vec{u} = 2\pi\vec{b}\times\vec{c}/V\)

\(\vec{v} = 2\pi\vec{b}\times\vec{c}/V\)

\(\vec{w} = 2\pi\vec{b}\times\vec{c}/V\)

Public Functions

void tolerance_to_alpha_ks(const double tolerance, const Configuration &config, double *alpha, int *kxmax, int *kymax, int *kzmax)

Recommend an alpha parameter for Ewald as described and implemented in LAMMPS https://docs.lammps.org/kspace_style.html https://doi.org/10.1080/08927029208049126 https://doi.org/10.1063/1.470043

param tolerance:

Estimated tolerance for the energy. This is relative to the energy exerted by two unit charges separated by a distance of one unit.

param alpha:

Return the alpha

param kxmax:

Return the maximum number of Fourier vectors in the x dimension.

param kymax:

as above but in y

param kzmax:

as above but in z

void update_wave_vectors(const Configuration &config, const double kmax_squared, std::vector<double> *wave_prefactor, std::vector<int> *wave_num, double *ux, double *uy, double *uz, double *vy, double *vz, double *wz) const

Precompute the wave vectors within cutoff, coefficients, and also resize the structure factors.

void update_struct_fact_eik(const Select &selection, const Configuration &config, const std::vector<double> &wave_prefactor, const std::vector<int> &wave_num, const double ux, const double uy, const double uz, const double vy, const double vz, const double wz, std::vector<double> *struct_fact_real, std::vector<double> *struct_fact_imag, std::vector<std::vector<std::vector<double>>> *eik_new) const

Compute new eiks and update the given structure factor.

void precompute(Configuration *config)

Process tolerance arguments and initialize wave vectors.

void compute(ModelOneBody *model, const ModelParams &model_params, Configuration *config, const int group_index = 0)

Compute interactions of entire group in configuration from scratch. This is not optimized for smaller perturbations to the configuration.

void compute(ModelOneBody *model, const ModelParams &model_params, const Select &selection, Configuration *config, const int group_index)

Compute by selection is optimized for perturbations to the system. Contrary to short range interactions, in this case the energy of the entire system is computed from the structure factor. Thus, to get an energy contribution of a selection, one must take a difference between the configuration with and without the selection’s contributions to the structure factor.

For MonteCarlo trials, this means that the type of trial must be taken into consideration, which is provided via Select::trial_state().

For 0, “old”, first prepare to return the (previously computed) energy of the entire system, and then remove the structure factor contributions of the selection, without updating the selection’s eik. This assumes that there will be a follow up of exactly one state 1, “move” which will update the eik. Ewald is not currently supported with more than one step in TrialStage.

For 1, “move”, eik are updated and their contributions are added to the structure factor. Return the new system energy.

For 2, “remove”, remove the eik contributions to the structure factor without updating the eik. Return the old minus new energy.

For 3, “add”, eik are updated and their contributions are added to the structure factor. Return the new minus old energy.

For 4, “volume”, use the compute without selection.

void change_volume(const double delta_volume, const int dimension, Configuration *config)

Change the volume.

int wave_num(const int vector_index, const int dim) const

Return the wave vector number for a given dimension.

const std::vector<double> &wave_prefactor() const

Return the prefactor in the fourier-space term for each vector.

int num_kx() const

Return the number of vectors in the x dimension.

int num_ky() const

Return the number of vectors in the y dimension.

int num_kz() const

Return the number of vectors in the z dimension.

int kxmax() const

Return the cutoff of the wave vector in the x dimension.

int kymax() const

Return the cutoff of the wave vector in the y dimension.

int kzmax() const

Return the cutoff of the wave vector in the z dimension.

double kmax_squared() const

Return the spherical cutoff of the wave vectors.

const std::vector<std::vector<std::vector<double>>> &eik() const

Same as above, but for all particles, sites and wave vectors.

double struct_fact_real(const int vector_index) const

Return the real part of the structure factor for a given vector index corresponding with wave_prefactor and wave_num.

double struct_fact_imag(const int vector_index) const

Return the imaginary part of the structure factor for a given vector index corresponding with wave_prefactor and wave_num.

const std::vector<double> &struct_fact_imag() const

Return the imaginary part of the structure factor.

double net_charge(const Configuration &config) const

Return the net charge of the configuration.

Arguments

  • tolerance: determine the alpha parameter and number of wave vectors by specifying the accuracy relative to the energy of two unit charges separated by a distance of one unit. As described and implemented in LAMMPS, see: https://docs.lammps.org/kspace_style.html, https://doi.org/10.1080/08927029208049126, https://doi.org/10.1063/1.470043.

  • tolerance_num_sites: for setting parameters with the tolerance, optionally set the number of sites to be used for the parameter calculation rather than the currently existing number of sites (default).

  • alpha: optionally specify the alpha parameter in units of inverse length.

  • kxmax: optionally specify the maximum wave vectors in the first dimension.

  • kymax: same as above, but in the second dimension.

  • kzmax: same as above, but in the third dimension.

  • kmax_squared: optionally set the squared maximum integer wave vector for cubic domains only, which also sets kxmax, etc.