[WIP] Magnetic space group identification

moyo/Cargo.toml

@@ -24,7 +24,8 @@ thiserror = "1.0"
 union-find = "0.4"
 strum = "0.25"
 strum_macros = "0.25"
-kiddo = "4.2"
+kiddo = "5.0.3"
+once_cell = "1.20.2"
 
 [dev-dependencies]
 rand = "0.8"

moyo/src/base/error.rs

@@ -27,6 +27,10 @@ pub enum MoyoError {
     ArithmeticCrystalClassIdentificationError,
     #[error("Space group type identification failed")]
     SpaceGroupTypeIdentificationError,
+    #[error("Construct type identification failed")]
+    ConstructTypeIdentificationError,
+    #[error("Magnetic space group type identification failed")]
+    MagneticSpaceGroupTypeIdentificationError,
     #[error("Standardization failed")]
     StandardizationError,
     #[error("Wyckoff position assignment failed")]

moyo/src/base/transformation.rs

@@ -1,10 +1,12 @@
+use std::ops::Mul;
+
 use itertools::iproduct;
 use nalgebra::base::{Matrix3, Vector3};
 
 use super::cell::Cell;
 use super::lattice::Lattice;
 use super::magnetic_cell::{MagneticCell, MagneticMoment};
-use super::operation::Operation;
+use super::operation::{MagneticOperation, Operation};
 use crate::math::SNF;
 
 pub type UnimodularLinear = Matrix3<i32>;
@@ -73,6 +75,24 @@ impl UnimodularTransformation {
             .collect()
     }
 
+    pub fn transform_magnetic_operation(
+        &self,
+        magnetic_operation: &MagneticOperation,
+    ) -> MagneticOperation {
+        let new_operation = self.transform_operation(&magnetic_operation.operation);
+        MagneticOperation::from_operation(new_operation, magnetic_operation.time_reversal)
+    }
+
+    pub fn transform_magnetic_operations(
+        &self,
+        magnetic_operations: &[MagneticOperation],
+    ) -> Vec<MagneticOperation> {
+        magnetic_operations
+            .iter()
+            .map(|mops| self.transform_magnetic_operation(mops))
+            .collect()
+    }
+
     pub fn transform_cell(&self, cell: &Cell) -> Cell {
         let new_lattice = self.transform_lattice(&cell.lattice);
         // (P, p)^-1 x = P^-1 (x - p)
@@ -98,6 +118,17 @@ impl UnimodularTransformation {
     }
 }
 
+impl Mul for UnimodularTransformation {
+    type Output = Self;
+
+    // (P_lhs, p_lhs) * (P_rhs, p_rhs) = (P_lhs * P_rhs, P_lhs * p_rhs + p_lhs)
+    fn mul(self, rhs: Self) -> Self::Output {
+        let new_linear = self.linear * rhs.linear;
+        let new_origin_shift = self.linear.map(|e| e as f64) * rhs.origin_shift + self.origin_shift;
+        Self::new(new_linear, new_origin_shift)
+    }
+}
+
 /// Represent change of origin and basis for an affine space
 #[derive(Debug, Clone)]
 pub struct Transformation {
@@ -188,6 +219,48 @@ impl Transformation {
             .collect()
     }
 
+    pub fn transform_magnetic_operation(
+        &self,
+        magnetic_operation: &MagneticOperation,
+    ) -> Option<MagneticOperation> {
+        let new_operation = self.transform_operation(&magnetic_operation.operation)?;
+        Some(MagneticOperation::from_operation(
+            new_operation,
+            magnetic_operation.time_reversal,
+        ))
+    }
+
+    pub fn transform_magnetic_operations(
+        &self,
+        magnetic_operations: &[MagneticOperation],
+    ) -> Vec<MagneticOperation> {
+        magnetic_operations
+            .iter()
+            .filter_map(|mops| self.transform_magnetic_operation(mops))
+            .collect()
+    }
+
+    pub fn inverse_transform_magnetic_operation(
+        &self,
+        magnetic_operation: &MagneticOperation,
+    ) -> Option<MagneticOperation> {
+        let new_operation = self.inverse_transform_operation(&magnetic_operation.operation)?;
+        Some(MagneticOperation::from_operation(
+            new_operation,
+            magnetic_operation.time_reversal,
+        ))
+    }
+
+    pub fn inverse_transform_magnetic_operations(
+        &self,
+        magnetic_operations: &[MagneticOperation],
+    ) -> Vec<MagneticOperation> {
+        magnetic_operations
+            .iter()
+            .filter_map(|mops| self.inverse_transform_magnetic_operation(mops))
+            .collect()
+    }
+
     // The transformation may increase the number of atoms in the cell.
     // Return the transformed cell and mapping from sites in the transformed cell to sites in the original cell.
     pub fn transform_cell(&self, cell: &Cell) -> (Cell, Vec<usize>) {

moyo/src/data.rs

@@ -11,13 +11,18 @@ mod wyckoff;
 
 pub use arithmetic_crystal_class::ArithmeticNumber;
 pub use centering::Centering;
-pub use hall_symbol::HallSymbol;
+pub use hall_symbol::{HallSymbol, MagneticHallSymbol};
 pub use hall_symbol_database::{hall_symbol_entry, HallNumber, HallSymbolEntry, Number};
+pub use magnetic_hall_symbol_database::{magnetic_hall_symbol_entry, MagneticHallSymbolEntry};
+pub use magnetic_space_group::{
+    get_magnetic_space_group_type, ConstructType, UNINumber, NUM_MAGNETIC_SPACE_GROUP_TYPES,
+};
 pub use setting::Setting;
 
 pub(super) use arithmetic_crystal_class::{
     arithmetic_crystal_class_entry, iter_arithmetic_crystal_entry,
 };
 pub(super) use classification::{CrystalSystem, GeometricCrystalClass, LatticeSystem};
+pub(super) use magnetic_space_group::uni_number_range;
 pub(super) use point_group::PointGroupRepresentative;
 pub(super) use wyckoff::{iter_wyckoff_positions, WyckoffPosition, WyckoffPositionSpace};

moyo/src/data/hall_symbol.rs

@@ -6,6 +6,8 @@ use nalgebra::{matrix, Vector3};
 
 use super::centering::Centering;
 use super::hall_symbol_database::{hall_symbol_entry, HallNumber};
+use super::magnetic_hall_symbol_database::magnetic_hall_symbol_entry;
+use super::magnetic_space_group::UNINumber;
 use crate::base::{
     MagneticOperation, MagneticOperations, Operation, Operations, OriginShift, Rotation,
     TimeReversal, Transformation, Translation, EPS,
@@ -102,6 +104,12 @@ impl HallSymbol {
         operations
     }
 
+    pub fn primitive_traverse(&self) -> Operations {
+        let operations = self.traverse();
+        Transformation::from_linear(self.centering.linear())
+            .inverse_transform_operations(&operations)
+    }
+
     pub fn from_hall_number(hall_number: HallNumber) -> Option<Self> {
         if let Some(entry) = hall_symbol_entry(hall_number) {
             Self::new(entry.hall_symbol)
@@ -206,6 +214,25 @@ impl MagneticHallSymbol {
 
         operations
     }
+
+    pub fn primitive_traverse(&self) -> MagneticOperations {
+        let operations = self.traverse();
+        Transformation::from_linear(self.centering.linear())
+            .inverse_transform_magnetic_operations(&operations)
+    }
+
+    pub fn from_uni_number(uni_number: UNINumber) -> Option<Self> {
+        if let Some(entry) = magnetic_hall_symbol_entry(uni_number) {
+            Self::new(entry.magnetic_hall_symbol)
+        } else {
+            None
+        }
+    }
+
+    pub fn primitive_generators(&self) -> MagneticOperations {
+        Transformation::from_linear(self.centering.linear())
+            .inverse_transform_magnetic_operations(&self.generators)
+    }
 }
 
 /// Tokenize string by whitespaces

moyo/src/data/magnetic_hall_symbol_database.rs

@@ -1,4 +1,8 @@
-use super::magnetic_space_group::UNINumber;
+use super::hall_symbol_database::HallNumber;
+use super::magnetic_space_group::{
+    get_magnetic_space_group_type, ConstructType, UNINumber, NUM_MAGNETIC_SPACE_GROUP_TYPES,
+};
+use super::setting::Setting;
 
 #[derive(Debug, Clone)]
 pub struct MagneticHallSymbolEntry {
@@ -13,23 +17,29 @@ impl MagneticHallSymbolEntry {
             uni_number,
         }
     }
-}
 
-// TODO:
-// smallest and largest UNI numbers for each Hall number
-// let hall_number_to_uni_numbers: HashMap<HallNumber, (UNINumber, UNINumber)>;
+    pub fn construct_type(&self) -> ConstructType {
+        get_magnetic_space_group_type(self.uni_number)
+            .unwrap()
+            .construct_type
+    }
 
-// TODO:
-// smallest and largest Hall numbers for XSG of each UNI number
-// let uni_number_to_xsg_hall_numbers: HashMap<UNINumber, (HallNumber, HallNumber)>;
+    pub fn reference_hall_number(&self) -> HallNumber {
+        let number = get_magnetic_space_group_type(self.uni_number)
+            .unwrap()
+            .number;
+        Setting::Standard.hall_number(number).unwrap()
+    }
+}
 
 pub fn magnetic_hall_symbol_entry(uni_number: UNINumber) -> Option<MagneticHallSymbolEntry> {
     MAGNETIC_HALL_SYMBOL_DATABASE
         .get((uni_number - 1) as usize)
         .cloned()
 }
 
-const MAGNETIC_HALL_SYMBOL_DATABASE: [MagneticHallSymbolEntry; 1651] = [
+// For type-IV MSG, we assume their Magnetic Hall symbol explicitly contains the anti-translation operation.
+const MAGNETIC_HALL_SYMBOL_DATABASE: [MagneticHallSymbolEntry; NUM_MAGNETIC_SPACE_GROUP_TYPES] = [
     MagneticHallSymbolEntry::new("P 1", 1),
     MagneticHallSymbolEntry::new("P 1 1'", 2),
     MagneticHallSymbolEntry::new("P 1 1c'", 3),
@@ -1686,7 +1696,7 @@ const MAGNETIC_HALL_SYMBOL_DATABASE: [MagneticHallSymbolEntry; 1651] = [
 #[cfg(test)]
 mod tests {
     use super::*;
-    use crate::data::hall_symbol::MagneticHallSymbol;
+    use crate::data::MagneticHallSymbol;
 
     fn iter_magnetic_hall_symbol_entry() -> impl Iterator<Item = &'static MagneticHallSymbolEntry> {
         MAGNETIC_HALL_SYMBOL_DATABASE.iter()