materialsproject / pymatgen

Python Materials Genomics (pymatgen) is a robust materials analysis code that defines classes for structures and molecules with support for many electronic structure codes. It powers the Materials Project.
https://pymatgen.org
Other
1.48k stars 855 forks source link

Periodic boundary condition is not considered in the interpolator of VolumetricData #3787

Open goodwilling opened 4 months ago

goodwilling commented 4 months ago

Python version

Python 3.11.4

Pymatgen version

pymatgen 2024.3.1

Operating system version

Windows 10

Current behavior

When I use the interpolation method implemented in pymatgen.io.common.VolumetricData, it is found that the results for valence charge density are slightly different compared to those by other methods like c2x. Moreover, the results become considerably different for all-electron charge density and the expected symmetry is broken. In the source codes, it seems that periodic boundary condition (PBC) is not considered in the interpolator of VolumetricData as in the following, since in PBC the points of data are not from 0.0 to 1.0 in fractional coordinates.

https://github.com/materialsproject/pymatgen/blob/master/pymatgen/io/common.py#L77

  self.xpoints = np.linspace(0.0, 1.0, num=self.dim[0])
    self.ypoints = np.linspace(0.0, 1.0, num=self.dim[1])
    self.zpoints = np.linspace(0.0, 1.0, num=self.dim[2])
    self.interpolator = RegularGridInterpolator(
        (self.xpoints, self.ypoints, self.zpoints),
        self.data["total"],
        bounds_error=True,
    )

Expected Behavior

The interpolation results of electronic density charge should be the same as expected, in particular that the symmetry should not be broken.

Minimal example

No response

Relevant files to reproduce this bug

No response

goodwilling commented 4 months ago

Using the following revision, expected interpolation results have been obtained.

    self.xpoints = np.linspace(0.0, (self.dim[0]-1)/self.dim[0], num=self.dim[0])
    self.ypoints = np.linspace(0.0, (self.dim[1]-1)/self.dim[1], num=self.dim[1])
    self.zpoints = np.linspace(0.0, (self.dim[2]-1)/self.dim[2], num=self.dim[2])