Lattice, Periodic Boundary Conditions
*************************************


The "lattice" attribute
=======================

The "lattice" property of the ChemicalSystem contains the lattice
vector information through the "Lattice" class.

The ChemicalSystem can handle periodic systems with arbitrary numbers
of periodic boundary conditions (i.e. 0, 1, 2, 3).

Here are some examples:

   from scm.base import ChemicalSystem
   from scm.base import Lattice

   # Create an empty chemical system
   chain = ChemicalSystem()

   # Add a helium atom at the origin
   chain.add_atom("He", [0, 0, 0])

   # Make it a 1D periodic chain by adding a lattice vector along X (units: Bohr)
   chain.lattice = Lattice([[2, 0, 0]])
   print(chain)

   # You can also initialize it from a System Block:
   graphene = ChemicalSystem("""
   System
        Atoms
             C    0.0000000000000000  0.0000000000000000  0.0000000000000000
             C    1.2300000000000000  0.7101408311032394  0.0000000000000000
        End
        Lattice
              2.4600000000000000   0.0000000000000000   0.0000000000000000
             -1.2300000000000000   2.1304224933097191   0.0000000000000000
        End
   End
   """)

   print(graphene)


Lattice class
=============

These are the methods in the Lattice class.

class Lattice
class Lattice(vectors: ArrayLike, unit: str = 'angstrom')

   A class representing the lattice of a ChemicalSystem.

   cartesian_to_fractional(cartesian_coords: ArrayLike, unit: str = 'angstrom') -> ndarray[Any, dtype[float64]]

      Convert Cartesian coordinates to fractional coordinates for a
      set of points.

      The shape of the output array will be the same as of the input
      array "cartesian_coords". The rows of "cartesian_coords"
      represent the independent points to convert the coordinates for.
      For <3D periodic systems, the number of columns should be ">=
      num_vectors" and "<= 3". No conversion will take place in non-
      periodic directions: the coordinate values are just copied to
      the output array.

   copy() -> Lattice

      Creates a deep copy of the Lattice.

   fractional_to_cartesian(fractional_coords: ArrayLike, unit: str = 'angstrom') -> ndarray[Any, dtype[float64]]

      Convert fractional coordinates to Cartesian coordinates for a
      set of points.

      The shape of the output array will be the same as of the input
      array "fractional_coords". The rows of "fractional_coords"
      represent the independent points to convert the coordinates for.
      For <3D periodic systems, the number of columns should be ">=
      num_vectors" and "<= 3". No conversion will take place in non-
      periodic directions: the coordinate values are just copied to
      the output array.

   classmethod from_kf(kf: KFFile, section: str) -> Lattice

      Constructs and returns a new Lattice from a section on a KF
      file.

   get_3x3_lattice_vectors(unit: str = 'angstrom') -> ndarray[Any, dtype[float64]]

      Get the lattice vectors as a 3x3 matrix padded out with zeros.

   get_angles(unit: str = 'rad') -> ndarray[Any, dtype[float64]]

      Get the angles between the lattice vectors (note: always returns
      3 numbers).

   get_cell_depth(unit: str = 'angstrom') -> ndarray[Any, dtype[float64]]

      Return unit cell depth, or thickness, in each dimension.

      The cell depth is the distance between the two opposite bounding
      planes of the parallelepiped spanned by the lattice vectors. For
      three lattice vectors "v_i", "v_j", and "v_k", the "i"-th
      component of the cell depth is the height of lattice vector i
      over the plane spanned by "v_j" and "v_k".

      Always returns a 3-component vector. For < 3D periodic systems,
      the remaining components for the non-periodic directions at the
      end are set to 0.

   get_lattice_translation(d: ArrayLike, unit: str = 'angstrom') -> ndarray[Any, dtype[float64]]

      Returns the overall translation given a multiple of the lattice
      vectors.

      The size of "d" shall "num_vectors <= len(d) <= 3". If "len(d) >
      num_vectors", the extra values are ignored.

   get_lengths(unit: str = 'angstrom') -> ndarray[Any, dtype[float64]]

      Get the length of the lattice vectors (note: always returns 3
      numbers).

   get_reciprocal_lattice_vectors(unit: str = 'angstrom-1') -> ndarray[Any, dtype[float64]]

      Get the reciprocal lattice vectors as a numpy array.

   get_volume(unit: str = 'angstrom') -> float

      Get the volume of the unit cell (for 2D: area, for 1D length)

   is_close(other: Lattice, tol: float = 0.001, unit: str = 'angstrom') -> bool
   is_close(other: Lattice, length_tol: float, angle_tol: float, length_unit: str = 'angstrom', angle_unit: str = 'rad') -> bool

      Checks if two sets of lattice vectors are close to each other.

   is_orthogonal() -> bool

      Is the lattice orthogonal?

   map_vector_to_central_cell(vec: ArrayLike, unit: str = 'angstrom') -> ndarray[Any, dtype[float64]]

      Maps a vector into the central cell, such that its fractional
      coordinates are in range [-0.5,+0.5). No mapping will take place
      in non-periodic directions.

   map_vectors_to_central_cell(vecs: ArrayLike, unit: str = 'angstrom') -> ndarray[Any, dtype[float64]]

      Maps a set of vectors into the central cell, such that their
      fractional coordinates are in range [-0.5,+0.5). No mapping will
      take place in non-periodic directions.

   write_kf(kf: KFFile, section: str) -> None

      Writes a Lattice to a section on a KF file.

   property num_vectors: int

      Returns the number of lattice vectors (i.e. the number of PBC).

   property vectors: AngstromCoordsArray

      The lattice vectors in angstrom.


ChemicalSystem methods for periodic systems
===========================================

Here are methods from the "ChemicalSystem" that relate to periodic
systems:

property ChemicalSystem.lattice: Lattice

   The lattice of the system.

ChemicalSystem.has_lattice() -> bool

   Checks if the system has a lattice.

   Note that the "lattice" attribute of a ChemicalSystem is never
   None. A system is considered to not have a lattice if the number of
   lattice vectors is 0. Using this method is the same a checking
   "lattice.num_vectors == 0."

ChemicalSystem.set_num_lattice_vectors(nvec: int) -> None

   Sets the number of lattice vectors.

   Can only be used to reduce the number of periodic directions. Note
   that this may change the number of bonds in the system: bonds that
   have lattice displacements in the no longer periodic directions
   will be removed.

ChemicalSystem.set_fractional_coordinates(frac_coords: ArrayLike, unit: str = 'angstrom') -> None

   Sets the fractional coordinates of all atoms. In non-periodic
   directions the plain coordinates should be passed in the specified
   unit.

ChemicalSystem.get_fractional_coordinates(unit: str = 'angstrom') -> ndarray[Any, dtype[float64]]

   Returns the fractional coordinates of all atoms with respect to the
   cell. In non-periodic directions the plain coordinates are returned
   in the specified unit.

ChemicalSystem.map_atoms(start_range: float | ArrayLike) -> Tuple[bool, ndarray[Any, dtype[int64]]]

   Maps all atoms into a unit cell from [start_range:start_range+1] in
   fractional coordinates.

   No mapping will take place in non-periodic directions. The returned
   tuple signifies if any mapping has taken place and the shifts in
   units of the lattice vectors that would move the mapped atoms back
   to their original position.

ChemicalSystem.map_atoms_around_atom(atom: int | Atom) -> Tuple[bool, ndarray[Any, dtype[int64]]]

   Map all atoms around the chosen atom that will be in the center of
   the new unit cell.

   No mapping will take place in non-periodic directions. The returned
   tuple signifies if any mapping has taken place and the shifts in
   units of the lattice vectors that would move the mapped atoms back
   to their original position.

   Note that the atom around which we are mapping does not move in the
   process. If you want the central atom to be in a special position
   (e.g. the origin), you will have to shift all atoms separately
   using the "translate" method.

ChemicalSystem.map_atoms_continuous() -> bool

   Map all atoms that are bonded across the periodic boundary back.

   Imagine a box full of molecules, with all atoms in the range [0,1]
   in fractional coordinates. Molecules that cross between cells will
   have bonds with lattice displacements in this case. Calling this
   method will remove all lattice displacements and move atoms outside
   of the range [0,1] in fractional coordinates. After calling it, all
   molecules are "continuous" and extend outside of the original cell.

   No mapping will take place in non-periodic directions.

   The return value indicates whether the method was able to
   completely remove all lattice displacements. In cases where a
   structure is connected to itself on both ends, it is not possible
   to make it continuous. An example of this would be a polymer chain
   going through the cell and connecting to itself across two cell
   boundaries at the same time. In this case the return value is
   "False", but disconnected parts of the system would still have been
   made continuous.

ChemicalSystem.apply_strain(strain_matrix: ArrayLike) -> None

   Apply a strain deformation to a periodic system.

   The atoms in the unit cell will be strained accordingly, keeping
   the fractional atomic coordinates constant.

   The strain_matrix argument should be (or be convertable to) an NxN
   numpy array, where N is the number of lattice vectors of the
   system.

ChemicalSystem.apply_strain_voigt(strain_voigt: Sequence[float]) -> None

   Apply a strain deformation to a periodic system.

   The atoms in the unit cell will be strained accordingly, keeping
   the fractional atomic coordinates constant.

   The strain_voigt argument should be passed in Voigt form:
      * for 3D periodic systems: [e_xx, e_yy, e_zz, gamma_yz,
        gamma_xz, gamma_xy]

      * for 2D periodic systems: [e_xx, e_yy, gamma_xy]

      * for 1D periodic systems: [e_xx]

   with e.g. e_xy = gamma_xy/2.

ChemicalSystem.perturb_lattice(noise_level: float) -> None

   Perturb the lattice vectors by applying random strain with matrix
   elements between [-noise_level,noise_level].

   This can be useful if you want to deviate from an ideal symmetric
   geometry, for example if you look for a phase change due to high
   pressure.

ChemicalSystem.make_supercell(supercell: ArrayLike) -> ChemicalSystem

   Create a supercell by scaling the lattice vectors.

   Returns a new system where copies of atoms will have the same
   properties as in the initial unit cell. Bonds are replicated
   between copies and across the new unit cell boundary.

ChemicalSystem.make_supercell_trafo(supercell: ArrayLike) -> ChemicalSystem

   Create a supercell by creating the matrix product of the lattice
   vectors and the supercell matrix.

   Returns a new system where copies of atoms will have the same
   properties as in the initial unit cell. Bonds are replicated
   between copies and across the new unit cell boundary.

ChemicalSystem.make_slab_thickness(miller: ArrayLike, top: float, bottom: float, translate: float = 0.0, unit: str = 'angstrom', ref_atom: int | None = None) -> ChemicalSystem

   Create a new system which is a 2D slab, by slicing a 3D system and
   specifying the upper and lower bound of the resulting slab.

   The "miller" argument accepts either three Miller indices "(h, k,
   l)" or four Miller-Bravais indices "(h, k, i, l)" for hexagonal
   lattices, in which case "h + k + i" must equal 0.

   The optional "ref_atom" specifies which atom should be used as the
   reference atom for the slice. It must be an integer in the range "0
   <= ref_atom < num_atoms". If it is not specified, slabs are
   generated for candidate reference atoms from a hard-coded search
   over at most the first 50 atoms, bonds are guessed for each
   candidate, and the slab with the highest sum of the minimum per-
   element coordination numbers is chosen. Ties are broken by
   preferring the smallest spread in those per-element minima.

ChemicalSystem.make_slab_layers(miller: ArrayLike, num_layers: int, translate: float = 0.0, unit: str = 'angstrom', ref_atom: int | None = None) -> ChemicalSystem

   Create a new system which is a 2D slab, by slicing a 3D system and
   specifying the number of layers in the resulting slab.

   The "miller" argument accepts either three Miller indices "(h, k,
   l)" or four Miller-Bravais indices "(h, k, i, l)" for hexagonal
   lattices, in which case "h + k + i" must equal 0.

   The optional "ref_atom" specifies which atom should be used as the
   reference atom for the slice. It must be an integer in the range "0
   <= ref_atom < num_atoms". If it is not specified, slabs are
   generated for candidate reference atoms from a hard-coded search
   over at most the first 50 atoms, bonds are guessed for each
   candidate, and the slab with the highest sum of the minimum per-
   element coordination numbers is chosen. Ties are broken by
   preferring the smallest spread in those per-element minima.

ChemicalSystem.make_primitive_cell(precision: float = 0.1) -> ChemicalSystem

   Create a new system by converting a 2D or 3D cell to its primitive
   representation.

ChemicalSystem.make_conventional_cell(precision: float = 0.1) -> ChemicalSystem

   Create a new system by converting a 2D or 3D cell to its
   conventional representation.

ChemicalSystem.density(unit: str = 'dalton/angstrom3') -> float

   Returns the density of the system in the specified unit. Only valid
   for 3D periodic systems.

ChemicalSystem.set_density(target_density: float, unit: str = 'dalton/angstrom3') -> None

   Applies a uniform strain to match the specified target density.
   Only valid for 3D periodic systems.
