qsym2/target/determinant/
mod.rs

1//! Slater determinants.
2
3use std::fmt;
4use std::iter::Sum;
5
6use anyhow::{self, format_err};
7use approx;
8use derive_builder::Builder;
9use itertools::Itertools;
10use log;
11use ndarray::{s, Array1, Array2, Ix2};
12use ndarray_einsum_beta::*;
13use ndarray_linalg::types::Lapack;
14use num::ToPrimitive;
15use num_complex::{Complex, ComplexFloat};
16use num_traits::float::{Float, FloatConst};
17
18use crate::angmom::spinor_rotation_3d::{SpinConstraint, SpinOrbitCoupled, StructureConstraint};
19use crate::auxiliary::molecule::Molecule;
20use crate::basis::ao::BasisAngularOrder;
21use crate::target::density::{DensitiesOwned, Density};
22use crate::target::orbital::MolecularOrbital;
23
24#[cfg(test)]
25mod determinant_tests;
26
27pub mod determinant_analysis;
28mod determinant_transformation;
29
30// ==================
31// Struct definitions
32// ==================
33
34/// Structure to manage single-determinantal wavefunctions.
35#[derive(Builder, Clone)]
36#[builder(build_fn(validate = "Self::validate"))]
37pub struct SlaterDeterminant<'a, T, SC>
38where
39    T: ComplexFloat + Lapack,
40    SC: StructureConstraint + fmt::Display,
41{
42    /// The structure constraint associated with the coefficients describing this determinant.
43    structure_constraint: SC,
44
45    /// The angular order of the basis functions with respect to which the coefficients are
46    /// expressed.
47    bao: &'a BasisAngularOrder<'a>,
48
49    /// A boolean indicating if inner products involving this determinant should be the
50    /// complex-symmetric bilinear form, rather than the conventional Hermitian sesquilinear form.
51    complex_symmetric: bool,
52
53    /// A boolean indicating if the determinant has been acted on by an antiunitary operation. This
54    /// is so that the correct metric can be used during overlap evaluation.
55    #[builder(default = "false")]
56    complex_conjugated: bool,
57
58    /// The associated molecule.
59    mol: &'a Molecule,
60
61    /// The coefficients describing this determinant.
62    #[builder(setter(custom))]
63    coefficients: Vec<Array2<T>>,
64
65    /// The occupation patterns of the molecular orbitals in [`Self::coefficients`].
66    #[builder(setter(custom))]
67    occupations: Vec<Array1<<T as ComplexFloat>::Real>>,
68
69    /// The energies of the molecular orbitals in [`Self::coefficients`].
70    #[builder(default = "None")]
71    mo_energies: Option<Vec<Array1<T>>>,
72
73    /// The energy of this determinant.
74    #[builder(default = "Err(\"Determinant energy not yet set.\".to_string())")]
75    energy: Result<T, String>,
76
77    /// The threshold for comparing determinants.
78    threshold: <T as ComplexFloat>::Real,
79}
80
81impl<'a, T, SC> SlaterDeterminantBuilder<'a, T, SC>
82where
83    T: ComplexFloat + Lapack,
84    SC: StructureConstraint + fmt::Display,
85{
86    pub fn coefficients(&mut self, cs: &[Array2<T>]) -> &mut Self {
87        self.coefficients = Some(cs.to_vec());
88        self
89    }
90
91    pub fn occupations(&mut self, occs: &[Array1<<T as ComplexFloat>::Real>]) -> &mut Self {
92        self.occupations = Some(occs.to_vec());
93        self
94    }
95
96    fn validate(&self) -> Result<(), String> {
97        let bao = self
98            .bao
99            .ok_or("No `BasisAngularOrder` found.".to_string())?;
100        let nbas = bao.n_funcs();
101        let coefficients = self
102            .coefficients
103            .as_ref()
104            .ok_or("No coefficients found.".to_string())?;
105        let structcons = self
106            .structure_constraint
107            .as_ref()
108            .ok_or("No structure constraints found.".to_string())?;
109        let coefficients_length_check = if coefficients.len() != structcons.n_coefficient_matrices()
110        {
111            log::error!(
112                    "Unexpected number of coefficient matrices: {} found, but {} expected for the structure constraint {}.",
113                    coefficients.len(),
114                    structcons.n_coefficient_matrices(),
115                    structcons
116                );
117            false
118        } else {
119            true
120        };
121        let coefficients_shape_check = {
122            let nrows = nbas * structcons.n_explicit_comps_per_coefficient_matrix();
123            if !coefficients.iter().all(|c| c.shape()[0] == nrows) {
124                log::error!(
125                    "Unexpected shapes of coefficient matrices: {} {} expected for all coefficient matrices, but {} found.",
126                    nrows,
127                    if nrows == 1 {"row"} else {"rows"},
128                    coefficients.iter().map(|c| c.shape()[0].to_string()).join(", ")
129                );
130                false
131            } else {
132                true
133            }
134        };
135
136        let occupations = self
137            .occupations
138            .as_ref()
139            .ok_or("No occupations found.".to_string())?;
140        let occupations_length_check = if occupations.len() != structcons.n_coefficient_matrices() {
141            log::error!(
142                "Unexpected number of occupation vectors: {} found, but {} expected for the structure constraint {}.",
143                occupations.len(),
144                structcons.n_coefficient_matrices(),
145                structcons
146            );
147            false
148        } else {
149            true
150        };
151        let occupations_shape_check = if !occupations
152            .iter()
153            .zip(coefficients.iter())
154            .all(|(o, c)| o.len() == c.shape()[1])
155        {
156            log::error!(
157                "Mismatched occupations and numbers of orbitals: {}",
158                occupations
159                    .iter()
160                    .zip(coefficients.iter())
161                    .map(|(o, c)| format!("{} vs. {}", o.len(), c.shape()[1]))
162                    .join(", ")
163            );
164            false
165        } else {
166            true
167        };
168
169        let mol = self.mol.ok_or("No molecule found.".to_string())?;
170        let natoms_check = mol.atoms.len() == bao.n_atoms();
171        if !natoms_check {
172            log::error!("The number of atoms in the molecule does not match the number of local sites in the basis.");
173        }
174
175        if coefficients_length_check
176            && coefficients_shape_check
177            && occupations_length_check
178            && occupations_shape_check
179            && natoms_check
180        {
181            Ok(())
182        } else {
183            Err("Slater determinant validation failed.".to_string())
184        }
185    }
186}
187
188impl<'a, T, SC> SlaterDeterminant<'a, T, SC>
189where
190    T: ComplexFloat + Clone + Lapack,
191    SC: StructureConstraint + Clone + fmt::Display,
192{
193    /// Returns a builder to construct a new [`SlaterDeterminant`].
194    pub fn builder() -> SlaterDeterminantBuilder<'a, T, SC> {
195        SlaterDeterminantBuilder::default()
196    }
197}
198
199impl<'a, T, SC> SlaterDeterminant<'a, T, SC>
200where
201    T: ComplexFloat + Clone + Lapack,
202    SC: StructureConstraint + Clone + fmt::Display,
203{
204    /// Extracts the molecular orbitals in this Slater determinant.
205    ///
206    /// # Returns
207    ///
208    /// A vector of the molecular orbitals constituting this Slater determinant. In the restricted
209    /// spin constraint, the identical molecular orbitals across different spin spaces are only
210    /// given once. Each molecular orbital does contain an index of the spin space it is in.
211    pub fn to_orbitals(&self) -> Vec<Vec<MolecularOrbital<'a, T, SC>>> {
212        self.coefficients
213            .iter()
214            .enumerate()
215            .map(|(spini, cs_spini)| {
216                cs_spini
217                    .columns()
218                    .into_iter()
219                    .enumerate()
220                    .map(move |(i, c)| {
221                        MolecularOrbital::builder()
222                            .coefficients(c.to_owned())
223                            .energy(self.mo_energies.as_ref().map(|moes| moes[spini][i]))
224                            .bao(self.bao)
225                            .mol(self.mol)
226                            .structure_constraint(self.structure_constraint.clone())
227                            .component_index(spini)
228                            .complex_symmetric(self.complex_symmetric)
229                            .threshold(self.threshold)
230                            .build()
231                            .expect("Unable to construct a molecular orbital.")
232                    })
233                    .collect::<Vec<_>>()
234            })
235            .collect::<Vec<_>>()
236    }
237}
238
239impl<'a, T, SC> SlaterDeterminant<'a, T, SC>
240where
241    T: ComplexFloat + Clone + Lapack,
242    SC: StructureConstraint + fmt::Display,
243{
244    /// Returns the complex-symmetric flag of the determinant.
245    pub fn complex_symmetric(&self) -> bool {
246        self.complex_symmetric
247    }
248
249    /// Returns the complex-conjugated flag of the determinant.
250    pub(crate) fn complex_conjugated(&self) -> bool {
251        self.complex_conjugated
252    }
253
254    /// Returns the constraint imposed on the coefficients.
255    pub fn structure_constraint(&self) -> &SC {
256        &self.structure_constraint
257    }
258
259    /// Returns the basis angular order information of the basis set in which the coefficients are
260    /// expressed.
261    pub fn bao(&self) -> &BasisAngularOrder {
262        self.bao
263    }
264
265    /// Returns the molecule associated with this Slater determinant.
266    pub fn mol(&self) -> &Molecule {
267        self.mol
268    }
269
270    /// Returns the determinantal energy.
271    pub fn energy(&self) -> Result<&T, &String> {
272        self.energy.as_ref()
273    }
274
275    /// Returns the molecular-orbital energies.
276    pub fn mo_energies(&self) -> Option<&Vec<Array1<T>>> {
277        self.mo_energies.as_ref()
278    }
279
280    /// Returns the occupation patterns of the molecular orbitals.
281    pub fn occupations(&self) -> &Vec<Array1<<T as ComplexFloat>::Real>> {
282        &self.occupations
283    }
284
285    /// Returns a shared reference to a vector of coefficient arrays.
286    pub fn coefficients(&self) -> &Vec<Array2<T>> {
287        &self.coefficients
288    }
289
290    /// Returns the threshold with which determinants are compared.
291    pub fn threshold(&self) -> <T as ComplexFloat>::Real {
292        self.threshold
293    }
294
295    /// Returns the total number of electrons in the determinant.
296    pub fn nelectrons(&self) -> <T as ComplexFloat>::Real
297    where
298        <T as ComplexFloat>::Real: Sum + From<u16>,
299    {
300        let implicit_factor = self
301            .structure_constraint
302            .implicit_factor()
303            .expect("Unable to retrieve the implicit factor from the structure constraint.");
304        <T as ComplexFloat>::Real::from(
305            implicit_factor
306                .to_u16()
307                .expect("Unable to convert the implicit factor to `u16`."),
308        ) * self
309            .occupations
310            .iter()
311            .map(|occ| occ.iter().copied().sum())
312            .sum()
313    }
314}
315
316impl<'a, T> SlaterDeterminant<'a, T, SpinConstraint>
317where
318    T: ComplexFloat + Clone + Lapack,
319{
320    /// Augments the encoding of coefficients in this Slater determinant to that in the
321    /// corresponding generalised spin constraint.
322    ///
323    /// # Returns
324    ///
325    /// The equivalent Slater determinant with the coefficients encoded in the generalised spin
326    /// constraint.
327    pub fn to_generalised(&self) -> Self {
328        match self.structure_constraint {
329            SpinConstraint::Restricted(n) => {
330                log::debug!(
331                    "Restricted Slater determinant will be augmented to generalised Slater determinant."
332                );
333                let nbas = self.bao.n_funcs();
334
335                let cr = &self.coefficients[0];
336                let occr = &self.occupations[0];
337                let norb = cr.ncols();
338                let mut cg = Array2::<T>::zeros((nbas * usize::from(n), norb * usize::from(n)));
339                let mut occg = Array1::<<T as ComplexFloat>::Real>::zeros((norb * usize::from(n),));
340                (0..usize::from(n)).for_each(|i| {
341                    let row_start = nbas * i;
342                    let row_end = nbas * (i + 1);
343                    let col_start = norb * i;
344                    let col_end = norb * (i + 1);
345                    cg.slice_mut(s![row_start..row_end, col_start..col_end])
346                        .assign(cr);
347                    occg.slice_mut(s![col_start..col_end]).assign(occr);
348                });
349                let moeg_opt = self.mo_energies.as_ref().map(|moer| {
350                    let mut moeg = Array1::<T>::zeros((norb * usize::from(n),));
351                    (0..usize::from(n)).for_each(|i| {
352                        let col_start = norb * i;
353                        let col_end = norb * (i + 1);
354                        moeg.slice_mut(s![col_start..col_end]).assign(&moer[0]);
355                    });
356                    vec![moeg]
357                });
358                Self::builder()
359                    .coefficients(&[cg])
360                    .occupations(&[occg])
361                    .mo_energies(moeg_opt)
362                    .energy(self.energy.clone())
363                    .bao(self.bao)
364                    .mol(self.mol)
365                    .structure_constraint(SpinConstraint::Generalised(n, false))
366                    .complex_symmetric(self.complex_symmetric)
367                    .threshold(self.threshold)
368                    .build()
369                    .expect("Unable to spin-generalise a `SlaterDeterminant`.")
370            }
371            SpinConstraint::Unrestricted(n, increasingm) => {
372                log::debug!(
373                    "Unrestricted Slater determinant will be augmented to generalised Slater determinant."
374                );
375                let nbas = self.bao.n_funcs();
376                let norb_tot = self.coefficients.iter().map(|c| c.ncols()).sum();
377                let mut cg = Array2::<T>::zeros((nbas * usize::from(n), norb_tot));
378                let mut occg = Array1::<<T as ComplexFloat>::Real>::zeros((norb_tot,));
379
380                let col_boundary_indices = (0..usize::from(n))
381                    .scan(0, |acc, ispin| {
382                        let start_index = *acc;
383                        *acc += self.coefficients[ispin].shape()[1];
384                        Some((start_index, *acc))
385                    })
386                    .collect::<Vec<_>>();
387                (0..usize::from(n)).for_each(|i| {
388                    let row_start = nbas * i;
389                    let row_end = nbas * (i + 1);
390                    let (col_start, col_end) = col_boundary_indices[i];
391                    cg.slice_mut(s![row_start..row_end, col_start..col_end])
392                        .assign(&self.coefficients[i]);
393                    occg.slice_mut(s![col_start..col_end])
394                        .assign(&self.occupations[i]);
395                });
396
397                let moeg_opt = self.mo_energies.as_ref().map(|moer| {
398                    let mut moeg = Array1::<T>::zeros((norb_tot,));
399                    (0..usize::from(n)).for_each(|i| {
400                        let (col_start, col_end) = col_boundary_indices[i];
401                        moeg.slice_mut(s![col_start..col_end]).assign(&moer[i]);
402                    });
403                    vec![moeg]
404                });
405
406                Self::builder()
407                    .coefficients(&[cg])
408                    .occupations(&[occg])
409                    .mo_energies(moeg_opt)
410                    .energy(self.energy.clone())
411                    .bao(self.bao)
412                    .mol(self.mol)
413                    .structure_constraint(SpinConstraint::Generalised(n, increasingm))
414                    .complex_symmetric(self.complex_symmetric)
415                    .threshold(self.threshold)
416                    .build()
417                    .expect("Unable to spin-generalise a `SlaterDeterminant`.")
418            }
419            SpinConstraint::Generalised(_, _) => self.clone(),
420        }
421    }
422}
423
424impl<'a> SlaterDeterminant<'a, f64, SpinConstraint> {
425    /// Constructs a vector of real densities, one for each spin space in a Slater determinant.
426    ///
427    /// For restricted and unrestricted spin constraints, spin spaces are well-defined. For
428    /// generalised spin constraints, each spin-space density is constructed from the corresponding
429    /// diagonal block of the overall density matrix.
430    ///
431    /// Occupation numbers are also incorporated in the formation of density matrices.
432    ///
433    /// # Returns
434    ///
435    /// A vector of real densities, one for each spin space.
436    pub fn to_densities(
437        &'a self,
438    ) -> Result<DensitiesOwned<'a, f64, SpinConstraint>, anyhow::Error> {
439        let densities = match self.structure_constraint {
440            SpinConstraint::Restricted(_) | SpinConstraint::Unrestricted(_, _) => {
441                self.coefficients().iter().zip(self.occupations().iter()).map(|(c, o)| {
442                    let denmat = einsum(
443                        "i,mi,ni->mn",
444                        &[&o.view(), &c.view(), &c.view()]
445                    )
446                    .expect("Unable to construct a density matrix from a determinant coefficient matrix.")
447                    .into_dimensionality::<Ix2>()
448                    .expect("Unable to convert the resultant density matrix to two dimensions.");
449                    Density::<f64>::builder()
450                        .density_matrix(denmat)
451                        .bao(self.bao())
452                        .mol(self.mol())
453                        .complex_symmetric(self.complex_symmetric())
454                        .threshold(self.threshold())
455                        .build()
456                    })
457                    .collect::<Result<Vec<_>, _>>()?
458            }
459            SpinConstraint::Generalised(nspins, _) => {
460                let denmat = einsum(
461                    "i,mi,ni->mn",
462                    &[&self.occupations[0].view(), &self.coefficients[0].view(), &self.coefficients[0].view()]
463                )
464                .expect("Unable to construct a density matrix from a determinant coefficient matrix.")
465                .into_dimensionality::<Ix2>()
466                .expect("Unable to convert the resultant density matrix to two dimensions.");
467                let nspatial = self.bao.n_funcs();
468                (0..usize::from(nspins)).map(|ispin| {
469                    let ispin_denmat = denmat.slice(
470                        s![ispin*nspatial..(ispin + 1)*nspatial, ispin*nspatial..(ispin + 1)*nspatial]
471                    ).to_owned();
472                    Density::<f64>::builder()
473                        .density_matrix(ispin_denmat)
474                        .bao(self.bao())
475                        .mol(self.mol())
476                        .complex_symmetric(self.complex_symmetric())
477                        .threshold(self.threshold())
478                        .build()
479                }).collect::<Result<Vec<_>, _>>()?
480            }
481        };
482        DensitiesOwned::builder()
483            .structure_constraint(self.structure_constraint.clone())
484            .densities(densities)
485            .build()
486            .map_err(|err| format_err!(err))
487    }
488}
489
490impl<'a, T> SlaterDeterminant<'a, Complex<T>, SpinConstraint>
491where
492    T: Float + FloatConst + Lapack,
493    Complex<T>: Lapack,
494{
495    /// Constructs a vector of complex densities, one from each coefficient matrix in a Slater
496    /// determinant.
497    ///
498    /// For restricted and unrestricted spin constraints, spin spaces are well-defined. For
499    /// generalised spin constraints, each spin-space density is constructed from the corresponding
500    /// diagonal block of the overall density matrix.
501    ///
502    /// Occupation numbers are also incorporated in the formation of density matrices.
503    ///
504    /// # Arguments
505    ///
506    /// * `sd` - A Slater determinant.
507    ///
508    /// # Returns
509    ///
510    /// A vector of complex densities, one for each spin space.
511    pub fn to_densities(
512        &'a self,
513    ) -> Result<DensitiesOwned<'a, Complex<T>, SpinConstraint>, anyhow::Error> {
514        let densities = match self.structure_constraint {
515            SpinConstraint::Restricted(_) | SpinConstraint::Unrestricted(_, _) => {
516                self.coefficients().iter().zip(self.occupations().iter()).map(|(c, o)| {
517                    let denmat = einsum(
518                        "i,mi,ni->mn",
519                        &[&o.map(Complex::<T>::from).view(), &c.view(), &c.map(Complex::conj).view()]
520                    )
521                    .expect("Unable to construct a density matrix from a determinant coefficient matrix.")
522                    .into_dimensionality::<Ix2>()
523                    .expect("Unable to convert the resultant density matrix to two dimensions.");
524                    Density::<Complex<T>>::builder()
525                        .density_matrix(denmat)
526                        .bao(self.bao())
527                        .mol(self.mol())
528                        .complex_symmetric(self.complex_symmetric())
529                        .threshold(self.threshold())
530                        .build()
531                    })
532                    .collect::<Result<Vec<_>, _>>()?
533            }
534            SpinConstraint::Generalised(nspins, _) => {
535                let denmat = einsum(
536                    "i,mi,ni->mn",
537                    &[
538                        &self.occupations[0].map(Complex::<T>::from).view(),
539                        &self.coefficients[0].view(),
540                        &self.coefficients[0].map(Complex::conj).view()
541                    ]
542                )
543                .expect("Unable to construct a density matrix from a determinant coefficient matrix.")
544                .into_dimensionality::<Ix2>()
545                .expect("Unable to convert the resultant density matrix to two dimensions.");
546                let nspatial = self.bao.n_funcs();
547                (0..usize::from(nspins)).map(|ispin| {
548                    let ispin_denmat = denmat.slice(
549                        s![ispin*nspatial..(ispin + 1)*nspatial, ispin*nspatial..(ispin + 1)*nspatial]
550                    ).to_owned();
551                    Density::<Complex<T>>::builder()
552                        .density_matrix(ispin_denmat)
553                        .bao(self.bao())
554                        .mol(self.mol())
555                        .complex_symmetric(self.complex_symmetric())
556                        .threshold(self.threshold())
557                        .build()
558                }).collect::<Result<Vec<_>, _>>()?
559            }
560        };
561        DensitiesOwned::builder()
562            .structure_constraint(self.structure_constraint.clone())
563            .densities(densities)
564            .build()
565            .map_err(|err| format_err!(err))
566    }
567}
568
569impl<'a, T> SlaterDeterminant<'a, Complex<T>, SpinOrbitCoupled>
570where
571    T: Float + FloatConst + Lapack,
572    Complex<T>: Lapack,
573{
574    /// Constructs a vector of complex densities, one for each component in the multi-component
575    /// j-adapted basis.
576    ///
577    /// Occupation numbers are also incorporated in the formation of density matrices.
578    ///
579    /// # Arguments
580    ///
581    /// * `sd` - A Slater determinant.
582    ///
583    /// # Returns
584    ///
585    /// A vector of complex densities.
586    pub fn to_densities(
587        &'a self,
588    ) -> Result<DensitiesOwned<'a, Complex<T>, SpinOrbitCoupled>, anyhow::Error> {
589        let densities = match self.structure_constraint {
590            SpinOrbitCoupled::JAdapted(ncomps) => {
591                let denmat = einsum(
592                    "i,mi,ni->mn",
593                    &[
594                        &self.occupations[0].map(Complex::<T>::from).view(),
595                        &self.coefficients[0].view(),
596                        &self.coefficients[0].map(Complex::conj).view()
597                    ]
598                )
599                .expect("Unable to construct a density matrix from a determinant coefficient matrix.")
600                .into_dimensionality::<Ix2>()
601                .expect("Unable to convert the resultant density matrix to two dimensions.");
602                let nspatial = self.bao.n_funcs();
603                (0..usize::from(ncomps)).map(|icomp| {
604                    let icomp_denmat = denmat.slice(
605                        s![icomp*nspatial..(icomp + 1)*nspatial, icomp*nspatial..(icomp + 1)*nspatial]
606                    ).to_owned();
607                    Density::<Complex<T>>::builder()
608                        .density_matrix(icomp_denmat)
609                        .bao(self.bao())
610                        .mol(self.mol())
611                        .complex_symmetric(self.complex_symmetric())
612                        .threshold(self.threshold())
613                        .build()
614                }).collect::<Result<Vec<_>, _>>()?
615            }
616        };
617        DensitiesOwned::builder()
618            .structure_constraint(self.structure_constraint.clone())
619            .densities(densities)
620            .build()
621            .map_err(|err| format_err!(err))
622    }
623}
624
625// =====================
626// Trait implementations
627// =====================
628
629// ----
630// From
631// ----
632impl<'a, T, SC> From<SlaterDeterminant<'a, T, SC>> for SlaterDeterminant<'a, Complex<T>, SC>
633where
634    T: Float + FloatConst + Lapack,
635    Complex<T>: Lapack,
636    SC: StructureConstraint + Clone + fmt::Display,
637{
638    fn from(value: SlaterDeterminant<'a, T, SC>) -> Self {
639        SlaterDeterminant::<'a, Complex<T>, SC>::builder()
640            .coefficients(
641                &value
642                    .coefficients
643                    .into_iter()
644                    .map(|coeffs| coeffs.map(Complex::from))
645                    .collect::<Vec<_>>(),
646            )
647            .occupations(&value.occupations)
648            .mo_energies(value.mo_energies.map(|moes| {
649                moes.iter()
650                    .map(|moe| moe.map(Complex::from))
651                    .collect::<Vec<_>>()
652            }))
653            .bao(value.bao)
654            .mol(value.mol)
655            .structure_constraint(value.structure_constraint)
656            .complex_symmetric(value.complex_symmetric)
657            .threshold(value.threshold)
658            .build()
659            .expect("Unable to complexify a `SlaterDeterminant`.")
660    }
661}
662
663// ---------
664// PartialEq
665// ---------
666impl<'a, T, SC> PartialEq for SlaterDeterminant<'a, T, SC>
667where
668    T: ComplexFloat<Real = f64> + Lapack,
669    SC: StructureConstraint + PartialEq + fmt::Display,
670{
671    fn eq(&self, other: &Self) -> bool {
672        let thresh = (self.threshold * other.threshold).sqrt();
673        let coefficients_eq =
674            self.coefficients.len() == other.coefficients.len()
675                && self.coefficients.iter().zip(other.coefficients.iter()).all(
676                    |(scoeffs, ocoeffs)| {
677                        approx::relative_eq!(
678                            (scoeffs - ocoeffs)
679                                .map(|x| ComplexFloat::abs(*x).powi(2))
680                                .sum()
681                                .sqrt(),
682                            0.0,
683                            epsilon = thresh,
684                            max_relative = thresh,
685                        )
686                    },
687                );
688        let occs_eq = self.occupations.len() == other.occupations.len()
689            && self
690                .occupations
691                .iter()
692                .zip(other.occupations.iter())
693                .all(|(soccs, ooccs)| {
694                    approx::relative_eq!(
695                        (soccs - ooccs).map(|x| x.abs().powi(2)).sum().sqrt(),
696                        0.0,
697                        epsilon = thresh,
698                        max_relative = thresh,
699                    )
700                });
701        self.structure_constraint == other.structure_constraint
702            && self.bao == other.bao
703            && self.mol == other.mol
704            && coefficients_eq
705            && occs_eq
706    }
707}
708
709// --
710// Eq
711// --
712impl<'a, T, SC> Eq for SlaterDeterminant<'a, T, SC>
713where
714    T: ComplexFloat<Real = f64> + Lapack,
715    SC: StructureConstraint + Eq + fmt::Display,
716{
717}
718
719// -----
720// Debug
721// -----
722impl<'a, T, SC> fmt::Debug for SlaterDeterminant<'a, T, SC>
723where
724    T: fmt::Debug + ComplexFloat + Lapack,
725    <T as ComplexFloat>::Real: Sum + From<u16> + fmt::Debug,
726    SC: StructureConstraint + fmt::Debug + fmt::Display,
727{
728    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
729        write!(
730            f,
731            "SlaterDeterminant[{:?}: {:?} electrons, {} coefficient {} of dimensions {}]",
732            self.structure_constraint,
733            self.nelectrons(),
734            self.coefficients.len(),
735            if self.coefficients.len() != 1 {
736                "matrices"
737            } else {
738                "matrix"
739            },
740            self.coefficients
741                .iter()
742                .map(|coeff| coeff
743                    .shape()
744                    .iter()
745                    .map(|x| x.to_string())
746                    .collect::<Vec<_>>()
747                    .join("×"))
748                .collect::<Vec<_>>()
749                .join(", ")
750        )?;
751        Ok(())
752    }
753}
754
755// -------
756// Display
757// -------
758impl<'a, T, SC> fmt::Display for SlaterDeterminant<'a, T, SC>
759where
760    T: fmt::Display + ComplexFloat + Lapack,
761    <T as ComplexFloat>::Real: Sum + From<u16> + fmt::Display,
762    SC: StructureConstraint + fmt::Display,
763{
764    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
765        write!(
766            f,
767            "SlaterDeterminant[{}: {} electrons, {} coefficient {} of dimensions {}]",
768            self.structure_constraint,
769            self.nelectrons(),
770            self.coefficients.len(),
771            if self.coefficients.len() != 1 {
772                "matrices"
773            } else {
774                "matrix"
775            },
776            self.coefficients
777                .iter()
778                .map(|coeff| coeff
779                    .shape()
780                    .iter()
781                    .map(|x| x.to_string())
782                    .collect::<Vec<_>>()
783                    .join("×"))
784                .collect::<Vec<_>>()
785                .join(", ")
786        )?;
787        Ok(())
788    }
789}