Remove own system in favor of metatensor.torch.atomistic.System

docs/src/references/lib/index.rst

@@ -9,4 +9,3 @@ are used for the meshLODE calculators.
 
    fourier_convolution
    mesh_interpolator
-   system

examples/library-tutorial.py

@@ -1,7 +1,6 @@
 """
 Basic Tutorial for Library functions
 ====================================
-
 This examples provides an illustration of the functioning of the underlaying library
 functions of ``meshlode`` and the construction LODE descriptors (`Grisafi 2019
 <https://doi.org/10.1063/1.5128375>`__, `Grisafi 2021
@@ -18,6 +17,7 @@
 import matplotlib.pyplot as plt
 import numpy as np
 import torch
+from metatensor.torch.atomistic import System
 
 import meshlode
 
@@ -48,22 +48,19 @@ def sliceplot(mesh, sz=12, cmap="viridis", vmin=None, vmax=None):
 # %%
 # Builds the structure
 # --------------------
-#
 # Builds a CsCl structure by replicating the primitive cell using ase and convert it to
-# a :py:class:`list` of :py:class:`meshlode.System`. We add a bit of noise to make
-# it less boring!
+# a :py:class:`list` of :py:class:`metatensor.torch.atomistic.System`. We add a bit of
+# noise to make it less boring!
 #
 
 positions = torch.tensor([[0, 0, 0], [0.5, 0.5, 0.5]]) * 4
-atomic_types = torch.tensor([55, 17])  # Cs and Cl
+types = torch.tensor([55, 17])  # Cs and Cl
 cell = torch.eye(3) * 4
-ase_frame = ase.Atoms(positions=positions, cell=cell, numbers=atomic_types).repeat(
-    [2, 2, 2]
-)
+ase_frame = ase.Atoms(positions=positions, cell=cell, numbers=types).repeat([2, 2, 2])
 ase_frame.positions[:] += np.random.normal(size=ase_frame.positions.shape) * 0.1
 charges = torch.tensor([1.0, -1.0] * 8)
-frame = meshlode.System(
-    species=torch.tensor(ase_frame.numbers),
+system = System(
+    types=torch.tensor(ase_frame.numbers),
     positions=torch.tensor(np.array(ase_frame.positions)),
     cell=torch.tensor(ase_frame.cell),
 )
@@ -85,7 +82,6 @@ def sliceplot(mesh, sz=12, cmap="viridis", vmin=None, vmax=None):
 # %%
 # MeshInterpolator
 # ----------------
-#
 # ``MeshInterpolator`` serves as a utility class to compute a mesh
 # representation of points, and/or to project a function defined on the
 # mesh on a set of points. Computing the mesh representation is a two-step
@@ -95,15 +91,15 @@ def sliceplot(mesh, sz=12, cmap="viridis", vmin=None, vmax=None):
 #
 
 interpol = meshlode.lib.mesh_interpolator.MeshInterpolator(
-    frame.cell, torch.tensor([16, 16, 16]), interpolation_order=3
+    system.cell, torch.tensor([16, 16, 16]), interpolation_order=3
 )
 
-interpol.compute_interpolation_weights(frame.positions)
+interpol.compute_interpolation_weights(system.positions)
 
 
 # %%
 # We use two sets of weights: ones (giving the atom density irrespective
-# of the species) and charges (giving a smooth representation of the point
+# of the types) and charges (giving a smooth representation of the point
 # charges).
 #
 
@@ -126,7 +122,6 @@ def sliceplot(mesh, sz=12, cmap="viridis", vmin=None, vmax=None):
 # %%
 # Fourier filter
 # --------------
-#
 # This module computes a Fourier-domain filter, that can be used e.g. to
 # smear the density and/or compute a 1/r^p potential field. This can also
 # be easily extended to compute an arbitrary filter
@@ -137,15 +132,15 @@ def sliceplot(mesh, sz=12, cmap="viridis", vmin=None, vmax=None):
 # %%
 # plain atomic_smearing
 rho_mesh = fsc.compute(
-    mesh_values=mesh, cell=frame.cell, potential_exponent=0, atomic_smearing=1
+    mesh_values=mesh, cell=system.cell, potential_exponent=0, atomic_smearing=1
 )
 
 sliceplot(rho_mesh[0, :, :, :5])
 
 # %%
 # coulomb-like potential, no atomic_smearing
 coulomb_mesh = fsc.compute(
-    mesh_values=mesh, cell=frame.cell, potential_exponent=1, atomic_smearing=0
+    mesh_values=mesh, cell=system.cell, potential_exponent=1, atomic_smearing=0
 )
 
 sliceplot(coulomb_mesh[1, :, :, :5], cmap="seismic")
@@ -154,7 +149,6 @@ def sliceplot(mesh, sz=12, cmap="viridis", vmin=None, vmax=None):
 # %%
 # Back-interpolation (on the same points)
 # ---------------------------------------
-#
 # The same ``MeshInterpolator`` object can be used to compute a field on
 # the same points used initially to generate the atom density
 #
@@ -167,21 +161,20 @@ def sliceplot(mesh, sz=12, cmap="viridis", vmin=None, vmax=None):
 # %%
 # Back-interpolation (on different points)
 # ----------------------------------------
-#
 # In order to compute the field on a different set of points, it is
 # sufficient to build another ``MeshInterpolator`` object and to compute
 # it with the desired field. One can also use a different
 # ``interpolation_order``, if wanted.
 #
 
 interpol_slice = meshlode.lib.mesh_interpolator.MeshInterpolator(
-    frame.cell, torch.tensor([16, 16, 16]), interpolation_order=4
+    system.cell, torch.tensor([16, 16, 16]), interpolation_order=4
 )
 
 # Compute a denser grid on a 2D slice
 n_points = 50
-x = torch.linspace(0, frame.cell[0, 0], n_points + 1)[:n_points]
-y = torch.linspace(0, frame.cell[1, 1], n_points + 1)[:n_points]
+x = torch.linspace(0, system.cell[0, 0], n_points + 1)[:n_points]
+y = torch.linspace(0, system.cell[1, 1], n_points + 1)[:n_points]
 xx, yy = torch.meshgrid(x, y, indexing="ij")
 
 # Flatten xx and yy, and concatenate with a zero column for the z-coordinate

examples/madelung.py

@@ -1,7 +1,6 @@
 """
 Compute Madelung Constants
 ==========================
-
 In this tutorial we show how to calculate the Madelung constants and total electrostatic
 energy of atomic structures using the :py:class:`meshlode.MeshPotential` and
 :py:class:`meshlode.metatensor.MeshPotential` calculator.
@@ -11,24 +10,24 @@
 import math
 
 import torch
+from metatensor.torch.atomistic import System
 
 import meshlode
 
 
 # %%
 # Define simple example structure having the CsCl structure and compute the reference
-# values. MeshPotential by default outputs the species sorted according to the atomic
+# values. MeshPotential by default outputs the types sorted according to the atomic
 # number. Thus, we input the compound "CsCl" and "ClCs" since Cl and Cs have atomic
 # numbers 17 and 55, respectively.
-atomic_types = torch.tensor([17, 55])  # Cl and Cs
+types = torch.tensor([17, 55])  # Cl and Cs
+positions = torch.tensor([[0, 0, 0], [0.5, 0.5, 0.5]])
 charges = torch.tensor([-1.0, 1.0])
 cell = torch.eye(3)
-positions = torch.tensor([[0, 0, 0], [0.5, 0.5, 0.5]])
-frame = meshlode.System(species=atomic_types, positions=positions, cell=torch.eye(3))
 
 # %%
 # Define the expected values of the energy
-n_atoms = len(positions)
+n_atoms = len(types)
 madelung = 2 * 1.7626 / math.sqrt(3)
 energies_ref = -madelung * torch.ones((n_atoms, 1))
 
@@ -43,7 +42,6 @@
 # %%
 # Computation using ``meshlode``
 # ------------------------------
-#
 # Compute features using
 
 MP = meshlode.MeshPotential(
@@ -52,7 +50,7 @@
     interpolation_order=interpolation_order,
     subtract_self=True,
 )
-potentials_torch = MP.compute(frame)
+potentials_torch = MP.compute(types=types, positions=positions, cell=cell)
 
 # %%
 # The "potentials" that have been computed so far are not the actual electrostatic
@@ -66,7 +64,7 @@
     for idx_n in range(n_atoms):
         # The coulomb potential between atoms i and j is charge_i * charge_j / d_ij
         # The features are simply computing a pure 1/r potential with no prefactors.
-        # Thus, to compute the energy between atoms of species i and j, we need to
+        # Thus, to compute the energy between atoms of types i and j, we need to
         # multiply by the charges of i and j.
         print(charges[idx_c] * charges[idx_n], potentials_torch[idx_n, idx_c])
         atomic_energies_torch[idx_c] += (
@@ -88,33 +86,36 @@
 # %%
 # Computation using ``meshlode.metatensor``
 # -----------------------------------------
-#
-# We now compute the same constants using the metatensor based calculator
+# We now compute the same constants using the metatensor based calculator. To achieve
+# this we first store our system parameters like the ``types``, ``positions`` and the
+# ``cell`` defined above into a :py:class:`metatensor.torch.atomistic.System` class.
+
+system = System(types=types, positions=positions, cell=cell)
 
 MP = meshlode.metatensor.MeshPotential(
     atomic_smearing=atomic_smearing,
     mesh_spacing=mesh_spacing,
     interpolation_order=interpolation_order,
     subtract_self=True,
 )
-potential_metatensor = MP.compute(frame)
+potential_metatensor = MP.compute(system)
 
 
 # %%
 # To get the Madelung constant, we again need to take a linear combination
 # of the "potentials" weighted by the charges of the atoms.
 
 atomic_energies_metatensor = torch.zeros((n_atoms, 1))
-for idx_c, c in enumerate(atomic_types):
-    for idx_n, n in enumerate(atomic_types):
-        # Take the coefficients with the correct center atom and neighbor atom species
+for idx_c, c in enumerate(types):
+    for idx_n, n in enumerate(types):
+        # Take the coefficients with the correct center atom and neighbor atom types
         block = potential_metatensor.block(
             {"center_type": int(c), "neighbor_type": int(n)}
         )
 
         # The coulomb potential between atoms i and j is charge_i * charge_j / d_ij
         # The features are simply computing a pure 1/r potential with no prefactors.
-        # Thus, to compute the energy between atoms of species i and j, we need to
+        # Thus, to compute the energy between atoms of types i and j, we need to
         # multiply by the charges of i and j.
         print(c, n, charges[idx_c] * charges[idx_n], block.values[0, 0])
         atomic_energies_metatensor[idx_c] += (

pyproject.toml

@@ -34,7 +34,7 @@ keywords = [
     "Atomistic Simulations",
 ]
 dependencies = [
-    "torch >= 1.11",
+    "torch >=1.11",
 ]
 dynamic = ["version"]
 
@@ -44,7 +44,7 @@ examples = [
     "matplotlib",
 ]
 metatensor = [
-    "metatensor[torch]",
+    "metatensor-torch >=0.3",
 ]
 
 [project.urls]

src/meshlode/__init__.py

@@ -1,11 +1,10 @@
 from .calculators.meshpotential import MeshPotential
-from .lib.system import System
 
 try:
     from . import metatensor  # noqa
 except ImportError:
     pass
 
 
-__all__ = ["MeshPotential", "System"]
+__all__ = ["MeshPotential"]
 __version__ = "0.0.0-dev"