Sunday, November 27, 2016

Interesting Links

1. Should I Get a PhD (H/T John D. Cook)

2. Google and Lip Reading (the Verge)

3. Magnus Carlsen and Chess (Vice)

4. The Theranos Whistleblower (WSJ)

Thursday, November 17, 2016

A Visual Interpretation of IVPs

Let us  try to explore the structure of an IVP graphically by considering a specific example:\[\dfrac{dy}{dt} = y - y^2, \quad \quad \quad y(0) = 0.1.\] Suppose, we are interested in the solution \(y(t)\) over the domain \(t \in [0, 10]\).

A jupyter notebook, which accompanies this blog, is available here on github.

click to enlarge
Let us consider the 2D domain (y versus t) on which the solution to the IVP is shown as a thick blue line. This is the solution y(t), which satisfies the initial condition, and the differential equation.

We can look at any point (t,y) on this domain, and ask "what is f(y,t) here?"

Note \(dy/dt = f(y,t)\) is the "slope"; it can be positive, negative or zero. The only restriction  is that it cannot point "backwards" in the direction of negative t.

To visualize the slope at each grid point, we can set the horizontal component \(u = 1\), and the vertical component \(v = f(y,t)\), and normalize by dividing each component by \(\sqrt{u^2 + v^2}\).

Therefore, the function f(y,t) completely defines the field of arrows. The streamlines, if you will.

The initial condition \((t_0, y_0)\) tells us where to start in this field; where to "drop the feather" in the river, to be guided and carried away by the streamlines.

The jupyter notebook linked to above, lets you play around with different problems, and different initial conditions. Take a spin.

Friday, November 11, 2016

Quotes from Twitter

Thursday, November 10, 2016

Virial Pressure and PBCs

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!


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.

Tuesday, November 8, 2016

Roger Cotes of Newton-Cotes Fame

We have been discussing Newton-Cotes integration rules in our class this week.

Everyone knows Newton. Not everyone knows the other guy: Robert Cotes.

Robert Cotes died young: he was only 33, when he died of a violent fever. Nevertheless, he made three important contributions:

  • helped Newton edit the second edition of the Principia
  • helped develop the Newton-Cotes forumulae
  • discovered a version of the famour Euler Formula: \[e^{i\theta} = \cos \theta+i \sin \theta.\]
They are beautifully summarized in this blog.

Wednesday, November 2, 2016

Kipnis and Title IX

Recently, we have a Title IX training seminar in our department. A colleague brought to our attention the fascinating case of Laura Kipnis at Northwestern University. This video summarizes the contours of the case:

Here's the rough outline. She wrote a piece in the Chronicle entitled "Sexual Paranoia Strikes Academe" (pdf), which prompted student outrage. She became the subject of a Title IX investigation, which she chronicled in another article (pdf), describing the opaque process.

Eventually, the charges against her were dropped, but the amount of time, money, and effort wasted on a stupid bureaucratic process raised important questions.
University President Morton Schapiro said many people have questioned the University’s decision to investigate the complaints at all. 
“The idea that a student shouldn’t be able to bring a Title IX complaint against a faculty member because of the faculty member’s protection under the First Amendment — that’s not my decision. That’s not Northwestern’s decision,” Schapiro told The Daily. “That’s federal law.”
The unfortunate outcome of such frivolous charges are that it helps de-legitimize/de-sensitize people from really outrageous violations. It is sad that the scope of the original legislation has been elastically expanded to become a political weapon.