ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merge pull request #4312 from navlalli/qtpie

Browse Source 
Add fix qtpie/reaxff
This commit is contained in:
Axel Kohlmeyer
2024-11-13 07:42:07 -05:00
committed by GitHub
parent 69526c87d8 2d2bc029a4
commit 5928389c5e
20 changed files with 2142 additions and 19 deletions
Download Patch File Download Diff File Expand all files Collapse all files

1
doc/src/Commands_fix.rst
View File

@ -187,6 +187,7 @@ OPT.
* :doc:`qeq/slater <fix_qeq>` * :doc:`qeq/slater <fix_qeq>`
* :doc:`qmmm <fix_qmmm>` * :doc:`qmmm <fix_qmmm>`
* :doc:`qtb <fix_qtb>` * :doc:`qtb <fix_qtb>`
* :doc:`qtpie/reaxff <fix_qtpie_reaxff>`
* :doc:`rattle <fix_shake>` * :doc:`rattle <fix_shake>`
* :doc:`reaxff/bonds (k) <fix_reaxff_bonds>` * :doc:`reaxff/bonds (k) <fix_reaxff_bonds>`
* :doc:`reaxff/species (k) <fix_reaxff_species>` * :doc:`reaxff/species (k) <fix_reaxff_species>`

1
doc/src/fix.rst
View File

@ -366,6 +366,7 @@ accelerated styles exist.
* :doc:`qeq/slater <fix_qeq>` - charge equilibration via Slater method * :doc:`qeq/slater <fix_qeq>` - charge equilibration via Slater method
* :doc:`qmmm <fix_qmmm>` - functionality to enable a quantum mechanics/molecular mechanics coupling * :doc:`qmmm <fix_qmmm>` - functionality to enable a quantum mechanics/molecular mechanics coupling
* :doc:`qtb <fix_qtb>` - implement quantum thermal bath scheme * :doc:`qtb <fix_qtb>` - implement quantum thermal bath scheme
* :doc:`qtpie/reaxff <fix_qtpie_reaxff>` - apply QTPIE charge equilibration
* :doc:`rattle <fix_shake>` - RATTLE constraints on bonds and/or angles * :doc:`rattle <fix_shake>` - RATTLE constraints on bonds and/or angles
* :doc:`reaxff/bonds <fix_reaxff_bonds>` - write out ReaxFF bond information * :doc:`reaxff/bonds <fix_reaxff_bonds>` - write out ReaxFF bond information
* :doc:`reaxff/species <fix_reaxff_species>` - write out ReaxFF molecule information * :doc:`reaxff/species <fix_reaxff_species>` - write out ReaxFF molecule information

6
doc/src/fix_acks2_reaxff.rst
View File

@ -111,7 +111,8 @@ LAMMPS was built with that package. See the :doc:`Build package
This fix does not correctly handle interactions involving multiple This fix does not correctly handle interactions involving multiple
periodic images of the same atom. Hence, it should not be used for periodic images of the same atom. Hence, it should not be used for
periodic cell dimensions less than :math:`10~\AA`. periodic cell dimensions smaller than the non-bonded cutoff radius,
which is typically :math:`10~\AA` for ReaxFF simulations.
This fix may be used in combination with :doc:`fix efield <fix_efield>` This fix may be used in combination with :doc:`fix efield <fix_efield>`
and will apply the external electric field during charge equilibration, and will apply the external electric field during charge equilibration,
@ -122,7 +123,8 @@ components in non-periodic directions.
Related commands Related commands
"""""""""""""""" """"""""""""""""
:doc:`pair_style reaxff <pair_reaxff>`, :doc:`fix qeq/reaxff <fix_qeq_reaxff>` :doc:`pair_style reaxff <pair_reaxff>`, :doc:`fix qeq/reaxff <fix_qeq_reaxff>`,
:doc:`fix qtpi/reaxff <fix_qtpie_reaxff>`
Default Default
""""""" """""""

6
doc/src/fix_qeq_reaxff.rst
View File

@ -124,7 +124,8 @@ LAMMPS was built with that package. See the :doc:`Build package
This fix does not correctly handle interactions involving multiple This fix does not correctly handle interactions involving multiple
periodic images of the same atom. Hence, it should not be used for periodic images of the same atom. Hence, it should not be used for
periodic cell dimensions less than 10 Angstroms. periodic cell dimensions smaller than the non-bonded cutoff radius,
which is typically :math:`10~\AA` for ReaxFF simulations.
This fix may be used in combination with :doc:`fix efield <fix_efield>` This fix may be used in combination with :doc:`fix efield <fix_efield>`
and will apply the external electric field during charge equilibration, and will apply the external electric field during charge equilibration,
@ -138,7 +139,8 @@ as an atom-style variable using the *potential* keyword for `fix efield`.
Related commands Related commands
"""""""""""""""" """"""""""""""""
:doc:`pair_style reaxff <pair_reaxff>`, :doc:`fix qeq/shielded <fix_qeq>` :doc:`pair_style reaxff <pair_reaxff>`, :doc:`fix qeq/shielded <fix_qeq>`,
:doc:`fix acks2/reaxff <fix_acks2_reaxff>`, :doc:`fix qtpie/reaxff <fix_qtpie_reaxff>`
Default Default
""""""" """""""

198
doc/src/fix_qtpie_reaxff.rst Normal file
View File

@ -0,0 +1,198 @@
.. index:: fix qtpie/reaxff
fix qtpie/reaxff command
========================
Syntax
""""""
.. code-block:: LAMMPS
fix ID group-ID qtpie/reaxff Nevery cutlo cuthi tolerance params gfile args
* ID, group-ID are documented in :doc:`fix <fix>` command
* qtpie/reaxff = style name of this fix command
* Nevery = perform QTPIE every this many steps
* cutlo,cuthi = lo and hi cutoff for Taper radius
* tolerance = precision to which charges will be equilibrated
* params = reaxff or a filename
* gfile = the name of a file containing Gaussian orbital exponents
* one or more keywords or keyword/value pairs may be appended
.. parsed-literal::
keyword = *maxiter*
*maxiter* N = limit the number of iterations to *N*
Examples
""""""""
.. code-block:: LAMMPS
fix 1 all qtpie/reaxff 1 0.0 10.0 1.0e-6 reaxff exp.qtpie
fix 1 all qtpie/reaxff 1 0.0 10.0 1.0e-6 params.qtpie exp.qtpie maxiter 500
Description
"""""""""""
The QTPIE charge equilibration method is an extension of the QEq charge
equilibration method. With QTPIE, the partial charges on individual atoms
are computed by minimizing the electrostatic energy of the system in the
same way as the QEq method but where the absolute electronegativity,
:math:`\chi_i`, of each atom in the QEq charge equilibration scheme
:ref:`(Rappe and Goddard) <Rappe3>` is replaced with an effective
electronegativity given by :ref:`(Chen) <qtpie-Chen>`
.. math::
\chi_{\mathrm{eff},i} = \frac{\sum_{j=1}^{N} (\chi_i - \chi_j) S_{ij}}
{\sum_{m=1}^{N}S_{im}},
which acts to penalize long-range charge transfer seen with the QEq charge
equilibration scheme. In this equation, :math:`N` is the number of atoms in
the system and :math:`S_{ij}` is the overlap integral between atom :math:`i`
and atom :math:`j`.
The effect of an external electric field can be incorporated into the QTPIE
method by modifying the absolute or effective electronegativities of each
atom :ref:`(Chen) <qtpie-Chen>`. This fix models the effect of an external
electric field by using the effective electronegativity given in
:ref:`(Gergs) <Gergs>`:
.. math::
\chi_{\mathrm{eff},i} = \frac{\sum_{j=1}^{N} (\chi_i - \chi_j + \phi_i - \phi_j) S_{ij}}
{\sum_{m=1}^{N}S_{im}},
where :math:`\phi_i` and :math:`\phi_j` are the electric
potentials at the positions of atom :math:`i` and :math:`j`
due to the external electric field.
This fix is typically used in conjunction with the ReaxFF force
field model as implemented in the :doc:`pair_style reaxff <pair_reaxff>`
command, but it can be used with any potential in LAMMPS, so long as it
defines and uses charges on each atom. For more technical details about the
charge equilibration performed by `fix qtpie/reaxff`, which is the same as in
:doc:`fix qeq/reaxff <fix_qeq_reaxff>` except for the use of
:math:`\chi_{\mathrm{eff},i}`, please refer to :ref:`(Aktulga) <qeq-Aktulga2>`.
To be explicit, this fix replaces :math:`\chi_k` of eq. 3 in
:ref:`(Aktulga) <qeq-Aktulga2>` with :math:`\chi_{\mathrm{eff},k}`.
This fix requires the absolute electronegativity, :math:`\chi`, in eV, the
self-Coulomb potential, :math:`\eta`, in eV, and the shielded Coulomb
constant, :math:`\gamma`, in :math:`\AA^{-1}`. If the *params* setting above
is the word "reaxff", then these are extracted from the
:doc:`pair_style reaxff <pair_reaxff>` command and the ReaxFF force field
file it reads in. If a file name is specified for *params*, then the
parameters are taken from the specified file and the file must contain
one line for each atom type. The latter form must be used when performing
QTPIE with a non-ReaxFF potential. Each line should be formatted as follows,
ensuring that the parameters are given in units of eV, eV, and :math:`\AA^{-1}`,
respectively:
.. parsed-literal::
itype chi eta gamma
where *itype* is the atom type from 1 to Ntypes. Note that eta is
defined here as twice the eta value in the ReaxFF file.
The overlap integrals in the equation for :math:`\chi_{\mathrm{eff},i}`
are computed by using normalized 1s Gaussian type orbitals. The Gaussian
orbital exponents, :math:`\alpha`, that are needed to compute the overlap
integrals are taken from the file given by *gfile*.
This file must contain one line for each atom type and provide the Gaussian
orbital exponent for each atom type in units of inverse square Bohr radius.
Each line should be formatted as follows:
.. parsed-literal::
itype alpha
Empty lines or any text following the pound sign (#) are ignored. An example
*gfile* for a system with two atom types is
.. parsed-literal::
# An example gfile. Exponents are taken from Table 2.2 of Chen, J. (2009).
# Theory and applications of fluctuating-charge models.
# The units of the exponents are 1 / (Bohr radius)^2 .
1 0.2240 # O
2 0.5434 # H
The optional *maxiter* keyword allows changing the max number
of iterations in the linear solver. The default value is 200.
.. 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.
Restart, fix_modify, output, run start/stop, minimize info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
No information about this fix is written to :doc:`binary restart files
<restart>`. This fix computes a global scalar (the number of
iterations) and a per-atom vector (the effective electronegativity), which
can be accessed by various :doc:`output commands <Howto_output>`.
No parameter of this fix can be used with the *start/stop* keywords of
the :doc:`run <run>` command.
This fix is invoked during :doc:`energy minimization <minimize>`.
Restrictions
""""""""""""
This fix is part of the REAXFF package. It is only enabled if
LAMMPS was built with that package. See the :doc:`Build package
<Build_package>` page for more info.
This fix does not correctly handle interactions involving multiple
periodic images of the same atom. Hence, it should not be used for
periodic cell dimensions smaller than the non-bonded cutoff radius,
which is typically :math:`10~\AA` for ReaxFF simulations.
This fix may be used in combination with :doc:`fix efield <fix_efield>`
and will apply the external electric field during charge equilibration,
but there may be only one fix efield instance used and the electric field
must be applied to all atoms in the system. Consequently, `fix efield` must
be used with *group-ID* all and must not be used with the keyword *region*.
Equal-style variables can be used for electric field vector
components without any further settings. Atom-style variables can be used
for spatially-varying electric field vector components, but the resulting
electric potential must be specified as an atom-style variable using
the *potential* keyword for `fix efield`.
Related commands
""""""""""""""""
:doc:`pair_style reaxff <pair_reaxff>`, :doc:`fix qeq/reaxff <fix_qeq_reaxff>`,
:doc:`fix acks2/reaxff <fix_acks2_reaxff>`
Default
"""""""
maxiter 200
----------
.. _Rappe3:
**(Rappe)** Rappe and Goddard III, Journal of Physical Chemistry, 95,
3358-3363 (1991).
.. _qtpie-Chen:
**(Chen)** Chen, Jiahao. Theory and applications of fluctuating-charge models.
University of Illinois at Urbana-Champaign, 2009.
.. _Gergs:
**(Gergs)** Gergs, Dirkmann and Mussenbrock.
Journal of Applied Physics 123.24 (2018).
.. _qeq-Aktulga2:
**(Aktulga)** Aktulga, Fogarty, Pandit, Grama, Parallel Computing, 38,
245-259 (2012).

26
doc/src/pair_reaxff.rst
View File

@ -20,7 +20,7 @@ Syntax
.. parsed-literal:: .. parsed-literal::
keyword = *checkqeq* or *lgvdw* or *safezone* or *mincap* or *minhbonds* or *tabulate* or *list/blocking* keyword = *checkqeq* or *lgvdw* or *safezone* or *mincap* or *minhbonds* or *tabulate* or *list/blocking*
*checkqeq* value = *yes* or *no* = whether or not to require qeq/reaxff or acks2/reaxff fix *checkqeq* value = *yes* or *no* = whether or not to require one of fix qeq/reaxff, fix acks2/reaxff or fix qtpie/reaxff
*enobonds* value = *yes* or *no* = whether or not to tally energy of atoms with no bonds *enobonds* value = *yes* or *no* = whether or not to tally energy of atoms with no bonds
*lgvdw* value = *yes* or *no* = whether or not to use a low gradient vdW correction *lgvdw* value = *yes* or *no* = whether or not to use a low gradient vdW correction
*safezone* = factor used for array allocation *safezone* = factor used for array allocation
@ -120,20 +120,22 @@ up that process.
The ReaxFF parameter files provided were created using a charge The ReaxFF parameter files provided were created using a charge
equilibration (QEq) model for handling the electrostatic interactions. equilibration (QEq) model for handling the electrostatic interactions.
Therefore, by default, LAMMPS requires that either the Therefore, by default, LAMMPS requires that
:doc:`fix qeq/reaxff <fix_qeq_reaxff>` or the :doc:`fix qeq/reaxff <fix_qeq_reaxff>` or :doc:`fix qeq/shielded <fix_qeq>`
:doc:`fix qeq/shielded <fix_qeq>` or :doc:`fix acks2/reaxff <fix_acks2_reaxff>` or :doc:`fix acks2/reaxff <fix_acks2_reaxff>`
command be used with or :doc:`fix qtpie/reaxff <fix_qtpie_reaxff>`
*pair_style reaxff* when simulating a ReaxFF model, to equilibrate is used with *pair_style reaxff* when simulating a ReaxFF model,
the charges each timestep. to equilibrate the charges at each timestep.
See the :doc:`fix qeq/reaxff <fix_qeq_reaxff>` or :doc:`fix qeq/shielded <fix_qeq>`
or :doc:`fix acks2/reaxff <fix_acks2_reaxff>`
or :doc:`fix qtpie/reaxff <fix_qtpie_reaxff>`
command documentation for more details.
Using the keyword *checkqeq* with the value *no* turns off the check Using the keyword *checkqeq* with the value *no* turns off the check
for the QEq fixes, allowing a simulation to be run without charge for the QEq fixes, allowing a simulation to be run without charge
equilibration. In this case, the static charges you assign to each equilibration. In this case, the static charges you assign to each
atom will be used for computing the electrostatic interactions in atom will be used for computing the electrostatic interactions in
the system. See the :doc:`fix qeq/reaxff <fix_qeq_reaxff>` or the system.
:doc:`fix qeq/shielded <fix_qeq>` or :doc:`fix acks2/reaxff <fix_acks2_reaxff>`
command documentation for more details.
Using the optional keyword *lgvdw* with the value *yes* turns on the Using the optional keyword *lgvdw* with the value *yes* turns on the
low-gradient correction of ReaxFF for long-range London Dispersion, low-gradient correction of ReaxFF for long-range London Dispersion,
@ -372,8 +374,8 @@ Related commands
"""""""""""""""" """"""""""""""""
:doc:`pair_coeff <pair_coeff>`, :doc:`fix qeq/reaxff <fix_qeq_reaxff>`, :doc:`pair_coeff <pair_coeff>`, :doc:`fix qeq/reaxff <fix_qeq_reaxff>`,
:doc:`fix acks2/reaxff <fix_acks2_reaxff>`, :doc:`fix reaxff/bonds <fix_reaxff_bonds>`, :doc:`fix acks2/reaxff <fix_acks2_reaxff>`, :doc:`fix qtpie/reaxff <fix_qtpie_reaxff>`,
:doc:`fix reaxff/species <fix_reaxff_species>`, :doc:`fix reaxff/bonds <fix_reaxff_bonds>`, :doc:`fix reaxff/species <fix_reaxff_species>`,
:doc:`compute reaxff/atom <compute_reaxff_atom>` :doc:`compute reaxff/atom <compute_reaxff_atom>`
Default Default

8
doc/utils/sphinx-config/false_positives.txt
View File

@ -815,6 +815,7 @@ dipoleflag
dir dir
Direc Direc
dirname dirname
Dirkmann
discoverable discoverable
discretization discretization
discretized discretized
@ -976,6 +977,7 @@ elaplong
elastance elastance
Electroneg Electroneg
electronegative electronegative
electronegativities
electronegativity electronegativity
electroneutral electroneutral
electroneutrality electroneutrality
@ -1293,6 +1295,7 @@ Geocomputing
georg georg
Georg Georg
Geotechnica Geotechnica
Gergs
germain germain
Germann Germann
Germano Germano
@ -1305,6 +1308,7 @@ gettimeofday
geturl geturl
gewald gewald
Gezelter Gezelter
gfile
Gflop Gflop
gfortran gfortran
ghostneigh ghostneigh
@ -1712,6 +1716,7 @@ Jewett
jgissing jgissing
ji ji
Jiang Jiang
Jiahao
Jiao Jiao
jik jik
JIK JIK
@ -2366,6 +2371,7 @@ mui
Mukherjee Mukherjee
Mulders Mulders
Müller Müller
Mulliken
mult mult
multi multi
multibody multibody
@ -2390,6 +2396,7 @@ Murdick
Murtola Murtola
Murty Murty
Muser Muser
Mussenbrock
mutexes mutexes
Muto Muto
muVT muVT
@ -3073,6 +3080,7 @@ qqr
qqrd qqrd
Qsb Qsb
qtb qtb
qtpie
quadratically quadratically
quadrupolar quadrupolar
quadrupole quadrupole

5
examples/reaxff/water/gauss_exp.txt Normal file
View File

@ -0,0 +1,5 @@
# Gaussian orbital exponents (required for fix qtpie/reaxff) taken from Table 2.2
# of Chen, J. (2009). Theory and applications of fluctuating-charge models.
# The units of the exponents are 1 / (Bohr radius)^2 .
1 0.2240 # O
2 0.5434 # H

29
examples/reaxff/water/in.water.qtpie Normal file
View File

@ -0,0 +1,29 @@
# QTPIE Water
boundary p p p
units real
atom_style charge
read_data data.water
variable x index 1
variable y index 1
variable z index 1
replicate $x $y $z
pair_style reaxff NULL safezone 3.0 mincap 150
pair_coeff * * qeq_ff.water O H
neighbor 0.5 bin
neigh_modify every 1 delay 0 check yes
velocity all create 300.0 4928459 rot yes dist gaussian
fix 1 all qtpie/reaxff 1 0.0 10.0 1.0e-6 reaxff gauss_exp.txt
fix 2 all nvt temp 300 300 50.0
timestep 0.5
thermo 10
thermo_style custom step temp press density vol
run 20

30
examples/reaxff/water/in.water.qtpie.field Normal file
View File

@ -0,0 +1,30 @@
# QTPIE Water
boundary p p p
units real
atom_style charge
read_data data.water
variable x index 1
variable y index 1
variable z index 1
replicate $x $y $z
pair_style reaxff NULL safezone 3.0 mincap 150
pair_coeff * * qeq_ff.water O H
neighbor 0.5 bin
neigh_modify every 1 delay 0 check yes
velocity all create 300.0 4928459 rot yes dist gaussian
fix 1 all qtpie/reaxff 1 0.0 10.0 1.0e-6 reaxff gauss_exp.txt
fix 2 all nvt temp 300 300 50.0
fix 3 all efield 0.0 0.0 0.05
timestep 0.5
thermo 10
thermo_style custom step temp press density vol
run 20

127
examples/reaxff/water/log.29Aug24.reaxff.water-qtpie-field.g++.1 Normal file
View File

@ -0,0 +1,127 @@
LAMMPS (29 Aug 2024 - Development - patch_29Aug2024-76-g3f232caf9b)
using 1 OpenMP thread(s) per MPI task
# QTPIE Water
boundary p p p
units real
atom_style charge
read_data data.water
Reading data file ...
orthogonal box = (0 0 0) to (31.043046 31.043046 31.043046)
1 by 1 by 1 MPI processor grid
reading atoms ...
3000 atoms
read_data CPU = 0.056 seconds
variable x index 1
variable y index 1
variable z index 1
replicate $x $y $z
replicate 1 $y $z
replicate 1 1 $z
replicate 1 1 1
Replication is creating a 1x1x1 = 1 times larger system...
orthogonal box = (0 0 0) to (31.043046 31.043046 31.043046)
1 by 1 by 1 MPI processor grid
3000 atoms
replicate CPU = 0.001 seconds
pair_style reaxff NULL safezone 3.0 mincap 150
pair_coeff * * qeq_ff.water O H
WARNING: Changed valency_val to valency_boc for X (src/REAXFF/reaxff_ffield.cpp:294)
neighbor 0.5 bin
neigh_modify every 1 delay 0 check yes
velocity all create 300.0 4928459 rot yes dist gaussian
fix 1 all qtpie/reaxff 1 0.0 10.0 1.0e-6 reaxff gauss_exp.txt
fix 2 all nvt temp 300 300 50.0
fix 3 all efield 0.0 0.0 0.05
timestep 0.5
thermo 10
thermo_style custom step temp press density vol
run 20
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- Type Label Framework: https://doi.org/10.1021/acs.jpcb.3c08419
@Article{Gissinger24,
author = {Jacob R. Gissinger, Ilia Nikiforov, Yaser Afshar, Brendon Waters, Moon-ki Choi, Daniel S. Karls, Alexander Stukowski, Wonpil Im, Hendrik Heinz, Axel Kohlmeyer, and Ellad B. Tadmor},
title = {Type Label Framework for Bonded Force Fields in LAMMPS},
journal = {J. Phys. Chem. B},
year = 2024,
volume = 128,
number = 13,
pages = {3282–-3297}
}
- pair reaxff command: doi:10.1016/j.parco.2011.08.005
@Article{Aktulga12,
author = {H. M. Aktulga and J. C. Fogarty and S. A. Pandit and A. Y. Grama},
title = {Parallel Reactive Molecular Dynamics: {N}umerical Methods and Algorithmic Techniques},
journal = {Parallel Computing},
year = 2012,
volume = 38,
number = {4--5},
pages = {245--259}
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 10.5
ghost atom cutoff = 10.5
binsize = 5.25, bins = 6 6 6
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair reaxff, perpetual
attributes: half, newton off, ghost
pair build: half/bin/ghost/newtoff
stencil: full/ghost/bin/3d
bin: standard
(2) fix qtpie/reaxff, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 539.2 | 539.2 | 539.2 Mbytes
Step Temp Press Density Volume
0 300 10137.041 1 29915.273
10 296.09128 3564.7969 1 29915.273
20 293.04308 10299.201 1 29915.273
Loop time of 10.7863 on 1 procs for 20 steps with 3000 atoms
Performance: 0.080 ns/day, 299.620 hours/ns, 1.854 timesteps/s, 5.563 katom-step/s
100.0% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 4.7275 | 4.7275 | 4.7275 | 0.0 | 43.83
Neigh | 0.17533 | 0.17533 | 0.17533 | 0.0 | 1.63
Comm | 0.0017376 | 0.0017376 | 0.0017376 | 0.0 | 0.02
Output | 8.2065e-05 | 8.2065e-05 | 8.2065e-05 | 0.0 | 0.00
Modify | 5.8812 | 5.8812 | 5.8812 | 0.0 | 54.52
Other | | 0.0005226 | | | 0.00
Nlocal: 3000 ave 3000 max 3000 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 11077 ave 11077 max 11077 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 971775 ave 971775 max 971775 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 971775
Ave neighs/atom = 323.925
Neighbor list builds = 2
Dangerous builds = 0
Total wall time: 0:00:12

127
examples/reaxff/water/log.29Aug24.reaxff.water-qtpie-field.g++.4 Normal file
View File

@ -0,0 +1,127 @@
LAMMPS (29 Aug 2024 - Development - patch_29Aug2024-76-g3f232caf9b)
using 1 OpenMP thread(s) per MPI task
# QTPIE Water
boundary p p p
units real
atom_style charge
read_data data.water
Reading data file ...
orthogonal box = (0 0 0) to (31.043046 31.043046 31.043046)
1 by 2 by 2 MPI processor grid
reading atoms ...
3000 atoms
read_data CPU = 0.053 seconds
variable x index 1
variable y index 1
variable z index 1
replicate $x $y $z
replicate 1 $y $z
replicate 1 1 $z
replicate 1 1 1
Replication is creating a 1x1x1 = 1 times larger system...
orthogonal box = (0 0 0) to (31.043046 31.043046 31.043046)
1 by 2 by 2 MPI processor grid
3000 atoms
replicate CPU = 0.002 seconds
pair_style reaxff NULL safezone 3.0 mincap 150
pair_coeff * * qeq_ff.water O H
WARNING: Changed valency_val to valency_boc for X (src/REAXFF/reaxff_ffield.cpp:294)
neighbor 0.5 bin
neigh_modify every 1 delay 0 check yes
velocity all create 300.0 4928459 rot yes dist gaussian
fix 1 all qtpie/reaxff 1 0.0 10.0 1.0e-6 reaxff gauss_exp.txt
fix 2 all nvt temp 300 300 50.0
fix 3 all efield 0.0 0.0 0.05
timestep 0.5
thermo 10
thermo_style custom step temp press density vol
run 20
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- Type Label Framework: https://doi.org/10.1021/acs.jpcb.3c08419
@Article{Gissinger24,
author = {Jacob R. Gissinger, Ilia Nikiforov, Yaser Afshar, Brendon Waters, Moon-ki Choi, Daniel S. Karls, Alexander Stukowski, Wonpil Im, Hendrik Heinz, Axel Kohlmeyer, and Ellad B. Tadmor},
title = {Type Label Framework for Bonded Force Fields in LAMMPS},
journal = {J. Phys. Chem. B},
year = 2024,
volume = 128,
number = 13,
pages = {3282–-3297}
}
- pair reaxff command: doi:10.1016/j.parco.2011.08.005
@Article{Aktulga12,
author = {H. M. Aktulga and J. C. Fogarty and S. A. Pandit and A. Y. Grama},
title = {Parallel Reactive Molecular Dynamics: {N}umerical Methods and Algorithmic Techniques},
journal = {Parallel Computing},
year = 2012,
volume = 38,
number = {4--5},
pages = {245--259}
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 10.5
ghost atom cutoff = 10.5
binsize = 5.25, bins = 6 6 6
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair reaxff, perpetual
attributes: half, newton off, ghost
pair build: half/bin/ghost/newtoff
stencil: full/ghost/bin/3d
bin: standard
(2) fix qtpie/reaxff, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 260.5 | 262.2 | 263.6 Mbytes
Step Temp Press Density Volume
0 300 10137.041 1 29915.273
10 296.09128 3564.7969 1 29915.273
20 293.04308 10299.201 1 29915.273
Loop time of 3.14492 on 4 procs for 20 steps with 3000 atoms
Performance: 0.275 ns/day, 87.359 hours/ns, 6.359 timesteps/s, 19.078 katom-step/s
99.6% CPU use with 4 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 1.6557 | 1.6847 | 1.7281 | 2.1 | 53.57
Neigh | 0.086503 | 0.086968 | 0.087627 | 0.2 | 2.77
Comm | 0.003309 | 0.046699 | 0.075729 | 12.4 | 1.48
Output | 5.0156e-05 | 5.483e-05 | 6.8111e-05 | 0.0 | 0.00
Modify | 1.3254 | 1.3261 | 1.3266 | 0.0 | 42.16
Other | | 0.0004552 | | | 0.01
Nlocal: 750 ave 760 max 735 min
Histogram: 1 0 0 0 1 0 0 0 0 2
Nghost: 6230.5 ave 6253 max 6193 min
Histogram: 1 0 0 0 0 0 1 0 1 1
Neighs: 276995 ave 280886 max 271360 min
Histogram: 1 0 0 0 1 0 0 0 1 1
Total # of neighbors = 1107981
Ave neighs/atom = 369.327
Neighbor list builds = 2
Dangerous builds = 0
Total wall time: 0:00:03

126
examples/reaxff/water/log.29Aug24.reaxff.water-qtpie.g++.1 Normal file
View File

@ -0,0 +1,126 @@
LAMMPS (29 Aug 2024 - Development - patch_29Aug2024-76-g3f232caf9b)
using 1 OpenMP thread(s) per MPI task
# QTPIE Water
boundary p p p
units real
atom_style charge
read_data data.water
Reading data file ...
orthogonal box = (0 0 0) to (31.043046 31.043046 31.043046)
1 by 1 by 1 MPI processor grid
reading atoms ...
3000 atoms
read_data CPU = 0.055 seconds
variable x index 1
variable y index 1
variable z index 1
replicate $x $y $z
replicate 1 $y $z
replicate 1 1 $z
replicate 1 1 1
Replication is creating a 1x1x1 = 1 times larger system...
orthogonal box = (0 0 0) to (31.043046 31.043046 31.043046)
1 by 1 by 1 MPI processor grid
3000 atoms
replicate CPU = 0.001 seconds
pair_style reaxff NULL safezone 3.0 mincap 150
pair_coeff * * qeq_ff.water O H
WARNING: Changed valency_val to valency_boc for X (src/REAXFF/reaxff_ffield.cpp:294)
neighbor 0.5 bin
neigh_modify every 1 delay 0 check yes
velocity all create 300.0 4928459 rot yes dist gaussian
fix 1 all qtpie/reaxff 1 0.0 10.0 1.0e-6 reaxff gauss_exp.txt
fix 2 all nvt temp 300 300 50.0
timestep 0.5
thermo 10
thermo_style custom step temp press density vol
run 20
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- Type Label Framework: https://doi.org/10.1021/acs.jpcb.3c08419
@Article{Gissinger24,
author = {Jacob R. Gissinger, Ilia Nikiforov, Yaser Afshar, Brendon Waters, Moon-ki Choi, Daniel S. Karls, Alexander Stukowski, Wonpil Im, Hendrik Heinz, Axel Kohlmeyer, and Ellad B. Tadmor},
title = {Type Label Framework for Bonded Force Fields in LAMMPS},
journal = {J. Phys. Chem. B},
year = 2024,
volume = 128,
number = 13,
pages = {3282–-3297}
}
- pair reaxff command: doi:10.1016/j.parco.2011.08.005
@Article{Aktulga12,
author = {H. M. Aktulga and J. C. Fogarty and S. A. Pandit and A. Y. Grama},
title = {Parallel Reactive Molecular Dynamics: {N}umerical Methods and Algorithmic Techniques},
journal = {Parallel Computing},
year = 2012,
volume = 38,
number = {4--5},
pages = {245--259}
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 10.5
ghost atom cutoff = 10.5
binsize = 5.25, bins = 6 6 6
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair reaxff, perpetual
attributes: half, newton off, ghost
pair build: half/bin/ghost/newtoff
stencil: full/ghost/bin/3d
bin: standard
(2) fix qtpie/reaxff, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 539.2 | 539.2 | 539.2 Mbytes
Step Temp Press Density Volume
0 300 10138.375 1 29915.273
10 295.97879 3575.2769 1 29915.273
20 292.76583 10309.128 1 29915.273
Loop time of 10.8138 on 1 procs for 20 steps with 3000 atoms
Performance: 0.080 ns/day, 300.383 hours/ns, 1.849 timesteps/s, 5.548 katom-step/s
99.9% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 4.7177 | 4.7177 | 4.7177 | 0.0 | 43.63
Neigh | 0.17607 | 0.17607 | 0.17607 | 0.0 | 1.63
Comm | 0.0017295 | 0.0017295 | 0.0017295 | 0.0 | 0.02
Output | 8.5431e-05 | 8.5431e-05 | 8.5431e-05 | 0.0 | 0.00
Modify | 5.9177 | 5.9177 | 5.9177 | 0.0 | 54.72
Other | | 0.0004911 | | | 0.00
Nlocal: 3000 ave 3000 max 3000 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 11077 ave 11077 max 11077 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 971830 ave 971830 max 971830 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 971830
Ave neighs/atom = 323.94333
Neighbor list builds = 2
Dangerous builds = 0
Total wall time: 0:00:12

126
examples/reaxff/water/log.29Aug24.reaxff.water-qtpie.g++.4 Normal file
View File

@ -0,0 +1,126 @@
LAMMPS (29 Aug 2024 - Development - patch_29Aug2024-76-g3f232caf9b)
using 1 OpenMP thread(s) per MPI task
# QTPIE Water
boundary p p p
units real
atom_style charge
read_data data.water
Reading data file ...
orthogonal box = (0 0 0) to (31.043046 31.043046 31.043046)
1 by 2 by 2 MPI processor grid
reading atoms ...
3000 atoms
read_data CPU = 0.053 seconds
variable x index 1
variable y index 1
variable z index 1
replicate $x $y $z
replicate 1 $y $z
replicate 1 1 $z
replicate 1 1 1
Replication is creating a 1x1x1 = 1 times larger system...
orthogonal box = (0 0 0) to (31.043046 31.043046 31.043046)
1 by 2 by 2 MPI processor grid
3000 atoms
replicate CPU = 0.002 seconds
pair_style reaxff NULL safezone 3.0 mincap 150
pair_coeff * * qeq_ff.water O H
WARNING: Changed valency_val to valency_boc for X (src/REAXFF/reaxff_ffield.cpp:294)
neighbor 0.5 bin
neigh_modify every 1 delay 0 check yes
velocity all create 300.0 4928459 rot yes dist gaussian
fix 1 all qtpie/reaxff 1 0.0 10.0 1.0e-6 reaxff gauss_exp.txt
fix 2 all nvt temp 300 300 50.0
timestep 0.5
thermo 10
thermo_style custom step temp press density vol
run 20
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- Type Label Framework: https://doi.org/10.1021/acs.jpcb.3c08419
@Article{Gissinger24,
author = {Jacob R. Gissinger, Ilia Nikiforov, Yaser Afshar, Brendon Waters, Moon-ki Choi, Daniel S. Karls, Alexander Stukowski, Wonpil Im, Hendrik Heinz, Axel Kohlmeyer, and Ellad B. Tadmor},
title = {Type Label Framework for Bonded Force Fields in LAMMPS},
journal = {J. Phys. Chem. B},
year = 2024,
volume = 128,
number = 13,
pages = {3282–-3297}
}
- pair reaxff command: doi:10.1016/j.parco.2011.08.005
@Article{Aktulga12,
author = {H. M. Aktulga and J. C. Fogarty and S. A. Pandit and A. Y. Grama},
title = {Parallel Reactive Molecular Dynamics: {N}umerical Methods and Algorithmic Techniques},
journal = {Parallel Computing},
year = 2012,
volume = 38,
number = {4--5},
pages = {245--259}
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 10.5
ghost atom cutoff = 10.5
binsize = 5.25, bins = 6 6 6
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair reaxff, perpetual
attributes: half, newton off, ghost
pair build: half/bin/ghost/newtoff
stencil: full/ghost/bin/3d
bin: standard
(2) fix qtpie/reaxff, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 260.5 | 262.2 | 263.6 Mbytes
Step Temp Press Density Volume
0 300 10138.375 1 29915.273
10 295.97879 3575.2769 1 29915.273
20 292.76583 10309.128 1 29915.273
Loop time of 3.13598 on 4 procs for 20 steps with 3000 atoms
Performance: 0.276 ns/day, 87.111 hours/ns, 6.378 timesteps/s, 19.133 katom-step/s
99.6% CPU use with 4 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 1.6622 | 1.695 | 1.7252 | 2.2 | 54.05
Neigh | 0.086543 | 0.087117 | 0.087848 | 0.2 | 2.78
Comm | 0.0048192 | 0.035002 | 0.067754 | 15.4 | 1.12
Output | 4.8033e-05 | 5.3375e-05 | 6.6893e-05 | 0.0 | 0.00
Modify | 1.3176 | 1.3183 | 1.3189 | 0.0 | 42.04
Other | | 0.0004753 | | | 0.02
Nlocal: 750 ave 760 max 735 min
Histogram: 1 0 0 0 1 0 0 0 0 2
Nghost: 6229.5 ave 6253 max 6191 min
Histogram: 1 0 0 0 0 0 1 0 1 1
Neighs: 277011 ave 280900 max 271380 min
Histogram: 1 0 0 0 1 0 0 0 1 1
Total # of neighbors = 1108044
Ave neighs/atom = 369.348
Neighbor list builds = 2
Dangerous builds = 0
Total wall time: 0:00:03

2
src/.gitignore vendored
View File

@ -988,6 +988,8 @@
/fix_qeq_reaxff.h /fix_qeq_reaxff.h
/fix_qmmm.cpp /fix_qmmm.cpp
/fix_qmmm.h /fix_qmmm.h
/fix_qtpie_reaxff.cpp
/fix_qtpie_reaxff.h
/fix_reaxff.cpp /fix_reaxff.cpp
/fix_reaxff.h /fix_reaxff.h
/fix_reaxff_bonds.cpp /fix_reaxff_bonds.cpp

2
src/REAXFF/fix_qeq_reaxff.cpp
View File

@ -1158,7 +1158,7 @@ void FixQEqReaxFF::get_chi_field()
for (int i = 0; i < nlocal; i++) { for (int i = 0; i < nlocal; i++) {
if (mask[i] & efgroupbit) { if (mask[i] & efgroupbit) {
if (region && !region->match(x[i][0],x[i][1],x[i][2])) continue; if (region && !region->match(x[i][0],x[i][1],x[i][2])) continue;
chi_field[i] = -efield->efield[i][3]; chi_field[i] = efield->efield[i][3];
} }
} }
} }

1193
src/REAXFF/fix_qtpie_reaxff.cpp Normal file
View File

File diff suppressed because it is too large Load Diff

141
src/REAXFF/fix_qtpie_reaxff.h Normal file
View File

@ -0,0 +1,141 @@
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
https://www.lammps.org/, Sandia National Laboratories
LAMMPS development team: developers@lammps.org
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#ifdef FIX_CLASS
// clang-format off
FixStyle(qtpie/reaxff,FixQtpieReaxFF);
// clang-format on
#else
#ifndef LMP_FIX_QTPIE_REAXFF_H
#define LMP_FIX_QTPIE_REAXFF_H
#include "fix.h"
namespace LAMMPS_NS {
class FixQtpieReaxFF : public Fix {
public:
FixQtpieReaxFF(class LAMMPS *, int, char **);
~FixQtpieReaxFF() override;
int setmask() override;
void post_constructor() override;
void init() override;
void init_list(int, class NeighList *) override;
void init_storage();
void setup_pre_force(int) override;
void pre_force(int) override;
void setup_pre_force_respa(int, int) override;
void pre_force_respa(int, int, int) override;
void min_setup_pre_force(int);
void min_pre_force(int) override;
double compute_scalar() override;
protected:
int nevery, reaxflag;
int matvecs;
int nn, nt, m_fill;
int n_cap, nmax, m_cap;
int pack_flag;
int nlevels_respa;
class NeighList *list;
class PairReaxFF *reaxff;
class FixEfield *efield;
int *ilist, *jlist, *numneigh, **firstneigh;
double swa, swb; // lower/upper Taper cutoff radius
double Tap[8]; // Taper function
double tolerance; // tolerance for the norm of the rel residual in CG
double *chi, *eta, *gamma; // qeq parameters
double **shld;
// fictitious charges
double *s, *t;
double **s_hist, **t_hist;
int nprev;
typedef struct {
int n, m;
int *firstnbr;
int *numnbrs;
int *jlist;
double *val;
} sparse_matrix;
sparse_matrix H;
double *Hdia_inv;
double *b_s, *b_t;
double *b_prc, *b_prm;
double *chi_eff; // array of effective electronegativities
//CG storage
double *p, *q, *r, *d;
int imax, maxwarn;
char *pertype_option; // argument to determine how per-type info is obtained
char *gauss_file; // input file for gaussian orbital exponents
double *gauss_exp; // array of gaussian orbital exponents for each atom type
double dist_cutoff; // separation distance beyond which to neglect overlap integrals
void pertype_parameters(char *);
void init_shielding();
void init_taper();
void allocate_storage();
void deallocate_storage();
void reallocate_storage();
void allocate_matrix();
void deallocate_matrix();
void reallocate_matrix();
void init_matvec();
void init_H();
void compute_H();
double calculate_H(double, double);
void calculate_Q();
int CG(double *, double *);
void sparse_matvec(sparse_matrix *, double *, double *);
int pack_forward_comm(int, int *, double *, int, int *) override;
void unpack_forward_comm(int, int, double *) override;
int pack_reverse_comm(int, int, double *) override;
void unpack_reverse_comm(int, int *, double *) override;
double memory_usage() override;
void grow_arrays(int) override;
void copy_arrays(int, int, int) override;
int pack_exchange(int, double *) override;
int unpack_exchange(int, double *) override;
double parallel_norm(double *, int);
double parallel_dot(double *, double *, int);
double parallel_vector_acc(double *, int);
void vector_sum(double *, double, double *, double, double *, int);
void vector_add(double *, double, double *, int);
void calc_chi_eff();
double find_min_exp(const double*, const int);
double distance(const double*, const double*);
int matvecs_s, matvecs_t; // Iteration count for each system
};
} // namespace LAMMPS_NS
#endif
#endif

6
src/REAXFF/pair_reaxff.cpp
View File

@ -340,11 +340,13 @@ void PairReaxFF::init_style()
auto acks2_fixes = modify->get_fix_by_style("^acks2/reax"); auto acks2_fixes = modify->get_fix_by_style("^acks2/reax");
int have_qeq = modify->get_fix_by_style("^qeq/reax").size() int have_qeq = modify->get_fix_by_style("^qeq/reax").size()
+ modify->get_fix_by_style("^qeq/shielded").size() + acks2_fixes.size(); + modify->get_fix_by_style("^qeq/shielded").size() + acks2_fixes.size()
+ modify->get_fix_by_style("^qtpie/reax").size();
if (qeqflag && (have_qeq != 1)) if (qeqflag && (have_qeq != 1))
error->all(FLERR,"Pair style reaxff requires use of exactly one of the " error->all(FLERR,"Pair style reaxff requires use of exactly one of the "
"fix qeq/reaxff or fix qeq/shielded or fix acks2/reaxff commands"); "fix qeq/reaxff or fix qeq/shielded or fix acks2/reaxff or "
"fix qtpie/reaxff commands");
api->system->acks2_flag = acks2_fixes.size(); api->system->acks2_flag = acks2_fixes.size();
if (api->system->acks2_flag) if (api->system->acks2_flag)

1
src/fix_efield.h
View File

@ -26,6 +26,7 @@ namespace LAMMPS_NS {
class FixEfield : public Fix { class FixEfield : public Fix {
friend class FixQEqReaxFF; friend class FixQEqReaxFF;
friend class FixQtpieReaxFF;
public: public:
FixEfield(class LAMMPS *, int, char **); FixEfield(class LAMMPS *, int, char **);