Merge pull request #1891 from charlessievers/fix_numerical_differentiation

   * :doc:`npt/eff <fix_nh_eff>`

   * :doc:`npt/sphere (o) <fix_npt_sphere>`

   * :doc:`npt/uef <fix_nh_uef>`

   * :doc:`numdiff <fix_numdiff>`

   * :doc:`nve (iko) <fix_nve>`

   * :doc:`nve/asphere (i) <fix_nve_asphere>`

   * :doc:`nve/asphere/noforce <fix_nve_asphere_noforce>`

*Compute ID for fix ave/time does not exist*

   Self-explanatory.

*Compute ID for fix numdiff does not exist*

   Self-explanatory.

*Compute ID for fix store/state does not exist*

   Self-explanatory.

   Self-explanatory.  The change in the box tilt is too extreme

   on a short timescale.

*Fix numdiff requires an atom map, see atom_modify*

   Self-explanatory. Efficient loop over all atoms for numerical

   difference requires an atom map.

*Fix numdiff requires consecutive atom IDs*

   Self-explanatory. Efficient loop over all atoms for numerical

   difference requires consecutive atom IDs.

*Fix nve/asphere requires extended particles*

   This fix can only be used for particles with a shape setting.

Related commands

""""""""""""""""

:doc:`fix phonon <fix_phonon>`

:doc:`fix phonon <fix_phonon>`, :doc:`fix numdiff <fix_numdiff>`, 

:doc:`compute hma <compute_hma>` uses an analytic formulation of the

Hessian provided by a pair_style's Pair::single_hessian() function,

* :doc:`npt/eff <fix_nh_eff>` - NPT for  nuclei and electrons in the electron force field model

* :doc:`npt/sphere <fix_npt_sphere>` - NPT for spherical particles

* :doc:`npt/uef <fix_nh_uef>` - NPT style time integration with diagonal flow

* :doc:`numdiff <fix_numdiff>` - compute derivatives of per-atom data from finite differences

* :doc:`nve <fix_nve>` - constant NVE time integration

* :doc:`nve/asphere <fix_nve_asphere>` - NVE for aspherical particles

* :doc:`nve/asphere/noforce <fix_nve_asphere_noforce>` - NVE for aspherical particles without forces

.. index:: fix numdiff

fix numdiff command

====================

Syntax

""""""

.. parsed-literal::

   fix ID group-ID numdiff Nevery delta

* ID, group-ID are documented in :doc:`fix <fix>` command

* numdiff = style name of this fix command

* Nevery = calculate force by finite difference every this many timesteps

* delta = finite difference displacement length (distance units)

Examples

""""""""

.. code-block:: LAMMPS

   fix 1 all numdiff 1 0.0001

   fix 1 all numdiff 10 1e-6

   fix 1 all numdiff 100 0.01

Description

"""""""""""

Calculate forces through finite difference calculations of energy

versus position.  These forces can be compared to analytic forces

computed by pair styles, bond styles, etc.  E.g. for debugging

purposes.

The group specified with the command means only atoms within the group

have their averages computed.  Results are set to 0.0 for atoms not in

the group.

This fix performs a loop over all atoms (in the group).  For each atom

and each component of force it adds *delta* to the position, and

computes the new energy of the entire system.  It then subtracts

*delta* from the original position and again computes the new energy

of the system.  It then restores the original position.  That

component of force is calculated as the difference in energy divided

by two times *delta*.

.. note::

   It is important to choose a suitable value for delta, the magnitude of

   atom displacements that are used to generate finite difference

   approximations to the exact forces.  For typical systems, a value in

   the range of 1 part in 1e4 to 1e5 of the typical separation distance

   between atoms in the liquid or solid state will be sufficient.

   However, the best value will depend on a multitude of factors

   including the stiffness of the interatomic potential, the thermodynamic

   state of the material being probed, and so on. The only way to be sure

   that you have made a good choice is to do a sensitivity study on a

   representative atomic configuration, sweeping over a wide range of

   values of delta.  If delta is too small, the output forces will vary

   erratically due to truncation effects. If delta is increased beyond a

   certain point, the output forces will start to vary smoothly with

   delta, due to growing contributions from higher order derivatives. In

   between these two limits, the numerical force values should be largely

   independent of delta.

.. note::

   The cost of each energy evaluation is essentially the cost of an MD

   timestep.  This invoking this fix once has a cost of 2N timesteps,

   where N is the total number of atoms in the system (assuming all atoms

   are included in the group).  So this fix can be very expensive to use

   for large systems.

----------

The *Nevery* argument specifies on what timesteps the force will

be used calculated by finite difference.

The *delta* argument specifies the positional displacement each

atom will undergo.

----------

**Restart, fix_modify, output, run start/stop, minimize info:**

No information about this fix is written to :doc:`binary restart files

<restart>`.  None of the :doc:`fix_modify <fix_modify>` options are

relevant to this fix.

This fix produces a per-atom array which can be accessed by various

:doc:`output commands <Howto_output>`, which stores the components of

the force on each atom as calculated by finite difference.  The

per-atom values can only be accessed on timesteps that are multiples

of *Nevery* since that is when the finite difference forces are

calculated.

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

""""""""""""

 none

Related commands

""""""""""""""""

:doc:`dynamical_matrix <dynamical_matrix>`,

**Default:** none