Axel Kohlmeyer
@ -186,6 +186,7 @@ OPT.
|
||||
* :doc:`qeq/fire <fix_qeq>`
|
||||
* :doc:`qeq/point <fix_qeq>`
|
||||
* :doc:`qeq/reaxff (ko) <fix_qeq_reaxff>`
|
||||
* :doc:`qeq/rel/reaxff <fix_qeq_rel_reaxff>`
|
||||
* :doc:`qeq/shielded <fix_qeq>`
|
||||
* :doc:`qeq/slater <fix_qeq>`
|
||||
* :doc:`qmmm <fix_qmmm>`
|
||||
|
||||
@ -365,6 +365,7 @@ accelerated styles exist.
|
||||
* :doc:`qeq/fire <fix_qeq>` - charge equilibration via FIRE minimizer
|
||||
* :doc:`qeq/point <fix_qeq>` - charge equilibration via point method
|
||||
* :doc:`qeq/reaxff <fix_qeq_reaxff>` - charge equilibration for ReaxFF potential
|
||||
* :doc:`qeq/rel/reaxff <fix_qeq_rel_reaxff>` - charge equilibration for ReaxFF potential with alternate efield implementation
|
||||
* :doc:`qeq/shielded <fix_qeq>` - charge equilibration via shielded method
|
||||
* :doc:`qeq/slater <fix_qeq>` - charge equilibration via Slater method
|
||||
* :doc:`qmmm <fix_qmmm>` - functionality to enable a quantum mechanics/molecular mechanics coupling
|
||||
|
||||
@ -123,8 +123,10 @@ components in non-periodic directions.
|
||||
Related commands
|
||||
""""""""""""""""
|
||||
|
||||
:doc:`pair_style reaxff <pair_reaxff>`, :doc:`fix qeq/reaxff <fix_qeq_reaxff>`,
|
||||
:doc:`fix qtpi/reaxff <fix_qtpie_reaxff>`
|
||||
:doc:`pair_style reaxff <pair_reaxff>`,
|
||||
:doc:`fix qeq/reaxff <fix_qeq_reaxff>`,
|
||||
:doc:`fix qtpie/reaxff <fix_qtpie_reaxff>`,
|
||||
:doc:`fix qeq/rel/reaxff <fix_qeq_rel_reaxff>`
|
||||
|
||||
Default
|
||||
"""""""
|
||||
|
||||
@ -59,7 +59,7 @@ 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 QeQ with a non-ReaxFF potential.
|
||||
form must be used when performing QEq with a non-ReaxFF potential.
|
||||
Each line should be formatted as follows:
|
||||
|
||||
.. parsed-literal::
|
||||
@ -140,7 +140,8 @@ Related commands
|
||||
""""""""""""""""
|
||||
|
||||
: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>`
|
||||
:doc:`fix acks2/reaxff <fix_acks2_reaxff>`, :doc:`fix qtpie/reaxff <fix_qtpie_reaxff>`,
|
||||
:doc:`fix qeq/rel/reaxff <fix_qeq_rel_reaxff>`
|
||||
|
||||
Default
|
||||
"""""""
|
||||
|
||||
195
doc/src/fix_qeq_rel_reaxff.rst
Normal file
195
doc/src/fix_qeq_rel_reaxff.rst
Normal file
@ -0,0 +1,195 @@
|
||||
.. index:: fix qeq/rel/reaxff
|
||||
|
||||
fix qeq/rel/reaxff command
|
||||
==========================
|
||||
|
||||
Syntax
|
||||
""""""
|
||||
|
||||
.. code-block:: LAMMPS
|
||||
|
||||
fix ID group-ID qeq/rel/reaxff Nevery cutlo cuthi tolerance params gfile args
|
||||
|
||||
* ID, group-ID are documented in :doc:`fix <fix>` command
|
||||
* qeq/rel/reaxff = style name of this fix command
|
||||
* Nevery = perform QEqR 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 = *scale* or *maxiter* or *nowarn*
|
||||
*scale* beta = set value of scaling factor *beta* (determines strength of electric polarization)
|
||||
*maxiter* N = limit the number of iterations to *N*
|
||||
*nowarn* = do not print a warning message if the maximum number of iterations is reached
|
||||
|
||||
Examples
|
||||
""""""""
|
||||
|
||||
.. code-block:: LAMMPS
|
||||
|
||||
fix 1 all qeq/rel/reaxff 1 0.0 10.0 1.0e-6 reaxff exp.qeqr
|
||||
fix 1 all qeq/rel/reaxff 1 0.0 10.0 1.0e-6 params.qeqr exp.qeqr scale 1.5 maxiter 500 nowarn
|
||||
|
||||
Description
|
||||
"""""""""""
|
||||
|
||||
.. versionadded:: 19Nov2024
|
||||
|
||||
This fix implements the QEqR method for charge equilibration, which
|
||||
differs from the QEq charge equilibration method :ref:`(Rappe and
|
||||
Goddard) <Rappe4>` only in how external electric fields are accounted
|
||||
for. This fix therefore raises a warning when used without :doc:`fix
|
||||
efield <fix_efield>` since :doc:`fix qeq/reaxff <fix_qeq_reaxff>` should
|
||||
be used without an external electric field. Charges are computed with
|
||||
the QEqR method 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 method is replaced with an
|
||||
effective electronegativity given by
|
||||
|
||||
.. math::
|
||||
\chi_{\mathrm{r}i} = \chi_i + \frac{\sum_{j=1}^{N} \beta(\phi_i - \phi_j) S_{ij}}
|
||||
{\sum_{m=1}^{N}S_{im}},
|
||||
|
||||
where :math:`N` is the number of atoms in the system, :math:`\beta` is a
|
||||
scaling factor, :math:`\phi_i` and :math:`\phi_j` are the electric
|
||||
potentials at the positions of atoms :math:`i` and :math:`j` due to the
|
||||
external electric field and :math:`S_{ij}` is the overlap integral
|
||||
between atoms :math:`i` and :math:`j`. This formulation is advantageous
|
||||
over the method used by :doc:`fix qeq/reaxff <fix_qeq_reaxff>` to
|
||||
account for an external electric field in that it permits periodic
|
||||
boundaries in the direction of an external electric field and in that it
|
||||
does not worsen long-range charge transfer seen with QEq.
|
||||
|
||||
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 qeq/rel/reaxff*, which is the
|
||||
same as in :doc:`fix qeq/reaxff <fix_qeq_reaxff>` except for the use of
|
||||
:math:`\chi_{\mathrm{r}i}`, please refer to :ref:`(Aktulga)
|
||||
<qeq-Aktulga3>`. To be explicit, *fix qeq/rel/reaxff* replaces
|
||||
:math:`\chi_k` of eq. 3 in :ref:`(Aktulga) <qeq-Aktulga3>` with
|
||||
:math:`\chi_{\mathrm{r}k}` when an external electric field is applied.
|
||||
|
||||
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
|
||||
using this fix 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 :math:`S_{ij}` 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 *scale* keyword sets the value of :math:`\beta` in the
|
||||
equation for :math:`\chi_{\mathrm{r}i}`. The default value is 1.0.
|
||||
|
||||
The optional *maxiter* keyword allows changing the max number of
|
||||
iterations in the linear solver. The default value is 200.
|
||||
|
||||
The optional *nowarn* keyword silences the warning message printed when
|
||||
the maximum number of iterations is reached. This can be useful for
|
||||
comparing serial and parallel results where having the same fixed number
|
||||
of iterations is desired, which can be achieved by using a very small
|
||||
tolerance and setting *maxiter* to the desired number of iterations.
|
||||
|
||||
.. 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>`, :doc:`fix qtpie/reaxff <fix_qtpie_reaxff>`
|
||||
|
||||
Default
|
||||
"""""""
|
||||
|
||||
scale = 1.0 and maxiter = 200
|
||||
|
||||
----------
|
||||
|
||||
.. _Rappe4:
|
||||
|
||||
**(Rappe)** Rappe and Goddard III, Journal of Physical Chemistry, 95,
|
||||
3358-3363 (1991).
|
||||
|
||||
.. _qeq-Aktulga3:
|
||||
|
||||
**(Aktulga)** Aktulga, Fogarty, Pandit, Grama, Parallel Computing, 38,
|
||||
245-259 (2012).
|
||||
@ -21,8 +21,10 @@ Syntax
|
||||
|
||||
.. parsed-literal::
|
||||
|
||||
keyword = *maxiter*
|
||||
keyword = *scale* or *maxiter* or *nowarn*
|
||||
*scale* beta = set value of scaling factor *beta* (determines strength of electric polarization)
|
||||
*maxiter* N = limit the number of iterations to *N*
|
||||
*nowarn* = do not print a warning message if the maximum number of iterations is reached
|
||||
|
||||
Examples
|
||||
""""""""
|
||||
@ -30,7 +32,7 @@ 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
|
||||
fix 1 all qtpie/reaxff 1 0.0 10.0 1.0e-6 params.qtpie exp.qtpie scale 1.5 maxiter 500 nowarn
|
||||
|
||||
Description
|
||||
"""""""""""
|
||||
@ -46,7 +48,7 @@ same way as the QEq method but where the absolute electronegativity,
|
||||
electronegativity given by :ref:`(Chen) <qtpie-Chen>`
|
||||
|
||||
.. math::
|
||||
\chi_{\mathrm{eff},i} = \frac{\sum_{j=1}^{N} (\chi_i - \chi_j) S_{ij}}
|
||||
\tilde{\chi}_{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
|
||||
@ -61,11 +63,11 @@ 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}}
|
||||
\tilde{\chi}_{\mathrm{r}i} = \frac{\sum_{j=1}^{N} (\chi_i - \chi_j + \beta(\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`
|
||||
where :math:`\beta` is a scaling factor and :math:`\phi_i` and :math:`\phi_j`
|
||||
are the electric potentials at the positions of atoms :math:`i` and :math:`j`
|
||||
due to the external electric field.
|
||||
|
||||
This fix is typically used in conjunction with the ReaxFF force
|
||||
@ -74,9 +76,12 @@ 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>`.
|
||||
:math:`\tilde{\chi}_{i}` or :math:`\tilde{\chi}_{\mathrm{r}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}`.
|
||||
:ref:`(Aktulga) <qeq-Aktulga2>` with :math:`\tilde{\chi}_{k}` when no external
|
||||
electric field is applied and with :math:`\tilde{\chi}_{\mathrm{r}k}` when an
|
||||
external electric field is applied.
|
||||
|
||||
This fix requires the absolute electronegativity, :math:`\chi`, in eV, the
|
||||
self-Coulomb potential, :math:`\eta`, in eV, and the shielded Coulomb
|
||||
@ -97,7 +102,7 @@ respectively:
|
||||
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}`
|
||||
The overlap integrals :math:`S_{ij}`
|
||||
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*.
|
||||
@ -120,15 +125,26 @@ Empty lines or any text following the pound sign (#) are ignored. An example
|
||||
1 0.2240 # O
|
||||
2 0.5434 # H
|
||||
|
||||
The optional *scale* keyword sets the value of :math:`\beta` in the equation for
|
||||
:math:`\tilde{\chi}_{\mathrm{r}i}`. This keyword only affects the computed charges
|
||||
when :doc:`fix efield <fix_efield>` is used. The default value is 1.0.
|
||||
|
||||
The optional *maxiter* keyword allows changing the max number
|
||||
of iterations in the linear solver. The default value is 200.
|
||||
|
||||
The optional *nowarn* keyword silences the warning message printed
|
||||
when the maximum number of iterations is reached. This can be
|
||||
useful for comparing serial and parallel results where having the
|
||||
same fixed number of iterations is desired, which can be achieved
|
||||
by using a very small tolerance and setting *maxiter* to the desired
|
||||
number of iterations.
|
||||
|
||||
.. 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
|
||||
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
|
||||
@ -170,12 +186,13 @@ Related commands
|
||||
""""""""""""""""
|
||||
|
||||
:doc:`pair_style reaxff <pair_reaxff>`, :doc:`fix qeq/reaxff <fix_qeq_reaxff>`,
|
||||
:doc:`fix acks2/reaxff <fix_acks2_reaxff>`
|
||||
:doc:`fix acks2/reaxff <fix_acks2_reaxff>`,
|
||||
:doc:`fix qeq/rel/reaxff <fix_qeq_rel_reaxff>`
|
||||
|
||||
Default
|
||||
"""""""
|
||||
|
||||
maxiter 200
|
||||
scale = 1.0 and maxiter = 200
|
||||
|
||||
----------
|
||||
|
||||
|
||||
@ -725,6 +725,7 @@ dashpot
|
||||
dat
|
||||
datafile
|
||||
datatype
|
||||
dataset
|
||||
datums
|
||||
Davidchack
|
||||
Daw
|
||||
@ -3120,9 +3121,11 @@ qE
|
||||
qeff
|
||||
qelectron
|
||||
qeq
|
||||
qeqr
|
||||
Qamar
|
||||
QeQ
|
||||
QEq
|
||||
QEqR
|
||||
qfactor
|
||||
qfile
|
||||
qi
|
||||
|
||||
29
examples/reaxff/water/in.water.qeqr
Normal file
29
examples/reaxff/water/in.water.qeqr
Normal file
@ -0,0 +1,29 @@
|
||||
# Water with QEqR
|
||||
|
||||
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 qeqr/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.qeqr.field
Normal file
30
examples/reaxff/water/in.water.qeqr.field
Normal file
@ -0,0 +1,30 @@
|
||||
# Water with QEqR
|
||||
|
||||
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 qeqr/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
|
||||
115
examples/reaxff/water/log.20Mar25.reaxff.water-qeqr-field.g++.1
Normal file
115
examples/reaxff/water/log.20Mar25.reaxff.water-qeqr-field.g++.1
Normal file
@ -0,0 +1,115 @@
|
||||
LAMMPS (4 Feb 2025 - Development - patch_4Feb2025-444-gbb8b6590d5-modified)
|
||||
using 1 OpenMP thread(s) per MPI task
|
||||
# Water with QEqR
|
||||
|
||||
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.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 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:300)
|
||||
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 qeqr/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:
|
||||
|
||||
- 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 qeqr/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 778.75601 1 29915.273
|
||||
10 301.42845 5423.6612 1 29915.273
|
||||
20 298.24707 1549.2257 1 29915.273
|
||||
Loop time of 10.6859 on 1 procs for 20 steps with 3000 atoms
|
||||
|
||||
Performance: 0.081 ns/day, 296.830 hours/ns, 1.872 timesteps/s, 5.615 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.7595 | 4.7595 | 4.7595 | 0.0 | 44.54
|
||||
Neigh | 0.17605 | 0.17605 | 0.17605 | 0.0 | 1.65
|
||||
Comm | 0.0017511 | 0.0017511 | 0.0017511 | 0.0 | 0.02
|
||||
Output | 8.3809e-05 | 8.3809e-05 | 8.3809e-05 | 0.0 | 0.00
|
||||
Modify | 5.748 | 5.748 | 5.748 | 0.0 | 53.79
|
||||
Other | | 0.0005279 | | | 0.00
|
||||
|
||||
Nlocal: 3000 ave 3000 max 3000 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
Nghost: 11075 ave 11075 max 11075 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
Neighs: 971785 ave 971785 max 971785 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
|
||||
Total # of neighbors = 971785
|
||||
Ave neighs/atom = 323.92833
|
||||
Neighbor list builds = 2
|
||||
Dangerous builds = 0
|
||||
Total wall time: 0:00:12
|
||||
115
examples/reaxff/water/log.20Mar25.reaxff.water-qeqr-field.g++.4
Normal file
115
examples/reaxff/water/log.20Mar25.reaxff.water-qeqr-field.g++.4
Normal file
@ -0,0 +1,115 @@
|
||||
LAMMPS (4 Feb 2025 - Development - patch_4Feb2025-444-gbb8b6590d5-modified)
|
||||
using 1 OpenMP thread(s) per MPI task
|
||||
# Water with QEqR
|
||||
|
||||
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:300)
|
||||
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 qeqr/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:
|
||||
|
||||
- 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 qeqr/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 778.75601 1 29915.273
|
||||
10 301.42845 5423.6623 1 29915.273
|
||||
20 298.24708 1549.2264 1 29915.273
|
||||
Loop time of 3.10467 on 4 procs for 20 steps with 3000 atoms
|
||||
|
||||
Performance: 0.278 ns/day, 86.241 hours/ns, 6.442 timesteps/s, 19.326 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.668 | 1.6843 | 1.7266 | 1.9 | 54.25
|
||||
Neigh | 0.08549 | 0.086004 | 0.086638 | 0.2 | 2.77
|
||||
Comm | 0.0135 | 0.055821 | 0.072105 | 10.4 | 1.80
|
||||
Output | 4.9632e-05 | 5.4515e-05 | 6.8384e-05 | 0.0 | 0.00
|
||||
Modify | 1.2774 | 1.2781 | 1.2786 | 0.0 | 41.17
|
||||
Other | | 0.000458 | | | 0.01
|
||||
|
||||
Nlocal: 750 ave 760 max 735 min
|
||||
Histogram: 1 0 0 0 1 0 0 0 0 2
|
||||
Nghost: 6230.75 ave 6255 max 6191 min
|
||||
Histogram: 1 0 0 0 0 1 0 0 1 1
|
||||
Neighs: 276996 ave 280553 max 271385 min
|
||||
Histogram: 1 0 0 0 0 1 0 0 0 2
|
||||
|
||||
Total # of neighbors = 1107985
|
||||
Ave neighs/atom = 369.32833
|
||||
Neighbor list builds = 2
|
||||
Dangerous builds = 0
|
||||
Total wall time: 0:00:03
|
||||
116
examples/reaxff/water/log.20Mar25.reaxff.water-qeqr.g++.1
Normal file
116
examples/reaxff/water/log.20Mar25.reaxff.water-qeqr.g++.1
Normal file
@ -0,0 +1,116 @@
|
||||
LAMMPS (4 Feb 2025 - Development - patch_4Feb2025-444-gbb8b6590d5-modified)
|
||||
using 1 OpenMP thread(s) per MPI task
|
||||
# Water with QEqR
|
||||
|
||||
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:300)
|
||||
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 qeqr/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:
|
||||
|
||||
- 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
|
||||
|
||||
WARNING: Use fix qeq/reaxff instead of fix qeqr/reaxff when not using fix efield
|
||||
(src/REAXFF/fix_qtpie_reaxff.cpp:493)
|
||||
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 qeqr/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 780.33989 1 29915.273
|
||||
10 301.29205 5433.7414 1 29915.273
|
||||
20 297.90652 1572.6111 1 29915.273
|
||||
Loop time of 6.87447 on 1 procs for 20 steps with 3000 atoms
|
||||
|
||||
Performance: 0.126 ns/day, 190.957 hours/ns, 2.909 timesteps/s, 8.728 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.8461 | 4.8461 | 4.8461 | 0.0 | 70.49
|
||||
Neigh | 0.17595 | 0.17595 | 0.17595 | 0.0 | 2.56
|
||||
Comm | 0.001787 | 0.001787 | 0.001787 | 0.0 | 0.03
|
||||
Output | 8.5794e-05 | 8.5794e-05 | 8.5794e-05 | 0.0 | 0.00
|
||||
Modify | 1.8501 | 1.8501 | 1.8501 | 0.0 | 26.91
|
||||
Other | | 0.0004811 | | | 0.01
|
||||
|
||||
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: 971826 ave 971826 max 971826 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 0
|
||||
|
||||
Total # of neighbors = 971826
|
||||
Ave neighs/atom = 323.942
|
||||
Neighbor list builds = 2
|
||||
Dangerous builds = 0
|
||||
Total wall time: 0:00:07
|
||||
116
examples/reaxff/water/log.20Mar25.reaxff.water-qeqr.g++.4
Normal file
116
examples/reaxff/water/log.20Mar25.reaxff.water-qeqr.g++.4
Normal file
@ -0,0 +1,116 @@
|
||||
LAMMPS (4 Feb 2025 - Development - patch_4Feb2025-444-gbb8b6590d5-modified)
|
||||
using 1 OpenMP thread(s) per MPI task
|
||||
# Water with QEqR
|
||||
|
||||
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.082 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:300)
|
||||
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 qeqr/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:
|
||||
|
||||
- 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
|
||||
|
||||
WARNING: Use fix qeq/reaxff instead of fix qeqr/reaxff when not using fix efield
|
||||
(src/REAXFF/fix_qtpie_reaxff.cpp:493)
|
||||
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 qeqr/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 780.34006 1 29915.273
|
||||
10 301.29205 5433.7414 1 29915.273
|
||||
20 297.90652 1572.6112 1 29915.273
|
||||
Loop time of 2.52349 on 4 procs for 20 steps with 3000 atoms
|
||||
|
||||
Performance: 0.342 ns/day, 70.097 hours/ns, 7.926 timesteps/s, 23.777 katom-step/s
|
||||
99.0% 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.7081 | 1.7518 | 1.7812 | 2.3 | 69.42
|
||||
Neigh | 0.10017 | 0.10116 | 0.10315 | 0.4 | 4.01
|
||||
Comm | 0.014848 | 0.044256 | 0.087941 | 14.7 | 1.75
|
||||
Output | 5.1199e-05 | 5.663e-05 | 7.1837e-05 | 0.0 | 0.00
|
||||
Modify | 0.62379 | 0.62575 | 0.62671 | 0.1 | 24.80
|
||||
Other | | 0.000504 | | | 0.02
|
||||
|
||||
Nlocal: 750 ave 759 max 735 min
|
||||
Histogram: 1 0 0 0 0 1 0 0 0 2
|
||||
Nghost: 6230.5 ave 6256 max 6190 min
|
||||
Histogram: 1 0 0 0 0 1 0 0 1 1
|
||||
Neighs: 277008 ave 280943 max 271394 min
|
||||
Histogram: 1 0 0 0 0 1 0 0 1 1
|
||||
|
||||
Total # of neighbors = 1108032
|
||||
Ave neighs/atom = 369.344
|
||||
Neighbor list builds = 2
|
||||
Dangerous builds = 0
|
||||
Total wall time: 0:00:02
|
||||
2
src/.gitignore
vendored
2
src/.gitignore
vendored
@ -998,6 +998,8 @@
|
||||
/fix_qeq_fire.h
|
||||
/fix_qeq_reaxff.cpp
|
||||
/fix_qeq_reaxff.h
|
||||
/fix_qeq_rel_reaxff.cpp
|
||||
/fix_qeq_rel_reaxff.h
|
||||
/fix_qmmm.cpp
|
||||
/fix_qmmm.h
|
||||
/fix_qtpie_reaxff.cpp
|
||||
|
||||
@ -106,13 +106,15 @@ void PairReaxFFOMP::init_style()
|
||||
|
||||
auto acks2_fixes = modify->get_fix_by_style("^acks2/reax");
|
||||
int have_qeq = modify->get_fix_by_style("^qeq/reax").size()
|
||||
+ modify->get_fix_by_style("^qeq/rel/reax").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))
|
||||
error->all(FLERR,"Pair style reaxff/omp requires use of exactly one of the "
|
||||
"fix qeq/reaxff or fix qeq/shielded or fix acks2/reaxff commands");
|
||||
error->all(FLERR,"Pair style reaxff requires use of exactly one of the "
|
||||
"fix qeq/reaxff or fix qeq/shielded or fix acks2/reaxff or "
|
||||
"fix qtpie/reaxff or fix qeq/rel/reaxff commands");
|
||||
|
||||
api->system->acks2_flag = acks2_fixes.size();
|
||||
if (api->system->acks2_flag)
|
||||
|
||||
109
src/REAXFF/fix_qeq_rel_reaxff.cpp
Normal file
109
src/REAXFF/fix_qeq_rel_reaxff.cpp
Normal file
@ -0,0 +1,109 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
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.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing authors:
|
||||
Navraj S Lalli, Imperial College London (navrajsinghlalli@gmail.com)
|
||||
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include "fix_qeq_rel_reaxff.h"
|
||||
|
||||
#include "atom.h"
|
||||
#include "comm.h"
|
||||
#include "error.h"
|
||||
#include "fix_efield.h"
|
||||
#include "force.h"
|
||||
#include "neighbor.h"
|
||||
|
||||
#include "pair_reaxff.h"
|
||||
#include "reaxff_api.h"
|
||||
|
||||
#include <cmath>
|
||||
#include <cstring>
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
using namespace FixConst;
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
FixQEqRelReaxFF::FixQEqRelReaxFF(LAMMPS *lmp, int narg, char **arg) : FixQtpieReaxFF(lmp, narg, arg)
|
||||
{
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixQEqRelReaxFF::calc_chi_eff()
|
||||
{
|
||||
memset(&chi_eff[0], 0, atom->nmax * sizeof(double));
|
||||
|
||||
const auto x = (const double *const *) atom->x;
|
||||
const int *type = atom->type;
|
||||
|
||||
double dx, dy, dz, dist_sq, overlap, sum_n, sum_d, chia, phia, phib;
|
||||
int i, j;
|
||||
|
||||
// check ghost atoms are stored up to the distance cutoff for overlap integrals
|
||||
const double comm_cutoff = MAX(neighbor->cutneighmax, comm->cutghostuser);
|
||||
if (comm_cutoff*comm_cutoff < dist_cutoff_sq) {
|
||||
error->all(FLERR, Error::NOLASTLINE,
|
||||
"Comm cutoff {} is smaller than distance cutoff {} for overlap integrals in fix {}. "
|
||||
"Increase accordingly using comm_modify cutoff",
|
||||
comm_cutoff, sqrt(dist_cutoff_sq), style);
|
||||
}
|
||||
|
||||
// efield energy is in real units of kcal/mol, factor needed for conversion to eV
|
||||
const double qe2f = force->qe2f;
|
||||
const double factor = 1.0 / qe2f;
|
||||
|
||||
if (efield) {
|
||||
if (efield->varflag != FixEfield::CONSTANT) efield->update_efield_variables();
|
||||
|
||||
// compute chi_eff for each local atom
|
||||
for (i = 0; i < nn; i++) {
|
||||
chia = chi[type[i]];
|
||||
if (efield->varflag != FixEfield::ATOM) {
|
||||
phia = -factor * (x[i][0] * efield->ex + x[i][1] * efield->ey + x[i][2] * efield->ez);
|
||||
} else { // atom-style potential from FixEfield
|
||||
phia = efield->efield[i][3];
|
||||
}
|
||||
|
||||
sum_n = 0.0;
|
||||
sum_d = 0.0;
|
||||
|
||||
for (j = 0; j < nt; j++) {
|
||||
dx = x[i][0] - x[j][0];
|
||||
dy = x[i][1] - x[j][1];
|
||||
dz = x[i][2] - x[j][2];
|
||||
dist_sq = (dx*dx + dy*dy + dz*dz);
|
||||
|
||||
if (dist_sq < dist_cutoff_sq) {
|
||||
|
||||
// overlap integral of two normalised 1s Gaussian type orbitals
|
||||
overlap = prefactor[type[i]][type[j]] * exp(-expfactor[type[i]][type[j]] * dist_sq);
|
||||
|
||||
if (efield->varflag != FixEfield::ATOM) {
|
||||
phib = -factor * (x[j][0] * efield->ex + x[j][1] * efield->ey + x[j][2] * efield->ez);
|
||||
} else { // atom-style potential from FixEfield
|
||||
phib = efield->efield[j][3];
|
||||
}
|
||||
sum_n += phib * overlap;
|
||||
sum_d += overlap;
|
||||
}
|
||||
}
|
||||
chi_eff[i] = chia + scale * (phia - sum_n / sum_d);
|
||||
}
|
||||
} else {
|
||||
for (i = 0; i < nn; i++) { chi_eff[i] = chi[type[i]]; }
|
||||
}
|
||||
}
|
||||
38
src/REAXFF/fix_qeq_rel_reaxff.h
Normal file
38
src/REAXFF/fix_qeq_rel_reaxff.h
Normal file
@ -0,0 +1,38 @@
|
||||
/* -*- 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(qeq/rel/reaxff,FixQEqRelReaxFF);
|
||||
// clang-format on
|
||||
#else
|
||||
|
||||
#ifndef LMP_FIX_QEQ_REL_REAXFF_H
|
||||
#define LMP_FIX_QEQ_REL_REAXFF_H
|
||||
|
||||
#include "fix_qtpie_reaxff.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class FixQEqRelReaxFF : public FixQtpieReaxFF {
|
||||
public:
|
||||
FixQEqRelReaxFF(class LAMMPS *, int, char **);
|
||||
|
||||
protected:
|
||||
void calc_chi_eff() override;
|
||||
};
|
||||
|
||||
} // namespace LAMMPS_NS
|
||||
|
||||
#endif
|
||||
#endif
|
||||
@ -16,7 +16,7 @@
|
||||
Contributing authors:
|
||||
Efstratios M Kritikos, California Institute of Technology
|
||||
(Implemented original version in LAMMMPS Aug 2019)
|
||||
Navraj S Lalli, Imperial College London
|
||||
Navraj S Lalli, Imperial College London (navrajsinghlalli@gmail.com)
|
||||
(Reimplemented QTPIE as a new fix in LAMMPS Aug 2024 and extended functionality)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
@ -51,7 +51,7 @@ using namespace FixConst;
|
||||
|
||||
static constexpr double CONV_TO_EV = 14.4;
|
||||
static constexpr double QSUMSMALL = 0.00001;
|
||||
static constexpr double ANGSTROM_TO_BOHRRADIUS = 1.8897261259;
|
||||
static constexpr double ANGSTROM_TO_BOHRRADIUS_SQ = 3.571064831;
|
||||
|
||||
static const char cite_fix_qtpie_reax[] =
|
||||
"fix qtpie/reaxff command: \n\n"
|
||||
@ -81,8 +81,9 @@ FixQtpieReaxFF::FixQtpieReaxFF(LAMMPS *lmp, int narg, char **arg) :
|
||||
|
||||
imax = 200;
|
||||
maxwarn = 1;
|
||||
scale = 1.0;
|
||||
|
||||
if ((narg < 9) || (narg > 12)) error->all(FLERR,"Illegal fix {} command", style);
|
||||
if ((narg < 9) || (narg > 14)) error->all(FLERR,"Illegal fix {} command", style);
|
||||
|
||||
nevery = utils::inumeric(FLERR,arg[3],false,lmp);
|
||||
if (nevery <= 0) error->all(FLERR,"Illegal fix {} command", style);
|
||||
@ -101,10 +102,17 @@ FixQtpieReaxFF::FixQtpieReaxFF(LAMMPS *lmp, int narg, char **arg) :
|
||||
error->all(FLERR,"Illegal fix {} command", style);
|
||||
imax = utils::numeric(FLERR,arg[iarg+1],false,lmp);
|
||||
iarg++;
|
||||
} else if (strcmp(arg[iarg],"scale") == 0) {
|
||||
if (iarg+1 > narg-1)
|
||||
error->all(FLERR,"Illegal fix {} command", style);
|
||||
scale = utils::numeric(FLERR,arg[iarg+1],false,lmp);
|
||||
iarg++;
|
||||
} else error->all(FLERR,"Illegal fix {} command", style);
|
||||
iarg++;
|
||||
}
|
||||
shld = nullptr;
|
||||
prefactor = nullptr;
|
||||
expfactor = nullptr;
|
||||
|
||||
nn = nt = n_cap = 0;
|
||||
nmax = 0;
|
||||
@ -166,8 +174,10 @@ FixQtpieReaxFF::~FixQtpieReaxFF()
|
||||
FixQtpieReaxFF::deallocate_matrix();
|
||||
|
||||
memory->destroy(shld);
|
||||
|
||||
memory->destroy(gauss_exp);
|
||||
memory->destroy(prefactor);
|
||||
memory->destroy(expfactor);
|
||||
|
||||
if (!reaxflag) {
|
||||
memory->destroy(chi);
|
||||
memory->destroy(eta);
|
||||
@ -179,7 +189,8 @@ FixQtpieReaxFF::~FixQtpieReaxFF()
|
||||
|
||||
void FixQtpieReaxFF::post_constructor()
|
||||
{
|
||||
if (lmp->citeme) lmp->citeme->add(cite_fix_qtpie_reax);
|
||||
if (utils::strmatch(style,"^qtpie/reax"))
|
||||
if (lmp->citeme) lmp->citeme->add(cite_fix_qtpie_reax);
|
||||
|
||||
grow_arrays(atom->nmax);
|
||||
for (int i = 0; i < atom->nmax; i++)
|
||||
@ -237,7 +248,7 @@ void FixQtpieReaxFF::pertype_parameters(char *arg)
|
||||
if (exp < 0)
|
||||
throw TokenizerException("Fix qtpie/reaxff: Invalid orbital exponent in gauss file",
|
||||
std::to_string(exp));
|
||||
gauss_exp[itype] = exp;
|
||||
gauss_exp[itype] = exp * ANGSTROM_TO_BOHRRADIUS_SQ;
|
||||
}
|
||||
fclose(fp);
|
||||
} catch (std::exception &e) {
|
||||
@ -247,11 +258,11 @@ void FixQtpieReaxFF::pertype_parameters(char *arg)
|
||||
|
||||
MPI_Bcast(gauss_exp,ntypes+1,MPI_DOUBLE,0,world);
|
||||
|
||||
// define a cutoff distance (in atomic units) beyond which overlap integrals are neglected
|
||||
// in calc_chi_eff()
|
||||
const double exp_min = find_min_exp(gauss_exp,ntypes+1);
|
||||
const int olap_cut = 10; // overlap integrals are neglected if less than pow(10,-olap_cut)
|
||||
dist_cutoff = sqrt(2*olap_cut/exp_min*log(10.0));
|
||||
// calculate a cutoff distance to neglect overlap integrals in calc_chi_eff()
|
||||
// when less than pow(10, -olap_cut)
|
||||
const double exp_min = find_min_exp(gauss_exp, ntypes+1);
|
||||
const int olap_cut = 10;
|
||||
dist_cutoff_sq = 2 * olap_cut * log(10.0) / exp_min;
|
||||
|
||||
// read chi, eta and gamma
|
||||
|
||||
@ -481,8 +492,13 @@ void FixQtpieReaxFF::init()
|
||||
if (efield->varflag == FixEfield::ATOM && efield->pstyle != FixEfield::ATOM)
|
||||
error->all(FLERR,"Atom-style external electric field requires atom-style "
|
||||
"potential variable when used with fix {}", style);
|
||||
} else {
|
||||
if (utils::strmatch(style,"^qeqr/reax") && comm->me == 0)
|
||||
error->warning(FLERR, "Use fix qeq/reaxff instead of fix {} when not using fix efield\n",
|
||||
style);
|
||||
}
|
||||
|
||||
|
||||
// we need a half neighbor list w/ Newton off
|
||||
// built whenever re-neighboring occurs
|
||||
|
||||
@ -490,6 +506,7 @@ void FixQtpieReaxFF::init()
|
||||
|
||||
init_shielding();
|
||||
init_taper();
|
||||
init_olap();
|
||||
|
||||
if (utils::strmatch(update->integrate_style,"^respa"))
|
||||
nlevels_respa = (dynamic_cast<Respa *>(update->integrate))->nlevels;
|
||||
@ -557,6 +574,31 @@ void FixQtpieReaxFF::init_taper()
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixQtpieReaxFF::init_olap()
|
||||
{
|
||||
int i,j;
|
||||
int ntypes;
|
||||
double expa,expb,expsum,expnorm;
|
||||
|
||||
ntypes = atom->ntypes;
|
||||
if (prefactor == nullptr)
|
||||
memory->create(prefactor,ntypes+1,ntypes+1,"qtpie:overlap_prefactor");
|
||||
if (expfactor == nullptr)
|
||||
memory->create(expfactor,ntypes+1,ntypes+1,"qtpie:overlap_expfactor");
|
||||
|
||||
for (i = 1; i <= ntypes; ++i)
|
||||
for (j = 1; j <= ntypes; ++j) {
|
||||
expa = gauss_exp[i];
|
||||
expb = gauss_exp[j];
|
||||
expsum = expa + expb;
|
||||
expnorm = expa * expb / expsum;
|
||||
prefactor[i][j] = pow((4.0 * expnorm / expsum), 0.75);
|
||||
expfactor[i][j] = expnorm;
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixQtpieReaxFF::setup_pre_force(int vflag)
|
||||
{
|
||||
if (reaxff) {
|
||||
@ -1119,15 +1161,16 @@ void FixQtpieReaxFF::calc_chi_eff()
|
||||
const auto x = (const double * const *)atom->x;
|
||||
const int *type = atom->type;
|
||||
|
||||
double dist,overlap,sum_n,sum_d,expa,expb,chia,chib,phia,phib,p,m;
|
||||
double dx,dy,dz,dist_sq,overlap,sum_n,sum_d,chia,chib,phia,phib;
|
||||
int i,j;
|
||||
|
||||
// check ghost atoms are stored up to the distance cutoff for overlap integrals
|
||||
const double comm_cutoff = MAX(neighbor->cutneighmax,comm->cutghostuser);
|
||||
if(comm_cutoff < dist_cutoff/ANGSTROM_TO_BOHRRADIUS) {
|
||||
error->all(FLERR,"comm cutoff = {} Angstrom is smaller than distance cutoff = {} Angstrom "
|
||||
"for overlap integrals in {}. Increase comm cutoff with comm_modify",
|
||||
comm_cutoff, dist_cutoff/ANGSTROM_TO_BOHRRADIUS, style);
|
||||
if(comm_cutoff*comm_cutoff < dist_cutoff_sq) {
|
||||
error->all(FLERR, Error::NOLASTLINE,
|
||||
"Comm cutoff {} is smaller than distance cutoff {} for overlap integrals in fix {}. "
|
||||
"Increase accordingly using comm_modify cutoff",
|
||||
comm_cutoff, sqrt(dist_cutoff_sq), style);
|
||||
}
|
||||
|
||||
// efield energy is in real units of kcal/mol, factor needed for conversion to eV
|
||||
@ -1141,7 +1184,6 @@ void FixQtpieReaxFF::calc_chi_eff()
|
||||
|
||||
// compute chi_eff for each local atom
|
||||
for (i = 0; i < nn; i++) {
|
||||
expa = gauss_exp[type[i]];
|
||||
chia = chi[type[i]];
|
||||
if (efield) {
|
||||
if (efield->varflag != FixEfield::ATOM) {
|
||||
@ -1155,16 +1197,16 @@ void FixQtpieReaxFF::calc_chi_eff()
|
||||
sum_d = 0.0;
|
||||
|
||||
for (j = 0; j < nt; j++) {
|
||||
dist = distance(x[i],x[j])*ANGSTROM_TO_BOHRRADIUS; // in atomic units
|
||||
dx = x[i][0] - x[j][0];
|
||||
dy = x[i][1] - x[j][1];
|
||||
dz = x[i][2] - x[j][2];
|
||||
dist_sq = (dx*dx + dy*dy + dz*dz);
|
||||
|
||||
if (dist < dist_cutoff) {
|
||||
expb = gauss_exp[type[j]];
|
||||
if (dist_sq < dist_cutoff_sq) {
|
||||
chib = chi[type[j]];
|
||||
|
||||
// overlap integral of two normalised 1s Gaussian type orbitals
|
||||
p = expa + expb;
|
||||
m = expa * expb / p;
|
||||
overlap = pow((4.0*m/p),0.75) * exp(-m*dist*dist);
|
||||
overlap = prefactor[type[i]][type[j]] * exp(-expfactor[type[i]][type[j]] * dist_sq);
|
||||
|
||||
if (efield) {
|
||||
if (efield->varflag != FixEfield::ATOM) {
|
||||
@ -1172,7 +1214,7 @@ void FixQtpieReaxFF::calc_chi_eff()
|
||||
} else { // atom-style potential from FixEfield
|
||||
phib = efield->efield[j][3];
|
||||
}
|
||||
sum_n += (chia - chib + phia - phib) * overlap;
|
||||
sum_n += (chia - chib + scale * (phia - phib)) * overlap;
|
||||
} else {
|
||||
sum_n += (chia - chib) * overlap;
|
||||
}
|
||||
@ -1197,14 +1239,3 @@ double FixQtpieReaxFF::find_min_exp(const double *array, const int array_length)
|
||||
}
|
||||
return exp_min;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
double FixQtpieReaxFF::distance(const double *posa, const double *posb)
|
||||
{
|
||||
double dx, dy, dz;
|
||||
dx = posb[0] - posa[0];
|
||||
dy = posb[1] - posa[1];
|
||||
dz = posb[2] - posa[2];
|
||||
return sqrt(dx*dx + dy*dy + dz*dz);
|
||||
}
|
||||
|
||||
@ -90,11 +90,15 @@ class FixQtpieReaxFF : public Fix {
|
||||
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
|
||||
double **prefactor; // factor used in computation of overlap integrals
|
||||
double **expfactor; // factor used in exponential term of overlap integrals
|
||||
double dist_cutoff_sq; // separation distance squared beyond which overlap integrals are neglected
|
||||
double scale; // scaling factor for electric polarization effects
|
||||
|
||||
void pertype_parameters(char *);
|
||||
void init_shielding();
|
||||
void init_taper();
|
||||
void init_olap();
|
||||
void allocate_storage();
|
||||
void deallocate_storage();
|
||||
void reallocate_storage();
|
||||
@ -128,9 +132,8 @@ class FixQtpieReaxFF : public Fix {
|
||||
void vector_sum(double *, double, double *, double, double *, int);
|
||||
void vector_add(double *, double, double *, int);
|
||||
|
||||
void calc_chi_eff();
|
||||
virtual 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
|
||||
};
|
||||
|
||||
@ -342,12 +342,13 @@ void PairReaxFF::init_style()
|
||||
auto acks2_fixes = modify->get_fix_by_style("^acks2/reax");
|
||||
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/rel/reax").size()
|
||||
+ modify->get_fix_by_style("^qtpie/reax").size();
|
||||
|
||||
if (qeqflag && (have_qeq != 1))
|
||||
error->all(FLERR,"Pair style reaxff requires use of exactly one of the "
|
||||
"fix qeq/reaxff or fix qeq/shielded or fix acks2/reaxff or "
|
||||
"fix qtpie/reaxff commands");
|
||||
"fix qtpie/reaxff or fix qeq/rel/reaxff commands");
|
||||
|
||||
api->system->acks2_flag = acks2_fixes.size();
|
||||
if (api->system->acks2_flag)
|
||||
|
||||
@ -26,6 +26,7 @@ namespace LAMMPS_NS {
|
||||
|
||||
class FixEfield : public Fix {
|
||||
friend class FixQEqReaxFF;
|
||||
friend class FixQEqRelReaxFF;
|
||||
friend class FixQtpieReaxFF;
|
||||
|
||||
public:
|
||||
|
||||
Reference in New Issue
Block a user