SourceForge.net Logo
MCCCS Towhee (Algorithms: Chemical Potential)

 

 

  • WARNING The treatment of normalization constants for rigid bonds and rigid angles changed yet again in version 6.2.10. Recomputing the chemical potential for such molecules is required when comparing between versions.
  • WARNING Prior to the 6.0.4 release the treatment of some normalization constants for rigid bond and angles was different. Chemical potentials for molecules with rigid bonds or angles computed in earlier versions are going to give a different result in newer versions. This is usually just a constant shift so a quick recalculation of the chemical potential is advised before using old values in a new Grand Canonical simulation.
  • WARNING Prior to the 4.5.2 release all chemical potentials (except the Gibbs total chemical potential) were defined differently than documented here.
Chemical Potential Definitions
    There are several different types of chemical potential output by the Towhee program. The various chemical potentials reported by the code are listed and defined in this essay. Consider the following equation as a starting point.
      μtotal(i) = μresidual(i) + μdensity(i) + μisolation(i)
    μtotal(i) is the total chemical potential of molecule i and is the appropriate value to consider when looking for total free energies of a system. This total chemical potential is then commonly broken down into several subcomponents either for computational convenience, or because the subcomponents relate to certain thermodynamic quantities. Towhee defines, and calculates, these subcomponents as follows.
  • The total Gibbs chemical potential is computed using the insertion Rosenbluth weight (containing all interactions - intermolecular and intramolecular) and is computed each time a multi-box swap move is attempted, and also chempotperstep times at the end of every Monte Carlo cycle.
      μtotal(i) = - kB T ln[ < W * V / ( [N(i)+1] * Λ3(i) ) > ]
    where kB is Boltzmann's constant, T is the temperature, V is the volume of the simulation box, N(i) is the number of molecules of type i in the simulation box, and Λ(i) is the thermal de Broglie wavelength for molecule type i. The < > brackets indicate an ensemble average of the quantity inside the brackets.
  • The ideal density dependent portion of the chemical potential is
      μdensity(i) = - kB T ln[ < V / N(i) > / Λ3(i)) ]
  • The isolated molecule portion of the chemical potential is the non-density dependent chemical potential of an isolated molecule in an ideal gas (which means it does not interact with any other molecule, but does interact with itself). This quantity is not normally computed by Towhee. However, this is computed in the special case where the total number of molecules (nmolectyp) of one of the species is set to zero. In that case, the isolated molecule chemical potential is computed chempotperstep times per cycle. It uses the full Rosenbluth weight (W) computed for growing the molecule in an extremely large box that has the Ewald sum and the van der Waals tail corrections disabled. This uses the standard formula for computing the chemical potential via configurational-bias.
      μisolation(i) = -kB T ln[ < Wisolation > ].
  • The NVT insertion chemical potential is computed in a very similar manner as the Gibbs total chemical potential. The only difference is the lack of any number density terms in the average. This chemical potential calculation is only formally correct for a single box canonical ensemble. However, experience has shown that this chemical potential is often better when the average number of molecules of the type of interest drops below 1. It turns out that systematic errors in the number density term from using the Gibbs total chemical potential (essentially the use of N+1) are significant in this case. This is discussed further in Martin and Siepmann 1998 (TCA).
      μNVT Insertion(i) = - kB T ln[ < W > ]
  • The NpT insertion chemical potential is computed in a very similar manner as the NVT insertion chemical potential. The only difference is the inclusion of a volume term in the average, that is then removed later by dividing out the average volume. This chemical potential calculation is only formally correct for a single box isobaric-isothermal ensemble. However, experience has shown that this chemical potential is often better when the average number of molecules of the type of interest drops below 1. It turns out that systematic errors in the number density term from using the Gibbs total chemical potential (essentially the use of N+1) are significant in this case. This is discussed further in Martin and Siepmann 1998 (TCA).
      μNpT Insertion(i) = - kB T ln[ < V W > / < V > ]
  • The residual chemical potential is not computed directly in Towhee. Instead it is inferred by computing the full insertion chemical potential and then substracting the isolated molecule chemical potential. This chemical potential is only computed in the special case when we are interested in Henry's Law coefficients and nmolectyp is set to 0 for at least one of the components. The residual chemical potential is computed in two slightly different ways depening on the ensemble. If we are using the canonical ensemble then the residual chemical potential is
      μresidual(i) = μNVT Insertion(i) - μIsolation(i)
    and if we are using the isobaric-isothermal ensemble (the typical ensemble for computing Henry's law) then the residual chemical potential is
      μresidual(i) = μNpT Insertion(i) - μIsolation(i).
Return to the Towhee algorithm page

 


Send comments to: Marcus G. Martin

Last updated: December 17, 2014