41 |
|
applies an external pressure to the facets comprising the convex |
42 |
|
hull surrounding the objects in the system. Additionally, a Langevin |
43 |
|
thermostat is applied to facets of the hull to mimic contact with an |
44 |
< |
external heat bath. This new method, the ``Langevin Hull'', |
45 |
< |
performs better than traditional affine transform methods for |
46 |
< |
systems containing heterogeneous mixtures of materials with |
47 |
< |
different compressibilities. It does not suffer from the edge |
48 |
< |
effects of boundary potential methods, and allows realistic |
49 |
< |
treatment of both external pressure and thermal conductivity to an |
50 |
< |
implicit solvents. We apply this method to several different |
51 |
< |
systems including bare nano-particles, nano-particles in explicit |
52 |
< |
solvent, as well as clusters of liquid water and ice. The predicted |
53 |
< |
mechanical and thermal properties of these systems are in good |
54 |
< |
agreement with experimental data. |
44 |
> |
external heat bath. This new method, the ``Langevin Hull'', performs |
45 |
> |
better than traditional affine transform methods for systems |
46 |
> |
containing heterogeneous mixtures of materials with different |
47 |
> |
compressibilities. It does not suffer from the edge effects of |
48 |
> |
boundary potential methods, and allows realistic treatment of both |
49 |
> |
external pressure and thermal conductivity to an implicit solvent. |
50 |
> |
We apply this method to several different systems including bare |
51 |
> |
nanoparticles, nanoparticles in an explicit solvent, as well as |
52 |
> |
clusters of liquid water and ice. The predicted mechanical and |
53 |
> |
thermal properties of these systems are in good agreement with |
54 |
> |
experimental data. |
55 |
|
\end{abstract} |
56 |
|
|
57 |
|
\newpage |
76 |
|
well as the scaled particle positions (but not the sizes of the |
77 |
|
particles). The most common constant pressure methods, including the |
78 |
|
Melchionna modification\cite{Melchionna1993} to the |
79 |
< |
Nos\'e-Hoover-Andersen equations of motion, the Berendsen pressure |
80 |
< |
bath, and the Langevin Piston, all utilize coordinate transformation |
81 |
< |
to adjust the box volume. |
79 |
> |
Nos\'e-Hoover-Andersen equations of |
80 |
> |
motion,\cite{Hoover85,ANDERSEN:1980vn,Sturgeon:2000kx} the Berendsen |
81 |
> |
pressure bath,\cite{ISI:A1984TQ73500045} and the Langevin |
82 |
> |
Piston,\cite{FELLER:1995fk,Jakobsen:2005uq} all utilize coordinate |
83 |
> |
transformation to adjust the box volume. |
84 |
|
|
85 |
+ |
As long as the material in the simulation box is essentially a bulk |
86 |
+ |
liquid which has a relatively uniform compressibility, the standard |
87 |
+ |
approach provides an excellent way of adjusting the volume of the |
88 |
+ |
system and applying pressure directly via the interactions between |
89 |
+ |
atomic sites. |
90 |
+ |
|
91 |
+ |
The problem with these approaches becomes apparent when the material |
92 |
+ |
being simulated is an inhomogeneous mixture in which portions of the |
93 |
+ |
simulation box are incompressible relative to other portions. |
94 |
+ |
Examples include simulations of metallic nanoparticles in liquid |
95 |
+ |
environments, proteins at interfaces, as well as other multi-phase or |
96 |
+ |
interfacial environments. In these cases, the affine transform of |
97 |
+ |
atomic coordinates will either cause numerical instability when the |
98 |
+ |
sites in the incompressible medium collide with each other, or lead to |
99 |
+ |
inefficient sampling of system volumes if the barostat is set slow |
100 |
+ |
enough to avoid collisions in the incompressible region. |
101 |
+ |
|
102 |
|
\begin{figure} |
103 |
|
\includegraphics[width=\linewidth]{AffineScale2} |
104 |
|
\caption{Affine Scaling constant pressure methods use box-length |
111 |
|
\label{affineScale} |
112 |
|
\end{figure} |
113 |
|
|
114 |
+ |
Additionally, one may often wish to simulate explicitly non-periodic |
115 |
+ |
systems, and the constraint that a periodic box must be used to |
116 |
|
|
96 |
– |
Heterogeneous mixtures of materials with different compressibilities? |
97 |
– |
|
117 |
|
Explicitly non-periodic systems |
118 |
|
|
119 |
|
Elastic Bag |
122 |
|
|
123 |
|
\section{Methodology} |
124 |
|
|
125 |
< |
A new method which uses a constant pressure and temperature bath that |
126 |
< |
interacts with the objects that are currently at the edge of the |
127 |
< |
system. |
125 |
> |
We have developed a new method which uses a constant pressure and |
126 |
> |
temperature bath. This bath interacts only with the objects that are |
127 |
> |
currently at the edge of the system. Since the edge is determined |
128 |
> |
dynamically as the simulation progresses, no {\it a priori} geometry |
129 |
> |
is defined. The pressure and temperature bath interacts {\it |
130 |
> |
directly} with the atoms on the edge and not with atoms interior to |
131 |
> |
the simulation. This means that there are no affine transforms |
132 |
> |
required. There are also no fictitious particles or bounding |
133 |
> |
potentials used in this approach. |
134 |
|
|
135 |
< |
Novel features: No a priori geometry is defined, No affine transforms, |
136 |
< |
No fictitious particles, No bounding potentials. |
135 |
> |
The basics of the method are as follows. The simulation starts as a |
136 |
> |
collection of atomic locations in three dimensions (a point cloud). |
137 |
> |
Delaunay triangulation is used to find all facets between coplanar |
138 |
> |
neighbors. In highly symmetric point clouds, facets can contain many |
139 |
> |
atoms, but in all but the most symmetric of cases one might experience |
140 |
> |
in a molecular dynamics simulation, the facets are simple triangles in |
141 |
> |
3-space that contain exactly three atoms. |
142 |
|
|
143 |
< |
Simulation starts as a collection of atomic locations in 3D (a point |
144 |
< |
cloud). |
145 |
< |
|
146 |
< |
Delaunay triangulation finds all facets between coplanar neighbors. |
147 |
< |
|
148 |
< |
The Convex Hull is the set of facets that have no concave corners at a |
119 |
< |
vertex. |
120 |
< |
|
121 |
< |
Molecules on the convex hull are dynamic. As they re-enter the |
122 |
< |
cluster, all interactions with the external bath are removed.The |
123 |
< |
external bath applies pressure to the facets of the convex hull in |
124 |
< |
direct proportion to the area of the facet. Thermal coupling depends on |
125 |
< |
the solvent temperature, friction and the size and shape of each |
126 |
< |
facet. |
143 |
> |
The convex hull is the set of facets that have {\it no concave |
144 |
> |
corners} at an atomic site. This eliminates all facets on the |
145 |
> |
interior of the point cloud, leaving only those exposed to the |
146 |
> |
bath. Sites on the convex hull are dynamic. As molecules re-enter the |
147 |
> |
cluster, all interactions between atoms on that molecule and the |
148 |
> |
external bath are removed. |
149 |
|
|
150 |
+ |
For atomic sites in the interior of the point cloud, the equations of |
151 |
+ |
motion are simple Newtonian dynamics, |
152 |
|
\begin{equation} |
153 |
< |
m_i \dot{\mathbf v}_i(t)=-{\mathbf \nabla}_i U |
153 |
> |
m_i \dot{\mathbf v}_i(t)=-{\mathbf \nabla}_i U, |
154 |
> |
\label{eq:Newton} |
155 |
|
\end{equation} |
156 |
< |
|
156 |
> |
where $m_i$ is the mass of site $i$, ${\mathbf v}_i(t)$ is the |
157 |
> |
instantaneous velocity of site $i$ at time $t$, and $U$ is the total |
158 |
> |
potential energy. For atoms on the exterior of the cluster |
159 |
> |
(i.e. those that occupy one of the vertices of the convex hull), the |
160 |
> |
equation of motion is modified with an external force, ${\mathbf |
161 |
> |
F}_i^{\mathrm ext}$, |
162 |
|
\begin{equation} |
163 |
< |
m_i \dot{\mathbf v}_i(t)=-{\mathbf \nabla}_i U + {\mathbf F}_i^{\mathrm ext} |
163 |
> |
m_i \dot{\mathbf v}_i(t)=-{\mathbf \nabla}_i U + {\mathbf F}_i^{\mathrm ext}. |
164 |
|
\end{equation} |
165 |
|
|
166 |
+ |
The external bath interacts directly with the facets of the convex |
167 |
+ |
hull. Since each vertex (or atom) provides one corner of a triangular |
168 |
+ |
facet, the force on the facets are divided equally to each vertex. |
169 |
+ |
However, each vertex can participate in multiple facets, so the resultant |
170 |
+ |
force is a sum over all facets $f$ containing vertex $i$: |
171 |
|
\begin{equation} |
172 |
|
{\mathbf F}_{i}^{\mathrm ext} = \sum_{\begin{array}{c}\mathrm{facets\ |
173 |
|
} f \\ \mathrm{containing\ } i\end{array}} \frac{1}{3}\ {\mathbf |
174 |
|
F}_f^{\mathrm ext} |
175 |
|
\end{equation} |
176 |
|
|
177 |
+ |
The external pressure bath applies a force to the facets of the convex |
178 |
+ |
hull in direct proportion to the area of the facet, while the thermal |
179 |
+ |
coupling depends on the solvent temperature, friction and the size and |
180 |
+ |
shape of each facet. The thermal interactions are expressed as a |
181 |
+ |
typical Langevin description of the forces, |
182 |
|
\begin{equation} |
183 |
|
\begin{array}{rclclcl} |
184 |
|
{\mathbf F}_f^{\text{ext}} & = & \text{external pressure} & + & \text{drag force} & + & \text{random force} \\ |
185 |
|
& = & -\hat{n}_f P A_f & - & \Xi_f(t) {\mathbf v}_f(t) & + & {\mathbf R}_f(t) |
186 |
|
\end{array} |
187 |
|
\end{equation} |
188 |
< |
|
188 |
> |
Here, $P$ is the external pressure, $A_f$ and $\hat{n}_f$ are the area |
189 |
> |
and normal vectors for facet $f$, respectively. ${\mathbf v}_f(t)$ is |
190 |
> |
the velocity of the facet, |
191 |
> |
\begin{equation} |
192 |
> |
{\mathbf v}_f(t) = \frac{1}{3} \sum_{i=1}^{3} {\mathbf v}_i, |
193 |
> |
\end{equation} |
194 |
> |
and $\Xi_f(t)$ is a ($3 \times 3$) hydrodynamic tensor that depends on |
195 |
> |
the geometry and surface area of facet $f$ and the viscosity of the |
196 |
> |
fluid (See Appendix A). The hydrodynamic tensor is related to the |
197 |
> |
fluctuations of the random force, $\mathbf{R}(t)$, by the |
198 |
> |
fluctuation-dissipation theorem, |
199 |
|
\begin{eqnarray} |
150 |
– |
A_f & = & \text{area of facet}\ f \\ |
151 |
– |
\hat{n}_f & = & \text{facet normal} \\ |
152 |
– |
P & = & \text{external pressure} |
153 |
– |
\end{eqnarray} |
154 |
– |
|
155 |
– |
\begin{eqnarray} |
156 |
– |
{\mathbf v}_f(t) & = & \text{velocity of facet} \\ |
157 |
– |
& = & \frac{1}{3} \sum_{i=1}^{3} {\mathbf v}_i \\ |
158 |
– |
\Xi_f(t) & = & \text{is a hydrodynamic tensor that depends} \\ |
159 |
– |
& & \text{on the geometry and surface area of} \\ |
160 |
– |
& & \text{facet} \ f\ \text{and the viscosity of the fluid.} |
161 |
– |
\end{eqnarray} |
162 |
– |
|
163 |
– |
\begin{eqnarray} |
200 |
|
\left< {\mathbf R}_f(t) \right> & = & 0 \\ |
201 |
|
\left<{\mathbf R}_f(t) {\mathbf R}_f^T(t^\prime)\right> & = & 2 k_B T\ |
202 |
< |
\Xi_f(t)\delta(t-t^\prime) |
202 |
> |
\Xi_f(t)\delta(t-t^\prime). |
203 |
> |
\label{eq:randomForce} |
204 |
|
\end{eqnarray} |
205 |
|
|
206 |
< |
Implemented in OpenMD.\cite{Meineke2005,openmd} |
206 |
> |
Once the hydrodynamic tensor is known for a given facet (see Appendix |
207 |
> |
A) obtaining a stochastic vector that has the properties in |
208 |
> |
Eq. (\ref{eq:randomForce}) can be done efficiently by carrying out a |
209 |
> |
one-time Cholesky decomposition to obtain the square root matrix of |
210 |
> |
the resistance tensor, |
211 |
> |
\begin{equation} |
212 |
> |
\Xi_f = {\bf S} {\bf S}^{T}, |
213 |
> |
\label{eq:Cholesky} |
214 |
> |
\end{equation} |
215 |
> |
where ${\bf S}$ is a lower triangular matrix.\cite{Schlick2002} A |
216 |
> |
vector with the statistics required for the random force can then be |
217 |
> |
obtained by multiplying ${\bf S}$ onto a random 3-vector ${\bf Z}$ which |
218 |
> |
has elements chosen from a Gaussian distribution, such that: |
219 |
> |
\begin{equation} |
220 |
> |
\langle {\bf Z}_i \rangle = 0, \hspace{1in} \langle {\bf Z}_i \cdot |
221 |
> |
{\bf Z}_j \rangle = \frac{2 k_B T}{\delta t} \delta_{ij}, |
222 |
> |
\end{equation} |
223 |
> |
where $\delta t$ is the timestep in use during the simulation. The |
224 |
> |
random force, ${\bf R}_{f} = {\bf S} {\bf Z}$, can be shown to |
225 |
> |
have the correct properties required by Eq. (\ref{eq:randomForce}). |
226 |
|
|
227 |
+ |
We have implemented this method by extending the Langevin dynamics |
228 |
+ |
integrator in our group code, OpenMD.\cite{Meineke2005,openmd} |
229 |
+ |
|
230 |
|
\section{Tests \& Applications} |
231 |
|
|
232 |
|
\subsection{Bulk modulus of gold nanoparticles} |