\(\renewcommand{\AA}{\text{Å}}\)
fix qeq/point command
fix qeq/shielded command
fix qeq/slater command
fix qeq/ctip command
fix qeq/dynamic command
fix qeq/fire command
Syntax
fix ID group-ID style Nevery cutoff tolerance maxiter qfile keyword ...
ID, group-ID are documented in fix command
style = qeq/point or qeq/shielded or qeq/slater or qeq/ctip or qeq/dynamic or qeq/fire
Nevery = perform charge equilibration every this many steps
cutoff = global cutoff for charge-charge interactions (distance unit)
tolerance = precision to which charges will be equilibrated
maxiter = maximum iterations to perform charge equilibration
qfile = a filename with QEq parameters or coul/streitz or coul/ctip or reaxff
zero or more keyword/value pairs may be appended
keyword = alpha or cdamp or maxrepeat or qdamp or qstep or warn
alpha value = Slater type orbital exponent (qeq/slater only) cdamp value = damping parameter for Coulomb interactions (qeq/ctip only) maxrepeat value = number of equilibration cycles allowed to ensure no atoms cross charge bounds (qeq/ctip only) qdamp value = damping factor for damped dynamics charge solver (qeq/dynamic and qeq/fire only) qstep value = time step size for damped dynamics charge solver (qeq/dynamic and qeq/fire only) warn value = do (=yes) or do not (=no) print a warning when the maximum number of iterations is reached
Examples
fix 1 all qeq/point 1 10 1.0e-6 200 param.qeq1
fix 1 qeq qeq/shielded 1 8 1.0e-6 100 param.qeq2
fix 1 all qeq/slater 5 10 1.0e-6 100 params alpha 0.2
fix 1 all qeq/ctip 1 12 1.0e-8 100 coul/ctip cdamp 0.30 maxrepeat 10
fix 1 qeq qeq/dynamic 1 12 1.0e-3 100 my_qeq
fix 1 all qeq/fire 1 10 1.0e-3 100 my_qeq qdamp 0.2 qstep 0.1
Description
Perform the charge equilibration (QEq) method as described in (Rappe and Goddard) and formulated in (Nakano) (also known as the matrix inversion method) and in (Rick and Stuart) (also known as the extended Lagrangian method) based on the electronegativity equilization principle.
These fixes can be used with any pair style in LAMMPS, so long as per-atom charges are defined. The most typical use-case is in conjunction with a pair style that performs charge equilibration periodically (e.g. every timestep), such as the ReaxFF or Streitz-Mintmire potential. But these fixes can also be used with potentials that normally assume per-atom charges are fixed, e.g. a Buckingham or LJ/Coulombic potential.
Because the charge equilibration calculation is effectively independent of the pair style, these fixes can also be used to perform a one-time assignment of charges to atoms. For example, you could define the QEq fix, perform a zero-timestep run via the run command without any pair style defined which would set per-atom charges (based on the current atom configuration), then remove the fix via the unfix command before performing further dynamics.
Note
Computing and using charge values different from published values defined for a fixed-charge potential like Buckingham or CHARMM or AMBER, can have a strong effect on energies and forces, and produces a different model than the published versions.
Note
The fix qeq/comb command must still be used to perform charge equilibration with the COMB potential. The fix qeq/reaxff command can be used to perform charge equilibration with the ReaxFF force field, although fix qeq/shielded yields the same results as fix qeq/reaxff if Nevery, cutoff, and tolerance are the same. Eventually the fix qeq/reaxff command will be deprecated.
The QEq method minimizes the electrostatic energy of the system (or equalizes the derivative of energy with respect to charge of all the atoms) by adjusting the partial charge on individual atoms based on interactions with their neighbors within cutoff. It requires a few parameters in the appropriate units for each atom type which are read from a file specified by qfile. The file has the following format:
1 chi eta gamma zeta qcore
2 chi eta gamma zeta qcore
...
Ntype chi eta gamma zeta qcore
except for fix style qeq/ctip where the format is:
1 chi eta gamma zeta qcore qmin qmax omega
2 chi eta gamma zeta qcore qmin qmax omega
...
Ntype chi eta gamma zeta qcore qmin qmax omega
There have to be parameters given for every atom type. Wildcard entries are possible using the same type range syntax as for “coeff” commands (i.e., n*m, n*, *m, *). Later entries will overwrite previous ones. Empty lines or any text following the pound sign (#) are ignored. Each line starts with the atom type followed by eight parameters. Only a subset of the parameters is used by each QEq style as described below, thus the others can be set to 0.0 if desired, but all eight entries per line are required.
chi = electronegativity in energy units
eta = self-Coulomb potential in energy units
gamma = shielded Coulomb constant defined by ReaxFF force field in distance units
zeta = Slater type orbital exponent defined by the Streitz-Mintmire potential in reverse distance units
qcore = charge of the nucleus defined by the Streitz-Mintmire potential potential in charge units
qmin = lower bound on the allowed charge defined by the CTIP potential in charge units
qmax = upper bound on the allowed charge defined by the CTIP potential in charge units
omega = penalty parameter used to enforce charge bounds defined by the CTIP potential in energy units
The fix qeq styles will print a warning if the charges are not equilibrated within tolerance by maxiter steps, unless the warn keyword is used with “no” as argument. This latter option may be useful for testing and benchmarking purposes, as it allows to use a fixed number of QEq iterations when tolerance is set to a small enough value to always reach the maxiter limit. Turning off warnings will avoid the excessive output in that case.
The qeq/point style describes partial charges on atoms as point charges. Interaction between a pair of charged particles is 1/r, which is the simplest description of the interaction between charges. Only the chi and eta parameters from the qfile file are used. Note that Coulomb catastrophe can occur if repulsion between the pair of charged particles is too weak. This style solves partial charges on atoms via the matrix inversion method. A tolerance of 1.0e-6 is usually a good number.
The qeq/shielded style describes partial charges on atoms also as point charges, but uses a shielded Coulomb potential to describe the interaction between a pair of charged particles. Interaction through the shielded Coulomb is given by equation (13) of the ReaxFF force field paper. The shielding accounts for charge overlap between charged particles at small separation. This style is the same as fix qeq/reaxff, and can be used with pair_style reaxff. Only the chi, eta, and gamma parameters from the qfile file are used. When using the string reaxff as filename, these parameters are extracted directly from an active reaxff pair style. This style solves partial charges on atoms via the matrix inversion method. A tolerance of 1.0e-6 is usually a good number.
The qeq/slater style describes partial charges on atoms as spherical charge densities centered around atoms via the Slater 1s orbital, so that the interaction between a pair of charged particles is the product of two Slater 1s orbitals. The expression for the Slater 1s orbital is given under equation (6) of the Streitz-Mintmire paper. Only the chi, eta, zeta, and qcore parameters from the qfile file are used. When using the string coul/streitz as filename, these parameters are extracted directly from an active coul/streitz pair style. This style solves partial charges on atoms via the matrix inversion method. A tolerance of 1.0e-6 is usually a good number. Keyword alpha can be used to change the Slater type orbital exponent.
Added in version 19Nov2024.
The qeq/ctip style describes partial charges on atoms in the same way as style qeq/shielded but also enables the definition of charge bounds. Only the chi, eta, gamma, qmin, qmax, and omega parameters from the qfile file are used. When using the string coul/ctip as filename, these parameters are extracted directly from an active coul/ctip pair style. This style solves partial charges on atoms via the matrix inversion method. Keyword cdamp can be used to change the damping parameter used to calculate Coulomb interactions. Keyword maxrepeat can be used to adjust the number of equilibration cycles allowed to ensure no atoms have crossed the charge bounds. A value of 10 is usually a good choice. A tolerance between 1.0e-6 and 1.0e-8 is usually a good choice but should be checked in conjunction with the timestep for adequate energy conservation during dynamic runs.
The qeq/dynamic style describes partial charges on atoms as point charges that interact through 1/r, but the extended Lagrangian method is used to solve partial charges on atoms. Only the chi and eta parameters from the qfile file are used. Note that Coulomb catastrophe can occur if repulsion between the pair of charged particles is too weak. A tolerance of 1.0e-3 is usually a good number. Keyword qdamp can be used to change the damping factor, while keyword qstep can be used to change the time step size.
The *qeq/fire* style describes the same charge model and charge solver as the qeq/dynamic style, but employs a FIRE minimization algorithm to solve for equilibrium charges. Keyword qdamp can be used to change the damping factor, while keyword qstep can be used to change the time step size.
Note that qeq/point, qeq/shielded, qeq/slater, and qeq/ctip describe different charge models, whereas the matrix inversion method and the extended Lagrangian method (qeq/dynamic and qeq/fire) are different solvers.
Note that qeq/point, qeq/dynamic and qeq/fire styles all describe charges as point charges that interact through 1/r relationship, but solve partial charges on atoms using different solvers. These three styles should yield comparable results if the QEq parameters and Nevery, cutoff, and tolerance are the same. Style qeq/point is typically faster, qeq/dynamic scales better on larger sizes, and qeq/fire is faster than qeq/dynamic.
Note
In order to solve the self-consistent equations for electronegativity equalization, LAMMPS imposes the additional constraint that all the charges in the fix group must add up to zero. The initial charge assignments should also satisfy this constraint. LAMMPS will print a warning if that is not the case.
Note
Developing QEq parameters (chi, eta, gamma, zeta, and qcore) is non-trivial. Charges on atoms are not guaranteed to equilibrate with arbitrary choices of these parameters. We do not develop these QEq parameters. See the examples/qeq directory for some examples.
Restart, fix_modify, output, run start/stop, minimize info
No information about these fixes is written to binary restart files. No global scalar or vector or per-atom quantities are stored by these fixes for access by various output commands. No parameter of these fixes can be used with the start/stop keywords of the run command.
Thexe fixes are invoked during energy minimization.
Restrictions
These fixes are part of the QEQ package. They are only enabled if LAMMPS was built with that package. See the Build package page for more info.
These qeq fixes will ignore electric field contributions from fix efield.
Default
warn yes
(Rappe and Goddard) A. K. Rappe and W. A. Goddard III, J Physical Chemistry, 95, 3358-3363 (1991).
(Nakano) A. Nakano, Computer Physics Communications, 104, 59-69 (1997).
(Rick and Stuart) S. W. Rick, S. J. Stuart, B. J. Berne, J Chemical Physics 101, 16141 (1994).
(Streitz-Mintmire) F. H. Streitz, J. W. Mintmire, Physical Review B, 50, 16, 11996 (1994)
(CTIP) G. Plummer, J. P. Tavenner, M. I. Mendelev, Z. Wu, J. W. Lawson, in preparation
(ReaxFF) A. C. T. van Duin, S. Dasgupta, F. Lorant, W. A. Goddard III, J Physical Chemistry, 105, 9396-9049 (2001)
(QEq/Fire) T.-R. Shan, A. P. Thompson, S. J. Plimpton, in preparation