\(\renewcommand{\AA}{\text{Å}}\)
compute pace command
Syntax
compute ID group-ID pace ace_potential_filename ... keyword values ...
ID, group-ID are documented in compute command
pace = style name of this compute command
ace_potential_filename = file name (in the .yace or .ace format from pace pair_style) including ACE hyper-parameters, bonds, and generalized coupling coefficients
keyword = bikflag or dgradflag
bikflag value = 0 or 1 0 = descriptors are summed over atoms of each type 1 = descriptors are listed separately for each atom dgradflag value = 0 or 1 0 = descriptor gradients are summed over atoms of each type 1 = descriptor gradients are listed separately for each atom pair
Examples
compute pace all pace coupling_coefficients.yace
compute pace all pace coupling_coefficients.yace 0 1
compute pace all pace coupling_coefficients.yace 1 1
Description
Added in version 7Feb2024.
This compute calculates a set of quantities related to the atomic cluster expansion (ACE) descriptors of the atoms in a group. ACE descriptors are highly general atomic descriptors, encoding the radial and angular distribution of neighbor atoms, up to arbitrary bond order (rank). The detailed mathematical definition is given in the paper by (Drautz). These descriptors are used in the pace pair_style. Quantities obtained from compute pace are related to those used in pace pair_style to evaluate atomic energies, forces, and stresses for linear ACE models.
For example, the energy for a linear ACE model is calculated as: \(E=\sum_i^{N\_atoms} \sum_{\boldsymbol{\nu}} c_{\boldsymbol{\nu}} B_{i,\boldsymbol{\boldsymbol{\nu}}}\). The ACE descriptors for atom i \(B_{i,\boldsymbol{\nu}}\), and \(c_{\nu}\) are linear model parameters. The detailed definition and indexing convention for ACE descriptors is given in (Drautz). In short, body order \(N\), angular character, radial character, and chemical elements in the N-body descriptor are encoded by \(\nu\). In the pace pair_style, the linear model parameters and the ACE descriptors are combined for efficient evaluation of energies and forces. The details and benefits of this efficient implementation are given in (Lysogorskiy), but the combined descriptors and linear model parameters for the purposes of compute pace may be expressed in terms of the ACE descriptors mentioned above.
\(c_{\boldsymbol{\nu}} B_{i,\boldsymbol{\nu}}= \sum_{\boldsymbol{\nu}' \in \boldsymbol{\nu} } \big[ c_{\boldsymbol{\nu}} C(\boldsymbol{\nu}') \big] A_{i,\boldsymbol{\nu}'}\)
where the bracketed terms on the right-hand side are the combined functions with linear model parameters typically provided in the <name>.yace potential file for pace pair_style. When these bracketed terms are multiplied by the products of the atomic base from (Drautz), \(A_{i,\boldsymbol{\nu'}}\), the ACE descriptors are recovered but they are also scaled by linear model parameters. The generalized coupling coefficients, written in short-hand here as \(C(\boldsymbol{\nu}')\), are the generalized Clebsch-Gordan or generalized Wigner symbols. It may be desirable to reverse the combination of these descriptors and the linear model parameters so that the ACE descriptors themselves may be used. The ACE descriptors and their gradients are often used when training ACE models, performing custom data analysis, generalizing ACE model forms, and other tasks that involve direct computation of descriptors. The key utility of compute pace is that it can compute the ACE descriptors and gradients so that these tasks can be performed during a LAMMPS simulation or so that LAMMPS can be used as a driver for tasks like ACE model parameterization. To see how this command can be used within a Python workflow to train ACE potentials, see the examples in FitSNAP. Examples on using outputs from this compute to construct general ACE potential forms are demonstrated in (Goff). The various keywords and inputs to compute pace determine what ACE descriptors and related quantities are returned in a compute array.
The coefficient file, <name>.yace, ultimately defines the number of ACE descriptors to be computed, their maximum body-order, the degree of angular character they have, the degree of radial character they have, the chemical character (which element-element interactions are encoded by descriptors), and other hyper-parameters defined in (Drautz). These may be modeled after the potential files in pace pair_style, and have the same format. Details on how to generate the coefficient files to train ACE models may be found in FitSNAP.
The keyword bikflag determines whether or not to list the descriptors of each atom separately, or sum them together and list in a single row. If bikflag is set to 0 then a single descriptor row is used, which contains the per-atom ACE descriptors \(B_{i,\boldsymbol{\nu}}\) summed over all atoms i to produce \(B_{\boldsymbol{\nu}}\). If bikflag is set to 1 this is replaced by a separate per-atom ACE descriptor row for each atom. In this case, the entries in the final column for these rows are set to zero.
The keyword dgradflag determines whether to sum atom gradients or list them separately. If dgradflag is set to 0, the ACE descriptor gradients w.r.t. atom j are summed over all atoms i’ of, which may be useful when training linear ACE models on atomic forces. If dgradflag is set to 1, gradients are listed separately for each pair of atoms. Each row corresponds to a single term \(\frac{\partial {B_{i,\boldsymbol{\nu}}}}{\partial {r}^a_j}\) where \({r}^a_j\) is the a-th position coordinate of the atom with global index j. This also changes the number of columns to be equal to the number of ACE descriptors, with 3 additional columns representing the indices \(i\), \(j\), and \(a\), as explained more in the Output info section below. The option dgradflag=1 requires that bikflag=1.
Note
It is noted here that in contrast to pace pair_style, the .yace file for compute pace typically should not contain linear parameters for an ACE potential. If \(c_{\nu}\) are included, the value of the descriptor will not be returned in the compute array, but instead, the energy contribution from that descriptor will be returned. Do not do this unless it is the desired behavior. In short, you should not plug in a ‘.yace’ for a pace potential into this compute to evaluate descriptors.
Note
Generalized Clebsch-Gordan or Generalized Wigner symbols (with appropriate factors) must be used to evaluate ACE descriptors with this compute. There are multiple ways to define the generalized coupling coefficients. Because of this, this compute will not revert your potential file to a coupling coefficient file. Instead this compute allows the user to supply coupling coefficients that follow any convention.
Note
Using dgradflag = 1 produces a global array with \(N + 3N^2 + 1\) rows which becomes expensive for systems with more than 1000 atoms.
Note
If you have a bonded system, then the settings of special_bonds command can remove pairwise interactions between atoms in the same bond, angle, or dihedral. This is the default setting for the special_bonds command, and means those pairwise interactions do not appear in the neighbor list. Because this fix uses the neighbor list, it also means those pairs will not be included in the calculation. One way to get around this, is to write a dump file, and use the rerun command to compute the ACE descriptors for snapshots in the dump file. The rerun script can use a special_bonds command that includes all pairs in the neighbor list.
Output info
Compute pace evaluates a global array. The columns are arranged into ntypes blocks, listed in order of atom type I. Each block contains one column for each ACE descriptor, the same as for compute sna/atomin compute snap. A final column contains the corresponding energy, force component on an atom, or virial stress component. The rows of the array appear in the following order:
1 row: pace average descriptor values for all atoms of type I
3*n force rows: quantities, with derivatives w.r.t. x, y, and z coordinate of atom i appearing in consecutive rows. The atoms are sorted based on atom ID and run up to the total number of atoms, n.
6 rows: virial quantities summed for all atoms of type I
For example, if \(\# \; B_{i, \boldsymbol{\nu}}\) =30 and ntypes=1, the number of columns in the The number of columns in the global array generated by pace are 31, and 931, respectively, while the number of rows is 1+3*n+6, where n is the total number of atoms.
If the bik keyword is set to 1, the structure of the pace array is expanded. The first \(N\) rows of the pace array correspond to \(\# \; B_{i,\boldsymbol{\nu}}\) instead of a single row summed over atoms \(i\). In this case, the entries in the final column for these rows are set to zero. Also, each row contains only non-zero entries for the columns corresponding to the type of that atom. This is not true in the case of dgradflag keyword = 1 (see below).
If the dgradflag keyword is set to 1, this changes the structure of the global array completely. Here the per-atom quantities are replaced with rows corresponding to descriptor gradient components on single atoms:
where \({r}^a_j\) is the a-th position coordinate of the atom with global index j. The rows are organized in chunks, where each chunk corresponds to an atom with global index \(j\). The rows in an atom \(j\) chunk correspond to atoms with global index \(i\). The total number of rows for these descriptor gradients is therefore \(3N^2\). The number of columns is equal to the number of ACE descriptors, plus 3 additional left-most columns representing the global atom indices \(i\), \(j\), and Cartesian direction \(a\) (0, 1, 2, for x, y, z). The first 3 columns of the first \(N\) rows belong to the reference potential force components. The remaining K columns contain the \(B_{i,\boldsymbol{\nu}}\) per-atom descriptors corresponding to the non-zero entries obtained when bikflag = 1. The first column of the last row, after the first \(N + 3N^2\) rows, contains the reference potential energy. The virial components are not used with this option. The total number of rows is therefore \(N + 3N^2 + 1\) and the number of columns is \(K + 3\).
These values can be accessed by any command that uses global values from a compute as input. See the Howto output doc page for an overview of LAMMPS output options.
Restrictions
These computes are part of the ML-PACE package. They are only enabled if LAMMPS was built with that package. See the Build package page for more info.
Default
The optional keyword defaults are bikflag = 0, dgradflag = 0
(Drautz) Drautz, Phys Rev B, 99, 014104 (2019).
(Lysogorskiy) Lysogorskiy, van der Oord, Bochkarev, Menon, Rinaldi, Hammerschmidt, Mrovec, Thompson, Csanyi, Ortner, Drautz, npj Comp Mat, 7, 97 (2021).
(Goff) Goff, Zhang, Negre, Rohskopf, Niklasson, Journal of Chemical Theory and Computation 19, no. 13 (2023).