Source code for mag4.lattice

"""Lattice geometry helpers (matrix construction, fractional ↔ Cartesian, PBC distances)."""

from __future__ import annotations

import numpy as np


[docs] def lattice_matrix(a: float, b: float, c: float, alpha_deg: float, beta_deg: float, gamma_deg: float) -> np.ndarray: """Build a 3×3 lattice matrix (columns = a, b, c) from cell parameters. Parameters ---------- a, b, c : float Lattice parameters in Å. alpha_deg, beta_deg, gamma_deg : float Cell angles in degrees. Returns ------- ndarray, shape (3, 3) Lattice matrix with the three lattice vectors as columns. """ alpha, beta, gamma = np.radians([alpha_deg, beta_deg, gamma_deg]) cos_a, cos_b, cos_g = np.cos(alpha), np.cos(beta), np.cos(gamma) sin_g = np.sin(gamma) v1 = np.array([a, 0.0, 0.0]) v2 = np.array([b * cos_g, b * sin_g, 0.0]) v3 = np.array([ c * cos_b, c * (cos_a - cos_b * cos_g) / sin_g, c * np.sqrt(max(0.0, 1 + 2 * cos_a * cos_b * cos_g - cos_a ** 2 - cos_b ** 2 - cos_g ** 2)) / sin_g, ]) return np.column_stack([v1, v2, v3])
[docs] def frac_to_cart(frac: np.ndarray, lat_mat: np.ndarray) -> np.ndarray: """Convert fractional → Cartesian coordinates.""" return (lat_mat @ frac.T).T
[docs] def get_cartesian_coords(atoms: list, lat_mat: np.ndarray) -> np.ndarray: """Cartesian coordinates of a list of atoms in the legacy ``[label, x, y, z, …]`` format.""" frac = np.array([[a[1], a[2], a[3]] for a in atoms]) return frac_to_cart(frac, lat_mat)
[docs] def pbc_distance(c1: np.ndarray, c2: np.ndarray, lat_vecs: np.ndarray, inv_lat: np.ndarray) -> float: """Minimum-image distance between two Cartesian points.""" df = inv_lat @ (c1 - c2) df -= np.round(df) return float(np.linalg.norm(lat_vecs.T @ df))
[docs] def find_site(new_cart: np.ndarray, cart_coords: np.ndarray, lat_vecs: np.ndarray, inv_lat: np.ndarray, tol: float = 0.05) -> int: """Index of the site in ``cart_coords`` closest to ``new_cart`` under PBC, or -1.""" best, best_d = -1, 1e9 for j, cj in enumerate(cart_coords): d = pbc_distance(new_cart, cj, lat_vecs, inv_lat) if d < best_d: best_d, best = d, j return best if best_d < tol else -1