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.51k stars 864 forks source link

Incompatability with numpy <2 on Windows for `Lattice.find_points_in_spheres` #3990

Closed Andrew-S-Rosen closed 2 months ago

Andrew-S-Rosen commented 2 months ago

Python version

3.11

Pymatgen version

2024.8.9

Operating system version

Windows

Current behavior

As of pyamtgen 2024.8.8, there is an incompatibility with Lattice.get_points_in_sphere on Windows when using numpy < 2. Everything seems to run fine on numpy 2.0.1 and/or non-Windows machines.

pip install pymatgen==2024.8.9 numpy==1.26.4
from pymatgen.core import Structure
structure = Structure(
  lattice=[[0, 2.13, 2.13], [2.13, 0, 2.13], [2.13, 2.13, 0]],
  species=["Mg", "O"],
  coords=[[0, 0, 0], [0.5, 0.5, 0.5]],
  )
structure.lattice.get_points_in_sphere(structure.frac_coords,(0,0,0),1)

Traceback:

----> 1 structure.lattice.get_points_in_sphere(structure.frac_coords,(0,0,0),1)

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\core\lattice.py:1387, in Lattice.get_points_in_sph)
   1384 pbc = np.ascontiguousarray(self.pbc, dtype=int)
   1385 center_coords = np.ascontiguousarray([center], dtype=float)
-> 1387 _, indices, images, distances = find_points_in_spheres(
   1388     all_coords=cart_coords, center_coords=center_coords, r=float(r), pbc=pbc, lattice=latt_matrix, tol=1e-8
   1389 )
   1390 if len(indices) < 1:
   1391     # Return empty np.array (not list or tuple) to ensure consistent return type
   1392     # whether sphere contains points or not
   1393     return np.array([]) if zip_results else tuple(np.array([]) for _ in range(4))

File src\\pymatgen\\optimization\\neighbors.pyx:48, in pymatgen.optimization.neighbors.find_points_in_sp)

ValueError: Buffer dtype mismatch, expected 'const int64_t' but got 'long'

This influences many other functions in Pymatgen, such as the example below:

from pymatgen.core.surface import generate_all_slabs
structure = Structure(
  lattice=[[0, 2.13, 2.13], [2.13, 0, 2.13], [2.13, 2.13, 0]],
  species=["Mg", "O"],
  coords=[[0, 0, 0], [0.5, 0.5, 0.5]],
  )
generate_all_slabs(structure,1,10,10)

Traceback:

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\core\surface.py:1650, in generate_all_slabs(structure, max_index, min_slab_size, min_vacuum_size, bonds, tol, ftol, max_broken_bonds, lll_reduce, center_slab, primitive, max_normal_search, symmetrize, repair, include_reconstructions, in_unit_planes)
   1638 for miller in get_symmetrically_distinct_miller_indices(structure, max_index):
   1639     gen = SlabGenerator(
   1640         structure,
   1641         miller,
   (...)
   1648         in_unit_planes=in_unit_planes,
   1649     )
-> 1650     slabs = gen.get_slabs(
   1651         bonds=bonds,
   1652         tol=tol,
   1653         ftol=ftol,
   1654         symmetrize=symmetrize,
   1655         max_broken_bonds=max_broken_bonds,
   1656         repair=repair,
   1657     )
   1659     if len(slabs) > 0:
   1660         logger.debug(f"{miller} has {len(slabs)} slabs... ")

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\core\surface.py:1348, in SlabGenerator.get_slabs(self, bonds, ftol, tol, max_broken_bonds, symmetrize, repair, ztol)
   1345 matcher = StructureMatcher(ltol=tol, stol=tol, primitive_cell=False, scale=False)
   1347 final_slabs: list[Slab] = []
-> 1348 for group in matcher.group_structures(slabs):
   1349     # For each unique slab, symmetrize the
   1350     # surfaces by removing sites from the bottom
   1351     if symmetrize:
   1352         sym_slabs = self.nonstoichiometric_symmetrized_slab(group[0])

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\analysis\structure_matcher.py:795, in StructureMatcher.group_structures(self, s_list, anonymous)
    793 s_list = self._process_species(s_list)
    794 # Prepare reduced structures beforehand
--> 795 s_list = [self._get_reduced_structure(s, self._primitive_cell, niggli=True) for s in s_list]
    797 # Use structure hash to pre-group structures
    798 if anonymous:

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\analysis\structure_matcher.py:795, in <listcomp>(.0)
    793 s_list = self._process_species(s_list)
    794 # Prepare reduced structures beforehand
--> 795 s_list = [self._get_reduced_structure(s, self._primitive_cell, niggli=True) for s in s_list]
    797 # Use structure hash to pre-group structures
    798 if anonymous:

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\analysis\structure_matcher.py:947, in StructureMatcher._get_reduced_structure(cls, struct, primitive_cell, niggli)
    945 reduced = struct.copy()
    946 if niggli:
--> 947     reduced = reduced.get_reduced_structure(reduction_algo="niggli")
    948 if primitive_cell:
    949     reduced = reduced.get_primitive_structure()

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\core\structure.py:2252, in IStructure.get_reduced_structure(self, reduction_algo)
   2242 """Get a reduced structure.
   2243
   2244 Args:
   (...)
   2249     Structure: Niggli- or LLL-reduced structure.
   2250 """
   2251 if reduction_algo == "niggli":
-> 2252     reduced_latt = self._lattice.get_niggli_reduced_lattice()
   2253 elif reduction_algo == "LLL":
   2254     reduced_latt = self._lattice.get_lll_reduced_lattice()

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\core\lattice.py:1227, in Lattice.get_niggli_reduced_lattice(self, tol)
   1224 gamma = math.acos(Y / 2 / a / b) / math.pi * 180
   1225 lattice = type(self).from_parameters(a, b, c, alpha, beta, gamma)
-> 1227 mapped = self.find_mapping(lattice, e, skip_rotation_matrix=True)
   1228 if mapped is not None:
   1229     if np.linalg.det(mapped[0].matrix) > 0:

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\core\lattice.py:1010, in Lattice.find_mapping(self, other_lattice, ltol, atol, skip_rotation_matrix)
    974 def find_mapping(
    975     self,
    976     other_lattice: Self,
   (...)
    979     skip_rotation_matrix: bool = False,
    980 ) -> tuple[Lattice, np.ndarray | None, np.ndarray] | None:
    981     """Find a mapping between current lattice and another lattice. There
    982     are an infinite number of choices of basis vectors for two entirely
    983     equivalent lattices. This method returns a mapping that maps
   (...)
   1008         None is returned if no matches are found.
   1009     """
-> 1010     return next(self.find_all_mappings(other_lattice, ltol, atol, skip_rotation_matrix), None)

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\core\lattice.py:947, in Lattice.find_all_mappings(self, other_lattice, ltol, atol, skip_rotation_matrix)
    944 lengths = other_lattice.lengths
    945 alpha, beta, gamma = other_lattice.angles
--> 947 frac, dist, _, _ = self.get_points_in_sphere(
    948     [[0, 0, 0]], [0, 0, 0], max(lengths) * (1 + ltol), zip_results=False
    949 )
    950 cart = self.get_cartesian_coords(frac)
    951 # This can't be broadcast because they're different lengths

File ~\miniconda\envs\test\Lib\site-packages\pymatgen\core\lattice.py:1387, in Lattice.get_points_in_sphere(self, frac_points, center, r, zip_results)
   1384 pbc = np.ascontiguousarray(self.pbc, dtype=int)
   1385 center_coords = np.ascontiguousarray([center], dtype=float)
-> 1387 _, indices, images, distances = find_points_in_spheres(
   1388     all_coords=cart_coords, center_coords=center_coords, r=float(r), pbc=pbc, lattice=latt_matrix, tol=1e-8
   1389 )
   1390 if len(indices) < 1:
   1391     # Return empty np.array (not list or tuple) to ensure consistent return type
   1392     # whether sphere contains points or not
   1393     return np.array([]) if zip_results else tuple(np.array([]) for _ in range(4))

File src\\pymatgen\\optimization\\neighbors.pyx:48, in pymatgen.optimization.neighbors.find_points_in_spheres()

ValueError: Buffer dtype mismatch, expected 'const int64_t' but got 'long
DanielYang59 commented 2 months ago

Thanks a lot for reporting. The unit test for chgnet is failing for the same reason, following our recent migration of build system to NumPy 2 #3894, as in NumPy 2 the default int type has been changed from int32 to int64 in Windows 64-bit system.

I would look into this ASAP.

DanielYang59 commented 2 months ago

@Andrew-S-Rosen I believe you just saved the world :) It turns out there're way more uncaught errors https://github.com/materialsproject/pymatgen/actions/runs/10330245550/job/28599007337