qsym2/bindings/python/
molecule_symmetrisation.rs

1//! Python bindings for QSym² molecule symmetrisation by bootstrapping.
2//!
3//! See [`crate::drivers::molecule_symmetrisation_bootstrap`] for more information.
4
5use std::path::PathBuf;
6
7use anyhow::format_err;
8use pyo3::exceptions::PyRuntimeError;
9use pyo3::prelude::*;
10
11use crate::auxiliary::atom::AtomKind;
12use crate::auxiliary::molecule::Molecule;
13use crate::bindings::python::symmetry_group_detection::PyMolecule;
14use crate::drivers::QSym2Driver;
15use crate::drivers::molecule_symmetrisation_bootstrap::{
16    MoleculeSymmetrisationBootstrapDriver, MoleculeSymmetrisationBootstrapParams,
17};
18#[allow(unused_imports)]
19use crate::io::QSym2FileType;
20
21/// Python-exposed function to perform molecule symmetrisation by bootstrapping and log the result
22/// via the `qsym2-output` logger at the `INFO` level.
23///
24/// See [`crate::drivers::molecule_symmetrisation_bootstrap`] for more information.
25///
26/// # Arguments
27///
28/// * `inp_xyz` - An optional string providing the path to an XYZ file containing the molecule to
29/// be symmetrised. Only one of `inp_xyz` or `inp_mol` can be specified.
30/// * `inp_mol` - An optional `PyMolecule` structure containing the molecule to be symmetrised. Only
31/// one of `inp_xyz` or `inp_mol` can be specified.
32/// * `out_target_sym` - An optional path for a [`QSym2FileType::Sym`] file to be saved that
33/// contains the symmetry-group detection results of the symmetrised molecule at the target
34/// thresholds.
35/// * `loose_moi_threshold` - The loose MoI threshold.
36/// * `loose_distance_threshold` - The loose distance threshold.
37/// * `target_moi_threshold` - The target (tight) MoI threshold.
38/// * `target_distance_threshold` - The target (tight) distance threshold.
39/// * `use_magnetic_group` - A boolean indicating if the magnetic group (*i.e.* the group including
40/// time-reversed operations) is to be used for the symmetrisation.
41/// * `reorientate_molecule` - A boolean indicating if the molecule is also reoriented to align its
42/// principal axes with the Cartesian axes.
43/// * `max_iterations` - The maximum number of iterations for the symmetrisation process.
44/// * `consistent_target_symmetry_iterations` - The number of consecutive iterations during which
45/// the symmetry group at the target level of threshold must be consistently found for convergence
46/// to be reached, if this group cannot become identical to the symmetry group at the loose level
47/// of threshold.
48/// * `verbose` - The print-out level.
49/// * `infinite_order_to_finite` - The finite order with which infinite-order generators are to be
50/// interpreted to form a finite subgroup of the prevailing infinite group. This finite subgroup
51/// will be used for the symmetrisation.
52///
53/// # Returns
54///
55/// The symmetrised molecule.
56///
57/// # Errors
58///
59/// Errors if any intermediate step in the symmetrisation procedure fails.
60#[allow(clippy::too_many_arguments)]
61#[pyfunction]
62#[pyo3(signature = (
63    inp_xyz,
64    inp_mol,
65    out_target_sym=None,
66    loose_moi_threshold=1e-2,
67    loose_distance_threshold=1e-2,
68    target_moi_threshold=1e-7,
69    target_distance_threshold=1e-7,
70    use_magnetic_group=false,
71    reorientate_molecule=true,
72    max_iterations=50,
73    consistent_target_symmetry_iterations=10,
74    verbose=0,
75    infinite_order_to_finite=None
76))]
77pub fn symmetrise_molecule(
78    py: Python<'_>,
79    inp_xyz: Option<PathBuf>,
80    inp_mol: Option<PyMolecule>,
81    out_target_sym: Option<PathBuf>,
82    loose_moi_threshold: f64,
83    loose_distance_threshold: f64,
84    target_moi_threshold: f64,
85    target_distance_threshold: f64,
86    use_magnetic_group: bool,
87    reorientate_molecule: bool,
88    max_iterations: usize,
89    consistent_target_symmetry_iterations: usize,
90    verbose: u8,
91    infinite_order_to_finite: Option<u32>,
92) -> PyResult<PyMolecule> {
93    py.detach(|| {
94        let mol = match (inp_xyz, inp_mol) {
95            (Some(xyz_path), None) => Molecule::from_xyz(xyz_path, 1e-7),
96            (None, Some(pymol)) => Molecule::from(pymol),
97            _ => {
98                return Err(PyRuntimeError::new_err(
99                    "One and only one of `inp_xyz` or `inp_mol` must be specified.",
100                ));
101            }
102        };
103
104        let msb_params = MoleculeSymmetrisationBootstrapParams::builder()
105            .reorientate_molecule(reorientate_molecule)
106            .use_magnetic_group(use_magnetic_group)
107            .loose_moi_threshold(loose_moi_threshold)
108            .loose_distance_threshold(loose_distance_threshold)
109            .target_moi_threshold(target_moi_threshold)
110            .target_distance_threshold(target_distance_threshold)
111            .infinite_order_to_finite(infinite_order_to_finite)
112            .max_iterations(max_iterations)
113            .consistent_target_symmetry_iterations(consistent_target_symmetry_iterations)
114            .verbose(verbose)
115            .symmetrised_result_save_name(out_target_sym)
116            .build()
117            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
118
119        let mut msb_driver = MoleculeSymmetrisationBootstrapDriver::builder()
120            .parameters(&msb_params)
121            .molecule(&mol)
122            .build()
123            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
124        msb_driver
125            .run()
126            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
127
128        let symmol = &msb_driver
129            .result()
130            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?
131            .symmetrised_molecule;
132
133        // Note that the magnetic field here will not have the original magnitude as it is back-deduded
134        // from the magnetic atoms which have been added using the normalised version of the original
135        // magnetic field.
136        let magnetic_field = symmol
137            .magnetic_atoms
138            .as_ref()
139            .map(|mag_atoms| {
140                if mag_atoms.len() != 2 {
141                    Err(format_err!("Only a uniform magnetic field is supported."))
142                } else {
143                    match (&mag_atoms[0].kind, &mag_atoms[1].kind) {
144                        (AtomKind::Magnetic(true), AtomKind::Magnetic(false)) => {
145                            let bvec = mag_atoms[0].coordinates - mag_atoms[1].coordinates;
146                            Ok([bvec[0], bvec[1], bvec[2]])
147                        }
148                        (AtomKind::Magnetic(false), AtomKind::Magnetic(true)) => {
149                            let bvec = mag_atoms[1].coordinates - mag_atoms[0].coordinates;
150                            Ok([bvec[0], bvec[1], bvec[2]])
151                        }
152                        _ => Err(format_err!("Invalid fictitious magnetic atoms detected.")),
153                    }
154                }
155            })
156            .transpose()
157            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
158
159        // Note that the electric field here will not have the original magnitude as it is back-deduded
160        // from the electric atoms which have been added using the normalised version of the original
161        // electric field.
162        let electric_field = symmol
163            .electric_atoms
164            .as_ref()
165            .map(|elec_atoms| {
166                if elec_atoms.len() != 1 {
167                    Err(format_err!("Only a uniform electric field is supported."))
168                } else {
169                    match &elec_atoms[0].kind {
170                        AtomKind::Electric(pos) => {
171                            let evec = if *pos {
172                                elec_atoms[0].coordinates
173                            } else {
174                                -elec_atoms[0].coordinates
175                            };
176                            Ok([evec[0], evec[1], evec[2]])
177                        }
178                        _ => Err(format_err!("Invalid fictitious electric atoms detected.")),
179                    }
180                }
181            })
182            .transpose()
183            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
184
185        let pymol = PyMolecule::new(
186            symmol
187                .atoms
188                .iter()
189                .map(|atom| {
190                    (
191                        atom.atomic_symbol.clone(),
192                        [
193                            atom.coordinates[0],
194                            atom.coordinates[1],
195                            atom.coordinates[2],
196                        ],
197                    )
198                })
199                .collect::<Vec<_>>(),
200            symmol.threshold,
201            magnetic_field,
202            electric_field,
203        );
204        Ok(pymol)
205    })
206}