| 1604 |
|
V_{\text{bond}}(b) = D_{ij} \left[ 1 - e^{-\beta_{ij} (b - b_{ij}^0)} \right]^2 |
| 1605 |
|
\end{equation} |
| 1606 |
|
|
| 1607 |
+ |
\begin{figure}[h] |
| 1608 |
+ |
\centering |
| 1609 |
+ |
\includegraphics[width=2.5in]{bond.pdf} |
| 1610 |
+ |
\caption[Bond coordinates]{The coordinate describing a |
| 1611 |
+ |
a bond between atoms $i$ and $j$ is $|r_{ij}|$, the length of the |
| 1612 |
+ |
$\vec{r}_{ij}$ vector. } |
| 1613 |
+ |
\label{fig:bond} |
| 1614 |
+ |
\end{figure} |
| 1615 |
+ |
|
| 1616 |
|
OpenMD can also simulate some less common types of bond potentials, |
| 1617 |
|
including {\tt Fixed} bonds (which are constrained to be at a |
| 1618 |
|
specified bond length), |
| 1684 |
|
The bending potential energy functions used in most force fields are |
| 1685 |
|
often simple functions of the angle between two bonds, |
| 1686 |
|
\begin{equation} |
| 1687 |
< |
\theta_{ijk} = \cos^{-1} \left(\frac{\vec{r}_{ij} \cdot |
| 1688 |
< |
\vec{r}_{jk}}{\left| \vec{r}_{ij} \right| \left| \vec{r}_{ij} |
| 1687 |
> |
\theta_{ijk} = \cos^{-1} \left(\frac{\vec{r}_{ji} \cdot |
| 1688 |
> |
\vec{r}_{jk}}{\left| \vec{r}_{ji} \right| \left| \vec{r}_{jk} |
| 1689 |
|
\right|} \right) |
| 1690 |
|
\end{equation} |
| 1691 |
|
Here atom $j$ is the central atom that is bonded to two partners $i$ |
| 1692 |
|
and $k$. |
| 1693 |
|
|
| 1694 |
+ |
\begin{figure}[h] |
| 1695 |
+ |
\centering |
| 1696 |
+ |
\includegraphics[width=3.5in]{bend.pdf} |
| 1697 |
+ |
\caption[Bend angle coordinates]{The coordinate describing a bend |
| 1698 |
+ |
between atoms $i$, $j$, and $k$ is the angle $\theta_{ijk} = |
| 1699 |
+ |
\cos^{-1} \left(\hat{r}_{ji} \cdot \hat{r}_{jk}\right)$ where $\hat{r}_{ji}$ is |
| 1700 |
+ |
the unit vector between atoms $j$ and $i$. } |
| 1701 |
+ |
\label{fig:bend} |
| 1702 |
+ |
\end{figure} |
| 1703 |
+ |
|
| 1704 |
+ |
|
| 1705 |
|
All BendTypes must specify three AtomType names ($i$, $j$ and $k$) |
| 1706 |
|
that describe when that bend potential should be applied, as well as |
| 1707 |
|
an equilibrium bending angle, $\theta_{ijk}^0$, in units of |
| 1823 |
|
\label{eq:torsPhi} |
| 1824 |
|
\end{equation} |
| 1825 |
|
Here, $\hat{\mathbf{r}}_{\alpha\beta}$ are the set of unit bond |
| 1826 |
< |
vectors between atoms $i$, $j$, $k$, and $l$. |
| 1826 |
> |
vectors between atoms $i$, $j$, $k$, and $l$. Note that some force |
| 1827 |
> |
fields define the zero of the $\phi_{ijkl}$ angle when atoms $i$ and |
| 1828 |
> |
$l$ are in the {\em trans} configuration, while most define the zero |
| 1829 |
> |
angle for when $i$ and $l$ are in the fully eclipsed orientation. The |
| 1830 |
> |
behavior of the torsion parser can be altered with the {\tt |
| 1831 |
> |
TorsionAngleConvention} keyword in the Options block. The default |
| 1832 |
> |
behavior is {\tt "180\_is\_trans"} while the opposite behavior can be |
| 1833 |
> |
invoked by setting this keyword to {\tt "0\_is\_trans"}. |
| 1834 |
> |
|
| 1835 |
> |
\begin{figure}[h] |
| 1836 |
> |
\centering |
| 1837 |
> |
\includegraphics[width=4.5in]{torsion.pdf} |
| 1838 |
> |
\caption[Torsion or dihedral angle coordinates]{The coordinate |
| 1839 |
> |
describing a torsion between atoms $i$, $j$, $k$, and $l$ is the |
| 1840 |
> |
dihedral angle $\phi_{ijkl}$ which measures the relative rotation of |
| 1841 |
> |
the two terminal atoms around the axis defined by the middle bond. } |
| 1842 |
> |
\label{fig:torsion} |
| 1843 |
> |
\end{figure} |
| 1844 |
|
|
| 1845 |
|
For computational efficiency, OpenMD recasts torsion potential in the |
| 1846 |
|
method of {\sc charmm},\cite{Brooks83} in which the angle series is |
| 2079 |
|
|
| 2080 |
|
\section{\label{section:electrostatics}Electrostatics} |
| 2081 |
|
|
| 2082 |
< |
To aid in performing simulations in more traditional force fields, we |
| 2083 |
< |
have added routines to carry out electrostatic interactions using a |
| 2084 |
< |
number of different electrostatic summation methods. These methods |
| 2085 |
< |
are extended from the damped and cutoff-neutralized Coulombic sum |
| 2086 |
< |
originally proposed by Wolf, {\it et al.}\cite{Wolf99} One of these, |
| 2087 |
< |
the damped shifted force method, shows a remarkable ability to |
| 2088 |
< |
reproduce the energetic and dynamic characteristics exhibited by |
| 2089 |
< |
simulations employing lattice summation techniques. The basic idea is |
| 2090 |
< |
to construct well-behaved real-space summation methods using two tricks: |
| 2082 |
> |
Because nearly all force fields involve electrostatic interactions in |
| 2083 |
> |
one form or another, OpenMD implements a number of different |
| 2084 |
> |
electrostatic summation methods. These methods are extended from the |
| 2085 |
> |
damped and cutoff-neutralized Coulombic sum originally proposed by |
| 2086 |
> |
Wolf, {\it et al.}\cite{Wolf99} One of these, the damped shifted force |
| 2087 |
> |
method, shows a remarkable ability to reproduce the energetic and |
| 2088 |
> |
dynamic characteristics exhibited by simulations employing lattice |
| 2089 |
> |
summation techniques. The basic idea is to construct well-behaved |
| 2090 |
> |
real-space summation methods using two tricks: |
| 2091 |
|
\begin{enumerate} |
| 2092 |
|
\item shifting through the use of image charges, and |
| 2093 |
|
\item damping the electrostatic interaction. |
| 2254 |
|
|
| 2255 |
|
\section{\label{section:cutoffGroups}Switching Functions} |
| 2256 |
|
|
| 2257 |
< |
If done poorly, calculating the the long-range interactions for $N$ |
| 2258 |
< |
atoms would involve $N(N-1)/2$ evaluations of atomic distances. To |
| 2259 |
< |
reduce the number of distance evaluations between pairs of atoms, {\sc |
| 2260 |
< |
OpenMD} allows the use of switched cutoffs with Verlet neighbor |
| 2261 |
< |
lists.\cite{Allen87} Neutral groups which contain charges will exhibit |
| 2262 |
< |
pathological forces unless the cutoff is applied to the neutral groups |
| 2263 |
< |
evenly instead of to the individual atoms.\cite{leach01:mm} {\sc |
| 2264 |
< |
OpenMD} allows users to specify cutoff groups which may contain an |
| 2265 |
< |
arbitrary number of atoms in the molecule. Atoms in a cutoff group |
| 2266 |
< |
are treated as a single unit for the evaluation of the switching |
| 2267 |
< |
function: |
| 2257 |
> |
Calculating the the long-range interactions for $N$ atoms involves |
| 2258 |
> |
$N(N-1)/2$ evaluations of atomic distances if it is done in a brute |
| 2259 |
> |
force manner. To reduce the number of distance evaluations between |
| 2260 |
> |
pairs of atoms, {\sc OpenMD} allows the use of hard or switched |
| 2261 |
> |
cutoffs with Verlet neighbor lists.\cite{Allen87} Neutral groups which |
| 2262 |
> |
contain charges can exhibit pathological forces unless the cutoff is |
| 2263 |
> |
applied to the neutral groups evenly instead of to the individual |
| 2264 |
> |
atoms.\cite{leach01:mm} {\sc OpenMD} allows users to specify cutoff |
| 2265 |
> |
groups which may contain an arbitrary number of atoms in the molecule. |
| 2266 |
> |
Atoms in a cutoff group are treated as a single unit for the |
| 2267 |
> |
evaluation of the switching function: |
| 2268 |
|
\begin{equation} |
| 2269 |
|
V_{\mathrm{long-range}} = \sum_{a} \sum_{b>a} s(r_{ab}) \sum_{i \in a} \sum_{j \in b} V_{ij}(r_{ij}), |
| 2270 |
|
\end{equation} |
| 2290 |
|
Here, $r_{\text{sw}}$ is the {\tt switchingRadius}, or the distance |
| 2291 |
|
beyond which interactions are reduced, and $r_{\text{cut}}$ is the |
| 2292 |
|
{\tt cutoffRadius}, or the distance at which interactions are |
| 2293 |
< |
truncated. |
| 2293 |
> |
truncated. |
| 2294 |
|
|
| 2295 |
|
Users of {\sc OpenMD} do not need to specify the {\tt cutoffRadius} or |
| 2296 |
< |
{\tt switchingRadius}. In simulations containing only Lennard-Jones |
| 2297 |
< |
atoms, the cutoff radius has a default value of $2.5\sigma_{ii}$, |
| 2298 |
< |
where $\sigma_{ii}$ is the largest Lennard-Jones length parameter |
| 2299 |
< |
present in the simulation. In simulations containing charged or |
| 2300 |
< |
dipolar atoms, the default cutoff radius is $15 \mbox{\AA}$. |
| 2296 |
> |
{\tt switchingRadius}. |
| 2297 |
> |
If the {\tt cutoffRadius} was not explicitly set, OpenMD will attempt |
| 2298 |
> |
to guess an appropriate choice. If the system contains electrostatic |
| 2299 |
> |
atoms, the default cutoff radius is 12 \AA. In systems without |
| 2300 |
> |
electrostatic (charge or multipolar) atoms, the atom types present in the simulation will be |
| 2301 |
> |
polled for suggested cutoff values (e.g. $2.5 max(\left\{ \sigma |
| 2302 |
> |
\right\})$ for Lennard-Jones atoms. The largest suggested value |
| 2303 |
> |
that was found will be used. |
| 2304 |
> |
|
| 2305 |
> |
By default, OpenMD uses shifted force potentials to force the |
| 2306 |
> |
potential energy and forces to smoothly approach zero at the cutoff |
| 2307 |
> |
radius. If the user would like to use another cutoff method |
| 2308 |
> |
they may do so by setting the {\tt cutoffMethod} parameter to: |
| 2309 |
> |
\begin{itemize} |
| 2310 |
> |
\item {\tt HARD} |
| 2311 |
> |
\item {\tt SWITCHED} |
| 2312 |
> |
\item {\tt SHIFTED\_FORCE} (default): |
| 2313 |
> |
\item {\tt TAYLOR\_SHIFTED} |
| 2314 |
> |
\item {\tt SHIFTED\_POTENTIAL} |
| 2315 |
> |
\end{itemize} |
| 2316 |
|
|
| 2317 |
|
The {\tt switchingRadius} is set to a default value of 95\% of the |
| 2318 |
|
{\tt cutoffRadius}. In the special case of a simulation containing |
| 2320 |
|
same value as the cutoff radius, and {\sc OpenMD} will use a shifted |
| 2321 |
|
potential to remove discontinuities in the potential at the cutoff. |
| 2322 |
|
Both radii may be specified in the meta-data file. |
| 2271 |
– |
|
| 2272 |
– |
|
| 2273 |
– |
\section{\label{section:WATER}The {\sc water} Force Field} |
| 2274 |
– |
|
| 2275 |
– |
In addition to the {\sc duff} force field's solvent description, a |
| 2276 |
– |
separate {\sc water} force field has been included for simulating most |
| 2277 |
– |
of the common rigid-body water models. This force field includes the |
| 2278 |
– |
simple and point-dipolar models (SSD, SSD1, SSD/E, SSD/RF, and DPD |
| 2279 |
– |
water), as well as the common charge-based models (SPC, SPC/E, TIP3P, |
| 2280 |
– |
TIP4P, and |
| 2281 |
– |
TIP5P).\cite{liu96:new_model,Ichiye03,fennell04,Marrink01,Berendsen81,Berendsen87,Jorgensen83,Mahoney00} |
| 2282 |
– |
In order to handle these models, charge-charge interactions were |
| 2283 |
– |
included in the force-loop: |
| 2284 |
– |
\begin{equation} |
| 2285 |
– |
V_{\text{charge}}(r_{ij}) = \sum_{ij}\frac{q_iq_je^2}{r_{ij}}, |
| 2286 |
– |
\end{equation} |
| 2287 |
– |
where $q$ represents the charge on particle $i$ or $j$, and $e$ is the |
| 2288 |
– |
charge of an electron in Coulombs. The charge-charge interaction |
| 2289 |
– |
support is rudimentary in the current version of {\sc OpenMD}. As with |
| 2290 |
– |
the other pair interactions, charges can be simulated with a pure |
| 2291 |
– |
cutoff or a reaction field. The various methods for performing the |
| 2292 |
– |
Ewald summation have not yet been included. The {\sc water} force |
| 2293 |
– |
field can be easily expanded through modification of the {\sc water} |
| 2294 |
– |
force field file ({\tt WATER.frc}). By adding atom types and inserting |
| 2295 |
– |
the appropriate parameters, it is possible to extend the force field |
| 2296 |
– |
to handle rigid molecules other than water. |
| 2323 |
|
|
| 2324 |
|
|
| 2299 |
– |
|
| 2325 |
|
\section{\label{section:pbc}Periodic Boundary Conditions} |
| 2326 |
|
|
| 2327 |
|
\newcommand{\roundme}{\operatorname{round}} |
