ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/tengDissertation/Introduction.tex
(Generate patch)

Comparing trunk/tengDissertation/Introduction.tex (file contents):
Revision 2700 by tim, Tue Apr 11 03:38:09 2006 UTC vs.
Revision 2702 by tim, Wed Apr 12 03:37:46 2006 UTC

# Line 27 | Line 27 | $F_ij$ be the force that particle $i$ exerts on partic
27   \end{equation}
28   A point mass interacting with other bodies moves with the
29   acceleration along the direction of the force acting on it. Let
30 < $F_ij$ be the force that particle $i$ exerts on particle $j$, and
31 < $F_ji$ be the force that particle $j$ exerts on particle $i$.
30 > $F_{ij}$ be the force that particle $i$ exerts on particle $j$, and
31 > $F_{ji}$ be the force that particle $j$ exerts on particle $i$.
32   Newton¡¯s third law states that
33   \begin{equation}
34 < F_ij = -F_ji
34 > F_{ij} = -F_{ji}
35   \label{introEquation:newtonThirdLaw}
36   \end{equation}
37  
# Line 396 | Line 396 | Liouville's theorem can be expresses in a variety of d
396   \label{introEquation:densityAndHamiltonian}
397   \end{equation}
398  
399 + \subsubsection{\label{introSection:phaseSpaceConservation}Conservation of Phase Space}
400 + Lets consider a region in the phase space,
401 + \begin{equation}
402 + \delta v = \int { \ldots \int {dq_1 } ...dq_f dp_1 } ..dp_f .
403 + \end{equation}
404 + If this region is small enough, the density $\rho$ can be regarded
405 + as uniform over the whole phase space. Thus, the number of phase
406 + points inside this region is given by,
407 + \begin{equation}
408 + \delta N = \rho \delta v = \rho \int { \ldots \int {dq_1 } ...dq_f
409 + dp_1 } ..dp_f.
410 + \end{equation}
411 +
412 + \begin{equation}
413 + \frac{{d(\delta N)}}{{dt}} = \frac{{d\rho }}{{dt}}\delta v + \rho
414 + \frac{d}{{dt}}(\delta v) = 0.
415 + \end{equation}
416 + With the help of stationary assumption
417 + (\ref{introEquation:stationary}), we obtain the principle of the
418 + \emph{conservation of extension in phase space},
419 + \begin{equation}
420 + \frac{d}{{dt}}(\delta v) = \frac{d}{{dt}}\int { \ldots \int {dq_1 }
421 + ...dq_f dp_1 } ..dp_f  = 0.
422 + \label{introEquation:volumePreserving}
423 + \end{equation}
424 +
425 + \subsubsection{\label{introSection:liouvilleInOtherForms}Liouville's Theorem in Other Forms}
426 +
427   Liouville's theorem can be expresses in a variety of different forms
428   which are convenient within different contexts. For any two function
429   $F$ and $G$ of the coordinates and momenta of a system, the Poisson
# Line 431 | Line 459 | expressed as
459   \label{introEquation:liouvilleTheoremInOperator}
460   \end{equation}
461  
434
462   \subsection{\label{introSection:ergodic}The Ergodic Hypothesis}
463  
464   Various thermodynamic properties can be calculated from Molecular
# Line 560 | Line 587 | H = \frac{1}{2}\left( {\frac{{\pi _1^2 }}{{I_1 }} + \f
587   }}{{I_2 }} + \frac{{\pi _3^2 }}{{I_3 }}} \right)
588   \end{equation}
589  
590 < \subsection{\label{introSection:geometricProperties}Geometric Properties}
590 > \subsection{\label{introSection:exactFlow}Exact Flow}
591 >
592   Let $x(t)$ be the exact solution of the ODE system,
593   \begin{equation}
594   \frac{{dx}}{{dt}} = f(x) \label{introEquation:ODE}
# Line 570 | Line 598 | space to itself. In most cases, it is not easy to find
598   x(t+\tau) =\varphi_\tau(x(t))
599   \]
600   where $\tau$ is a fixed time step and $\varphi$ is a map from phase
601 < space to itself. In most cases, it is not easy to find the exact
574 < flow $\varphi_\tau$. Instead, we use a approximate map, $\psi_\tau$,
575 < which is usually called integrator. The order of an integrator
576 < $\psi_\tau$ is $p$, if the Taylor series of $\psi_\tau$ agree to
577 < order $p$,
601 > space to itself. The flow has the continuous group property,
602   \begin{equation}
603 + \varphi _{\tau _1 }  \circ \varphi _{\tau _2 }  = \varphi _{\tau _1
604 + + \tau _2 } .
605 + \end{equation}
606 + In particular,
607 + \begin{equation}
608 + \varphi _\tau   \circ \varphi _{ - \tau }  = I
609 + \end{equation}
610 + Therefore, the exact flow is self-adjoint,
611 + \begin{equation}
612 + \varphi _\tau   = \varphi _{ - \tau }^{ - 1}.
613 + \end{equation}
614 + The exact flow can also be written in terms of the of an operator,
615 + \begin{equation}
616 + \varphi _\tau  (x) = e^{\tau \sum\limits_i {f_i (x)\frac{\partial
617 + }{{\partial x_i }}} } (x) \equiv \exp (\tau f)(x).
618 + \label{introEquation:exponentialOperator}
619 + \end{equation}
620 +
621 + In most cases, it is not easy to find the exact flow $\varphi_\tau$.
622 + Instead, we use a approximate map, $\psi_\tau$, which is usually
623 + called integrator. The order of an integrator $\psi_\tau$ is $p$, if
624 + the Taylor series of $\psi_\tau$ agree to order $p$,
625 + \begin{equation}
626   \psi_tau(x) = x + \tau f(x) + O(\tau^{p+1})
627   \end{equation}
628  
629 + \subsection{\label{introSection:geometricProperties}Geometric Properties}
630 +
631   The hidden geometric properties of ODE and its flow play important
632 < roles in numerical studies. Let $\varphi$ be the flow of Hamiltonian
633 < vector field, $\varphi$ is a \emph{symplectic} flow if it satisfies,
632 > roles in numerical studies. Many of them can be found in systems
633 > which occur naturally in applications.
634 >
635 > Let $\varphi$ be the flow of Hamiltonian vector field, $\varphi$ is
636 > a \emph{symplectic} flow if it satisfies,
637   \begin{equation}
638   '\varphi^T J '\varphi = J.
639   \end{equation}
# Line 593 | Line 645 | is the property must be preserved by the integrator. I
645   \begin{equation}
646   '\varphi ^T J '\varphi  = J \circ \varphi
647   \end{equation}
648 < is the property must be preserved by the integrator. It is possible
649 < to construct a \emph{volume-preserving} flow for a source free($
650 < \nabla \cdot f = 0 $) ODE, if the flow satisfies $ \det d\varphi  =
651 < 1$. Changing the variables $y = h(x)$ in a
652 < ODE\ref{introEquation:ODE} will result in a new system,
648 > is the property must be preserved by the integrator.
649 >
650 > It is possible to construct a \emph{volume-preserving} flow for a
651 > source free($ \nabla \cdot f = 0 $) ODE, if the flow satisfies $
652 > \det d\varphi  = 1$. One can show easily that a symplectic flow will
653 > be volume-preserving.
654 >
655 > Changing the variables $y = h(x)$ in a ODE\ref{introEquation:ODE}
656 > will result in a new system,
657   \[
658   \dot y = \tilde f(y) = ((dh \cdot f)h^{ - 1} )(y).
659   \]
660   The vector filed $f$ has reversing symmetry $h$ if $f = - \tilde f$.
661   In other words, the flow of this vector field is reversible if and
662 < only if $ h \circ \varphi ^{ - 1}  = \varphi  \circ h $. When
607 < designing any numerical methods, one should always try to preserve
608 < the structural properties of the original ODE and its flow.
662 > only if $ h \circ \varphi ^{ - 1}  = \varphi  \circ h $.
663  
664 + When designing any numerical methods, one should always try to
665 + preserve the structural properties of the original ODE and its flow.
666 +
667   \subsection{\label{introSection:constructionSymplectic}Construction of Symplectic Methods}
668   A lot of well established and very effective numerical methods have
669   been successful precisely because of their symplecticities even
# Line 622 | Line 679 | and difficult to use\cite{}. In dissipative systems, v
679   \end{enumerate}
680  
681   Generating function tends to lead to methods which are cumbersome
682 < and difficult to use\cite{}. In dissipative systems, variational
683 < methods can capture the decay of energy accurately\cite{}. Since
684 < their geometrically unstable nature against non-Hamiltonian
685 < perturbations, ordinary implicit Runge-Kutta methods are not
686 < suitable for Hamiltonian system. Recently, various high-order
687 < explicit Runge--Kutta methods have been developed to overcome this
682 > and difficult to use. In dissipative systems, variational methods
683 > can capture the decay of energy accurately. Since their
684 > geometrically unstable nature against non-Hamiltonian perturbations,
685 > ordinary implicit Runge-Kutta methods are not suitable for
686 > Hamiltonian system. Recently, various high-order explicit
687 > Runge--Kutta methods have been developed to overcome this
688   instability \cite{}. However, due to computational penalty involved
689   in implementing the Runge-Kutta methods, they do not attract too
690   much attention from Molecular Dynamics community. Instead, splitting
691   have been widely accepted since they exploit natural decompositions
692 < of the system\cite{Tuckerman92}. The main idea behind splitting
693 < methods is to decompose the discrete $\varphi_h$ as a composition of
694 < simpler flows,
692 > of the system\cite{Tuckerman92}.
693 >
694 > \subsubsection{\label{introSection:splittingMethod}Splitting Method}
695 >
696 > The main idea behind splitting methods is to decompose the discrete
697 > $\varphi_h$ as a composition of simpler flows,
698   \begin{equation}
699   \varphi _h  = \varphi _{h_1 }  \circ \varphi _{h_2 }  \ldots  \circ
700   \varphi _{h_n }
701   \label{introEquation:FlowDecomposition}
702   \end{equation}
703   where each of the sub-flow is chosen such that each represent a
704 < simpler integration of the system. Let $\phi$ and $\psi$ both be
705 < symplectic maps, it is easy to show that any composition of
706 < symplectic flows yields a symplectic map,
704 > simpler integration of the system.
705 >
706 > Suppose that a Hamiltonian system takes the form,
707 > \[
708 > H = H_1 + H_2.
709 > \]
710 > Here, $H_1$ and $H_2$ may represent different physical processes of
711 > the system. For instance, they may relate to kinetic and potential
712 > energy respectively, which is a natural decomposition of the
713 > problem. If $H_1$ and $H_2$ can be integrated using exact flows
714 > $\varphi_1(t)$ and $\varphi_2(t)$, respectively, a simple first
715 > order is then given by the Lie-Trotter formula
716   \begin{equation}
717 + \varphi _h  = \varphi _{1,h}  \circ \varphi _{2,h},
718 + \label{introEquation:firstOrderSplitting}
719 + \end{equation}
720 + where $\varphi _h$ is the result of applying the corresponding
721 + continuous $\varphi _i$ over a time $h$. By definition, as
722 + $\varphi_i(t)$ is the exact solution of a Hamiltonian system, it
723 + must follow that each operator $\varphi_i(t)$ is a symplectic map.
724 + It is easy to show that any composition of symplectic flows yields a
725 + symplectic map,
726 + \begin{equation}
727   (\varphi '\phi ')^T J\varphi '\phi ' = \phi '^T \varphi '^T J\varphi
728 < '\phi ' = \phi '^T J\phi ' = J.
728 > '\phi ' = \phi '^T J\phi ' = J,
729   \label{introEquation:SymplecticFlowComposition}
730   \end{equation}
731 < Suppose that a Hamiltonian system has a form with $H = T + V$
731 > where $\phi$ and $\psi$ both are symplectic maps. Thus operator
732 > splitting in this context automatically generates a symplectic map.
733  
734 + The Lie-Trotter splitting(\ref{introEquation:firstOrderSplitting})
735 + introduces local errors proportional to $h^2$, while Strang
736 + splitting gives a second-order decomposition,
737 + \begin{equation}
738 + \varphi _h  = \varphi _{1,h/2}  \circ \varphi _{2,h}  \circ \varphi
739 + _{1,h/2} ,
740 + \label{introEqaution:secondOrderSplitting}
741 + \end{equation}
742 + which has a local error proportional to $h^3$. Sprang splitting's
743 + popularity in molecular simulation community attribute to its
744 + symmetric property,
745 + \begin{equation}
746 + \varphi _h^{ - 1} = \varphi _{ - h}.
747 + \lable{introEquation:timeReversible}
748 + \end{equation}
749 +
750 + \subsubsection{\label{introSection:exampleSplittingMethod}Example of Splitting Method}
751 + The classical equation for a system consisting of interacting
752 + particles can be written in Hamiltonian form,
753 + \[
754 + H = T + V
755 + \]
756 + where $T$ is the kinetic energy and $V$ is the potential energy.
757 + Setting $H_1 = T, H_2 = V$ and applying Strang splitting, one
758 + obtains the following:
759 + \begin{align}
760 + q(\Delta t) &= q(0) + \dot{q}(0)\Delta t +
761 +    \frac{F[q(0)]}{m}\frac{\Delta t^2}{2}, %
762 + \label{introEquation:Lp10a} \\%
763 + %
764 + \dot{q}(\Delta t) &= \dot{q}(0) + \frac{\Delta t}{2m}
765 +    \biggl [F[q(0)] + F[q(\Delta t)] \biggr]. %
766 + \label{introEquation:Lp10b}
767 + \end{align}
768 + where $F(t)$ is the force at time $t$. This integration scheme is
769 + known as \emph{velocity verlet} which is
770 + symplectic(\ref{introEquation:SymplecticFlowComposition}),
771 + time-reversible(\ref{introEquation:timeReversible}) and
772 + volume-preserving (\ref{introEquation:volumePreserving}). These
773 + geometric properties attribute to its long-time stability and its
774 + popularity in the community. However, the most commonly used
775 + velocity verlet integration scheme is written as below,
776 + \begin{align}
777 + \dot{q}\biggl (\frac{\Delta t}{2}\biggr ) &=
778 +    \dot{q}(0) + \frac{\Delta t}{2m}\, F[q(0)], \label{introEquation:Lp9a}\\%
779 + %
780 + q(\Delta t) &= q(0) + \Delta t\, \dot{q}\biggl (\frac{\Delta t}{2}\biggr ),%
781 +    \label{introEquation:Lp9b}\\%
782 + %
783 + \dot{q}(\Delta t) &= \dot{q}\biggl (\frac{\Delta t}{2}\biggr ) +
784 +    \frac{\Delta t}{2m}\, F[q(0)]. \label{introEquation:Lp9c}
785 + \end{align}
786 + From the preceding splitting, one can see that the integration of
787 + the equations of motion would follow:
788 + \begin{enumerate}
789 + \item calculate the velocities at the half step, $\frac{\Delta t}{2}$, from the forces calculated at the initial position.
790 +
791 + \item Use the half step velocities to move positions one whole step, $\Delta t$.
792 +
793 + \item Evaluate the forces at the new positions, $\mathbf{r}(\Delta t)$, and use the new forces to complete the velocity move.
794 +
795 + \item Repeat from step 1 with the new position, velocities, and forces assuming the roles of the initial values.
796 + \end{enumerate}
797 +
798 + Simply switching the order of splitting and composing, a new
799 + integrator, the \emph{position verlet} integrator, can be generated,
800 + \begin{align}
801 + \dot q(\Delta t) &= \dot q(0) + \Delta tF(q(0))\left[ {q(0) +
802 + \frac{{\Delta t}}{{2m}}\dot q(0)} \right], %
803 + \label{introEquation:positionVerlet1} \\%
804 + %
805 + q(\Delta t) = q(0) + \frac{{\Delta t}}{2}\left[ {\dot q(0) + \dot
806 + q(\Delta t)} \right]. %
807 + \label{introEquation:positionVerlet1}
808 + \end{align}
809 +
810 + \subsubsection{\label{introSection:errorAnalysis}Error Analysis and Higher Order Methods}
811 +
812 + Baker-Campbell-Hausdorff formula can be used to determine the local
813 + error of splitting method in terms of commutator of the
814 + operators(\ref{introEquation:exponentialOperator}) associated with
815 + the sub-flow. For operators $hX$ and $hY$ which are associate to
816 + $\varphi_1(t)$ and $\varphi_2(t$ respectively , we have
817 + \begin{equation}
818 + \exp (hX + hY) = \exp (hZ)
819 + \end{equation}
820 + where
821 + \begin{equation}
822 + hZ = hX + hY + \frac{{h^2 }}{2}[X,Y] + \frac{{h^3 }}{2}\left(
823 + {[X,[X,Y]] + [Y,[Y,X]]} \right) +  \ldots .
824 + \end{equation}
825 + Here, $[X,Y]$ is the commutators of operator $X$ and $Y$ given by
826 + \[
827 + [X,Y] = XY - YX .
828 + \]
829 + Applying Baker-Campbell-Hausdorff formula to Sprang splitting, we
830 + can obtain
831 + \begin{eqnarray}
832 + \exp (h X/2)\exp (h Y)\exp (h X/2) & = & \exp (h X + h Y + h^2
833 + [X,Y]/4 + h^2 [Y,X]/4 \\ & & \mbox{} + h^2 [X,X]/8 + h^2 [Y,Y]/8 +
834 + h^3 [Y,[Y,X]]/12 - h^3 [X,[X,Y]]/24 +  \ldots )
835 + \end{eqnarray}
836 + Since \[ [X,Y] + [Y,X] = 0\] and \[ [X,X] = 0\], the dominant local
837 + error of Spring splitting is proportional to $h^3$. The same
838 + procedure can be applied to general splitting,  of the form
839 + \begin{equation}
840 + \varphi _{b_m h}^2  \circ \varphi _{a_m h}^1  \circ \varphi _{b_{m -
841 + 1} h}^2  \circ  \ldots  \circ \varphi _{a_1 h}^1 .
842 + \end{equation}
843 + Careful choice of coefficient $a_1 ,\ldot , b_m$ will lead to higher
844 + order method. Yoshida proposed an elegant way to compose higher
845 + order methods based on symmetric splitting. Given a symmetric second
846 + order base method $ \varphi _h^{(2)} $, a fourth-order symmetric
847 + method can be constructed by composing,
848 + \[
849 + \varphi _h^{(4)}  = \varphi _{\alpha h}^{(2)}  \circ \varphi _{\beta
850 + h}^{(2)}  \circ \varphi _{\alpha h}^{(2)}
851 + \]
852 + where $ \alpha  =  - \frac{{2^{1/3} }}{{2 - 2^{1/3} }}$ and $ \beta
853 + = \frac{{2^{1/3} }}{{2 - 2^{1/3} }}$. Moreover, a symmetric
854 + integrator $ \varphi _h^{(2n + 2)}$ can be composed by
855 + \begin{equation}
856 + \varphi _h^{(2n + 2)}  = \varphi _{\alpha h}^{(2n)}  \circ \varphi
857 + _{\beta h}^{(2n)}  \circ \varphi _{\alpha h}^{(2n)}
858 + \end{equation}
859 + , if the weights are chosen as
860 + \[
861 + \alpha  =  - \frac{{2^{1/(2n + 1)} }}{{2 - 2^{1/(2n + 1)} }},\beta =
862 + \frac{{2^{1/(2n + 1)} }}{{2 - 2^{1/(2n + 1)} }} .
863 + \]
864 +
865   \section{\label{introSection:molecularDynamics}Molecular Dynamics}
866  
867   As a special discipline of molecular modeling, Molecular dynamics
# Line 677 | Line 888 | Applications of dynamics of rigid bodies.
888  
889   \subsection{\label{introSection:otherRBMotionEquation}Other Formulations for Rigid Body Motion}
890  
680 %\subsection{\label{introSection:poissonBrackets}Poisson Brackets}
681
891   \section{\label{introSection:correlationFunctions}Correlation Functions}
892  
893   \section{\label{introSection:langevinDynamics}Langevin Dynamics}
# Line 729 | Line 938 | introEquation:motionHamiltonianMomentum},
938   \dot p &=  - \frac{{\partial H}}{{\partial x}}
939         &= m\ddot x
940         &= - \frac{{\partial W(x)}}{{\partial x}} - \sum\limits_{\alpha  = 1}^N {g_\alpha  \left( {x_\alpha   - \frac{{g_\alpha  }}{{m_\alpha  w_\alpha ^2 }}x} \right)}
941 < \label{introEq:Lp5}
941 > \label{introEquation:Lp5}
942   \end{align}
943   , and
944   \begin{align}

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines