I have been reading a recent paper. In it, the authors performed molecular dynamics (MD) simulations of parallel-plate supercapacitors, in which liquid resides between the parallel-plate electrodes. To simplify the situation, let us suppose that the liquid between the electrodes is argon liquid.

The system has a "slab" geometry, so the authors are only interested in variations of the liquid structure along the $z$ direction. Thus, the authors compute the particle number densities averaged over $x$ and $y$: $\bar{n}_\alpha(z)$, where $\alpha$ is a solvent species. (That is, in my simplified example, $\alpha$ is argon -- an argon atom.) $\bar{n} _\alpha(z)$ has dimensions of $\frac{\text{number}}{\text{length}^3}$ or simply $\text{length}^{-3}$, I think.

The $xy$-plane is given by the inequalities $-x_0 < x < x_0$ and $-y_0 < y < y_0$. The area $A_0$ of the $xy$-plane is thus given by $A_0 = 4x_0y_0$.

So, the authors define the particle number density averaged over $x$ and $y$ as follows: $$\bar{n}_\alpha(z) = A_0^{-1} \int_{-x_0}^{x_0} \int_{-y_0}^{y_0} dx^\prime dy^\prime n_\alpha(x^\prime, y^\prime, z)$$ where $A_0 = 4x_0y_0$ and $n_\alpha(x, y, z)$ is the local number density of $\alpha$ at $(x, y, z)$.

Thus, $\bar{n}_\alpha(z)$ is simply proportional to $n_\alpha$ integrated over $x$ and $y$. But, my question is, what is $n_\alpha(x, y, z)$? How is $n_\alpha(x, y, z)$ determined in practice?

As far as the computer is concerned, the argon atoms are point particles; they are modeled as having zero volume (although they interact by Lennard-Jones interactions). So how is it possible to define a number density?

Do we simply "cut" the "slab" in "slices" along $z$ and then assign the particles to these slices? There might be 5 particles in the first $z$ slice, 10 in the second, 7 in the third, and so on. If I then divide 5, 10, and 7 by the volume of the respective slice, then I have a sort of number density, with units of $\frac{\text{number}}{\text{length}^3}$ or simply $\text{length}^{-3}$. But how do I now integrate this $n_\alpha(x^\prime, y^\prime, z)$ over $x$ and $y$? Do I have to additionally perform binning in the $x$ and $y$ directions?


2 Answers 2


Without seeing the paper, it's hard to know for sure, but the actual particle density probably takes the form

$$n_\alpha(x,y,z) = \sum_{i\in\text{ particles}} \delta^{(3)}(x_i, y_i, z_i)$$

When you integrate this over $x$ and $y$ and some small range $\Delta z$, you get the number of particles in the region you integrated over. So a computer wouldn't actually have to do an integral, it would just count the number of particles in the region. In other words, the simulation probably works with $\bar{n}_\alpha$ directly, not $n_\alpha$.

  • $\begingroup$ Thanks. The paper is available here: http://pubs.rsc.org/en/Content/ArticleLanding/2012/FD/C1FD00086A. The particle number density is discussed in equation 1 on page 252. $\endgroup$
    – Andrew
    Jul 14, 2012 at 20:50
  • 1
    $\begingroup$ OK, thanks, I edited the link into your question. I'm not able to access the full text of the paper right now, but I'll come back and revisit this when I can. $\endgroup$
    – David Z
    Jul 14, 2012 at 21:27
  • $\begingroup$ This is right without seeing the paper--- it's not good to put caveats in, since it seeds doubt. $\endgroup$
    – Ron Maimon
    Jul 14, 2012 at 21:29

To calculate $\bar{n}_\alpha$ it's pretty much just what you said. You take the slice between, say $z=2.3$nm and $z=2.301$nm, and count the average number of atoms in it. Divide that number by the volume of the slice (cross-sectional area of the simulation box, multiplied by the slice thickness, i.e. 0.001nm). The answer you get is the number density at $z=2.3$nm

In practice: Each simulation snapshot, you write down the z-coordinate of each atom. As the simulation goes on, you get a larger and larger list of real numbers---all those z-coordinates. Now, plot those numbers in the form of a histogram. If you have a long enough simulation, you can make the bin size of your histogram very very small, so the histogram will look like a smooth curve. (Make sure you scale the histogram so that the integral under the curve is the total number of particles in the simulation divided by the cross-sectional area.)

You never have to explicitly bin or integrate over x and y, if all you need is $\bar{n}_\alpha$.

An alternate approach to calculating $\bar{n}_\alpha$---although it makes no sense to do it this way---is to calculate $n_\alpha$ first, then $\bar{n}_{\alpha}$ second. For the first step, you need to bin in the x,y,z directions---draw little cubes, count the average number of atoms in them, divide by volume. For the second step, you use the formula you cited to integrate $n_\alpha$ over x and y, then divide by cross-sectional area (or in simpler terms, take the mean value of $n_\alpha(x,y,z)$ as $x$ and $y$ vary but $z$ is fixed).

I think you may have gotten confused because the authors discuss the concept of averaging over $x$ and $y$, but you can and should calculate $\bar{n}_\alpha$ without actually explicitly doing that as a separate step.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.