Ewald sum matrix

Ewald sum matrix [1] can be viewed as a logical extension of the Coulomb matrix for periodic systems, as it models the interaction between atoms in a periodic crystal through electrostatic interaction. However, in order to calculate the Coulomb interaction between the atomic cores in a periodic crystal, the Ewald summation technique has to be used together with a uniform neutralizing background charge.

Notice that the terms of the Ewald sum matrix in DScribe differ slightly from the original article. We correct an issue in the energy term related to the self-energy and the charged-cell correction. In the original article [1] (equation 20) this correction for the off-diagonal elements was defined as:

\[\phi^{\text{self}}_{ij} + \phi^{\text{bg}}_{ij} = -\frac{\alpha}{\sqrt{\pi}}(Z_i^2 + Z_j^2) -\frac{\pi}{2V\alpha^2}(Z_i + Z_j)^2~\forall~i \neq j\]

This term does however not correspond to the correct Ewald energy, as seen in the case of two interacting atoms, one of which has no charge: \(Z_i = 0\), \(Z_j \neq 0\). In this case the interaction should be zero, but the given term is non-zero and additionally depends on the screening parameter \(\alpha\). DScribe instead uses the corrected term:

\[\phi^{\text{self}}_{ij} + \phi^{\text{bg}}_{ij} = -\frac{\pi}{2 V \alpha^2} Z_i Z_j~\forall~i \neq j\]

Setup

Instantiating the object that is used to create Ewald sum matrices can be done as follows:

The constructor takes the following parameters:

EwaldSumMatrix.__init__(n_atoms_max, permutation='sorted_l2', sigma=None, seed=None, flatten=True, sparse=False)

Parameters

• n_atoms_max (int) – The maximum nuber of atoms that any of the samples can have. This controls how much zeros need to be padded to the final result.

• permutation (string) –

  Defines the method for handling permutational invariance. Can be one of the following:

  • none: The matrix is returned in the order defined by the Atoms.

  • sorted_l2: The rows and columns are sorted by the L2 norm.

  • eigenspectrum: Only the eigenvalues are returned sorted by their absolute value in descending order.

  • random: The rows and columns are sorted by their L2 norm after applying Gaussian noise to the norms. The standard deviation of the noise is determined by the sigma-parameter.

• sigma (float) – Provide only when using the random-permutation option. Standard deviation of the gaussian distributed noise determining how much the rows and columns of the randomly sorted matrix are scrambled.

• seed (int) – Provide only when using the random-permutation option. A seed to use for drawing samples from a normal distribution.

• flatten (bool) – Whether the output of create() should be flattened to a 1D array.

• sparse (bool) – Whether the output should be a sparse matrix or a dense numpy array.

Creation

After the Ewald sum matrix has been set up, it may be used on periodic atomic structures with the create()-method.

The call syntax for the create-function is as follows:

EwaldSumMatrix.create(system, accuracy=1e-05, w=1, rcut=None, gcut=None, a=None, n_jobs=1, verbose=False)

Return the Coulomb matrix for the given systems.

Parameters

• system (ase.Atoms or list of ase.Atoms) – One or many atomic structures.

• accuracy (float) – The accuracy to which the sum is converged to. Corresponds to the variable \(A\) in https://doi.org/10.1080/08927022.2013.840898. Used only if gcut, rcut and a have not been specified. Provide either one value or a list of values for each system.

• w (float) – Weight parameter that represents the relative computational expense of calculating a term in real and reciprocal space. This has little effect on the total energy, but may influence speed of computation in large systems. Note that this parameter is used only when the cutoffs and a are set to None. Provide either one value or a list of values for each system.

• rcut (float) – Real space cutoff radius dictating how many terms are used in the real space sum. Provide either one value or a list of values for each system.

• gcut (float) – Reciprocal space cutoff radius. Provide either one value or a list of values for each system.

• a (float) – The screening parameter that controls the width of the Gaussians. If not provided, a default value of \(\alpha = \sqrt{\pi}\left(\frac{N}{V^2}\right)^{1/6}\) is used. Corresponds to the standard deviation of the Gaussians. Provide either one value or a list of values for each system.

• n_jobs (int) – Number of parallel jobs to instantiate. Parallellizes the calculation across samples. Defaults to serial calculation with n_jobs=1.

• verbose (bool) – Controls whether to print the progress of each job into to the console.

Returns

Ewald sum matrix for the given systems. The return type depends on the ‘sparse’ and ‘flatten’-attributes. For flattened output a single numpy array or sparse scipy.csr_matrix is returned. The first dimension is determined by the amount of systems.

Return type

np.ndarray | scipy.sparse.csr_matrix

Note that if you specify in n_atoms_max a lower number than atoms in your structure it will cause an error. The output will in this case be a flattened matrix, specifically a numpy array with size #atoms * #atoms. The number of features may be requested beforehand with the get_number_of_features()-method.

In the case of multiple samples, the creation can also be parallellized by using the n_jobs-parameter. This splits the list of structures into equally sized parts and spaws a separate process to handle each part.

Examples

The following examples demonstrate usage of the descriptor. These examples are also available in dscribe/examples/ewaldsummatrix.py.

Accuracy

Easiest way to control the accuracy of the Ewald summation is to use the accuracy-parameter. Lower values of this parameter correspond to tighter convergence criteria and better accuracy.

Another option is to directly provide the real- and reciprocal space cutoffs:

Total electrostatic energy

Let’s calculate the electrostatic energy of a crystal by using the information contained in the Ewald sum matrix.

We can compare the result against the Ewald summation implementation in pymatgen:

1(1,2)

Felix Faber, Alexander Lindmaa, O. Anatole von Lilienfeld, and Rickard Armiento. Crystal structure representations for machine learning models of formation energies. International Journal of Quantum Chemistry, 115(16):1094–1101, 8 2015. doi:10.1002/qua.24917.