I have a 3D charge density file, calculated from Quantum Espresso DFT.

If someone calculates the total charge (number of electrons) of an unit cell from that cube file (using a numerical approximation method, for example, a spherical integration), what should be the basic mathematical equation to express that?

If needed to see the charge density file: 3D cube file

For example, the DOS integration gives us the electron density at T=0 K using the following equation: enter image description here.

Now need to know the equation for total charge (using spherical integration) for the 3D cube file of charge density. Thanks in advance.


1 Answer 1


The cube file contains the density for each grid point.

The general formula for calculating the total number of electrons for that unit-cell would be: $$ n = \iiint \rho(x,y,z)dxdydz $$

Each grid point $\rho(x,y,z)$ contains the energy integrated density $$\rho(x,y,z)=\int_{-\infty}^{E_F} DOS(E,x,y,z) dE$$

In some cases, the cube files may contain LDOS quantities, that is charge in some energy range so that the 2nd equation has different bounds, but that may be code specific (I am not sure what QE can do here).

Numerically, say in Python, one would just do: np.sum(grid) * volume / grid.size.

  • $\begingroup$ Thanks @nickpaior. What can be the limit of integration. For example, the cube file given here has two Si atoms in a cubic unit cell. $\endgroup$
    – Sak
    Dec 9, 2022 at 4:06
  • 1
    $\begingroup$ What do you mean about limits? Could you see the edit, whether it makes more sense? $\endgroup$
    – nickpapior
    Dec 9, 2022 at 11:47
  • $\begingroup$ Thanks @nickpapior. It makes sense. I also see that the python code for the integration that you wrote. As one of the goals of the task is also to use a numerical approximation method (for example, any numerical approximation integration method), could you suggest any method and it's corresponding idea to implement in the code? Thanks in advance. $\endgroup$
    – Sak
    Dec 9, 2022 at 20:04
  • 1
    $\begingroup$ @Sak you should just add the values up. Although this seems poor at first glance, it is actually correct to extremely high-order in periodic boundary conditions. $\endgroup$ Dec 9, 2022 at 22:49
  • 1
    $\begingroup$ @Sak you would have to use the same broadening function and smearing width as the DFT code to get the exact answer $\endgroup$ Dec 11, 2022 at 2:16

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .