qsym2/bindings/python/projection/
density.rs

1//! Python bindings for QSym² symmetry projection of electron densities.
2
3use std::path::PathBuf;
4
5use indexmap::IndexMap;
6use itertools::Itertools;
7use num_complex::Complex;
8use pyo3::exceptions::{PyIOError, PyRuntimeError};
9use pyo3::{IntoPyObjectExt, prelude::*};
10
11use crate::bindings::python::integrals::PyBasisAngularOrder;
12use crate::bindings::python::projection::PyProjectionTarget;
13use crate::bindings::python::representation_analysis::density::{
14    PyDensity, PyDensityComplex, PyDensityReal,
15};
16use crate::drivers::QSym2Driver;
17use crate::drivers::projection::density::{DensityProjectionDriver, DensityProjectionParams};
18use crate::drivers::representation_analysis::{
19    CharacterTableDisplay, MagneticSymmetryAnalysisKind,
20};
21use crate::drivers::symmetry_group_detection::SymmetryGroupDetectionResult;
22use crate::io::format::qsym2_output;
23use crate::io::{QSym2FileType, read_qsym2_binary};
24use crate::symmetry::symmetry_group::UnitaryRepresentedSymmetryGroup;
25use crate::symmetry::symmetry_transformation::SymmetryTransformationKind;
26use crate::target::density::Density;
27
28type C128 = Complex<f64>;
29
30// =====================
31// Functions definitions
32// =====================
33
34/// Python-exposed function to perform symmetry projection for real and complex electron densities.
35///
36/// # Arguments
37///
38/// * `inp_sym` - A path to the [`QSym2FileType::Sym`] file containing the symmetry-group detection
39/// result for the system. This will be used to construct abstract groups and character tables for
40/// symmetry projection.
41/// * `pydens` - A sequence of Python-exposed electron densities whose density matrices are of type
42/// `float64` or `complex128`. Each density is accompanied by a description string.
43/// * `projection_targets` - A sequence of subspace labels for projection. Each label is either a
44/// symbolic string or a numerical index for the subspace in the character table of the prevailing
45/// group.
46/// * `pybao` - Python-exposed structure containing basis angular order information for the density
47/// matrices.
48/// * `use_magnetic_group` - An option indicating if the magnetic group is to be used for symmetry
49/// analysis, and if so, whether unitary representations or unitary-antiunitary corepresentations
50/// should be used.
51/// * `use_double_group` - A boolean indicating if the double group of the prevailing symmetry
52/// group is to be used for representation analysis instead.
53/// * `symmetry_transformation_kind` - An enumerated type indicating the type of symmetry
54/// transformations to be performed on the origin electron density to generate the orbit.
55/// * `write_character_table` - A boolean indicating if the character table of the prevailing
56/// symmetry group is to be printed out.
57/// * `infinite_order_to_finite` - The finite order with which infinite-order generators are to be
58/// interpreted to form a finite subgroup of the prevailing infinite group. This finite subgroup
59/// will be used for symmetry analysis.
60///
61/// # Returns
62///
63/// The result will be returned as a list of tuples, each of which contains the name/description of
64/// an original density and a dictionary in which the keys are the subspace labels and the values
65/// are the corresponding projected density.
66#[pyfunction]
67#[pyo3(signature = (
68    inp_sym,
69    pydens,
70    projection_targets,
71    pybao,
72    use_magnetic_group,
73    use_double_group,
74    symmetry_transformation_kind,
75    write_character_table=true,
76    infinite_order_to_finite=None,
77))]
78pub fn project_densities(
79    py: Python<'_>,
80    inp_sym: PathBuf,
81    pydens: Vec<(String, PyDensity)>,
82    projection_targets: Vec<PyProjectionTarget>,
83    pybao: &PyBasisAngularOrder,
84    use_magnetic_group: Option<MagneticSymmetryAnalysisKind>,
85    use_double_group: bool,
86    symmetry_transformation_kind: SymmetryTransformationKind,
87    write_character_table: bool,
88    infinite_order_to_finite: Option<u32>,
89) -> PyResult<Vec<(String, IndexMap<String, Py<PyAny>>)>> {
90    let pd_res: SymmetryGroupDetectionResult =
91        read_qsym2_binary(inp_sym.clone(), QSym2FileType::Sym)
92            .map_err(|err| PyIOError::new_err(err.to_string()))?;
93
94    let mut file_name = inp_sym.to_path_buf();
95    file_name.set_extension(QSym2FileType::Sym.ext());
96    qsym2_output!(
97        "Symmetry-group detection results read in from {}.",
98        file_name.display(),
99    );
100    qsym2_output!("");
101
102    let mol = &pd_res.pre_symmetry.recentred_molecule;
103    let bao = pybao
104        .to_qsym2(mol)
105        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
106
107    let symbolic_projection_targets = projection_targets
108        .iter()
109        .filter_map(|pt| match pt {
110            PyProjectionTarget::Symbolic(symbolic_pt) => Some(symbolic_pt.clone()),
111            PyProjectionTarget::Numeric(_) => None,
112        })
113        .collect_vec();
114    let numeric_projection_targets = projection_targets
115        .iter()
116        .filter_map(|pt| match pt {
117            PyProjectionTarget::Symbolic(_) => None,
118            PyProjectionTarget::Numeric(numeric_pt) => Some(*numeric_pt),
119        })
120        .collect_vec();
121
122    let dp_params = DensityProjectionParams::builder()
123        .symbolic_projection_targets(Some(symbolic_projection_targets))
124        .numeric_projection_targets(Some(numeric_projection_targets))
125        .use_magnetic_group(use_magnetic_group.clone())
126        .use_double_group(use_double_group)
127        .symmetry_transformation_kind(symmetry_transformation_kind)
128        .write_character_table(if write_character_table {
129            Some(CharacterTableDisplay::Symbolic)
130        } else {
131            None
132        })
133        .infinite_order_to_finite(infinite_order_to_finite)
134        .build()
135        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
136
137    let any_complex = pydens
138        .iter()
139        .any(|(_, pyden)| matches!(pyden, PyDensity::Complex(_)));
140
141    let projected_densities = if !any_complex {
142        // All density matrices are real.
143        let dens = pydens
144            .iter()
145            .map(|(_, pyden)| match pyden {
146                PyDensity::Real(pyden_r) => pyden_r.to_qsym2(&bao, mol),
147                PyDensity::Complex(_) => panic!("Unexpected complex density."),
148            })
149            .collect::<Result<Vec<_>, _>>()
150            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
151        let dens_ref = dens
152            .iter()
153            .zip(pydens.iter())
154            .map(|(den, (desc, _))| (desc.clone(), den))
155            .collect::<Vec<_>>();
156
157        let projected_densities = match &use_magnetic_group {
158            Some(MagneticSymmetryAnalysisKind::Corepresentation) => {
159                return Err(PyRuntimeError::new_err(
160                    "Projection using corepresentations is not yet supported.",
161                ));
162            }
163            Some(MagneticSymmetryAnalysisKind::Representation) | None => {
164                let mut dp_driver =
165                    DensityProjectionDriver::<UnitaryRepresentedSymmetryGroup, f64>::builder()
166                        .parameters(&dp_params)
167                        .densities(dens_ref)
168                        .symmetry_group(&pd_res)
169                        .build()
170                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
171                py.detach(|| {
172                    dp_driver
173                        .run()
174                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))
175                })?;
176                dp_driver
177                    .result()
178                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))?
179                    .projected_densities()
180                    .clone()
181            }
182        };
183        projected_densities
184            .into_iter()
185            .map(|(name, projected_dens)| {
186                (
187                    name,
188                    projected_dens
189                        .into_iter()
190                        .map(|(row, den_res)| {
191                            (
192                                row.to_string(),
193                                den_res
194                                    .and_then(|den| {
195                                        let pyden = PyDensityReal {
196                                            complex_symmetric: den.complex_symmetric(),
197                                            density_matrix: den.density_matrix().clone(),
198                                            threshold: den.threshold(),
199                                        };
200                                        pyden.into_py_any(py).map_err(|err| err.to_string())
201                                    })
202                                    .expect("Unable to convert a projected density into a Python object."),
203                            )
204                        })
205                        .collect::<IndexMap<_, _>>(),
206                )
207            })
208            .collect::<Vec<_>>()
209    } else {
210        // At least one of coefficients or sao_4c are not real.
211        let dens: Vec<Density<C128>> = pydens
212            .iter()
213            .map(|(_, pyden)| match pyden {
214                PyDensity::Real(pyden_r) => pyden_r.to_qsym2(&bao, mol).map(|den_r| den_r.into()),
215                PyDensity::Complex(pyden_c) => pyden_c.to_qsym2(&bao, mol),
216            })
217            .collect::<Result<Vec<_>, _>>()
218            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
219        let dens_ref = dens
220            .iter()
221            .zip(pydens.iter())
222            .map(|(den, (desc, _))| (desc.clone(), den))
223            .collect::<Vec<_>>();
224
225        let projected_densities = match &use_magnetic_group {
226            Some(MagneticSymmetryAnalysisKind::Corepresentation) => {
227                return Err(PyRuntimeError::new_err(
228                    "Projection using corepresentations is not yet supported.",
229                ));
230            }
231            Some(MagneticSymmetryAnalysisKind::Representation) | None => {
232                let mut dp_driver =
233                    DensityProjectionDriver::<UnitaryRepresentedSymmetryGroup, C128>::builder()
234                        .parameters(&dp_params)
235                        .densities(dens_ref)
236                        .symmetry_group(&pd_res)
237                        .build()
238                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
239                py.detach(|| {
240                    dp_driver
241                        .run()
242                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))
243                })?;
244                dp_driver
245                    .result()
246                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))?
247                    .projected_densities()
248                    .clone()
249            }
250        };
251        projected_densities
252            .into_iter()
253            .map(|(name, projected_dens)| {
254                (
255                    name,
256                    projected_dens
257                        .into_iter()
258                        .map(|(row, den_res)| {
259                            (
260                                row.to_string(),
261                                den_res
262                                    .and_then(|den| {
263                                        let pyden = PyDensityComplex {
264                                            complex_symmetric: den.complex_symmetric(),
265                                            density_matrix: den.density_matrix().clone(),
266                                            threshold: den.threshold(),
267                                        };
268                                        pyden.into_py_any(py).map_err(|err| err.to_string())
269                                    })
270                                    .expect("Unable to convert a projected density into a Python object."),
271                            )
272                        })
273                        .collect::<IndexMap<_, _>>(),
274                )
275            })
276            .collect::<Vec<_>>()
277    };
278    Ok(projected_densities)
279}
qsym2/bindings/python/projection/density.rs

qsym2/bindings/python/projection/
density.rs
