The total energy of a system can be decomposed into kinetic and potential parts:
\[E = K + U.\] Using the thermodynamic relation, \[p = - \left(\dfrac{dE}{dV}\right)_T.\] This translates to, \[P = \dfrac{NkT}{V} - \left\langle\dfrac{dU}{dV}\right\rangle_T,\] where the angle brackets denote a statistical average. I will drop the subscript \(T\) for brevity from here on.
Molecular simulations are typically carried out in boxes. Let us assume a cubic box of volume \(V = L^3\), with \(N\) particles. Let the coordinates of the particles be \(\mathbf{r}^N = \{\mathbf{r}_1, \mathbf{r}_2, ..., \mathbf{r}_N,\}\).
People often like to use scaled coordinates:\[\mathbf{r}_i = \mathbf{s}_i L.\] The potential energy describes interactions of the particles, and can be thought of \(U(\mathbf{r}^N, L)\). The explicit dependence on \(L\) will become clear shortly.
Now one can write:
\[\dfrac{dU}{dV} = \dfrac{dL}{dV}\left[\dfrac{\partial U}{\partial L} + \sum_{i=1}^N \dfrac{\partial U}{\partial \mathbf{r}_i} \cdot \dfrac{d\mathbf{r}_i}{dL} \right].\] Using \(V = L^3\) and \(d\mathbf{r}_i/dL = \mathbf{s}_i\), and identifying the force on particle \(i\) as \(f_i = -\partial U/\partial \mathbf{r}_i\), one can simplify this relation to obtain:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{3V} \sum_{i=1}^N \mathbf{r}_i \cdot \mathbf{f}_i - \dfrac{1}{3L^2} \dfrac{\partial U}{\partial L} \right\rangle.\]
Finite Non-periodic Systems
For finite, non-periodic systems, the potential energy depends only on atomic coordinates, \(U(\mathbf{r}^N)\). Hence, the partial with respect to \(L\) is equal to zero, and we obtain the familiar form:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{3V} \sum_{i=1}^N \mathbf{r}_i \cdot \mathbf{f}_i \right\rangle.\] It is easy to show that in \(d\) dimensions, \(V = L^d\), this equation can be generalized to: \[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{dV} \sum_{i=1}^N \mathbf{r}_i \cdot \mathbf{f}_i \right\rangle.\] For pairwise forces, \(\mathbf{f}_i = \sum_{j \neq i} \mathbf{f}_{ij}\), where \(\mathbf{f}_{ij} = - \mathbf{f}_{ji}\), this equation can be further simplified to:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{6V} \sum_{i=1}^N \sum_{j \neq i}^N \mathbf{r}_{ij} \cdot \mathbf{f}_{ij} \right\rangle.\]
Periodic Systems
Very often, periodic boundary conditions are used in molecular simulations.
Thus, a given particle interacts not only with particles in the given box, but all the periodic images. This includes, by the way, images of itself in other periodic boxes.
One now has to contend with the entire expression, including the partial with respect to box length:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{3V} \sum_{i=1}^N \mathbf{r}_i \cdot \mathbf{f}_i - \dfrac{1}{3L^2} \dfrac{\partial U}{\partial L} \right\rangle.\] Note that the dependence of \(U\) on \(L\) is real. If we expand the box, we have to choose whether we scale all the coordinates, or just expand the box, leaving the particle coordinates in the simulation box unchanged.
In either case, the interaction energy between particle \(i\) and periodic images is mediated by \(L\). Thus, ignoring this dependence makes virial pressure calculations for general (especially multibody) potentials incorrect!
However, for pairwise potentials (which constitute the bulk of simulations in the literature), it turns out that one can write the total potential energy as, \[U = \dfrac{1}{2} \sum_{i=1}^N \sum_{j \neq i}^{MN} u(r_{ij}),\] where \(M\) sums over the (potentially infinite) periodic boxes. Note that since we are summing over all the periodic boxes, we don't need an explicit dependence on \(L\) - the only dependence on \(L\) is implicit via scaling of \(\mathbf{r}_i\).
For affine deformation, in which particle coordinates are scaled, it turns out that we indeed recover the expression:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{6V} \sum_{i=1}^N \sum_{j \neq i}^N \mathbf{r}_{ij} \cdot \mathbf{f}_{ij} \right\rangle.\]
This is extraordinarily fortunate!
References
1. Manuel J. Louwerse, Evert Jan Baerends, ``Calculation of pressure in case of periodic boundary conditions", Chemical Physics Letters 421, 138-141, 2006.
This paper, while not the first, is the most recent reminder that PBCs don't play well with multibody potentials.
2. Aidan P. Thompson, Steven J. Plimpton, and William Mattson, ``General formulation of pressure and stress tensor for arbitrary many-body interaction potentials under periodic boundary conditions", Journal of Chemical Physics, 154107, 2009.
This paper offers a surprisingly straightforward solution for addressing PBCs and multibody potentials.
\[E = K + U.\] Using the thermodynamic relation, \[p = - \left(\dfrac{dE}{dV}\right)_T.\] This translates to, \[P = \dfrac{NkT}{V} - \left\langle\dfrac{dU}{dV}\right\rangle_T,\] where the angle brackets denote a statistical average. I will drop the subscript \(T\) for brevity from here on.
Molecular simulations are typically carried out in boxes. Let us assume a cubic box of volume \(V = L^3\), with \(N\) particles. Let the coordinates of the particles be \(\mathbf{r}^N = \{\mathbf{r}_1, \mathbf{r}_2, ..., \mathbf{r}_N,\}\).
Now one can write:
\[\dfrac{dU}{dV} = \dfrac{dL}{dV}\left[\dfrac{\partial U}{\partial L} + \sum_{i=1}^N \dfrac{\partial U}{\partial \mathbf{r}_i} \cdot \dfrac{d\mathbf{r}_i}{dL} \right].\] Using \(V = L^3\) and \(d\mathbf{r}_i/dL = \mathbf{s}_i\), and identifying the force on particle \(i\) as \(f_i = -\partial U/\partial \mathbf{r}_i\), one can simplify this relation to obtain:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{3V} \sum_{i=1}^N \mathbf{r}_i \cdot \mathbf{f}_i - \dfrac{1}{3L^2} \dfrac{\partial U}{\partial L} \right\rangle.\]
Finite Non-periodic Systems
For finite, non-periodic systems, the potential energy depends only on atomic coordinates, \(U(\mathbf{r}^N)\). Hence, the partial with respect to \(L\) is equal to zero, and we obtain the familiar form:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{3V} \sum_{i=1}^N \mathbf{r}_i \cdot \mathbf{f}_i \right\rangle.\] It is easy to show that in \(d\) dimensions, \(V = L^d\), this equation can be generalized to: \[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{dV} \sum_{i=1}^N \mathbf{r}_i \cdot \mathbf{f}_i \right\rangle.\] For pairwise forces, \(\mathbf{f}_i = \sum_{j \neq i} \mathbf{f}_{ij}\), where \(\mathbf{f}_{ij} = - \mathbf{f}_{ji}\), this equation can be further simplified to:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{6V} \sum_{i=1}^N \sum_{j \neq i}^N \mathbf{r}_{ij} \cdot \mathbf{f}_{ij} \right\rangle.\]
Periodic Systems
Very often, periodic boundary conditions are used in molecular simulations.
Thus, a given particle interacts not only with particles in the given box, but all the periodic images. This includes, by the way, images of itself in other periodic boxes.
One now has to contend with the entire expression, including the partial with respect to box length:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{3V} \sum_{i=1}^N \mathbf{r}_i \cdot \mathbf{f}_i - \dfrac{1}{3L^2} \dfrac{\partial U}{\partial L} \right\rangle.\] Note that the dependence of \(U\) on \(L\) is real. If we expand the box, we have to choose whether we scale all the coordinates, or just expand the box, leaving the particle coordinates in the simulation box unchanged.
In either case, the interaction energy between particle \(i\) and periodic images is mediated by \(L\). Thus, ignoring this dependence makes virial pressure calculations for general (especially multibody) potentials incorrect!
However, for pairwise potentials (which constitute the bulk of simulations in the literature), it turns out that one can write the total potential energy as, \[U = \dfrac{1}{2} \sum_{i=1}^N \sum_{j \neq i}^{MN} u(r_{ij}),\] where \(M\) sums over the (potentially infinite) periodic boxes. Note that since we are summing over all the periodic boxes, we don't need an explicit dependence on \(L\) - the only dependence on \(L\) is implicit via scaling of \(\mathbf{r}_i\).
For affine deformation, in which particle coordinates are scaled, it turns out that we indeed recover the expression:
\[P = \dfrac{NkT}{V} + \left\langle \dfrac{1}{6V} \sum_{i=1}^N \sum_{j \neq i}^N \mathbf{r}_{ij} \cdot \mathbf{f}_{ij} \right\rangle.\]
This is extraordinarily fortunate!
References
1. Manuel J. Louwerse, Evert Jan Baerends, ``Calculation of pressure in case of periodic boundary conditions", Chemical Physics Letters 421, 138-141, 2006.
This paper, while not the first, is the most recent reminder that PBCs don't play well with multibody potentials.
2. Aidan P. Thompson, Steven J. Plimpton, and William Mattson, ``General formulation of pressure and stress tensor for arbitrary many-body interaction potentials under periodic boundary conditions", Journal of Chemical Physics, 154107, 2009.
This paper offers a surprisingly straightforward solution for addressing PBCs and multibody potentials.
No comments:
Post a Comment