In condensed matter simulations, how is particle number density computed in practice? 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?
 A: 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$.
A: 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.
