Skip to main content

qsym2/bindings/python/representation_analysis/multideterminant/
multideterminant_fci_basis.rs

1//! Python bindings for QSym² symmetry analysis of multi-determinants with FCI bases.
2
3use std::path::PathBuf;
4
5use itertools::Itertools;
6use num_complex::Complex;
7use numpy::{PyArray1, PyArrayMethods};
8use pyo3::exceptions::{PyIOError, PyRuntimeError};
9use pyo3::prelude::*;
10
11use crate::analysis::EigenvalueComparisonMode;
12use crate::angmom::spinor_rotation_3d::{SpinConstraint, SpinOrbitCoupled};
13use crate::bindings::python::integrals::{PyBasisAngularOrder, PyStructureConstraint};
14use crate::bindings::python::representation_analysis::slater_determinant::PySlaterDeterminant;
15use crate::bindings::python::representation_analysis::{PyArray1RC, PyArray2RC};
16use crate::drivers::QSym2Driver;
17use crate::drivers::representation_analysis::angular_function::AngularFunctionRepAnalysisParams;
18use crate::drivers::representation_analysis::multideterminant::{
19    MultiDeterminantRepAnalysisDriver, MultiDeterminantRepAnalysisParams,
20};
21use crate::drivers::representation_analysis::{
22    CharacterTableDisplay, MagneticSymmetryAnalysisKind,
23};
24use crate::drivers::symmetry_group_detection::SymmetryGroupDetectionResult;
25use crate::io::format::qsym2_output;
26use crate::io::{QSym2FileType, read_qsym2_binary};
27use crate::symmetry::symmetry_group::{
28    MagneticRepresentedSymmetryGroup, UnitaryRepresentedSymmetryGroup,
29};
30use crate::symmetry::symmetry_transformation::SymmetryTransformationKind;
31use crate::target::determinant::SlaterDeterminant;
32use crate::target::noci::basis::FCIBasis;
33use crate::target::noci::multideterminant::MultiDeterminant;
34use crate::target::noci::multideterminants::MultiDeterminants;
35
36type C128 = Complex<f64>;
37
38/// Python-exposed function to perform representation symmetry analysis for real and complex
39/// multi-determinantal wavefunctions constructed from a FCI basis based on a reference Slater
40/// determinant and log the result via the `qsym2-output` logger at the `INFO` level.
41///
42/// If `symmetry_transformation_kind` includes spin transformation, the provided
43/// multi-determinantal wavefunctions will be augmented to generalised spin constraint
44/// automatically.
45///
46/// # Arguments
47///
48/// * `inp_sym` - A path to the [`QSym2FileType::Sym`] file containing the symmetry-group detection
49/// result for the system. This will be used to construct abstract groups and character tables for
50/// representation analysis.
51/// * `pydet` - The Python-exposed reference Slater determinants whose coefficients are of type
52/// `float64` or `complex128`. This determinant serves as the reference for the
53/// full-configuration-interaction basis and should contain the full set of occupied and virtual
54/// molecular orbitals.
55/// * `coefficients` - The coefficient matrix where each column gives the linear combination
56/// coefficients for one multi-determinantal wavefunction. The number of rows must match the number
57/// of elements in the full-configuration-interaction basis, which is also the number of elements in
58/// `occupation_patterns`. The elements are of type `float64` or `complex128`.
59/// * `energies` - The `float64` or `complex128` energies of the multi-determinantal wavefunctions.
60/// The number of terms must match the number of columns of `coefficients`.
61/// * `occupation_patterns` - The occupation patterns of the determinants in the
62/// full-configuration-interaction basis. Each element specifies one determinant in the basis.
63/// * `pybaos` - Python-exposed structures containing basis angular order information, one for each
64/// explicit component per coefficient matrix.
65/// * `integrality_threshold` - The threshold for verifying if subspace multiplicities are
66/// integral.
67/// * `linear_independence_threshold` - The threshold for determining the linear independence
68/// subspace via the non-zero eigenvalues of the orbit overlap matrix.
69/// * `orbit_overlap_symmetry_threshold` - The threshold for checking the symmetry of the orbit
70/// overlap matrix.
71/// * `use_magnetic_group` - An option indicating if the magnetic group is to be used for symmetry
72/// analysis, and if so, whether unitary representations or unitary-antiunitary corepresentations
73/// should be used.
74/// * `use_double_group` - A boolean indicating if the double group of the prevailing symmetry
75/// group is to be used for representation analysis instead.
76/// * `use_cayley_table` - A boolean indicating if the Cayley table for the group, if available,
77/// should be used to speed up the calculation of orbit overlap matrices.
78/// * `symmetry_transformation_kind` - An enumerated type indicating the type of symmetry
79/// transformations to be performed on the origin determinant to generate the orbit. If this
80/// contains spin transformation, the multi-determinant will be augmented to generalised spin
81/// constraint automatically.
82/// * `eigenvalue_comparison_mode` - An enumerated type indicating the mode of comparison of orbit
83/// overlap eigenvalues with the specified `linear_independence_threshold`.
84/// * `sao` - The atomic-orbital overlap matrix whose elements are of type `float64` or
85/// `complex128`.
86/// * `sao_h` - The optional complex-symmetric atomic-orbital overlap matrix whose elements
87/// are of type `float64` or `complex128`. This is required if antiunitary symmetry operations are
88/// involved.
89/// * `write_overlap_eigenvalues` - A boolean indicating if the eigenvalues of the determinant
90/// orbit overlap matrix are to be written to the output.
91/// * `write_character_table` - A boolean indicating if the character table of the prevailing
92/// symmetry group is to be printed out.
93/// * `infinite_order_to_finite` - The finite order with which infinite-order generators are to be
94/// interpreted to form a finite subgroup of the prevailing infinite group. This finite subgroup
95/// will be used for symmetry analysis.
96/// * `angular_function_integrality_threshold` - The threshold for verifying if subspace
97/// multiplicities are integral for the symmetry analysis of angular functions.
98/// * `angular_function_linear_independence_threshold` - The threshold for determining the linear
99/// independence subspace via the non-zero eigenvalues of the orbit overlap matrix for the symmetry
100/// analysis of angular functions.
101/// * `angular_function_orbit_overlap_symmetry_threshold` - The threshold for checking the symmetry
102/// of the orbit overlap matrices for the symmetry analysis of angular functions.
103/// * `angular_function_max_angular_momentum` - The maximum angular momentum order to be used in
104/// angular function symmetry analysis.
105#[allow(clippy::too_many_arguments)]
106#[pyfunction]
107#[pyo3(signature = (
108    inp_sym,
109    pydet,
110    coefficients,
111    energies,
112    occupation_patterns,
113    pybaos,
114    integrality_threshold,
115    linear_independence_threshold,
116    orbit_overlap_symmetry_threshold,
117    use_magnetic_group,
118    use_double_group,
119    use_cayley_table,
120    symmetry_transformation_kind,
121    eigenvalue_comparison_mode,
122    sao,
123    sao_h=None,
124    write_overlap_eigenvalues=true,
125    write_character_table=true,
126    infinite_order_to_finite=None,
127    angular_function_integrality_threshold=1e-7,
128    angular_function_linear_independence_threshold=1e-7,
129    angular_function_orbit_overlap_symmetry_threshold=1e-7,
130    angular_function_max_angular_momentum=2
131))]
132pub fn rep_analyse_multideterminants_fci_basis(
133    py: Python<'_>,
134    inp_sym: PathBuf,
135    pydet: PySlaterDeterminant,
136    coefficients: PyArray2RC,
137    energies: PyArray1RC,
138    occupation_patterns: Vec<Vec<Bound<'_, PyArray1<f64>>>>,
139    pybaos: Vec<PyBasisAngularOrder>,
140    integrality_threshold: f64,
141    linear_independence_threshold: f64,
142    orbit_overlap_symmetry_threshold: f64,
143    use_magnetic_group: Option<MagneticSymmetryAnalysisKind>,
144    use_double_group: bool,
145    use_cayley_table: bool,
146    symmetry_transformation_kind: SymmetryTransformationKind,
147    eigenvalue_comparison_mode: EigenvalueComparisonMode,
148    sao: PyArray2RC,
149    sao_h: Option<PyArray2RC>,
150    write_overlap_eigenvalues: bool,
151    write_character_table: bool,
152    infinite_order_to_finite: Option<u32>,
153    angular_function_integrality_threshold: f64,
154    angular_function_linear_independence_threshold: f64,
155    angular_function_orbit_overlap_symmetry_threshold: f64,
156    angular_function_max_angular_momentum: u32,
157) -> PyResult<()> {
158    // Read in point-group detection results
159    let pd_res: SymmetryGroupDetectionResult =
160        read_qsym2_binary(inp_sym.clone(), QSym2FileType::Sym)
161            .map_err(|err| PyIOError::new_err(err.to_string()))?;
162
163    let mut file_name = inp_sym.to_path_buf();
164    file_name.set_extension(QSym2FileType::Sym.ext());
165    qsym2_output!(
166        "Symmetry-group detection results read in from {}.",
167        file_name.display(),
168    );
169    qsym2_output!("");
170
171    // Set up basic parameters
172    let mol = &pd_res.pre_symmetry.recentred_molecule;
173
174    let baos = pybaos
175        .iter()
176        .map(|bao| {
177            bao.to_qsym2(mol)
178                .map_err(|err| PyRuntimeError::new_err(err.to_string()))
179        })
180        .collect::<Result<Vec<_>, _>>()?;
181    let baos_ref = baos.iter().collect::<Vec<_>>();
182    let augment_to_generalised = match symmetry_transformation_kind {
183        SymmetryTransformationKind::SpatialWithSpinTimeReversal
184        | SymmetryTransformationKind::Spin
185        | SymmetryTransformationKind::SpinSpatial => true,
186        SymmetryTransformationKind::Spatial => false,
187    };
188    let afa_params = AngularFunctionRepAnalysisParams::builder()
189        .integrality_threshold(angular_function_integrality_threshold)
190        .linear_independence_threshold(angular_function_linear_independence_threshold)
191        .orbit_overlap_symmetry_threshold(angular_function_orbit_overlap_symmetry_threshold)
192        .max_angular_momentum(angular_function_max_angular_momentum)
193        .build()
194        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
195    let mda_params = MultiDeterminantRepAnalysisParams::<f64>::builder()
196        .integrality_threshold(integrality_threshold)
197        .linear_independence_threshold(linear_independence_threshold)
198        .orbit_overlap_symmetry_threshold(orbit_overlap_symmetry_threshold)
199        .use_magnetic_group(use_magnetic_group.clone())
200        .use_double_group(use_double_group)
201        .use_cayley_table(use_cayley_table)
202        .symmetry_transformation_kind(symmetry_transformation_kind.clone())
203        .eigenvalue_comparison_mode(eigenvalue_comparison_mode)
204        .write_overlap_eigenvalues(write_overlap_eigenvalues)
205        .write_character_table(if write_character_table {
206            Some(CharacterTableDisplay::Symbolic)
207        } else {
208            None
209        })
210        .infinite_order_to_finite(infinite_order_to_finite)
211        .build()
212        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
213
214    let real = matches!(pydet, PySlaterDeterminant::Real(_));
215
216    let structure_constraint = match pydet {
217        PySlaterDeterminant::Real(ref pd) => pd.structure_constraint().clone(),
218        PySlaterDeterminant::Complex(ref pd) => pd.structure_constraint().clone(),
219    };
220
221    let occupation_patterns = occupation_patterns
222        .into_iter()
223        .map(|occs| {
224            occs.into_iter()
225                .map(|occ| occ.to_owned_array())
226                .collect_vec()
227        })
228        .collect_vec();
229
230    // Decision tree:
231    // - all real numerical data?
232    //   + yes:
233    //     - structure_constraint:
234    //       + SpinConstraint:
235    //         - use_magnetic_group:
236    //           + Some(Corepresentation)
237    //           + Some(Representation) | None
238    //       + SpinOrbitCoupled: not supported
239    //   - no:
240    //     - structure_constraint:
241    //       + SpinConstraint:
242    //         - use_magnetic_group:
243    //           + Some(Corepresentation)
244    //           + Some(Representation) | None
245    //       + SpinOrbitCoupled:
246    //         - use_magnetic_group:
247    //           + Some(Corepresentation)
248    //           + Some(Representation) | None
249    match (real, &coefficients, &energies, &sao) {
250        (
251            true,
252            PyArray2RC::Real(pycoefficients_r),
253            PyArray1RC::Real(pyenergies_r),
254            PyArray2RC::Real(pysao_r),
255        ) => {
256            // Real numeric data type
257            if matches!(
258                structure_constraint,
259                PyStructureConstraint::SpinOrbitCoupled(_)
260            ) {
261                return Err(PyRuntimeError::new_err(
262                    "Real determinants cannot support spin--orbit-coupled structure constraint.",
263                ));
264            }
265
266            // Preparation
267            let sao_r = pysao_r.to_owned_array();
268            let coefficients_r = pycoefficients_r.to_owned_array();
269            let energies_r = pyenergies_r.to_owned_array();
270            let det_r = if augment_to_generalised {
271                if let PySlaterDeterminant::Real(ref pydet_r) = pydet {
272                    pydet_r
273                        .to_qsym2(&baos_ref, mol)
274                        .map(|det_r_inner| det_r_inner.to_generalised())
275                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))
276                } else {
277                    Err(PyRuntimeError::new_err(
278                        "Unexpected complex type for a Slater determinant.".to_string(),
279                    ))
280                }
281            } else {
282                if let PySlaterDeterminant::Real(ref pydet_r) = pydet {
283                    pydet_r
284                        .to_qsym2(&baos_ref, mol)
285                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))
286                } else {
287                    Err(PyRuntimeError::new_err(
288                        "Unexpected complex type for a Slater determinant.".to_string(),
289                    ))
290                }
291            }?;
292
293            let n_energies = energies_r.len();
294            if coefficients_r.shape()[1] != n_energies {
295                return Err(PyRuntimeError::new_err(
296                    "Mismatched number of FCI energies and number of FCI states.",
297                ));
298            }
299
300            // Construct the FCI basis
301            let fci_basis = FCIBasis::builder()
302                .reference(det_r)
303                .occupation_patterns(occupation_patterns)
304                .build()
305                .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
306
307            // Construct the multi-determinantal wavefunctions
308            let multidets_vec = energies_r
309                .iter()
310                .zip(coefficients_r.columns())
311                .map(|(energy, coeffs)| {
312                    MultiDeterminant::builder()
313                        .basis(fci_basis.clone())
314                        .coefficients(coeffs.to_owned())
315                        .threshold(1e-7)
316                        .energy(Ok(*energy))
317                        .build()
318                })
319                .collect::<Result<Vec<_>, _>>()
320                .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
321            let multidets = MultiDeterminants::from_multideterminant_vec(
322                &multidets_vec.iter().collect::<Vec<_>>()
323            )
324            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
325
326            // Construct the driver
327            match &use_magnetic_group {
328                Some(MagneticSymmetryAnalysisKind::Corepresentation) => {
329                    let mut mda_driver = MultiDeterminantRepAnalysisDriver::<
330                        MagneticRepresentedSymmetryGroup,
331                        f64,
332                        _,
333                        SpinConstraint,
334                    >::builder()
335                    .parameters(&mda_params)
336                    .angular_function_parameters(&afa_params)
337                    .multidets(&multidets)
338                    .sao(&sao_r)
339                    .sao_h(None) // Real SAO.
340                    .symmetry_group(&pd_res)
341                    .build()
342                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
343
344                    // Run the driver
345                    py.detach(|| {
346                        mda_driver
347                            .run()
348                            .map_err(|err| PyRuntimeError::new_err(err.to_string()))
349                    })?
350                }
351                Some(MagneticSymmetryAnalysisKind::Representation) | None => {
352                    let mut mda_driver = MultiDeterminantRepAnalysisDriver::<
353                        UnitaryRepresentedSymmetryGroup,
354                        f64,
355                        _,
356                        SpinConstraint,
357                    >::builder()
358                    .parameters(&mda_params)
359                    .angular_function_parameters(&afa_params)
360                    .multidets(&multidets)
361                    .sao(&sao_r)
362                    .sao_h(None) // Real SAO.
363                    .symmetry_group(&pd_res)
364                    .build()
365                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
366
367                    // Run the driver
368                    py.detach(|| {
369                        mda_driver
370                            .run()
371                            .map_err(|err| PyRuntimeError::new_err(err.to_string()))
372                    })?
373                }
374            }
375        }
376        (_, _, _, _) => {
377            // Complex numeric data type
378
379            // Preparation
380            let sao_c = match sao {
381                PyArray2RC::Real(pysao_r) => pysao_r.to_owned_array().mapv(Complex::from),
382                PyArray2RC::Complex(pysao_c) => pysao_c.to_owned_array(),
383            };
384            let sao_h_c = sao_h.map(|pysao_h| match pysao_h {
385                // sao_spatial_h must have the same reality as sao_spatial.
386                PyArray2RC::Real(pysao_h_r) => pysao_h_r.to_owned_array().mapv(Complex::from),
387                PyArray2RC::Complex(pysao_h_c) => pysao_h_c.to_owned_array(),
388            });
389            let coefficients_c = match coefficients {
390                PyArray2RC::Real(pycoefficients_r) => {
391                    pycoefficients_r.to_owned_array().mapv(Complex::from)
392                }
393                PyArray2RC::Complex(pycoefficients_c) => pycoefficients_c.to_owned_array(),
394            };
395            let energies_c = match energies {
396                PyArray1RC::Real(pyenergies_r) => pyenergies_r.to_owned_array().mapv(Complex::from),
397                PyArray1RC::Complex(pyenergies_c) => pyenergies_c.to_owned_array(),
398            };
399
400            match structure_constraint {
401                PyStructureConstraint::SpinConstraint(_) => {
402                    let det_c = if augment_to_generalised {
403                        match pydet {
404                            PySlaterDeterminant::Real(ref pydet_r) => pydet_r
405                                .to_qsym2::<SpinConstraint>(&baos_ref, mol)
406                                .map(|det_r| {
407                                    SlaterDeterminant::<C128, SpinConstraint>::from(det_r)
408                                        .to_generalised()
409                                }),
410                            PySlaterDeterminant::Complex(ref pydet_c) => pydet_c
411                                .to_qsym2::<SpinConstraint>(&baos_ref, mol)
412                                .map(|det_c| det_c.to_generalised()),
413                        }
414                    } else {
415                        match pydet {
416                            PySlaterDeterminant::Real(ref pydet_r) => pydet_r
417                                .to_qsym2::<SpinConstraint>(&baos_ref, mol)
418                                .map(|det_r| {
419                                    SlaterDeterminant::<C128, SpinConstraint>::from(det_r)
420                                }),
421                            PySlaterDeterminant::Complex(ref pydet_c) => {
422                                pydet_c.to_qsym2::<SpinConstraint>(&baos_ref, mol)
423                            }
424                        }
425                    }
426                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
427
428                    let n_energies = energies_c.len();
429                    if coefficients_c.shape()[1] != n_energies {
430                        return Err(PyRuntimeError::new_err(
431                            "Mismatched number of FCI energies and number of FCI states.",
432                        ));
433                    }
434
435                    // Construct the FCI basis
436                    let fci_basis = FCIBasis::builder()
437                        .reference(det_c)
438                        .occupation_patterns(occupation_patterns)
439                        .build()
440                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
441
442                    // Construct the multi-determinantal wavefunctions
443                    let multidets_vec = energies_c
444                        .iter()
445                        .zip(coefficients_c.columns())
446                        .map(|(energy, coeffs)| {
447                            MultiDeterminant::builder()
448                                .basis(fci_basis.clone())
449                                .coefficients(coeffs.to_owned())
450                                .threshold(1e-7)
451                                .energy(Ok(*energy))
452                                .build()
453                        })
454                        .collect::<Result<Vec<_>, _>>()
455                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
456                    let multidets = MultiDeterminants::from_multideterminant_vec(
457                        &multidets_vec.iter().collect::<Vec<_>>()
458                    )
459                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
460
461                    // Construct the driver
462                    match &use_magnetic_group {
463                        Some(MagneticSymmetryAnalysisKind::Corepresentation) => {
464                            let mut mda_driver = MultiDeterminantRepAnalysisDriver::<
465                                MagneticRepresentedSymmetryGroup,
466                                C128,
467                                _,
468                                SpinConstraint,
469                            >::builder()
470                            .parameters(&mda_params)
471                            .angular_function_parameters(&afa_params)
472                            .multidets(&multidets)
473                            .sao(&sao_c)
474                            .sao_h(sao_h_c.as_ref())
475                            .symmetry_group(&pd_res)
476                            .build()
477                            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
478
479                            // Run the driver
480                            py.detach(|| {
481                                mda_driver
482                                    .run()
483                                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))
484                            })?
485                        }
486                        Some(MagneticSymmetryAnalysisKind::Representation) | None => {
487                            let mut mda_driver = MultiDeterminantRepAnalysisDriver::<
488                                UnitaryRepresentedSymmetryGroup,
489                                C128,
490                                _,
491                                SpinConstraint,
492                            >::builder()
493                            .parameters(&mda_params)
494                            .angular_function_parameters(&afa_params)
495                            .multidets(&multidets)
496                            .sao(&sao_c)
497                            .sao_h(sao_h_c.as_ref())
498                            .symmetry_group(&pd_res)
499                            .build()
500                            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
501
502                            // Run the driver
503                            py.detach(|| {
504                                mda_driver
505                                    .run()
506                                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))
507                            })?
508                        }
509                    }
510                }
511                PyStructureConstraint::SpinOrbitCoupled(_) => {
512                    let det_c = match pydet {
513                        PySlaterDeterminant::Real(ref pydet_r) => pydet_r
514                            .to_qsym2::<SpinOrbitCoupled>(&baos_ref, mol)
515                            .map(SlaterDeterminant::<C128, SpinOrbitCoupled>::from),
516                        PySlaterDeterminant::Complex(ref pydet_c) => {
517                            pydet_c.to_qsym2::<SpinOrbitCoupled>(&baos_ref, mol)
518                        }
519                    }
520                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
521
522                    let n_energies = energies_c.len();
523                    if coefficients_c.shape()[1] != n_energies {
524                        return Err(PyRuntimeError::new_err(
525                            "Mismatched number of FCI energies and number of FCI states.",
526                        ));
527                    }
528
529                    // Construct the FCI basis
530                    let fci_basis = FCIBasis::builder()
531                        .reference(det_c)
532                        .occupation_patterns(occupation_patterns)
533                        .build()
534                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
535
536                    // Construct the multi-determinantal wavefunctions
537                    let multidets_vec = energies_c
538                        .iter()
539                        .zip(coefficients_c.columns())
540                        .map(|(energy, coeffs)| {
541                            MultiDeterminant::builder()
542                                .basis(fci_basis.clone())
543                                .coefficients(coeffs.to_owned())
544                                .threshold(1e-7)
545                                .energy(Ok(*energy))
546                                .build()
547                        })
548                        .collect::<Result<Vec<_>, _>>()
549                        .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
550                    let multidets = MultiDeterminants::from_multideterminant_vec(
551                        &multidets_vec.iter().collect::<Vec<_>>()
552                    )
553                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
554
555                    // Construct the driver
556                    match &use_magnetic_group {
557                        Some(MagneticSymmetryAnalysisKind::Corepresentation) => {
558                            let mut mda_driver = MultiDeterminantRepAnalysisDriver::<
559                                MagneticRepresentedSymmetryGroup,
560                                C128,
561                                _,
562                                SpinOrbitCoupled,
563                            >::builder()
564                            .parameters(&mda_params)
565                            .angular_function_parameters(&afa_params)
566                            .multidets(&multidets)
567                            .sao(&sao_c)
568                            .sao_h(sao_h_c.as_ref())
569                            .symmetry_group(&pd_res)
570                            .build()
571                            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
572
573                            // Run the driver
574                            py.detach(|| {
575                                mda_driver
576                                    .run()
577                                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))
578                            })?
579                        }
580                        Some(MagneticSymmetryAnalysisKind::Representation) | None => {
581                            let mut mda_driver = MultiDeterminantRepAnalysisDriver::<
582                                UnitaryRepresentedSymmetryGroup,
583                                C128,
584                                _,
585                                SpinOrbitCoupled,
586                            >::builder()
587                            .parameters(&mda_params)
588                            .angular_function_parameters(&afa_params)
589                            .multidets(&multidets)
590                            .sao(&sao_c)
591                            .sao_h(sao_h_c.as_ref())
592                            .symmetry_group(&pd_res)
593                            .build()
594                            .map_err(|err| PyRuntimeError::new_err(err.to_string()))?;
595
596                            // Run the driver
597                            py.detach(|| {
598                                mda_driver
599                                    .run()
600                                    .map_err(|err| PyRuntimeError::new_err(err.to_string()))
601                            })?
602                        }
603                    }
604                }
605            }
606        }
607    }
608    Ok(())
609}