diff --git a/doc/src/fix_alchemy.rst b/doc/src/fix_alchemy.rst index 90c436139e..056cf60369 100644 --- a/doc/src/fix_alchemy.rst +++ b/doc/src/fix_alchemy.rst @@ -12,7 +12,7 @@ Syntax * ID, group-ID are documented in :doc:`fix ` command * alchemy = style name of this fix command -* v_name = variable with name that determines the :math:`\lambda_p` value +* v_name = variable with name that determines the :math:`\lambda_R` value Examples """""""" @@ -26,80 +26,121 @@ Description .. versionadded:: TBD -This fix command enables running an "alchemical transformation" MD -simulation between two topologies (i.e. the same number and positions of -atoms, but differences in atom parameters like type, charge, bonds, -angles and so on). For this a :ref:`multi-partition run ` is -required with exactly two partitions. During the MD run, the fix will -will determine a factor, :math:`\lambda_p`, for each partition *p* that -will be taken from an equal style or equivalent :doc:`variable -`. Typically, this variable would be chose to linearly ramp -*down* from 1.0 to 0.0 for the *first* partition (*p=0*) and linearly -ramp *up* from 0.0 to 1.0 for the *second* partition (*p=1*). The -forces used for the propagation of the atoms will be the sum of the -forces of the two systems combined and scaled with their respective -:math:`\lambda_p` factor. This allows to perform transformations that -are not easily possible with :doc:`pair style hybrid/scaled -`, :doc:`fix adapt ` or :doc:`fix adapt/fep -`. +This fix command enables an "alchemical transformation" to be performed +between two systems, whereby one system slowly transforms into the other +over the course of a molecular dynamics run. This is useful for +measuring thermodynamic differences between two different systems. It +also allows transformations that are not easily possible with the +:doc:`pair style hybrid/scaled `, :doc:`fix adapt +` or :doc:`fix adapt/fep ` commmands. -.. note:: +Example inputs are included in the ``examples/PACKAGES/alchemy`` +directory for (a) transforming a pure copper system into a +copper/aluminum bronze alloy and (b) transforming two water molecules +into a hydronium plus hydroxyl ion. - Since the definition of the variable to provide the :math:`\lambda_p` is - independent in the two partitions, no check is made that the two values - remain between 0.0 and 1.0 and that they add up to 1.0. So care needs to - be taken when defining those variables that this is the case. +The two systems must be defined in :doc:`separate replicas +` and run in separate partitions of processors using the +:doc:`-partition ` command-line switch. Exacty two +partitions must be specified and each partition must use the same number +of processors and the same domain decomposition. -Due to the specifics of the implementation, the initial geometry and -dimensions of the system must be exactly the same and the fix will -synchronize them during the run. It is thus not possible to initialize -the two partitions by reading different data files or creating different -systems from scratch, but rather they have to be started from the same -system and then the desired modifications need to be applied to the -system of the second partition. The commands :doc:`pair style -hybrid/scaled `, :doc:`fix adapt ` or :doc:`fix -adapt/fep ` could be used for simulations where the -requirements for fix alchemy are not given. +Because the forces applied to the atoms are the same mix of the forces +from each partition and the simulation starts with the same atom +positions across both partitions, they will generate the same trajectory +of coordinates for each atom, and the same simulation box size and +shape. The latter two conditions are *enforced* by this fix; it +exchanges coordinates and box information between the replicas. This is +not strictly required, but since MD simulations are an example of a +chaotic system, even the tiniest random difference will eventually grow +exponentially into an unwanted divergence. -The commands below demonstrate how the setup for the second partition -can be done for the example of transforming a pure copper system into a -copper/aluminum bronze. +Otherwise the properties of each atom (type, charge, bond and angle +partners, etc), as well as energy and forces between interacting atoms +(pair, bond, angle styles, etc) can be different in the two systems. -.. code-block:: LAMMPS +This can be initialized in the same input script by using commands which +only apply to one or the other replica. The example scripts use a +world-style :doc:`variable ` command along with +:doc:`if/then/else ` commmands for this purpose. The +:doc:`partition ` command can also be used. variable name world pure alloy +.. code-block:: LAMMPS + create_box 2 box create_atoms 1 box pair_style eam/alloy pair_coeff * * AlCu.eam.alloy Cu Al # replace 5% of copper with aluminum on the second partition only + variable name world pure alloy if "${name} == alloy" then & "set type 1 type/fraction 2 0.05 6745234" - # define ramp variable to combine the two different partitions - if "${name} == pure" then & - "variable ramp equal ramp(1.0,0.0)" & - else & - "variable ramp equal ramp(0.0,1.0)" +Both replicas must define an instance of this fix, but with a different +*v_name* variable. The named variable must be an equal-style or +equivalent :doc:`variable `. The two variables should be +defined so that one ramps *down* from 1.0 to 0.0 for the *first* replica +(*R=0*) and the other ramps *up* from 0.0 to 1.0 for the *second* +replica (*R=1*). A simple way is to do this is lineraly, which can be +done using the ramp() function of the :doc:`variable ` +command. You could also define a variable which returns a value between +0.0 and 1.0 as a non-linear function of the timestep. Here is a linear +example: +.. code-block:: LAMMPS + + partition yes 1 variable ramp equal ramp(1.0,0.0) + partition yes 2 variable ramp equal ramp(0.0,1.0) fix 2 all alchemy v_ramp +.. note:: -The ``examples/PACKAGES/alchemy`` folder contains complete example -inputs for this command. + For an alchemical transformation, the two variables should sum to + exactly 1.0 at any timestep. LAMMPS does *NOT* check that this is + the case. + +If you use the ``ramp()`` function to define the two variables, this fix +can easily be used across successive runs in the same input script by +ensuring each instance of the :doc:`run ` command specifies the +appropriate *start* or *stop* options. + +At each timestep of an MD run, the two instances of this fix evaluate +their respective variables as a :math:`\lambda_R` factor, where *R* = 0 +or 1 for each replica. The forces used by each system for the +propagation of theire atoms is set to the sum of the forces for the two +systems, each scaled by their respective :math:`\lambda_R` factor. Thus +during the MD run, the system will transform incrementally from from the +first system to the second system. + +.. note:: + + As mentioned above, the coordinates of the atoms and box size/shape + must be exactly the same in the two replicas. Thus it is generally + not a good idea to initialize the two replicas by reading different + data files or creating them from scratch. Rather, a single system + should be initialized and desired modifications applied to the system + of the second replica. If your input script somehow induces the two + systems to become different (e.g. by performing :doc:`atom_modify + sort ` differently, or by adding or depositing a + different number of atoms), then LAMMPS will detect the mismatchand + generate an error. This is done by ensuring that each step the + number and ordering of atoms is identical within each pair of + processors in the two replicas. ---------- Restart, fix_modify, output, run start/stop, minimize info """"""""""""""""""""""""""""""""""""""""""""""""""""""""""" -No information about this fix is written to :doc:`binary restart files `. -None of the :doc:`fix_modify ` options are relevant to this fix. +No information about this fix is written to :doc:`binary restart files +`. None of the :doc:`fix_modify ` options are +relevant to this fix. -This fix stores a global scalar (the current value of :math:`\lambda_p`) +This fix stores a global scalar (the current value of :math:`\lambda_R`) and a global vector of length 3 which contains the potential energy of the first partition, the second partition and the combined value, respectively. The global scalar is unitless and "intensive", the vector @@ -117,12 +158,9 @@ This fix is part of the REPLICA package. It is only enabled if LAMMPS was built with that package. See the :doc:`Build package ` page for more info. -There may be only one instance of this fix in use at any time. +There may be only one instance of this fix in use at a time within +each replica. -This fix requires to perform a :ref:`multi-partition run ` -with *exactly* two partitions. - -This fix is *not* compatible with :doc:`load balancing `. Related commands """"""""""""""""