git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@14913 f3b2605a-c512-4ea7-a41b-209d697bcdaa

This commit is contained in:
sjplimp
2016-04-30 18:08:26 +00:00
parent 6ff2012d88
commit 114742bc0d
9 changed files with 330 additions and 156 deletions

View File

@ -493,7 +493,7 @@ name links to a sub-section below with more details.</p>
<td>python</td>
<td>lib/python</td>
</tr>
<tr class="row-even"><td><a class="reference internal" href="#reax"><span>REAX</span></a></td>
<tr class="row-even"><td><a class="reference internal" href="Section_tools.html#reax"><span>REAX</span></a></td>
<td>ReaxFF potential</td>
<td>Aidan Thompson (Sandia)</td>
<td><a class="reference internal" href="pair_reax.html"><em>pair_style reax</em></a></td>

View File

@ -5,7 +5,7 @@
LAMMPS Documentation
====================
27 Apr 2016 version
28 Apr 2016 version
-------------------
Version info:

View File

@ -53,64 +53,70 @@ Description
Calculate the dynamical matrix from molecular dynamics simulations
based on fluctuation-dissipation theory for a group of atoms.
Consider a crystal with *N* unit cells in three dimensions labelled *l
= (l<sub>1</sub>,l<sub>2</sub>,l<sub>3</sub>)* where *l<sub>i</sub>*
Consider a crystal with :math:`N` unit cells in three dimensions labelled :math:`l = (l_1, l_2, l_3)` where :math:`l_i`
are integers. Each unit cell is defined by three linearly independent
vectors **a**<sub>1</sub>, **a**<sub>2</sub>, **a**<sub>3</sub> forming a
parallelipiped, containing *K* basis atoms labelled *k*.
vectors :math:`\mathbf{a}_1`, :math:`\mathbf{a}_2`, :math:`\mathbf{a}_3` forming a
parallelipiped, containing :math:`K` basis atoms labeled :math:`k`.
Based on fluctuation-dissipation theory, the force constant
coefficients of the system in reciprocal space are given by
(:ref:`Campa&ntilde;&aacute; <Campana>` , :ref:`Kong <Kong>`)
<center><b>&Phi;</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>) =
k<sub>B</sub>T
<b>G</b><sup>-1</sup><sub>k&alpha;,k'&beta;</sub>(<b>q</b>),</center>
(:ref:`Campana <Campana>` , :ref:`Kong <Kong>`)
where **G** is the Green's functions coefficients given by
.. raw:: html
.. math::
<center><b>G</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>) =
<<b>u</b><sub>k&alpha;</sub>(<b>q</b>)&#149;<b>u</b><sub>k'&beta;</sub><sup>*</sup>(<b>q</b>)>,</center>
\begin{equation}\mathbf{\Phi}_{k\alpha,k^\prime \beta}(\mathbf{q}) = k_B T \mathbf{G}^{-1}_{k\alpha,k^\prime \beta}(\mathbf{q})\end{equation}
where <...> denotes the ensemble average, and
<center>**u**<sub>k&alpha;</sub>(<b>q</b>) = &sum;<sub>l</sub>
<b>u</b><sub>lk&alpha;</sub> exp(i**qr**<sub>l</sub>)</center>
where :math:`\mathbf{G}` is the Green's functions coefficients given by
is the &alpha; component of the atomic displacement for the *k*th atom
in the unit cell in reciprocal space at **q**. In practice, the Green's
.. math::
\begin{equation}\mathbf{G}_{k\alpha,k^\prime \beta}(\mathbf{q}) = \left< \mathbf{u}_{k\alpha}(\mathbf{q}) \bullet \mathbf{u}_{k^\prime \beta}^*(\mathbf{q}) \right>\end{equation}
where :math:`\left< \ldots \right>` denotes the ensemble average, and
.. math::
\begin{equation}\mathbf{u}_{k\alpha}(\mathbf{q}) = \sum_l \mathbf{u}_{l k \alpha} \exp{(i\mathbf{qr}_l)}\end{equation}
is the :math:`\alpha` component of the atomic displacement for the :math:`k` th atom
in the unit cell in reciprocal space at :math:`\mathbf{q}`. In practice, the Green's
functions coefficients can also be measured according to the following
formula,
.. raw:: html
<center><b>G</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>) =
<<b>R</b><sub>k&alpha;</sub>(<b>q</b>)&#149;<b>R</b><sup>*</sup><sub>k'&beta;</sub>(<b>q</b>)>
- <<b>R</b>><sub>k&alpha;</sub>(<b>q</b>)&#149;<<b>R</b>><sup>*</sup><sub>k'&beta;</sub>(<b>q</b>),
</center>
.. math::
where **R** is the instantaneous positions of atoms, and <**R**> is the
\begin{equation}\mathbf{G}_{k\alpha,k^\prime \beta}(\mathbf{q}) =
\left< \mathbf{R}_{k \alpha}(\mathbf{q}) \bullet \mathbf{R}^*_{k^\prime \beta}(\mathbf{q}) \right>
- \left<\mathbf{R}\right>_{k \alpha}(\mathbf{q}) \bullet \left<\mathbf{R}\right>^*_{k^\prime \beta}(\mathbf{q})\end{equation}
where :math:`\mathbf{R}` is the instantaneous positions of atoms, and :math:`\left<\mathbf{R}\right>` is the
averaged atomic positions. It gives essentially the same results as
the displacement method and is easier to implement in an MD code.
Once the force constant matrix is known, the dynamical matrix **D** can
Once the force constant matrix is known, the dynamical matrix :math:`\mathbf{D}` can
then be obtained by
.. raw:: html
<center><b>D</b><sub>k&alpha;, k'&beta;</sub>(<b>q</b>) = (m<sub>k</sub>m<sub>k'</sub>)<sup>-1/2</sup> <b>&Phi;</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>)</center>
.. math::
whose eigenvalues are exactly the phonon frequencies at **q**.
\begin{equation}\mathbf{D}_{k\alpha, k^\prime\beta}(\mathbf{q}) =
(m_k m_{k^\prime})^{-\frac{1}{2}} \mathbf{\Phi}_{k \alpha, k^\prime \beta}(\mathbf{q})\end{equation}
whose eigenvalues are exactly the phonon frequencies at :math:`\mathbf{q}`.
This fix uses positions of atoms in the specified group and calculates
two-point correlations. To achieve this. the positions of the atoms
are examined every *Nevery* steps and are Fourier-transformed into
reciprocal space, where the averaging process and correlation
computation is then done. After every *Noutput* measurements, the
matrix **G**(**q**) is calculated and inverted to obtain the elastic
matrix :math:`\mathbf{G}(\mathbf{q})` is calculated and inverted to obtain the elastic
stiffness coefficients. The dynamical matrices are then constructed
and written to *prefix*.bin.timestep files in binary format and to the
file *prefix*.log for each wavevector **q**.
file *prefix*.log for each wavevector :math:`\mathbf{q}`.
A detailed description of this method can be found in
(:ref:`Kong2011 <Kong2011>`).
@ -122,7 +128,7 @@ model a 2D or 3D system, the phonon dispersion of a 1D atomic chain
can be computed using *sysdim* = 1.
The *nasr* keyword is optional. An iterative procedure is employed to
enforce the acoustic sum rule on &Phi; at &Gamma;, and the number
enforce the acoustic sum rule on :math:`\Phi` at :math:`\Gamma`, and the number
provided by keyword *nasr* gives the total number of iterations. For a
system whose unit cell has only one atom, *nasr* = 1 is sufficient;
for other systems, *nasr* = 10 is typically sufficient.
@ -202,19 +208,19 @@ the `dimension <dimension>`_ command, and nasr = 20.
**(Campa&ntilde;&aacute;)** C. Campa&ntilde;&aacute; and
M. H. M&uuml;ser, *Practical Green's function approach to the
**(Campana)** C. Campana and
M. H. Muser, *Practical Green's function approach to the
simulation of elastic semi-infinite solids*, `Phys. Rev. B [74], 075420 (2006) <http://dx.doi.org/10.1103/PhysRevB.74.075420>`_
.. _Kong:
**(Kong)** L.T. Kong, G. Bartels, C. Campa&ntilde;&aacute;,
C. Denniston, and Martin H. M&uuml;ser, *Implementation of Green's
**(Kong)** L.T. Kong, G. Bartels, C. Campana,
C. Denniston, and Martin H. Muser, *Implementation of Green's
function molecular dynamics: An extension to LAMMPS*, `Computer Physics Communications [180](6):1004-1010 (2009). <http://dx.doi.org/10.1016/j.cpc.2008.12.035>`_
L.T. Kong, C. Denniston, and Martin H. M&uuml;ser,
L.T. Kong, C. Denniston, and Martin H. Muser,
*An improved version of the Green's function molecular dynamics
method*, `Computer Physics Communications [182](2):540-541 (2011). <http://dx.doi.org/10.1016/j.cpc.2010.10.006>`_

View File

@ -12,11 +12,12 @@ Syntax
* ID, group-ID are documented in :doc:`fix <fix>` command
* property/atom = style name of this fix command
* vec1,vec2,... = *mol* or *q* or *i_name* or *d_name*
* vec1,vec2,... = *mol* or *q* or *rmass* or *i_name* or *d_name*
.. parsed-literal::
*mol* = molecule IDs
*q* = charge
*rmass* = per-atom mass
*i_name* = new integer vector referenced by name
*d_name* = new floating-point vector referenced by name
@ -49,25 +50,36 @@ properties, as explained on the :doc:`atom_style <atom_style>` and
:doc:`read_data <read_data>` doc pages. The latter command allows these
properties to be defined for each atom in the system when a data file
is read. This fix will augment the set of properties with new custom
ones.
ones. This can be useful in several scenarios.
This can be useful in at least two scenarios.
If the atom style does not define molecule IDs, per-atom charge,
or per-atom mass, they can be added using the *mol*, *q* or *rmass*
keywords. This can be useful, e.g, to define "molecules" to use as
rigid bodies with the :doc:`fix rigid <fix_rigid>` command, or just to
carry around an extra flag with the atoms (stored as a molecule ID)
that can be used to group atoms without having to use the group
command (which is limited to a total of 32 groups including *all*).
If the atom style does not define molecule IDs or per-atom charge,
they can be added using the *mol* or *q* keywords. This can be
useful, e.g, to define "molecules" to use as rigid bodies with the
:doc:`fix rigid <fix_rigid>` command, or just to carry around an extra
flag with the atoms (stored as a molecule ID). An alternative is to
use an atom style that does define molecule IDs or charge or to use a
hybrid atom style that combines two styles to allow for molecule IDs
or charge, but that has 2 practical drawbacks. First it typically
necessitates changing the format of the data file. And it may define
additional properties that aren't needed such as bond lists, which has
some overhead when there are no bonds.
Another application would be to use the *rmass* flag in order to have
per-atom masses instead of per-type masses, for example this can be
useful to study isotope effects with partial isotope substitution.
Please :ref:`see below <isotopes>` for an example of simulating a mixture
of light and heavy water with the TIP4P water potential.
An alternative to using fix *property/atom* in these ways is to
use an atom style that does define molecule IDs or charge or per-atom
mass (indirectly via diameter and density) or to use a hybrid atom
style that combines two or more atom styles
to provide the union of all atom properties. However, this has two
practical drawbacks: first, it typically necessitates changing the
format of the data file, which can be tedious for large systems;
and second, it may define additional properties that are not needed
such as bond lists, which has some overhead when there are no bonds.
In the future, we may add additional per-atom properties similar to
*mol* or *q*, which "turn-on" specific properties defined by some atom
styles, so they can be used by atom styles that don't define them.
*mol*, *q* or *rmass*, which "turn-on" specific properties defined
by some atom styles, so they can be used by atom styles that do not
define them.
More generally, the *i_name* and *d_name* vectors allow one or more
new custom per-atom properties to be defined. Each name must be
@ -92,9 +104,9 @@ new properties are also defined for the ghost atoms.
.. note::
If you use this command with the *mol* or *charge* vectors than
you most likely want to set *ghost* yes, since these properties are
stored with ghost atoms if you use an :doc:`atom_style <atom_style>`
If you use this command with the *mol*, *q* or *rmass* vectors,
then you most likely want to set *ghost* yes, since these properties
are stored with ghost atoms if you use an :doc:`atom_style <atom_style>`
that defines them, and many LAMMPS operations that use molecule IDs or
charge, such as neighbor lists and pair styles, will expect ghost
atoms to have these valuse. LAMMPS will issue a warning it you define
@ -208,6 +220,54 @@ on how the properties can be accessed from added classes.
----------
.. _isotopes:
Example for using per-atom masses with TIP4P water to study isotope
effects. When setting up simulations with the :ref:`TIP4P pair styles <howto_8>` for water, you have to provide
exactly one atom type each to identify the water oxygen and hydrogen
atoms. Since the atom mass is normally tied to the atom type, this
makes it impossible to study multiple isotopes in the same simulation.
With *fix property/atom rmass* however, the per-type masses are
replaced by per-atom masses. Asumming you have a working input deck
for regular TIP4P water, where water oxygen is atom type 1 and
water hydrogen is atom type 2, the following lines of input script
convert this to using per-atom masses:
.. parsed-literal::
fix Isotopes all property/atom rmass ghost yes
set type 1 mass 15.9994
set type 2 mass 1.008
When writing out the system data with the :doc:`write_data <write_data>`
command, there will be a new section named with the fix-ID
(i.e. *Isotopes* in this case). Alternatively, you can take an
existing data file and just add this *Isotopes* section with
one line per atom containing atom-ID and mass. Either way, the
extended data file can be read back with:
.. parsed-literal::
fix Isotopes all property/atom rmass ghost yes
read_data tip4p-isotopes.data fix Isotopes NULL Isotopes
Please note that the first *Isotopes* refers to the fix-ID
and the second to the name of the section. The following input
script code will now change the first 100 water molecules in this
example to heavy water:
.. parsed-literal::
group hwat id 2:300:3
group hwat id 3:300:3
set group hwat mass 2.0141018
----------
Restart, fix_modify, output, run start/stop, minimize info
""""""""""""""""""""""""""""""""""""""""""""""""""""""""""

View File

@ -171,44 +171,47 @@ fix 1 all phonon 10 5000 500000 GAMMA EAM0D nasr 100
<h2>Description<a class="headerlink" href="#description" title="Permalink to this headline"></a></h2>
<p>Calculate the dynamical matrix from molecular dynamics simulations
based on fluctuation-dissipation theory for a group of atoms.</p>
<p>Consider a crystal with <em>N</em> unit cells in three dimensions labelled <em>l
= (l&lt;sub&gt;1&lt;/sub&gt;,l&lt;sub&gt;2&lt;/sub&gt;,l&lt;sub&gt;3&lt;/sub&gt;)</em> where <em>l&lt;sub&gt;i&lt;/sub&gt;</em>
<p>Consider a crystal with <span class="math">\(N\)</span> unit cells in three dimensions labelled <span class="math">\(l = (l_1, l_2, l_3)\)</span> where <span class="math">\(l_i\)</span>
are integers. Each unit cell is defined by three linearly independent
vectors <a href="#id1"><span class="problematic" id="id2">**</span></a>a**&lt;sub&gt;1&lt;/sub&gt;, <a href="#id3"><span class="problematic" id="id4">**</span></a>a**&lt;sub&gt;2&lt;/sub&gt;, <a href="#id5"><span class="problematic" id="id6">**</span></a>a**&lt;sub&gt;3&lt;/sub&gt; forming a
parallelipiped, containing <em>K</em> basis atoms labelled <em>k</em>.</p>
vectors <span class="math">\(\mathbf{a}_1\)</span>, <span class="math">\(\mathbf{a}_2\)</span>, <span class="math">\(\mathbf{a}_3\)</span> forming a
parallelipiped, containing <span class="math">\(K\)</span> basis atoms labeled <span class="math">\(k\)</span>.</p>
<p>Based on fluctuation-dissipation theory, the force constant
coefficients of the system in reciprocal space are given by
(<a class="reference internal" href="#campana"><span>Campa&amp;ntilde;&amp;aacute;</span></a> , <a class="reference internal" href="#kong"><span>Kong</span></a>)
&lt;center&gt;&lt;b&gt;&amp;Phi;&lt;/b&gt;&lt;sub&gt;k&amp;alpha;,k&#8217;&amp;beta;&lt;/sub&gt;(&lt;b&gt;q&lt;/b&gt;) =
k&lt;sub&gt;B&lt;/sub&gt;T
&lt;b&gt;G&lt;/b&gt;&lt;sup&gt;-1&lt;/sup&gt;&lt;sub&gt;k&amp;alpha;,k&#8217;&amp;beta;&lt;/sub&gt;(&lt;b&gt;q&lt;/b&gt;),&lt;/center&gt;</p>
<p>where <strong>G</strong> is the Green&#8217;s functions coefficients given by</p>
<center><b>G</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>) =
<<b>u</b><sub>k&alpha;</sub>(<b>q</b>)&#149;<b>u</b><sub>k'&beta;</sub><sup>*</sup>(<b>q</b>)>,</center><p>where &lt;...&gt; denotes the ensemble average, and
&lt;center&gt;**u**&lt;sub&gt;k&amp;alpha;&lt;/sub&gt;(&lt;b&gt;q&lt;/b&gt;) = &amp;sum;&lt;sub&gt;l&lt;/sub&gt;
&lt;b&gt;u&lt;/b&gt;&lt;sub&gt;lk&amp;alpha;&lt;/sub&gt; exp(i**qr**&lt;sub&gt;l&lt;/sub&gt;)&lt;/center&gt;</p>
<p>is the &amp;alpha; component of the atomic displacement for the <em>k*th atom
in the unit cell in reciprocal space at **q*</em>. In practice, the Green&#8217;s
(<a class="reference internal" href="#campana"><span>Campana</span></a> , <a class="reference internal" href="#kong"><span>Kong</span></a>)</p>
<div class="math">
\[\begin{equation}\mathbf{\Phi}_{k\alpha,k^\prime \beta}(\mathbf{q}) = k_B T \mathbf{G}^{-1}_{k\alpha,k^\prime \beta}(\mathbf{q})\end{equation}\]</div>
<p>where <span class="math">\(\mathbf{G}\)</span> is the Green&#8217;s functions coefficients given by</p>
<div class="math">
\[\begin{split}\begin{equation}\mathbf{G}_{k\alpha,k^\prime \beta}(\mathbf{q}) = \left&lt; \mathbf{u}_{k\alpha}(\mathbf{q}) \bullet \mathbf{u}_{k^\prime \beta}^*(\mathbf{q}) \right&gt;\end{equation}\end{split}\]</div>
<p>where <span class="math">\(\left&lt; \ldots \right&gt;\)</span> denotes the ensemble average, and</p>
<div class="math">
\[\begin{equation}\mathbf{u}_{k\alpha}(\mathbf{q}) = \sum_l \mathbf{u}_{l k \alpha} \exp{(i\mathbf{qr}_l)}\end{equation}\]</div>
<p>is the <span class="math">\(\alpha\)</span> component of the atomic displacement for the <span class="math">\(k\)</span> th atom
in the unit cell in reciprocal space at <span class="math">\(\mathbf{q}\)</span>. In practice, the Green&#8217;s
functions coefficients can also be measured according to the following
formula,</p>
<center><b>G</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>) =
<<b>R</b><sub>k&alpha;</sub>(<b>q</b>)&#149;<b>R</b><sup>*</sup><sub>k'&beta;</sub>(<b>q</b>)>
- <<b>R</b>><sub>k&alpha;</sub>(<b>q</b>)&#149;<<b>R</b>><sup>*</sup><sub>k'&beta;</sub>(<b>q</b>),
</center><p>where <strong>R</strong> is the instantaneous positions of atoms, and &lt;<strong>R</strong>&gt; is the
<div class="math">
\[\begin{split}\begin{equation}\mathbf{G}_{k\alpha,k^\prime \beta}(\mathbf{q}) =
\left&lt; \mathbf{R}_{k \alpha}(\mathbf{q}) \bullet \mathbf{R}^*_{k^\prime \beta}(\mathbf{q}) \right&gt;
- \left&lt;\mathbf{R}\right&gt;_{k \alpha}(\mathbf{q}) \bullet \left&lt;\mathbf{R}\right&gt;^*_{k^\prime \beta}(\mathbf{q})\end{equation}\end{split}\]</div>
<p>where <span class="math">\(\mathbf{R}\)</span> is the instantaneous positions of atoms, and <span class="math">\(\left&lt;\mathbf{R}\right&gt;\)</span> is the
averaged atomic positions. It gives essentially the same results as
the displacement method and is easier to implement in an MD code.</p>
<p>Once the force constant matrix is known, the dynamical matrix <strong>D</strong> can
<p>Once the force constant matrix is known, the dynamical matrix <span class="math">\(\mathbf{D}\)</span> can
then be obtained by</p>
<center><b>D</b><sub>k&alpha;, k'&beta;</sub>(<b>q</b>) = (m<sub>k</sub>m<sub>k'</sub>)<sup>-1/2</sup> <b>&Phi;</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>)</center><p>whose eigenvalues are exactly the phonon frequencies at <strong>q</strong>.</p>
<div class="math">
\[\begin{equation}\mathbf{D}_{k\alpha, k^\prime\beta}(\mathbf{q}) =
(m_k m_{k^\prime})^{-\frac{1}{2}} \mathbf{\Phi}_{k \alpha, k^\prime \beta}(\mathbf{q})\end{equation}\]</div>
<p>whose eigenvalues are exactly the phonon frequencies at <span class="math">\(\mathbf{q}\)</span>.</p>
<p>This fix uses positions of atoms in the specified group and calculates
two-point correlations. To achieve this. the positions of the atoms
are examined every <em>Nevery</em> steps and are Fourier-transformed into
reciprocal space, where the averaging process and correlation
computation is then done. After every <em>Noutput</em> measurements, the
matrix <strong>G**(**q</strong>) is calculated and inverted to obtain the elastic
matrix <span class="math">\(\mathbf{G}(\mathbf{q})\)</span> is calculated and inverted to obtain the elastic
stiffness coefficients. The dynamical matrices are then constructed
and written to <em>prefix</em>.bin.timestep files in binary format and to the
file <em>prefix</em>.log for each wavevector <strong>q</strong>.</p>
file <em>prefix</em>.log for each wavevector <span class="math">\(\mathbf{q}\)</span>.</p>
<p>A detailed description of this method can be found in
(<a class="reference internal" href="#kong2011"><span>Kong2011</span></a>).</p>
<p>The <em>sysdim</em> keyword is optional. If specified with a value smaller
@ -217,7 +220,7 @@ for the dynamical matrix calculation. For example, using LAMMPS ot
model a 2D or 3D system, the phonon dispersion of a 1D atomic chain
can be computed using <em>sysdim</em> = 1.</p>
<p>The <em>nasr</em> keyword is optional. An iterative procedure is employed to
enforce the acoustic sum rule on &amp;Phi; at &amp;Gamma;, and the number
enforce the acoustic sum rule on <span class="math">\(\Phi\)</span> at <span class="math">\(\Gamma\)</span>, and the number
provided by keyword <em>nasr</em> gives the total number of iterations. For a
system whose unit cell has only one atom, <em>nasr</em> = 1 is sufficient;
for other systems, <em>nasr</em> = 10 is typically sufficient.</p>
@ -274,13 +277,13 @@ LAMMPS was built with that package. See the <a class="reference internal" href=
<p>The option defaults are sysdim = the same dimemsion as specified by
the <a class="reference external" href="dimension">dimension</a> command, and nasr = 20.</p>
<hr class="docutils" />
<p id="campana"><strong>(Campa&amp;ntilde;&amp;aacute;)</strong> C. Campa&amp;ntilde;&amp;aacute; and
M. H. M&amp;uuml;ser, <em>Practical Green&#8217;s function approach to the
<p id="campana"><strong>(Campana)</strong> C. Campana and
M. H. Muser, <em>Practical Green&#8217;s function approach to the
simulation of elastic semi-infinite solids</em>, <a class="reference external" href="http://dx.doi.org/10.1103/PhysRevB.74.075420">Phys. Rev. B [74], 075420 (2006)</a></p>
<p id="kong"><strong>(Kong)</strong> L.T. Kong, G. Bartels, C. Campa&amp;ntilde;&amp;aacute;,
C. Denniston, and Martin H. M&amp;uuml;ser, <em>Implementation of Green&#8217;s
<p id="kong"><strong>(Kong)</strong> L.T. Kong, G. Bartels, C. Campana,
C. Denniston, and Martin H. Muser, <em>Implementation of Green&#8217;s
function molecular dynamics: An extension to LAMMPS</em>, <a class="reference external" href="http://dx.doi.org/10.1016/j.cpc.2008.12.035">Computer Physics Communications [180](6):1004-1010 (2009).</a></p>
<p>L.T. Kong, C. Denniston, and Martin H. M&amp;uuml;ser,
<p>L.T. Kong, C. Denniston, and Martin H. Muser,
<em>An improved version of the Green&#8217;s function molecular dynamics
method</em>, <a class="reference external" href="http://dx.doi.org/10.1016/j.cpc.2010.10.006">Computer Physics Communications [182](2):540-541 (2011).</a></p>
<p id="kong2011"><strong>(Kong2011)</strong> L.T. Kong, <em>Phonon dispersion measured directly from

View File

@ -43,58 +43,67 @@ fix 1 all phonon 10 5000 500000 GAMMA EAM0D nasr 100 :pre
Calculate the dynamical matrix from molecular dynamics simulations
based on fluctuation-dissipation theory for a group of atoms.
Consider a crystal with {N} unit cells in three dimensions labelled {l
= (l<sub>1</sub>,l<sub>2</sub>,l<sub>3</sub>)} where {l<sub>i</sub>}
Consider a crystal with \(N\) unit cells in three dimensions labelled \(l = (l_1, l_2, l_3)\) where \(l_i\)
are integers. Each unit cell is defined by three linearly independent
vectors [a]<sub>1</sub>, [a]<sub>2</sub>, [a]<sub>3</sub> forming a
parallelipiped, containing {K} basis atoms labelled {k}.
vectors \(\mathbf\{a\}_1\), \(\mathbf\{a\}_2\), \(\mathbf\{a\}_3\) forming a
parallelipiped, containing \(K\) basis atoms labeled \(k\).
Based on fluctuation-dissipation theory, the force constant
coefficients of the system in reciprocal space are given by
("Campa&ntilde;&aacute;"_#Campana , "Kong"_#Kong)
<center><b>&Phi;</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>) =
k<sub>B</sub>T
<b>G</b><sup>-1</sup><sub>k&alpha;,k'&beta;</sub>(<b>q</b>),</center>
("Campana"_#Campana , "Kong"_#Kong)
where [G] is the Green's functions coefficients given by
\begin\{equation\}
\mathbf\{\Phi\}_\{k\alpha,k^\prime \beta\}(\mathbf\{q\}) = k_B T \mathbf\{G\}^\{-1\}_\{k\alpha,k^\prime \beta\}(\mathbf\{q\})
\end\{equation\}
<center><b>G</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>) =
<<b>u</b><sub>k&alpha;</sub>(<b>q</b>)&#149;<b>u</b><sub>k'&beta;</sub><sup>*</sup>(<b>q</b>)>,</center>
where \(\mathbf\{G\}\) is the Green's functions coefficients given by
where <...> denotes the ensemble average, and
<center>[u]<sub>k&alpha;</sub>(<b>q</b>) = &sum;<sub>l</sub>
<b>u</b><sub>lk&alpha;</sub> exp(i[qr]<sub>l</sub>)</center>
\begin\{equation\}
\mathbf\{G\}_\{k\alpha,k^\prime \beta\}(\mathbf\{q\}) = \left< \mathbf\{u\}_\{k\alpha\}(\mathbf\{q\}) \bullet \mathbf\{u\}_\{k^\prime \beta\}^*(\mathbf\{q\}) \right>
\end\{equation\}
is the &alpha; component of the atomic displacement for the {k}th atom
in the unit cell in reciprocal space at [q]. In practice, the Green's
where \(\left< \ldots \right>\) denotes the ensemble average, and
\begin\{equation\}
\mathbf\{u\}_\{k\alpha\}(\mathbf\{q\}) = \sum_l \mathbf\{u\}_\{l k \alpha\} \exp\{(i\mathbf\{qr\}_l)\}
\end\{equation\}
is the \(\alpha\) component of the atomic displacement for the \(k\) th atom
in the unit cell in reciprocal space at \(\mathbf\{q\}\). In practice, the Green's
functions coefficients can also be measured according to the following
formula,
<center><b>G</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>) =
<<b>R</b><sub>k&alpha;</sub>(<b>q</b>)&#149;<b>R</b><sup>*</sup><sub>k'&beta;</sub>(<b>q</b>)>
- <<b>R</b>><sub>k&alpha;</sub>(<b>q</b>)&#149;<<b>R</b>><sup>*</sup><sub>k'&beta;</sub>(<b>q</b>),
</center>
\begin\{equation\}
\mathbf\{G\}_\{k\alpha,k^\prime \beta\}(\mathbf\{q\}) =
\left< \mathbf\{R\}_\{k \alpha\}(\mathbf\{q\}) \bullet \mathbf\{R\}^*_\{k^\prime \beta\}(\mathbf\{q\}) \right>
- \left<\mathbf\{R\}\right>_\{k \alpha\}(\mathbf\{q\}) \bullet \left<\mathbf\{R\}\right>^*_\{k^\prime \beta\}(\mathbf\{q\})
\end\{equation\}
where [R] is the instantaneous positions of atoms, and <[R]> is the
where \(\mathbf\{R\}\) is the instantaneous positions of atoms, and \(\left<\mathbf\{R\}\right>\) is the
averaged atomic positions. It gives essentially the same results as
the displacement method and is easier to implement in an MD code.
Once the force constant matrix is known, the dynamical matrix [D] can
Once the force constant matrix is known, the dynamical matrix \(\mathbf\{D\}\) can
then be obtained by
<center><b>D</b><sub>k&alpha;, k'&beta;</sub>(<b>q</b>) = (m<sub>k</sub>m<sub>k'</sub>)<sup>-1/2</sup> <b>&Phi;</b><sub>k&alpha;,k'&beta;</sub>(<b>q</b>)</center>
\begin\{equation\}
\mathbf\{D\}_\{k\alpha, k^\prime\beta\}(\mathbf\{q\}) =
(m_k m_\{k^\prime\})^\{-\frac\{1\}\{2\}\} \mathbf\{\Phi\}_\{k \alpha, k^\prime \beta\}(\mathbf\{q\})
\end\{equation\}
whose eigenvalues are exactly the phonon frequencies at [q].
whose eigenvalues are exactly the phonon frequencies at \(\mathbf\{q\}\).
This fix uses positions of atoms in the specified group and calculates
two-point correlations. To achieve this. the positions of the atoms
are examined every {Nevery} steps and are Fourier-transformed into
reciprocal space, where the averaging process and correlation
computation is then done. After every {Noutput} measurements, the
matrix [G]([q]) is calculated and inverted to obtain the elastic
matrix \(\mathbf\{G\}(\mathbf\{q\})\) is calculated and inverted to obtain the elastic
stiffness coefficients. The dynamical matrices are then constructed
and written to {prefix}.bin.timestep files in binary format and to the
file {prefix}.log for each wavevector [q].
file {prefix}.log for each wavevector \(\mathbf\{q\}\).
A detailed description of this method can be found in
("Kong2011"_#Kong2011).
@ -106,7 +115,7 @@ model a 2D or 3D system, the phonon dispersion of a 1D atomic chain
can be computed using {sysdim} = 1.
The {nasr} keyword is optional. An iterative procedure is employed to
enforce the acoustic sum rule on &Phi; at &Gamma;, and the number
enforce the acoustic sum rule on \(\Phi\) at \(\Gamma\), and the number
provided by keyword {nasr} gives the total number of iterations. For a
system whose unit cell has only one atom, {nasr} = 1 is sufficient;
for other systems, {nasr} = 10 is typically sufficient.
@ -178,19 +187,19 @@ the "dimension"_dimension command, and nasr = 20.
:line
:link(Campana)
[(Campa&ntilde;&aacute;)] C. Campa&ntilde;&aacute; and
M. H. M&uuml;ser, {Practical Green's function approach to the
[(Campana)] C. Campana and
M. H. Muser, {Practical Green's function approach to the
simulation of elastic semi-infinite solids}, "Phys. Rev. B \[74\],
075420 (2006)"_http://dx.doi.org/10.1103/PhysRevB.74.075420
:link(Kong)
[(Kong)] L.T. Kong, G. Bartels, C. Campa&ntilde;&aacute;,
C. Denniston, and Martin H. M&uuml;ser, {Implementation of Green's
[(Kong)] L.T. Kong, G. Bartels, C. Campana,
C. Denniston, and Martin H. Muser, {Implementation of Green's
function molecular dynamics: An extension to LAMMPS}, "Computer
Physics Communications \[180\](6):1004-1010
(2009)."_http://dx.doi.org/10.1016/j.cpc.2008.12.035
L.T. Kong, C. Denniston, and Martin H. M&uuml;ser,
L.T. Kong, C. Denniston, and Martin H. Muser,
{An improved version of the Green's function molecular dynamics
method}, "Computer Physics Communications \[182\](2):540-541
(2011)."_http://dx.doi.org/10.1016/j.cpc.2010.10.006

View File

@ -134,11 +134,12 @@
<ul class="simple">
<li>ID, group-ID are documented in <a class="reference internal" href="fix.html"><em>fix</em></a> command</li>
<li>property/atom = style name of this fix command</li>
<li>vec1,vec2,... = <em>mol</em> or <em>q</em> or <em>i_name</em> or <em>d_name</em></li>
<li>vec1,vec2,... = <em>mol</em> or <em>q</em> or <em>rmass</em> or <em>i_name</em> or <em>d_name</em></li>
</ul>
<pre class="literal-block">
<em>mol</em> = molecule IDs
<em>q</em> = charge
<em>rmass</em> = per-atom mass
<em>i_name</em> = new integer vector referenced by name
<em>d_name</em> = new floating-point vector referenced by name
</pre>
@ -168,22 +169,32 @@ properties, as explained on the <a class="reference internal" href="atom_style.h
<a class="reference internal" href="read_data.html"><em>read_data</em></a> doc pages. The latter command allows these
properties to be defined for each atom in the system when a data file
is read. This fix will augment the set of properties with new custom
ones.</p>
<p>This can be useful in at least two scenarios.</p>
<p>If the atom style does not define molecule IDs or per-atom charge,
they can be added using the <em>mol</em> or <em>q</em> keywords. This can be
useful, e.g, to define &#8220;molecules&#8221; to use as rigid bodies with the
<a class="reference internal" href="fix_rigid.html"><em>fix rigid</em></a> command, or just to carry around an extra
flag with the atoms (stored as a molecule ID). An alternative is to
use an atom style that does define molecule IDs or charge or to use a
hybrid atom style that combines two styles to allow for molecule IDs
or charge, but that has 2 practical drawbacks. First it typically
necessitates changing the format of the data file. And it may define
additional properties that aren&#8217;t needed such as bond lists, which has
some overhead when there are no bonds.</p>
ones. This can be useful in several scenarios.</p>
<p>If the atom style does not define molecule IDs, per-atom charge,
or per-atom mass, they can be added using the <em>mol</em>, <em>q</em> or <em>rmass</em>
keywords. This can be useful, e.g, to define &#8220;molecules&#8221; to use as
rigid bodies with the <a class="reference internal" href="fix_rigid.html"><em>fix rigid</em></a> command, or just to
carry around an extra flag with the atoms (stored as a molecule ID)
that can be used to group atoms without having to use the group
command (which is limited to a total of 32 groups including <em>all</em>).</p>
<p>Another application would be to use the <em>rmass</em> flag in order to have
per-atom masses instead of per-type masses, for example this can be
useful to study isotope effects with partial isotope substitution.
Please <a class="reference internal" href="#isotopes"><span>see below</span></a> for an example of simulating a mixture
of light and heavy water with the TIP4P water potential.</p>
<p>An alternative to using fix <em>property/atom</em> in these ways is to
use an atom style that does define molecule IDs or charge or per-atom
mass (indirectly via diameter and density) or to use a hybrid atom
style that combines two or more atom styles
to provide the union of all atom properties. However, this has two
practical drawbacks: first, it typically necessitates changing the
format of the data file, which can be tedious for large systems;
and second, it may define additional properties that are not needed
such as bond lists, which has some overhead when there are no bonds.</p>
<p>In the future, we may add additional per-atom properties similar to
<em>mol</em> or <em>q</em>, which &#8220;turn-on&#8221; specific properties defined by some atom
styles, so they can be used by atom styles that don&#8217;t define them.</p>
<em>mol</em>, <em>q</em> or <em>rmass</em>, which &#8220;turn-on&#8221; specific properties defined
by some atom styles, so they can be used by atom styles that do not
define them.</p>
<p>More generally, the <em>i_name</em> and <em>d_name</em> vectors allow one or more
new custom per-atom properties to be defined. Each name must be
unique and can use alphanumeric or underscore characters. These
@ -204,9 +215,9 @@ when ghost atoms are created (at every re-neighboring) to insure the
new properties are also defined for the ghost atoms.</p>
<div class="admonition note">
<p class="first admonition-title">Note</p>
<p class="last">If you use this command with the <em>mol</em> or <em>charge</em> vectors than
you most likely want to set <em>ghost</em> yes, since these properties are
stored with ghost atoms if you use an <a class="reference internal" href="atom_style.html"><em>atom_style</em></a>
<p class="last">If you use this command with the <em>mol</em>, <em>q</em> or <em>rmass</em> vectors,
then you most likely want to set <em>ghost</em> yes, since these properties
are stored with ghost atoms if you use an <a class="reference internal" href="atom_style.html"><em>atom_style</em></a>
that defines them, and many LAMMPS operations that use molecule IDs or
charge, such as neighbor lists and pair styles, will expect ghost
atoms to have these valuse. LAMMPS will issue a warning it you define
@ -291,6 +302,41 @@ dump 1 all custom 100 tmp.dump id x y z c_1[1] c_1[2]
<a class="reference internal" href="fix.html"><em>fixes</em></a>, or <a class="reference internal" href="compute.html"><em>computes</em></a> that use the per-atom
properties defined by this fix, see <a class="reference internal" href="Section_modify.html#mod-1"><span>Section modify</span></a> of the manual which has some details
on how the properties can be accessed from added classes.</p>
<hr class="docutils" />
<p id="isotopes">Example for using per-atom masses with TIP4P water to study isotope
effects. When setting up simulations with the <a class="reference internal" href="Section_howto.html#howto-8"><span>TIP4P pair styles</span></a> for water, you have to provide
exactly one atom type each to identify the water oxygen and hydrogen
atoms. Since the atom mass is normally tied to the atom type, this
makes it impossible to study multiple isotopes in the same simulation.
With <em>fix property/atom rmass</em> however, the per-type masses are
replaced by per-atom masses. Asumming you have a working input deck
for regular TIP4P water, where water oxygen is atom type 1 and
water hydrogen is atom type 2, the following lines of input script
convert this to using per-atom masses:</p>
<div class="highlight-python"><div class="highlight"><pre>fix Isotopes all property/atom rmass ghost yes
set type 1 mass 15.9994
set type 2 mass 1.008
</pre></div>
</div>
<p>When writing out the system data with the <a class="reference internal" href="write_data.html"><em>write_data</em></a>
command, there will be a new section named with the fix-ID
(i.e. <em>Isotopes</em> in this case). Alternatively, you can take an
existing data file and just add this <em>Isotopes</em> section with
one line per atom containing atom-ID and mass. Either way, the
extended data file can be read back with:</p>
<div class="highlight-python"><div class="highlight"><pre>fix Isotopes all property/atom rmass ghost yes
read_data tip4p-isotopes.data fix Isotopes NULL Isotopes
</pre></div>
</div>
<p>Please note that the first <em>Isotopes</em> refers to the fix-ID
and the second to the name of the section. The following input
script code will now change the first 100 water molecules in this
example to heavy water:</p>
<div class="highlight-python"><div class="highlight"><pre>group hwat id 2:300:3
group hwat id 3:300:3
set group hwat mass 2.0141018
</pre></div>
</div>
</div>
<hr class="docutils" />
<div class="section" id="restart-fix-modify-output-run-start-stop-minimize-info">

View File

@ -14,9 +14,10 @@ fix ID group-ID property/atom vec1 vec2 ... keyword value ... :pre
ID, group-ID are documented in "fix"_fix.html command :ulb,l
property/atom = style name of this fix command :l
vec1,vec2,... = {mol} or {q} or {i_name} or {d_name} :l
vec1,vec2,... = {mol} or {q} or {rmass} or {i_name} or {d_name} :l
{mol} = molecule IDs
{q} = charge
{rmass} = per-atom mass
{i_name} = new integer vector referenced by name
{d_name} = new floating-point vector referenced by name :pre
zero of more keyword/value pairs may be appended :l
@ -41,25 +42,36 @@ properties, as explained on the "atom_style"_atom_style.html and
"read_data"_read_data.html doc pages. The latter command allows these
properties to be defined for each atom in the system when a data file
is read. This fix will augment the set of properties with new custom
ones.
ones. This can be useful in several scenarios.
This can be useful in at least two scenarios.
If the atom style does not define molecule IDs, per-atom charge,
or per-atom mass, they can be added using the {mol}, {q} or {rmass}
keywords. This can be useful, e.g, to define "molecules" to use as
rigid bodies with the "fix rigid"_fix_rigid.html command, or just to
carry around an extra flag with the atoms (stored as a molecule ID)
that can be used to group atoms without having to use the group
command (which is limited to a total of 32 groups including {all}).
If the atom style does not define molecule IDs or per-atom charge,
they can be added using the {mol} or {q} keywords. This can be
useful, e.g, to define "molecules" to use as rigid bodies with the
"fix rigid"_fix_rigid.html command, or just to carry around an extra
flag with the atoms (stored as a molecule ID). An alternative is to
use an atom style that does define molecule IDs or charge or to use a
hybrid atom style that combines two styles to allow for molecule IDs
or charge, but that has 2 practical drawbacks. First it typically
necessitates changing the format of the data file. And it may define
additional properties that aren't needed such as bond lists, which has
some overhead when there are no bonds.
Another application would be to use the {rmass} flag in order to have
per-atom masses instead of per-type masses, for example this can be
useful to study isotope effects with partial isotope substitution.
Please "see below"_#isotopes for an example of simulating a mixture
of light and heavy water with the TIP4P water potential.
An alternative to using fix {property/atom} in these ways is to
use an atom style that does define molecule IDs or charge or per-atom
mass (indirectly via diameter and density) or to use a hybrid atom
style that combines two or more atom styles
to provide the union of all atom properties. However, this has two
practical drawbacks: first, it typically necessitates changing the
format of the data file, which can be tedious for large systems;
and second, it may define additional properties that are not needed
such as bond lists, which has some overhead when there are no bonds.
In the future, we may add additional per-atom properties similar to
{mol} or {q}, which "turn-on" specific properties defined by some atom
styles, so they can be used by atom styles that don't define them.
{mol}, {q} or {rmass}, which "turn-on" specific properties defined
by some atom styles, so they can be used by atom styles that do not
define them.
More generally, the {i_name} and {d_name} vectors allow one or more
new custom per-atom properties to be defined. Each name must be
@ -82,9 +94,9 @@ keyword with a value of {yes}. This will invoke extra communication
when ghost atoms are created (at every re-neighboring) to insure the
new properties are also defined for the ghost atoms.
NOTE: If you use this command with the {mol} or {charge} vectors than
you most likely want to set {ghost} yes, since these properties are
stored with ghost atoms if you use an "atom_style"_atom_style.html
NOTE: If you use this command with the {mol}, {q} or {rmass} vectors,
then you most likely want to set {ghost} yes, since these properties
are stored with ghost atoms if you use an "atom_style"_atom_style.html
that defines them, and many LAMMPS operations that use molecule IDs or
charge, such as neighbor lists and pair styles, will expect ghost
atoms to have these valuse. LAMMPS will issue a warning it you define
@ -184,6 +196,44 @@ on how the properties can be accessed from added classes.
:line
:link(isotopes)
Example for using per-atom masses with TIP4P water to study isotope
effects. When setting up simulations with the "TIP4P pair
styles"_Section_howto.html#howto_8 for water, you have to provide
exactly one atom type each to identify the water oxygen and hydrogen
atoms. Since the atom mass is normally tied to the atom type, this
makes it impossible to study multiple isotopes in the same simulation.
With {fix property/atom rmass} however, the per-type masses are
replaced by per-atom masses. Asumming you have a working input deck
for regular TIP4P water, where water oxygen is atom type 1 and
water hydrogen is atom type 2, the following lines of input script
convert this to using per-atom masses:
fix Isotopes all property/atom rmass ghost yes
set type 1 mass 15.9994
set type 2 mass 1.008 :pre
When writing out the system data with the "write_data"_write_data.html
command, there will be a new section named with the fix-ID
(i.e. {Isotopes} in this case). Alternatively, you can take an
existing data file and just add this {Isotopes} section with
one line per atom containing atom-ID and mass. Either way, the
extended data file can be read back with:
fix Isotopes all property/atom rmass ghost yes
read_data tip4p-isotopes.data fix Isotopes NULL Isotopes :pre
Please note that the first {Isotopes} refers to the fix-ID
and the second to the name of the section. The following input
script code will now change the first 100 water molecules in this
example to heavy water:
group hwat id 2:300:3
group hwat id 3:300:3
set group hwat mass 2.0141018 :pre
:line
[Restart, fix_modify, output, run start/stop, minimize info:]
This fix writes the per-atom values it stores to "binary restart

File diff suppressed because one or more lines are too long