diff --git a/doc/src/fix_msst.rst b/doc/src/fix_msst.rst index 36f2bd2ae0..80bfaf24bd 100644 --- a/doc/src/fix_msst.rst +++ b/doc/src/fix_msst.rst @@ -35,7 +35,7 @@ Examples """""""" -.. parsed-literal:: +.. code-block:: LAMMPS fix 1 all msst y 100.0 q 1.0e5 mu 1.0e5 fix 2 all msst z 50.0 q 1.0e4 mu 1.0e4 v0 4.3419e+03 p0 3.7797e+03 e0 -9.72360e+02 tscale 0.01 @@ -101,16 +101,17 @@ timestep. To do this, the fix creates its own computes of style "temp" "pressure", and "pe", as if these commands had been issued: -.. parsed-literal:: +.. code-block:: LAMMPS compute fix-ID_MSST_temp all temp compute fix-ID_MSST_press all pressure fix-ID_MSST_temp compute fix-ID_MSST_pe all pe -See the :doc:`compute temp ` and :doc:`compute pressure ` commands for details. Note that the -IDs of the new computes are the fix-ID + "_MSST\_temp`or `_ -or "_MSST\_pe". The group for the new computes is "all". +See the :doc:`compute temp ` and :doc:`compute pressure +` commands for details. Note that the IDs of the +new computes are the fix-ID + "_MSST\_temp" or "MSST\_press" or +"_MSST\_pe". The group for the new computes is "all". ---------- @@ -132,9 +133,10 @@ timestepping. DFTB+ will communicate its info to LAMMPS via that fix. **Restart, fix\_modify, output, run start/stop, minimize info:** -This fix writes the state of all internal variables to :doc:`binary restart files `. See the :doc:`read_restart ` command -for info on how to re-specify a fix in an input script that reads a -restart file, so that the operation of the fix continues in an +This fix writes the state of all internal variables to :doc:`binary +restart files `. See the :doc:`read_restart ` +command for info on how to re-specify a fix in an input script that +reads a restart file, so that the operation of the fix continues in an uninterrupted fashion. The progress of the MSST can be monitored by printing the global @@ -161,7 +163,7 @@ To print these quantities to the log file with descriptive column headers, the following LAMMPS commands are suggested: -.. parsed-literal:: +.. code-block:: LAMMPS fix msst all msst z fix_modify msst energy yes @@ -172,15 +174,17 @@ headers, the following LAMMPS commands are suggested: thermo_style custom step temp ke pe lz pzz etotal v_dhug v_dray v_lgr_vel v_lgr_pos f_msst These fixes compute a global scalar and a global vector of 4 -quantities, which can be accessed by various :doc:`output commands `. The scalar values calculated by this fix -are "extensive"; the vector values are "intensive". +quantities, which can be accessed by various :doc:`output commands +`. The scalar values calculated by this fix are +"extensive"; the vector values are "intensive". Restrictions """""""""""" This fix style is part of the SHOCK package. It is only enabled if -LAMMPS was built with that package. See the :doc:`Build package ` doc page for more info. +LAMMPS was built with that package. See the :doc:`Build package +` doc page for more info. All cell dimensions must be periodic. This fix can not be used with a triclinic cell. The MSST fix has been tested only for the group-ID diff --git a/doc/src/fix_qbmsst.rst b/doc/src/fix_qbmsst.rst index 282fe8487b..fac772686c 100644 --- a/doc/src/fix_qbmsst.rst +++ b/doc/src/fix_qbmsst.rst @@ -26,7 +26,7 @@ Syntax *v0* value = initial simulation cell volume in the shock equations (distance\^3 units) *e0* value = initial total energy (energy units) *tscale* value = reduction in initial temperature (unitless fraction between 0.0 and 1.0) - *damp* value = damping parameter (time units) inverse of friction γ + *damp* value = damping parameter (time units) inverse of friction *gamma* *seed* value = random number seed (positive integer) *f_max* value = upper cutoff frequency of the vibration spectrum (1/time units) *N_f* value = number of frequency bins (positive integer) @@ -40,44 +40,47 @@ Examples """""""" -.. parsed-literal:: +.. code-block:: LAMMPS - fix 1 all qbmsst z 0.122 q 25 mu 0.9 tscale 0.01 damp 200 seed 35082 f_max 0.3 N_f 100 eta 1 beta 400 T_init 110 (liquid methane modeled with the REAX force field, real units) - fix 2 all qbmsst z 72 q 40 tscale 0.05 damp 1 seed 47508 f_max 120.0 N_f 100 eta 1.0 beta 500 T_init 300 (quartz modeled with the BKS force field, metal units) + # (liquid methane modeled with the REAX force field, real units) + fix 1 all qbmsst z 0.122 q 25 mu 0.9 tscale 0.01 damp 200 seed 35082 f_max 0.3 N_f 100 eta 1 beta 400 T_init 110 + # (quartz modeled with the BKS force field, metal units) + fix 2 all qbmsst z 72 q 40 tscale 0.05 damp 1 seed 47508 f_max 120.0 N_f 100 eta 1.0 beta 500 T_init 300 Two example input scripts are given, including shocked -:math:`\alpha-\mathrm{quartz}` and shocked liquid methane. -The input script first equilibrate an initial state with the quantum -thermal bath at the target temperature and then apply the qbmsst to -simulate shock compression with quantum nuclear correction. The -following two figures plot related quantities for shocked alpha quartz. +:math:`\alpha\textrm{-quartz}` and shocked liquid methane. The input +script first equilibrates an initial state with the quantum thermal +bath at the target temperature and then applies *fix qbmsst* to simulate +shock compression with quantum nuclear correction. The following two +figures plot relevant quantities for shocked +:math:`\alpha\textrm{-quartz}`. .. image:: JPG/qbmsst_init.jpg :align: center -Figure 1. Classical temperature Tcl = ∑ -mivi2/3NkB vs. time -for coupling the alpha quartz initial state with the quantum thermal -bath at target quantum temperature Tqm = 300 K. The -NpH ensemble is used for time integration while QTB provides the -colored random force. Tcl converges at the timescale -of *damp* which is set to be 1 ps. +Figure 1. Classical temperature +:math:`T_{cl} = \sum \frac{m_iv_i^2}{3Nk_B}` vs. time for coupling the +:math:`\alpha\textrm{-quartz}` initial state with the quantum thermal +bath at target quantum temperature :math:`T^{qm} = 300 K`. The NpH +ensemble is used for time integration while QTB provides the colored +random force. :math:`T^{cl}` converges at the timescale of *damp* +which is set to be 1 ps. .. image:: JPG/qbmsst_shock.jpg :align: center Figure 2. Quantum temperature and pressure vs. time for simulating -shocked alpha quartz with the QBMSST. The shock propagates along the z -direction. Restart of the QBMSST command is demonstrated in the -example input script. Thermodynamic quantities stay continuous before -and after the restart. +shocked :math:`\alpha\textrm{-quartz}` with *fix qbmsst*\. The shock +propagates along the z direction. Restart of the *fix qbmsst* command +is demonstrated in the example input script. Thermodynamic quantities +stay continuous before and after the restart. Description """"""""""" This command performs the Quantum-Bath coupled Multi-Scale Shock Technique (QBMSST) integration. See :ref:`(Qi) ` for a detailed -description of this method. The QBMSST provides description of the +description of this method. QBMSST provides description of the thermodynamics and kinetics of shock processes while incorporating quantum nuclear effects. The *shockvel* setting determines the steady shock velocity that will be simulated along direction *dir*\ . @@ -107,14 +110,14 @@ in the command :doc:`fix msst `. The values of *e0*\ , *p0*\ , or parameter of *damp*\ , *f\_max*, and *N\_f* are described in the command :doc:`fix qtb `. -The fix qbmsst command couples the shock system to a quantum thermal +The *fix qbmsst* command couples the shock system to a quantum thermal bath with a rate that is proportional to the change of the total -energy of the shock system, etot - etot0. -Here etot consists of both the system energy and a thermal -term, see :ref:`(Qi) `, and etot0 = *e0* is the +energy of the shock system, :math:`E^{tot} - E^{tot}_0`. +Here :math:`E^{etot}` consists of both the system energy and a thermal +term, see :ref:`(Qi) `, and :math:`E^{tot}_0 = e0` is the initial total energy. -The *eta* (η) parameter is a unitless coupling constant +The *eta* (:math:`\eta`) parameter is a unitless coupling constant between the shock system and the quantum thermal bath. A small *eta* value cannot adjust the quantum temperature fast enough during the temperature ramping period of shock compression while large *eta* @@ -123,21 +126,18 @@ leads to big temperature oscillation. A value of *eta* between 0.3 and compression. We observe that different values of *eta* lead to almost the same final thermodynamic state behind the shock, as expected. -The quantum temperature is updated every *beta* (β) steps +The quantum temperature is updated every *beta* (:math:`\beta`) steps with an integration time interval *beta* times longer than the -simulation time step. In that case, etot is taken as its +simulation time step. In that case, :math:`E^{tot}` is taken as its average over the past *beta* steps. The temperature of the quantum -thermal bath Tqm changes dynamically according to -the following equation where Δt is the MD time step and -γ is the friction constant which is equal to the inverse +thermal bath :math:`T^{qm}` changes dynamically according to +the following equation where :math:`\Delta_t` is the MD time step and +:math:`\gamma` is the friction constant which is equal to the inverse of the *damp* parameter. -.. raw:: html +.. math:: -
dTqm/dt = - γη∑βl = - 1[etot(t-lΔt)-etot0]/3βNkB -
+ \frac{dT^{qm}}{dt} = \gamma\eta\sum^\beta_{l=1}\frac{E^{tot}(t-l\Delta t) - E^{tot}_0}{3\beta N k_B} The parameter *T\_init* is the initial temperature of the quantum thermal bath and the system before shock loading. @@ -172,12 +172,12 @@ vector contains five values in this order: 2. *drayleigh* is the departure from the Rayleigh line (pressure units). 3. *lagrangian\_speed* is the laboratory-frame Lagrangian speed (particle velocity) of the computational cell (velocity units). 4. *lagrangian\_position* is the computational cell position in the reference frame moving at the shock speed. This is the distance of the computational cell behind the shock front. -5. *quantum\_temperature* is the temperature of the quantum thermal bath Tqm. +5. *quantum\_temperature* is the temperature of the quantum thermal bath :math:`T^{qm}`. To print these quantities to the log file with descriptive column headers, the following LAMMPS commands are suggested. Here the :doc:`fix_modify ` energy command is also enabled to allow -the thermo keyword *etotal* to print the quantity etot. See +the thermo keyword *etotal* to print the quantity :math:`E^{tot}`. See also the :doc:`thermo_style ` command. @@ -193,10 +193,11 @@ also the :doc:`thermo_style ` command. thermo_style custom step temp ke pe lz pzz etotal v_dhug v_dray v_lgr_vel v_lgr_pos v_T_qm f_fix_id The global scalar under the entry f\_fix\_id is the quantity of thermo -energy as an extra part of etot. This global scalar and the -vector of 5 quantities can be accessed by various :doc:`output commands `. It is worth noting that the temp keyword +energy as an extra part of :math:`E^{tot}`. This global scalar and the +vector of 5 quantities can be accessed by various :doc:`output commands `. +It is worth noting that the temp keyword under the :doc:`thermo_style ` command print the -instantaneous classical temperature Tcl as described +instantaneous classical temperature :math:`T^{cl}` as described in the command :doc:`fix qtb `. diff --git a/doc/src/fix_qtb.rst b/doc/src/fix_qtb.rst index 235c28ca19..7242d0739b 100644 --- a/doc/src/fix_qtb.rst +++ b/doc/src/fix_qtb.rst @@ -19,7 +19,7 @@ Syntax .. parsed-literal:: *temp* value = target quantum temperature (temperature units) - *damp* value = damping parameter (time units) inverse of friction γ + *damp* value = damping parameter (time units) inverse of friction *gamma* *seed* value = random number seed (positive integer) *f_max* value = upper cutoff frequency of the vibration spectrum (1/time units) *N_f* value = number of frequency bins (positive integer) @@ -30,12 +30,14 @@ Examples """""""" -.. parsed-literal:: +.. code-block:: LAMMPS + # (liquid methane modeled with the REAX force field, real units) fix 1 all nve - fix 1 all qtb temp 110 damp 200 seed 35082 f_max 0.3 N_f 100 (liquid methane modeled with the REAX force field, real units) + fix 1 all qtb temp 110 damp 200 seed 35082 f_max 0.3 N_f 100 + # (quartz modeled with the BKS force field, metal units) fix 2 all nph iso 1.01325 1.01325 1 - fix 2 all qtb temp 300 damp 1 seed 47508 f_max 120.0 N_f 100 (quartz modeled with the BKS force field, metal units) + fix 2 all qtb temp 300 damp 1 seed 47508 f_max 120.0 N_f 100 Description """"""""""" @@ -56,61 +58,54 @@ atoms and thus higher classical limits. The equation of motion implemented by this command follows a Langevin form: -.. raw:: html +.. math:: -
miai = fi - + Ri - - miγvi.
+ m_i a_i = f_i + R_i - m_i\gamma v_i -Here mi, ai, fi -, Ri, γ and vi -represent mass, acceleration, force exerted by all other atoms, random +Here :math:`m_i, a_i, f_i, R_i, \gamma, \textrm{and} v_i` +represent in this order mass, acceleration, force exerted by all other atoms, random force, frictional coefficient (the inverse of damping parameter damp), -and velocity. The random force Ri is "colored" so -that any vibrational mode with frequency ω will have a -temperature-sensitive energy θ(ω,T) which +and velocity. The random force :math:`R_i` is "colored" so +that any vibrational mode with frequency :math:`\omega` will have a +temperature-sensitive energy :math:`\theta(\omega,T)` which resembles the energy expectation for a quantum harmonic oscillator with the same natural frequency: -.. raw:: html -
θ(ω,T) = - ℏω/2 + - ℏω[exp(ℏω/kBT)-1]-1 -
+.. math:: + + \theta(\omega T) = \frac{\hbar}{2} + \hbar\omega \left[\exp(\frac{\hbar\omega}{k_B T})-1 \right]^{-1} To efficiently generate the random forces, we employ the method of :ref:`(Barrat) `, that circumvents the need to generate all random forces for all times before the simulation. The memory requirement of this approach is less demanding and independent -of the simulation duration. Since the total random force Rtot +of the simulation duration. Since the total random force :math:`R_{tot}` does not necessarily vanish for a finite number of atoms, -Ri is replaced by Ri - Rtot/Ntot +:math:`R_i` is replaced by :math:`R_i - \frac{R_{tot}}{N_{tot}}` to avoid collective motion of the system. The *temp* parameter sets the target quantum temperature. LAMMPS will still have an output temperature in its thermo style. That is the -instantaneous classical temperature Tcl derived from +instantaneous classical temperature :math:`T^{cl}` derived from the atom velocities at thermal equilibrium. A non-zero -Tcl will be present even when the quantum +:math:`T^{cl}` will be present even when the quantum temperature approaches zero. This is associated with zero-point energy at low temperatures. -.. raw:: html +.. math:: -
Tcl = ∑ - mivi2/3NkB -
+ T^{cl} = \sum \frac{m_i v_i^2}{3 N k_B} The *damp* parameter is specified in time units, and it equals the -inverse of the frictional coefficient γ. γ +inverse of the frictional coefficient :math:`\gamma`. :math:`\gamma` should be as small as possible but slightly larger than the timescale of anharmonic coupling in the system which is about 10 ps to 100 -ps. When γ is too large, it gives an energy spectrum that -differs from the desired Bose-Einstein spectrum. When γ +ps. When :math:`\gamma` is too large, it gives an energy spectrum that +differs from the desired Bose-Einstein spectrum. When :math:`\gamma` is too small, the quantum thermal bath coupling to the system will be less significant than anharmonic effects, reducing to a classical -limit. We find that setting γ between 5 THz and 1 THz +limit. We find that setting :math:`\gamma` between 5 THz and 1 THz could be appropriate depending on the system. The random number *seed* is a positive integer used to initiate a @@ -121,23 +116,22 @@ runs on different numbers of processors. The *f\_max* parameter truncate the noise frequency domain so that vibrational modes with frequencies higher than *f\_max* will not be -modulated. If we denote Δt as the time interval for the +modulated. If we denote :math:`\Delta t` as the time interval for the MD integration, *f\_max* is always reset by the code to make -α = (int)(2*f\_max*Δt)-1 a +:math:`\alpha = (int)(2` *f\_max* :math:`\Delta t)^{-1}` a positive integer and print out relative information. An appropriate -value for the cutoff frequency *f\_max* would be around 2~3 -fD, where fD is the Debye -frequency. +value for the cutoff frequency *f\_max* would be around 2~3 :math:`f_D`, +where :math:`f_D` is the Debye frequency. The *N\_f* parameter is the frequency grid size, the number of points from 0 to *f\_max* in the frequency domain that will be -sampled. 3×2 *N\_f* per-atom random numbers are required +sampled. 3*2\ *N\_f* per-atom random numbers are required in the random force generation and there could be as many atoms as in the whole simulation that can migrate into every individual processor. A larger *N\_f* provides a more accurate sampling of the spectrum while consumes more memory. With fixed *f\_max* and -γ, *N\_f* should be big enough to converge the classical -temperature Tcl as a function of target quantum bath +:math:`\gamma`, *N\_f* should be big enough to converge the classical +temperature :math:`T^{cl}` as a function of target quantum bath temperature. Memory usage per processor could be from 10 to 100 Mbytes. @@ -147,10 +141,12 @@ Mbytes. Nose/Hoover thermostatting AND time integration, this fix does NOT perform time integration. It only modifies forces to a colored thermostat. Thus you must use a separate time integration fix, like - :doc:`fix nve ` or :doc:`fix nph ` to actually update the - velocities and positions of atoms (as shown in the - examples). Likewise, this fix should not normally be used with other - fixes or commands that also specify system temperatures , e.g. :doc:`fix nvt ` and :doc:`fix temp/rescale `. + :doc:`fix nve ` or :doc:`fix nph ` to actually + update the velocities and positions of atoms (as shown in the + examples). Likewise, this fix should not normally be used with + other fixes or commands that also specify system temperatures , + e.g. :doc:`fix nvt ` and :doc:`fix temp/rescale + `. ---------- @@ -158,10 +154,11 @@ Mbytes. **Restart, fix\_modify, output, run start/stop, minimizie info:** -No information about this fix is written to :doc:`binary restart files `. Because the state of the random number generator -is not saved in restart files, this means you cannot do "exact" -restarts with this fix. However, in a statistical sense, a restarted -simulation should produce similar behaviors of the system. +No information about this fix is written to :doc:`binary restart files +`. Because the state of the random number generator is not +saved in restart files, this means you cannot do "exact" restarts with +this fix. However, in a statistical sense, a restarted simulation +should produce similar behaviors of the system. This fix is not invoked during :doc:`energy minimization `. @@ -174,7 +171,8 @@ Restrictions This fix style is part of the USER-QTB package. It is only enabled if -LAMMPS was built with that package. See the :doc:`Build package ` doc page for more info. +LAMMPS was built with that package. See the :doc:`Build package +` doc page for more info. ---------- @@ -183,7 +181,8 @@ LAMMPS was built with that package. See the :doc:`Build package ` Related commands """""""""""""""" -:doc:`fix nve `, :doc:`fix nph `, :doc:`fix langevin `, :doc:`fix qbmsst ` +:doc:`fix nve `, :doc:`fix nph `, +:doc:`fix langevin `, :doc:`fix qbmsst ` ---------- diff --git a/doc/src/pair_comb.rst b/doc/src/pair_comb.rst index f7fe3dc55a..b893f56502 100644 --- a/doc/src/pair_comb.rst +++ b/doc/src/pair_comb.rst @@ -45,8 +45,8 @@ Description Style *comb* computes the second-generation variable charge COMB (Charge-Optimized Many-Body) potential. Style *comb3* computes the third-generation COMB potential. These COMB potentials are described -in :ref:`(COMB) ` and :ref:`(COMB3) `. Briefly, the total energy -*ET* of a system of atoms is given by +in :ref:`(COMB) ` and :ref:`(COMB3) `. Briefly, the +total energy :math:`E_T` of a system of atoms is given by .. math::