| 498 |
|
{\tt minimizer} & string & Chooses a minimizer & Possible minimizers |
| 499 |
|
are SD and CG. Either {\tt ensemble} or {\tt minimizer} must be specified. \\ |
| 500 |
|
{\tt ensemble} & string & Sets the ensemble. & Possible ensembles are |
| 501 |
< |
NVE, NVT, NPTi, NPAT, NPTf, and NPTxyz. Either {\tt ensemble} |
| 501 |
> |
NVE, NVT, NPTi, NPAT, NPTf, NPTxyz, and LD. Either {\tt ensemble} |
| 502 |
|
or {\tt minimizer} must be specified. \\ |
| 503 |
|
{\tt dt} & fs & Sets the time step. & Selection of {\tt dt} should be |
| 504 |
|
small enough to sample the fastest motion of the simulation. ({\tt |
| 1905 |
|
& (with changes to box shape) & \\ |
| 1906 |
|
NPTxyz & approximate isobaric-isothermal & {\tt ensemble = NPTxyz;} \\ |
| 1907 |
|
& (with separate barostats on each box dimension) & \\ |
| 1908 |
+ |
LD & Langevin Dynamics & {\tt ensemble = LD;} \\ |
| 1909 |
+ |
& (approximates the effects of an implicit solvent) & \\ |
| 1910 |
|
\end{tabular} |
| 1911 |
|
\end{center} |
| 1912 |
|
|
| 2371 |
|
for use with liquid simulations, or in cases where there is |
| 2372 |
|
orientational anisotropy in the system (i.e. in lipid bilayer |
| 2373 |
|
simulations). |
| 2372 |
– |
|
| 2373 |
– |
|
| 2374 |
– |
\section{\label{sec:constraints}Constraint Methods} |
| 2374 |
|
|
| 2375 |
< |
\subsection{\label{oopseSec:rattle}The {\sc rattle} Method for Bond |
| 2377 |
< |
Constraints} |
| 2378 |
< |
|
| 2379 |
< |
In order to satisfy the constraints of fixed bond lengths within {\sc |
| 2380 |
< |
oopse}, we have implemented the {\sc rattle} algorithm of |
| 2381 |
< |
Andersen.\cite{andersen83} {\sc rattle} is a velocity-Verlet |
| 2382 |
< |
formulation of the {\sc shake} method\cite{ryckaert77} for iteratively |
| 2383 |
< |
solving the Lagrange multipliers which maintain the holonomic |
| 2384 |
< |
constraints. Both methods are covered in depth in the |
| 2385 |
< |
literature,\cite{leach01:mm,Allen87} and a detailed description of |
| 2386 |
< |
this method would be redundant. |
| 2387 |
< |
|
| 2388 |
< |
\subsection{\label{oopseSec:zcons}The Z-Constraint Method} |
| 2389 |
< |
|
| 2390 |
< |
A force auto-correlation method based on the fluctuation-dissipation |
| 2391 |
< |
theorem was developed by Roux and Karplus to investigate the dynamics |
| 2392 |
< |
of ions inside ion channels.\cite{Roux91} The time-dependent friction |
| 2393 |
< |
coefficient can be calculated from the deviation of the instantaneous |
| 2394 |
< |
force from its mean value: |
| 2395 |
< |
\begin{equation} |
| 2396 |
< |
\xi(z,t)=\langle\delta F(z,t)\delta F(z,0)\rangle/k_{B}T, |
| 2397 |
< |
\end{equation} |
| 2398 |
< |
where% |
| 2399 |
< |
\begin{equation} |
| 2400 |
< |
\delta F(z,t)=F(z,t)-\langle F(z,t)\rangle. |
| 2401 |
< |
\end{equation} |
| 2375 |
> |
\section{Langevin Dynamics (LD)\label{LDRB}} |
| 2376 |
|
|
| 2377 |
< |
If the time-dependent friction decays rapidly, the static friction |
| 2378 |
< |
coefficient can be approximated by |
| 2379 |
< |
\begin{equation} |
| 2380 |
< |
\xi_{\text{static}}(z)=\int_{0}^{\infty}\langle\delta F(z,t)\delta F(z,0)\rangle dt. |
| 2381 |
< |
\end{equation} |
| 2382 |
< |
|
| 2409 |
< |
This allows the diffusion constant to then be calculated through the |
| 2410 |
< |
Einstein relation:\cite{Marrink94} |
| 2411 |
< |
\begin{equation} |
| 2412 |
< |
D(z)=\frac{k_{B}T}{\xi_{\text{static}}(z)}=\frac{(k_{B}T)^{2}}{\int_{0}^{\infty |
| 2413 |
< |
}\langle\delta F(z,t)\delta F(z,0)\rangle dt}.% |
| 2414 |
< |
\end{equation} |
| 2415 |
< |
|
| 2416 |
< |
The Z-Constraint method, which fixes the $z$ coordinates of a few |
| 2417 |
< |
``tagged'' molecules with respect to the center of the mass of the |
| 2418 |
< |
system is a technique that was proposed to obtain the forces required |
| 2419 |
< |
for the force auto-correlation calculation.\cite{Marrink94} However, |
| 2420 |
< |
simply resetting the coordinate will move the center of the mass of |
| 2421 |
< |
the whole system. To avoid this problem, we have developed a new |
| 2422 |
< |
method that is utilized in {\sc oopse}. Instead of resetting the |
| 2423 |
< |
coordinates, we reset the forces of $z$-constrained molecules and |
| 2424 |
< |
subtract the total constraint forces from the rest of the system after |
| 2425 |
< |
the force calculation at each time step. |
| 2426 |
< |
|
| 2427 |
< |
After the force calculation, the total force on molecule $\alpha$ is: |
| 2428 |
< |
\begin{equation} |
| 2429 |
< |
G_{\alpha} = \sum_i F_{\alpha i}, |
| 2430 |
< |
\label{oopseEq:zc1} |
| 2431 |
< |
\end{equation} |
| 2432 |
< |
where $F_{\alpha i}$ is the force in the $z$ direction on atom $i$ in |
| 2433 |
< |
$z$-constrained molecule $\alpha$. The forces on the atoms in the |
| 2434 |
< |
$z$-constrained molecule are then adjusted to remove the total force |
| 2435 |
< |
on molecule $\alpha$: |
| 2436 |
< |
\begin{equation} |
| 2437 |
< |
F_{\alpha i} = F_{\alpha i} - |
| 2438 |
< |
\frac{m_{\alpha i} G_{\alpha}}{\sum_i m_{\alpha i}}. |
| 2439 |
< |
\end{equation} |
| 2440 |
< |
Here, $m_{\alpha i}$ is the mass of atom $i$ in the $z$-constrained |
| 2441 |
< |
molecule. After the forces have been adjusted, the velocities must |
| 2442 |
< |
also be modified to subtract out molecule $\alpha$'s center-of-mass |
| 2443 |
< |
velocity in the $z$ direction. |
| 2444 |
< |
\begin{equation} |
| 2445 |
< |
v_{\alpha i} = v_{\alpha i} - |
| 2446 |
< |
\frac{\sum_i m_{\alpha i} v_{\alpha i}}{\sum_i m_{\alpha i}}, |
| 2447 |
< |
\end{equation} |
| 2448 |
< |
where $v_{\alpha i}$ is the velocity of atom $i$ in the $z$ direction. |
| 2449 |
< |
Lastly, all of the accumulated constraint forces must be subtracted |
| 2450 |
< |
from the rest of the unconstrained system to keep the system center of |
| 2451 |
< |
mass of the entire system from drifting. |
| 2452 |
< |
\begin{equation} |
| 2453 |
< |
F_{\beta i} = F_{\beta i} - \frac{m_{\beta i} \sum_{\alpha} G_{\alpha}} |
| 2454 |
< |
{\sum_{\beta}\sum_i m_{\beta i}}, |
| 2455 |
< |
\end{equation} |
| 2456 |
< |
where $\beta$ denotes all {\it unconstrained} molecules in the |
| 2457 |
< |
system. Similarly, the velocities of the unconstrained molecules must |
| 2458 |
< |
also be scaled: |
| 2459 |
< |
\begin{equation} |
| 2460 |
< |
v_{\beta i} = v_{\beta i} + \sum_{\alpha} \frac{\sum_i m_{\alpha i} |
| 2461 |
< |
v_{\alpha i}}{\sum_i m_{\alpha i}}. |
| 2462 |
< |
\end{equation} |
| 2463 |
< |
|
| 2464 |
< |
This method will pin down the centers-of-mass of all of the |
| 2465 |
< |
$z$-constrained molecules, and will also keep the entire system fixed |
| 2466 |
< |
at the original system center-of-mass location. |
| 2377 |
> |
{\sc oopse} implements a Langevin integrator in order to perform |
| 2378 |
> |
molecular dynamics simulations in implicit solvent environments. This |
| 2379 |
> |
can result in substantial performance gains when the detailed dynamics |
| 2380 |
> |
of the solvent is not important. Since {\sc oopse} also handles rigid |
| 2381 |
> |
bodies of arbitrary composition and shape, the Langevin integrator is |
| 2382 |
> |
by necessity somewhat more complex than in other simulation packages. |
| 2383 |
|
|
| 2384 |
< |
At the very beginning of the simulation, the molecules may not be at |
| 2469 |
< |
their desired positions. To steer a $z$-constrained molecule to its |
| 2470 |
< |
specified position, a simple harmonic potential is used: |
| 2384 |
> |
Consider the Langevin equations of motion in generalized coordinates |
| 2385 |
|
\begin{equation} |
| 2472 |
– |
U(t)=\frac{1}{2}k_{\text{Harmonic}}(z(t)-z_{\text{cons}})^{2},% |
| 2473 |
– |
\end{equation} |
| 2474 |
– |
where $k_{\text{Harmonic}}$ is an harmonic force constant, $z(t)$ is |
| 2475 |
– |
the current $z$ coordinate of the center of mass of the constrained |
| 2476 |
– |
molecule, and $z_{\text{cons}}$ is the desired constraint |
| 2477 |
– |
position. The harmonic force operating on the $z$-constrained molecule |
| 2478 |
– |
at time $t$ can be calculated by |
| 2479 |
– |
\begin{equation} |
| 2480 |
– |
F_{z_{\text{Harmonic}}}(t)=-\frac{\partial U(t)}{\partial z(t)}= |
| 2481 |
– |
-k_{\text{Harmonic}}(z(t)-z_{\text{cons}}). |
| 2482 |
– |
\end{equation} |
| 2483 |
– |
|
| 2484 |
– |
The user may also specify the use of a constant velocity method |
| 2485 |
– |
(steered molecular dynamics) to move the molecules to their desired |
| 2486 |
– |
initial positions. Based on concepts from atomic force microscopy, |
| 2487 |
– |
{\sc smd} has been used to study many processes which occur via rare |
| 2488 |
– |
events on the time scale of a few hundreds of picoseconds. For |
| 2489 |
– |
example,{\sc smd} has been used to observe the dissociation of |
| 2490 |
– |
Streptavidin-biotin Complex.\cite{smd} |
| 2491 |
– |
|
| 2492 |
– |
To use of the $z$-constraint method in an {\sc oopse} simulation, the |
| 2493 |
– |
molecules must be specified using the {\tt nZconstraints} keyword in |
| 2494 |
– |
the meta-data file. The other parameters for modifying the behavior |
| 2495 |
– |
of the $z$-constraint method are listed in table~\ref{table:zconParams}. |
| 2496 |
– |
|
| 2497 |
– |
\begin{longtable}[c]{ABCD} |
| 2498 |
– |
\caption{Meta-data Keywords: Z-Constraint Parameters} |
| 2499 |
– |
\\ |
| 2500 |
– |
{\bf keyword} & {\bf units} & {\bf use} & {\bf remarks} \\ \hline |
| 2501 |
– |
\endhead |
| 2502 |
– |
\hline |
| 2503 |
– |
\endfoot |
| 2504 |
– |
{\tt zconsTime} & fs & Sets the frequency at which the {\tt .fz} file |
| 2505 |
– |
is written & \\ |
| 2506 |
– |
{\tt zconsForcePolicy} & string & The strategy for subtracting |
| 2507 |
– |
the $z$-constraint force from the {\it unconstrained} molecules & Possible |
| 2508 |
– |
strategies are {\tt BYMASS} and {\tt BYNUMBER}. The default |
| 2509 |
– |
strategy is {\tt BYMASS}\\ |
| 2510 |
– |
{\tt zconsGap} & $\mbox{\AA}$ & Sets the distance between two adjacent |
| 2511 |
– |
constraint positions&Used mainly to move molecules through a |
| 2512 |
– |
simulation to estimate potentials of mean force. \\ |
| 2513 |
– |
{\tt zconsFixtime} & fs & Sets the length of time the $z$-constraint |
| 2514 |
– |
molecule is held fixed & {\tt zconsFixtime} must be set if {\tt |
| 2515 |
– |
zconsGap} is set\\ |
| 2516 |
– |
{\tt zconsUsingSMD} & logical & Flag for using Steered Molecular |
| 2517 |
– |
Dynamics to move the molecules to the correct constrained positions & |
| 2518 |
– |
Harmonic Forces are used by default |
| 2519 |
– |
\label{table:zconParams} |
| 2520 |
– |
\end{longtable} |
| 2521 |
– |
|
| 2522 |
– |
\section{Langevin Dynamics for Rigid Particles of Arbitrary Shape (LANGEVINDYNAMICS)\label{LDRB}} |
| 2523 |
– |
|
| 2524 |
– |
{\sc oopse} implements langevin dynamics integrator to perform the |
| 2525 |
– |
molecular dynamics simulations at implicit solvents environment to |
| 2526 |
– |
speed up the simulation when the properties of solvents are not |
| 2527 |
– |
important. Consider the Langevin equations of motion in generalized |
| 2528 |
– |
coordinates |
| 2529 |
– |
\begin{equation} |
| 2386 |
|
{\bf M} \dot{{\bf V}}(t) = {\bf F}_{s}(t) + |
| 2387 |
|
{\bf F}_{f}(t) + {\bf F}_{r}(t) |
| 2388 |
|
\label{LDGeneralizedForm} |
| 2436 |
|
2 k_B T \Xi_R \delta(t - t'). \label{eq:randomForce} |
| 2437 |
|
\end{equation} |
| 2438 |
|
$\Xi_R$ is the $6\times6$ resistance tensor at the center of |
| 2439 |
< |
resistance. Once this tensor is known for a given rigid body (as |
| 2440 |
< |
described in the previous section) obtaining a stochastic vector that |
| 2441 |
< |
has the properties in Eq. (\ref{eq:randomForce}) can be done |
| 2442 |
< |
efficiently by carrying out a one-time Cholesky decomposition to |
| 2443 |
< |
obtain the square root matrix of the resistance tensor, |
| 2439 |
> |
resistance. |
| 2440 |
> |
|
| 2441 |
> |
For atoms and ellipsoids, there are good approximations for this |
| 2442 |
> |
tensor that are based on Stokes' law. For arbitrary rigid bodies, the |
| 2443 |
> |
resistance tensor must be pre-computed before Langevin dynamics can be |
| 2444 |
> |
used. The {\sc oopse} distribution contains a utitilty program called |
| 2445 |
> |
Hydro that performs this computation. |
| 2446 |
> |
|
| 2447 |
> |
Once this tensor is known for a given {\tt integrableObject}, |
| 2448 |
> |
obtaining a stochastic vector that has the properties in |
| 2449 |
> |
Eq. (\ref{eq:randomForce}) can be done efficiently by carrying out a |
| 2450 |
> |
one-time Cholesky decomposition to obtain the square root matrix of |
| 2451 |
> |
the resistance tensor, |
| 2452 |
|
\begin{equation} |
| 2453 |
|
\Xi_R = {\bf S} {\bf S}^{T}, |
| 2454 |
|
\label{eq:Cholesky} |
| 2486 |
|
\begin{equation} |
| 2487 |
|
\frac{\partial}{\partial t}{\bf j}(t) = \tau^{~b}(t) |
| 2488 |
|
\end{equation} |
| 2489 |
< |
Embedding the friction and random forces into the the total force and |
| 2490 |
< |
torque, one can integrate the Langevin equations of motion for a rigid |
| 2491 |
< |
body of arbitrary shape in a velocity-Verlet style 2-part algorithm, |
| 2492 |
< |
where $h = \delta t$: |
| 2489 |
> |
By embedding the friction and random forces into the the total force |
| 2490 |
> |
and torque, {\sc oopse} integrates the Langevin equations of motion |
| 2491 |
> |
for a rigid body of arbitrary shape in a velocity-Verlet style 2-part |
| 2492 |
> |
algorithm, where $h = \delta t$: |
| 2493 |
|
|
| 2494 |
|
{\tt move A:} |
| 2495 |
|
\begin{align*} |
| 2594 |
|
+ \frac{h}{2} {\bf \tau}^{~b}(t + h) . |
| 2595 |
|
\end{align*} |
| 2596 |
|
|
| 2597 |
< |
The viscosity of the implicit of solvents must be specified using {\tt |
| 2598 |
< |
viscosity} keywords in the meta-data file to use langevin dynamics |
| 2599 |
< |
integrator. For simple shaped particles (spheres and ellipsoids), no |
| 2600 |
< |
further parameters are necessary. However, there are no analytical |
| 2601 |
< |
solutions for composite shaped particles, the approximate methods have |
| 2602 |
< |
to be applied to get the resistance tensor. The file which contains |
| 2603 |
< |
the information about hydro properties is indicated by {\tt |
| 2604 |
< |
HydroPropFile} keyword in meta-data file. The {\tt HydroPropFile} is |
| 2605 |
< |
precalculated by {\tt Hydro}. |
| 2597 |
> |
The viscosity of the implicit solvent must be specified using the {\tt |
| 2598 |
> |
viscosity} keyword in the meta-data file if the Langevin integrator is |
| 2599 |
> |
selected. For simple particles (spheres and ellipsoids), no further |
| 2600 |
> |
parameters are necessary. Since there are no analytic solutions for |
| 2601 |
> |
the resistance tensors for composite rigid bodies, the approximate |
| 2602 |
> |
tensors for these objects must also be specified in order to use |
| 2603 |
> |
Langevin dynamics. The meta-data file must therefore point to another |
| 2604 |
> |
file which contains the information about the hydrodynamic properties |
| 2605 |
> |
of all complex rigid bodies being used during the simulation. The |
| 2606 |
> |
{\tt HydroPropFile} keyword is used to specify the name of this |
| 2607 |
> |
file. A {\tt HydroPropFile} should be precalculated using the Hydro |
| 2608 |
> |
program that is included in the {\sc oopse} distribution. |
| 2609 |
|
|
| 2610 |
+ |
\begin{longtable}[c]{ABG} |
| 2611 |
+ |
\caption{Meta-data Keywords: Required parameters for the Langevin integrator} |
| 2612 |
+ |
\\ |
| 2613 |
+ |
{\bf keyword} & {\bf units} & {\bf use} \\ \hline |
| 2614 |
+ |
\endhead |
| 2615 |
+ |
\hline |
| 2616 |
+ |
\endfoot |
| 2617 |
+ |
{\tt viscosity} & centipoise & Sets the value of viscosity of the implicit |
| 2618 |
+ |
solvent \\ |
| 2619 |
+ |
{\tt targetTemp} & K & Sets the target temperature of the system. |
| 2620 |
+ |
This parameter must be specified to use Langevin dynamics. \\ |
| 2621 |
+ |
{\tt HydroPropFile} & string & Specifies the name of the resistance |
| 2622 |
+ |
tensor (usually a {\tt .diff} file) which is precalculated by {\tt |
| 2623 |
+ |
Hydro}. This keyworkd is not necessary if the simulation contains only |
| 2624 |
+ |
simple bodies (spheres and ellipsoids). \\ |
| 2625 |
+ |
{\tt beadSize} & $\mbox{\AA}$ & Sets the diameter of the beads to use |
| 2626 |
+ |
when the {\tt RoughShell} model is used to approximate the resistance |
| 2627 |
+ |
tensor. |
| 2628 |
+ |
\label{table:ldParameters} |
| 2629 |
+ |
\end{longtable} |
| 2630 |
+ |
|
| 2631 |
+ |
\section{\label{sec:constraints}Constraint Methods} |
| 2632 |
+ |
|
| 2633 |
+ |
\subsection{\label{oopseSec:rattle}The {\sc rattle} Method for Bond |
| 2634 |
+ |
Constraints} |
| 2635 |
+ |
|
| 2636 |
+ |
In order to satisfy the constraints of fixed bond lengths within {\sc |
| 2637 |
+ |
oopse}, we have implemented the {\sc rattle} algorithm of |
| 2638 |
+ |
Andersen.\cite{andersen83} {\sc rattle} is a velocity-Verlet |
| 2639 |
+ |
formulation of the {\sc shake} method\cite{ryckaert77} for iteratively |
| 2640 |
+ |
solving the Lagrange multipliers which maintain the holonomic |
| 2641 |
+ |
constraints. Both methods are covered in depth in the |
| 2642 |
+ |
literature,\cite{leach01:mm,Allen87} and a detailed description of |
| 2643 |
+ |
this method would be redundant. |
| 2644 |
+ |
|
| 2645 |
+ |
\subsection{\label{oopseSec:zcons}The Z-Constraint Method} |
| 2646 |
+ |
|
| 2647 |
+ |
A force auto-correlation method based on the fluctuation-dissipation |
| 2648 |
+ |
theorem was developed by Roux and Karplus to investigate the dynamics |
| 2649 |
+ |
of ions inside ion channels.\cite{Roux91} The time-dependent friction |
| 2650 |
+ |
coefficient can be calculated from the deviation of the instantaneous |
| 2651 |
+ |
force from its mean value: |
| 2652 |
+ |
\begin{equation} |
| 2653 |
+ |
\xi(z,t)=\langle\delta F(z,t)\delta F(z,0)\rangle/k_{B}T, |
| 2654 |
+ |
\end{equation} |
| 2655 |
+ |
where% |
| 2656 |
+ |
\begin{equation} |
| 2657 |
+ |
\delta F(z,t)=F(z,t)-\langle F(z,t)\rangle. |
| 2658 |
+ |
\end{equation} |
| 2659 |
+ |
|
| 2660 |
+ |
If the time-dependent friction decays rapidly, the static friction |
| 2661 |
+ |
coefficient can be approximated by |
| 2662 |
+ |
\begin{equation} |
| 2663 |
+ |
\xi_{\text{static}}(z)=\int_{0}^{\infty}\langle\delta F(z,t)\delta F(z,0)\rangle dt. |
| 2664 |
+ |
\end{equation} |
| 2665 |
+ |
|
| 2666 |
+ |
This allows the diffusion constant to then be calculated through the |
| 2667 |
+ |
Einstein relation:\cite{Marrink94} |
| 2668 |
+ |
\begin{equation} |
| 2669 |
+ |
D(z)=\frac{k_{B}T}{\xi_{\text{static}}(z)}=\frac{(k_{B}T)^{2}}{\int_{0}^{\infty |
| 2670 |
+ |
}\langle\delta F(z,t)\delta F(z,0)\rangle dt}.% |
| 2671 |
+ |
\end{equation} |
| 2672 |
+ |
|
| 2673 |
+ |
The Z-Constraint method, which fixes the $z$ coordinates of a few |
| 2674 |
+ |
``tagged'' molecules with respect to the center of the mass of the |
| 2675 |
+ |
system is a technique that was proposed to obtain the forces required |
| 2676 |
+ |
for the force auto-correlation calculation.\cite{Marrink94} However, |
| 2677 |
+ |
simply resetting the coordinate will move the center of the mass of |
| 2678 |
+ |
the whole system. To avoid this problem, we have developed a new |
| 2679 |
+ |
method that is utilized in {\sc oopse}. Instead of resetting the |
| 2680 |
+ |
coordinates, we reset the forces of $z$-constrained molecules and |
| 2681 |
+ |
subtract the total constraint forces from the rest of the system after |
| 2682 |
+ |
the force calculation at each time step. |
| 2683 |
+ |
|
| 2684 |
+ |
After the force calculation, the total force on molecule $\alpha$ is: |
| 2685 |
+ |
\begin{equation} |
| 2686 |
+ |
G_{\alpha} = \sum_i F_{\alpha i}, |
| 2687 |
+ |
\label{oopseEq:zc1} |
| 2688 |
+ |
\end{equation} |
| 2689 |
+ |
where $F_{\alpha i}$ is the force in the $z$ direction on atom $i$ in |
| 2690 |
+ |
$z$-constrained molecule $\alpha$. The forces on the atoms in the |
| 2691 |
+ |
$z$-constrained molecule are then adjusted to remove the total force |
| 2692 |
+ |
on molecule $\alpha$: |
| 2693 |
+ |
\begin{equation} |
| 2694 |
+ |
F_{\alpha i} = F_{\alpha i} - |
| 2695 |
+ |
\frac{m_{\alpha i} G_{\alpha}}{\sum_i m_{\alpha i}}. |
| 2696 |
+ |
\end{equation} |
| 2697 |
+ |
Here, $m_{\alpha i}$ is the mass of atom $i$ in the $z$-constrained |
| 2698 |
+ |
molecule. After the forces have been adjusted, the velocities must |
| 2699 |
+ |
also be modified to subtract out molecule $\alpha$'s center-of-mass |
| 2700 |
+ |
velocity in the $z$ direction. |
| 2701 |
+ |
\begin{equation} |
| 2702 |
+ |
v_{\alpha i} = v_{\alpha i} - |
| 2703 |
+ |
\frac{\sum_i m_{\alpha i} v_{\alpha i}}{\sum_i m_{\alpha i}}, |
| 2704 |
+ |
\end{equation} |
| 2705 |
+ |
where $v_{\alpha i}$ is the velocity of atom $i$ in the $z$ direction. |
| 2706 |
+ |
Lastly, all of the accumulated constraint forces must be subtracted |
| 2707 |
+ |
from the rest of the unconstrained system to keep the system center of |
| 2708 |
+ |
mass of the entire system from drifting. |
| 2709 |
+ |
\begin{equation} |
| 2710 |
+ |
F_{\beta i} = F_{\beta i} - \frac{m_{\beta i} \sum_{\alpha} G_{\alpha}} |
| 2711 |
+ |
{\sum_{\beta}\sum_i m_{\beta i}}, |
| 2712 |
+ |
\end{equation} |
| 2713 |
+ |
where $\beta$ denotes all {\it unconstrained} molecules in the |
| 2714 |
+ |
system. Similarly, the velocities of the unconstrained molecules must |
| 2715 |
+ |
also be scaled: |
| 2716 |
+ |
\begin{equation} |
| 2717 |
+ |
v_{\beta i} = v_{\beta i} + \sum_{\alpha} \frac{\sum_i m_{\alpha i} |
| 2718 |
+ |
v_{\alpha i}}{\sum_i m_{\alpha i}}. |
| 2719 |
+ |
\end{equation} |
| 2720 |
+ |
|
| 2721 |
+ |
This method will pin down the centers-of-mass of all of the |
| 2722 |
+ |
$z$-constrained molecules, and will also keep the entire system fixed |
| 2723 |
+ |
at the original system center-of-mass location. |
| 2724 |
+ |
|
| 2725 |
+ |
At the very beginning of the simulation, the molecules may not be at |
| 2726 |
+ |
their desired positions. To steer a $z$-constrained molecule to its |
| 2727 |
+ |
specified position, a simple harmonic potential is used: |
| 2728 |
+ |
\begin{equation} |
| 2729 |
+ |
U(t)=\frac{1}{2}k_{\text{Harmonic}}(z(t)-z_{\text{cons}})^{2},% |
| 2730 |
+ |
\end{equation} |
| 2731 |
+ |
where $k_{\text{Harmonic}}$ is an harmonic force constant, $z(t)$ is |
| 2732 |
+ |
the current $z$ coordinate of the center of mass of the constrained |
| 2733 |
+ |
molecule, and $z_{\text{cons}}$ is the desired constraint |
| 2734 |
+ |
position. The harmonic force operating on the $z$-constrained molecule |
| 2735 |
+ |
at time $t$ can be calculated by |
| 2736 |
+ |
\begin{equation} |
| 2737 |
+ |
F_{z_{\text{Harmonic}}}(t)=-\frac{\partial U(t)}{\partial z(t)}= |
| 2738 |
+ |
-k_{\text{Harmonic}}(z(t)-z_{\text{cons}}). |
| 2739 |
+ |
\end{equation} |
| 2740 |
+ |
|
| 2741 |
+ |
The user may also specify the use of a constant velocity method |
| 2742 |
+ |
(steered molecular dynamics) to move the molecules to their desired |
| 2743 |
+ |
initial positions. Based on concepts from atomic force microscopy, |
| 2744 |
+ |
{\sc smd} has been used to study many processes which occur via rare |
| 2745 |
+ |
events on the time scale of a few hundreds of picoseconds. For |
| 2746 |
+ |
example,{\sc smd} has been used to observe the dissociation of |
| 2747 |
+ |
Streptavidin-biotin Complex.\cite{smd} |
| 2748 |
+ |
|
| 2749 |
+ |
To use of the $z$-constraint method in an {\sc oopse} simulation, the |
| 2750 |
+ |
molecules must be specified using the {\tt nZconstraints} keyword in |
| 2751 |
+ |
the meta-data file. The other parameters for modifying the behavior |
| 2752 |
+ |
of the $z$-constraint method are listed in table~\ref{table:zconParams}. |
| 2753 |
+ |
|
| 2754 |
|
\begin{longtable}[c]{ABCD} |
| 2755 |
< |
\caption{Meta-data Keywords: Langevin Dynamics Parameters} |
| 2755 |
> |
\caption{Meta-data Keywords: Z-Constraint Parameters} |
| 2756 |
|
\\ |
| 2757 |
|
{\bf keyword} & {\bf units} & {\bf use} & {\bf remarks} \\ \hline |
| 2758 |
|
\endhead |
| 2759 |
|
\hline |
| 2760 |
|
\endfoot |
| 2761 |
< |
{\tt viscosity} & centipoise & Sets the value of viscosity of implicit |
| 2762 |
< |
solvents & \\ {\tt HydroPropFile} & string & specifies the resistance |
| 2763 |
< |
tensor file & usually a {\tt .diff} file which is precalculated by |
| 2764 |
< |
{\sc Hydro}. Not necessory for simple shaped particles (spheres and |
| 2765 |
< |
ellipsoids) \\ |
| 2766 |
< |
{\tt beadSize} & $\mbox{\AA}$ & Sets the diameters of |
| 2767 |
< |
beads when {\sc Rough Shell Model} is used to generate the resistance |
| 2768 |
< |
tensor file. \\ |
| 2769 |
< |
\label{table:ldParameters} |
| 2761 |
> |
{\tt zconsTime} & fs & Sets the frequency at which the {\tt .fz} file |
| 2762 |
> |
is written & \\ |
| 2763 |
> |
{\tt zconsForcePolicy} & string & The strategy for subtracting |
| 2764 |
> |
the $z$-constraint force from the {\it unconstrained} molecules & Possible |
| 2765 |
> |
strategies are {\tt BYMASS} and {\tt BYNUMBER}. The default |
| 2766 |
> |
strategy is {\tt BYMASS}\\ |
| 2767 |
> |
{\tt zconsGap} & $\mbox{\AA}$ & Sets the distance between two adjacent |
| 2768 |
> |
constraint positions&Used mainly to move molecules through a |
| 2769 |
> |
simulation to estimate potentials of mean force. \\ |
| 2770 |
> |
{\tt zconsFixtime} & fs & Sets the length of time the $z$-constraint |
| 2771 |
> |
molecule is held fixed & {\tt zconsFixtime} must be set if {\tt |
| 2772 |
> |
zconsGap} is set\\ |
| 2773 |
> |
{\tt zconsUsingSMD} & logical & Flag for using Steered Molecular |
| 2774 |
> |
Dynamics to move the molecules to the correct constrained positions & |
| 2775 |
> |
Harmonic Forces are used by default |
| 2776 |
> |
\label{table:zconParams} |
| 2777 |
|
\end{longtable} |
| 2778 |
|
|
| 2779 |
|
\chapter{\label{oopseSec:thermInt}Thermodynamic Integration} |
| 3473 |
|
\end{longtable} |
| 3474 |
|
|
| 3475 |
|
\section{\label{oopseSec:Hydro}Hydro} |
| 3476 |
< |
{\tt Hydro} generates {\tt .diff} file which is required when a Langevin |
| 3477 |
< |
Dynamics simulation using approximate models (supports Bead Model and |
| 3478 |
< |
Rough Shell Model) is performed. To generate the {\tt }.diff file, the |
| 3479 |
< |
meta-data file is needed as the input file. The viscosity of the fluid |
| 3480 |
< |
flow (solvent) and the temperature of the system have to be defined in |
| 3481 |
< |
meta-data file. If the approximate model is {\tt Rough Shell Model}, |
| 3482 |
< |
the {\tt beadSize} which is the diameter of every beads must be |
| 3483 |
< |
specified in meta-data file. |
| 3476 |
> |
{\tt Hydro} generates resistance tensor ({\tt .diff}) files which are |
| 3477 |
> |
required when using the Langevin integrator using complex rigid |
| 3478 |
> |
bodies. {\tt Hydro} supports two approximate models: the {\tt |
| 3479 |
> |
BeadModel} and {\tt RoughShell}. Additionally, {\tt Hydro} can |
| 3480 |
> |
generate resistance tensor files using analytic solutions for simple |
| 3481 |
> |
shapes. To generate a {\tt }.diff file, a meta-data file is needed as |
| 3482 |
> |
the input file. Since the resistance tensor depends on these |
| 3483 |
> |
quantities, the {\tt viscosity} of the solvent and the temperature |
| 3484 |
> |
({\tt targetTemp}) of the system must be defined in meta-data file. If |
| 3485 |
> |
the approximate model in use is the {\tt RoughShell} model the {\tt |
| 3486 |
> |
beadSize} (the diameter of the small beads used to approximate the |
| 3487 |
> |
surface of the body) must also be specified. |
| 3488 |
|
|
| 3489 |
|
The options available for Hydro are as follows: |
| 3490 |
|
\begin{longtable}[c]{|EFG|} |
| 3498 |
|
-V& {\tt -{}-version} & Print version and exit\\ |
| 3499 |
|
-i& {\tt -{}-input=filename} & input MetaData (md) file\\ |
| 3500 |
|
-o& {\tt -{}-output=STRING} & Output file name\\ |
| 3501 |
< |
& {\tt -{}-model=STRING} & hydrodynamics model (supports |
| 3502 |
< |
RoughShell and BeadModel)\\ |
| 3501 |
> |
& {\tt -{}-model=STRING} & hydrodynamics model (supports both |
| 3502 |
> |
{\tt RoughShell} and {\tt BeadModel})\\ |
| 3503 |
|
-b& {\tt -{}-beads} & generate the beads only, |
| 3504 |
< |
hydrodynamics will be performed (default=off)\\ |
| 3504 |
> |
hydrodynamic calculations will not be performed (default=off)\\ |
| 3505 |
|
\end{longtable} |
| 3506 |
|
|
| 3507 |
|
|
