1use std::collections::HashSet;
4use std::fmt;
5use std::hash::Hash;
6use std::ops::Mul;
7
8use anyhow::{self, bail, format_err};
9use derive_builder::Builder;
10use duplicate::duplicate_item;
11use ndarray::{s, Array2};
12use ndarray_linalg::types::Lapack;
13use num_complex::{Complex, ComplexFloat};
14use num_traits::Float;
15use serde::{Deserialize, Serialize};
16
17use crate::analysis::{
18 log_overlap_eigenvalues, EigenvalueComparisonMode, Overlap, ProjectionDecomposition,
19 RepAnalysis,
20};
21use crate::angmom::spinor_rotation_3d::{SpinConstraint, StructureConstraint};
22use crate::chartab::chartab_group::CharacterProperties;
23use crate::chartab::SubspaceDecomposable;
24use crate::drivers::representation_analysis::angular_function::{
25 find_angular_function_representation, find_spinor_function_representation,
26 AngularFunctionRepAnalysisParams,
27};
28use crate::drivers::representation_analysis::{
29 fn_construct_magnetic_group, fn_construct_unitary_group, log_bao, log_cc_transversal,
30 CharacterTableDisplay, MagneticSymmetryAnalysisKind,
31};
32use crate::drivers::symmetry_group_detection::SymmetryGroupDetectionResult;
33use crate::drivers::QSym2Driver;
34use crate::group::{GroupProperties, MagneticRepresentedGroup, UnitaryRepresentedGroup};
35use crate::io::format::{
36 log_micsec_begin, log_micsec_end, log_subtitle, nice_bool, qsym2_output, write_subtitle,
37 write_title, QSym2Output,
38};
39use crate::symmetry::symmetry_group::{
40 MagneticRepresentedSymmetryGroup, SymmetryGroupProperties, UnitaryRepresentedSymmetryGroup,
41};
42use crate::symmetry::symmetry_transformation::SymmetryTransformationKind;
43use crate::target::determinant::SlaterDeterminant;
44use crate::target::noci::basis::{Basis, EagerBasis, OrbitBasis};
45use crate::target::noci::multideterminant::multideterminant_analysis::MultiDeterminantSymmetryOrbit;
46use crate::target::noci::multideterminant::MultiDeterminant;
47
48#[cfg(test)]
49#[path = "multideterminant_tests.rs"]
50mod multideterminant_tests;
51
52const fn default_true() -> bool {
61 true
62}
63const fn default_symbolic() -> Option<CharacterTableDisplay> {
64 Some(CharacterTableDisplay::Symbolic)
65}
66
67#[derive(Clone, Builder, Debug, Serialize, Deserialize)]
70pub struct MultiDeterminantRepAnalysisParams<T: From<f64>> {
71 pub integrality_threshold: T,
73
74 pub linear_independence_threshold: T,
76
77 #[builder(default = "None")]
80 #[serde(default)]
81 pub use_magnetic_group: Option<MagneticSymmetryAnalysisKind>,
82
83 #[builder(default = "false")]
85 #[serde(default)]
86 pub use_double_group: bool,
87
88 #[builder(default = "true")]
91 #[serde(default = "default_true")]
92 pub use_cayley_table: bool,
93
94 #[builder(default = "SymmetryTransformationKind::Spatial")]
97 #[serde(default)]
98 pub symmetry_transformation_kind: SymmetryTransformationKind,
99
100 #[builder(default = "Some(CharacterTableDisplay::Symbolic)")]
103 #[serde(default = "default_symbolic")]
104 pub write_character_table: Option<CharacterTableDisplay>,
105
106 #[builder(default = "true")]
108 #[serde(default = "default_true")]
109 pub write_overlap_eigenvalues: bool,
110
111 #[builder(default = "EigenvalueComparisonMode::Modulus")]
113 #[serde(default)]
114 pub eigenvalue_comparison_mode: EigenvalueComparisonMode,
115
116 #[builder(default = "None")]
119 #[serde(default)]
120 pub infinite_order_to_finite: Option<u32>,
121}
122
123impl<T> MultiDeterminantRepAnalysisParams<T>
124where
125 T: Float + From<f64>,
126{
127 pub fn builder() -> MultiDeterminantRepAnalysisParamsBuilder<T> {
129 MultiDeterminantRepAnalysisParamsBuilder::default()
130 }
131}
132
133impl Default for MultiDeterminantRepAnalysisParams<f64> {
134 fn default() -> Self {
135 Self::builder()
136 .integrality_threshold(1e-7)
137 .linear_independence_threshold(1e-7)
138 .build()
139 .expect("Unable to construct a default `MultiDeterminantRepAnalysisParams<f64>`.")
140 }
141}
142
143impl<T> fmt::Display for MultiDeterminantRepAnalysisParams<T>
144where
145 T: From<f64> + fmt::LowerExp + fmt::Debug,
146{
147 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
148 writeln!(
149 f,
150 "Integrality threshold: {:.3e}",
151 self.integrality_threshold
152 )?;
153 writeln!(
154 f,
155 "Linear independence threshold: {:.3e}",
156 self.linear_independence_threshold
157 )?;
158 writeln!(
159 f,
160 "Orbit eigenvalue comparison mode: {}",
161 self.eigenvalue_comparison_mode
162 )?;
163 writeln!(
164 f,
165 "Write overlap eigenvalues: {}",
166 nice_bool(self.write_overlap_eigenvalues)
167 )?;
168 writeln!(f)?;
169 writeln!(
170 f,
171 "Use magnetic group for analysis: {}",
172 match self.use_magnetic_group {
173 None => "no",
174 Some(MagneticSymmetryAnalysisKind::Representation) =>
175 "yes, using unitary representations",
176 Some(MagneticSymmetryAnalysisKind::Corepresentation) =>
177 "yes, using magnetic corepresentations",
178 }
179 )?;
180 writeln!(
181 f,
182 "Use double group for analysis: {}",
183 nice_bool(self.use_double_group)
184 )?;
185 writeln!(
186 f,
187 "Use Cayley table for orbit overlap matrices: {}",
188 nice_bool(self.use_cayley_table)
189 )?;
190 if let Some(finite_order) = self.infinite_order_to_finite {
191 writeln!(f, "Infinite order to finite: {finite_order}")?;
192 }
193 writeln!(
194 f,
195 "Symmetry transformation kind: {}",
196 self.symmetry_transformation_kind
197 )?;
198 writeln!(f)?;
199 writeln!(
200 f,
201 "Write character table: {}",
202 if let Some(chartab_display) = self.write_character_table.as_ref() {
203 format!("yes, {}", chartab_display.to_string().to_lowercase())
204 } else {
205 "no".to_string()
206 }
207 )?;
208
209 Ok(())
210 }
211}
212
213#[derive(Clone, Builder)]
219pub struct MultiDeterminantRepAnalysisResult<'a, G, T, B, SC>
220where
221 G: SymmetryGroupProperties + Clone,
222 G::CharTab: SubspaceDecomposable<T>,
223 T: ComplexFloat + Lapack,
224 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug,
225 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
226 SC: StructureConstraint + Hash + Eq + fmt::Display,
227{
228 parameters: &'a MultiDeterminantRepAnalysisParams<<T as ComplexFloat>::Real>,
231
232 multidets: Vec<&'a MultiDeterminant<'a, T, B, SC>>,
234
235 group: G,
237
238 multidet_symmetries:
240 Vec<Result<<G::CharTab as SubspaceDecomposable<T>>::Decomposition, String>>,
241
242 multidet_symmetries_thresholds: Vec<(Option<T>, Option<T>)>,
245}
246
247impl<'a, G, T, B, SC> MultiDeterminantRepAnalysisResult<'a, G, T, B, SC>
248where
249 G: SymmetryGroupProperties + Clone,
250 G::CharTab: SubspaceDecomposable<T>,
251 T: ComplexFloat + Lapack,
252 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug,
253 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
254 SC: StructureConstraint + Clone + Hash + Eq + fmt::Display,
255{
256 fn builder() -> MultiDeterminantRepAnalysisResultBuilder<'a, G, T, B, SC> {
259 MultiDeterminantRepAnalysisResultBuilder::default()
260 }
261
262 pub fn multidet_symmetries(
264 &self,
265 ) -> &Vec<Result<<G::CharTab as SubspaceDecomposable<T>>::Decomposition, String>> {
266 &self.multidet_symmetries
267 }
268}
269
270impl<'a, G, T, B, SC> fmt::Display for MultiDeterminantRepAnalysisResult<'a, G, T, B, SC>
271where
272 G: SymmetryGroupProperties + Clone,
273 G::CharTab: SubspaceDecomposable<T>,
274 T: ComplexFloat + Lapack,
275 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug + fmt::Display,
276 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
277 SC: StructureConstraint + Clone + Hash + Eq + fmt::Display,
278{
279 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
280 write_subtitle(f, "Orbit-based symmetry analysis results")?;
281 writeln!(f)?;
282 writeln!(
283 f,
284 "> Group: {} ({})",
285 self.group
286 .finite_subgroup_name()
287 .map(|subgroup_name| format!("{} > {}", self.group.name(), subgroup_name))
288 .unwrap_or(self.group.name()),
289 self.group.group_type().to_string().to_lowercase()
290 )?;
291 writeln!(f)?;
292
293 let multidet_index_length = usize::try_from(self.multidets.len().ilog10() + 1).unwrap_or(4);
294 let multidet_symmetry_length = self
295 .multidet_symmetries
296 .iter()
297 .map(|multidet_sym| {
298 multidet_sym
299 .as_ref()
300 .map(|sym| sym.to_string())
301 .unwrap_or_else(|err| err.clone())
302 .chars()
303 .count()
304 })
305 .max()
306 .unwrap_or(0)
307 .max(8);
308 let multidet_energy_length = self
309 .multidets
310 .iter()
311 .map(|multidet| {
312 multidet
313 .energy()
314 .map(|v| format!("{v:+.7}").chars().count())
315 .unwrap_or(2)
316 })
317 .max()
318 .unwrap_or(6)
319 .max(6);
320
321 let multidet_eig_above_length: usize = self
322 .multidet_symmetries_thresholds
323 .iter()
324 .map(|(above, _)| {
325 above
326 .as_ref()
327 .map(|eig| format!("{eig:+.3e}"))
328 .unwrap_or("--".to_string())
329 .chars()
330 .count()
331 })
332 .max()
333 .unwrap_or(10)
334 .max(10);
335
336 let multidet_eig_below_length: usize = self
337 .multidet_symmetries_thresholds
338 .iter()
339 .map(|(_, below)| {
340 below
341 .as_ref()
342 .map(|eig| format!("{eig:+.3e}"))
343 .unwrap_or("--".to_string())
344 .chars()
345 .count()
346 })
347 .max()
348 .unwrap_or(10)
349 .max(10);
350
351 let table_width = 10
352 + multidet_index_length
353 + multidet_energy_length
354 + multidet_symmetry_length
355 + multidet_eig_above_length
356 + multidet_eig_below_length;
357
358 writeln!(f, "> Multi-determinantal results")?;
359 writeln!(
360 f,
361 " Structure constraint: {}",
362 self.multidets
363 .get(0)
364 .map(|multidet_0| multidet_0.structure_constraint().to_string().to_lowercase())
365 .unwrap_or("--".to_string())
366 )?;
367 writeln!(f, "{}", "┈".repeat(table_width))?;
368 writeln!(
369 f,
370 " {:>multidet_index_length$} {:<multidet_energy_length$} {:<multidet_symmetry_length$} {:<multidet_eig_above_length$} Eig. below",
371 "#",
372 "Energy",
373 "Symmetry",
374 "Eig. above",
375 )?;
376 writeln!(f, "{}", "┈".repeat(table_width))?;
377
378 for (multidet_i, multidet) in self.multidets.iter().enumerate() {
379 let multidet_energy_str = multidet
380 .energy()
381 .map(|multidet_energy| format!("{multidet_energy:>+multidet_energy_length$.7}"))
382 .unwrap_or("--".to_string());
383 let multidet_sym_str = self
384 .multidet_symmetries
385 .get(multidet_i)
386 .ok_or_else(|| format!("Unable to retrieve the symmetry of multideterminantal wavefunction index `{multidet_i}`."))
387 .and_then(|sym_res| sym_res.as_ref().map(|sym| sym.to_string()).map_err(|err| err.to_string()))
388 .unwrap_or_else(|err| err);
389
390 let (eig_above_str, eig_below_str) = self
391 .multidet_symmetries_thresholds
392 .get(multidet_i)
393 .map(|(eig_above_opt, eig_below_opt)| {
394 (
395 eig_above_opt
396 .map(|eig_above| format!("{eig_above:>+.3e}"))
397 .unwrap_or("--".to_string()),
398 eig_below_opt
399 .map(|eig_below| format!("{eig_below:>+.3e}"))
400 .unwrap_or("--".to_string()),
401 )
402 })
403 .unwrap_or(("--".to_string(), "--".to_string()));
404 writeln!(
405 f,
406 " {multidet_i:>multidet_index_length$} \
407 {multidet_energy_str:<multidet_energy_length$} \
408 {multidet_sym_str:<multidet_symmetry_length$} \
409 {eig_above_str:<multidet_eig_above_length$} \
410 {eig_below_str}"
411 )?;
412 }
413
414 writeln!(f, "{}", "┈".repeat(table_width))?;
415 writeln!(f)?;
416
417 Ok(())
418 }
419}
420
421impl<'a, G, T, B, SC> fmt::Debug for MultiDeterminantRepAnalysisResult<'a, G, T, B, SC>
422where
423 G: SymmetryGroupProperties + Clone,
424 G::CharTab: SubspaceDecomposable<T>,
425 T: ComplexFloat + Lapack,
426 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug + fmt::Display,
427 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
428 SC: StructureConstraint + Clone + Hash + Eq + fmt::Display,
429{
430 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
431 writeln!(f, "{self}")
432 }
433}
434
435#[derive(Clone, Builder)]
445#[builder(build_fn(validate = "Self::validate"))]
446pub struct MultiDeterminantRepAnalysisDriver<'a, G, T, B, SC>
447where
448 G: SymmetryGroupProperties + Clone,
449 G::CharTab: SubspaceDecomposable<T>,
450 T: ComplexFloat + Lapack,
451 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug,
452 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
453 SC: StructureConstraint + Hash + Eq + fmt::Display,
454{
455 parameters: &'a MultiDeterminantRepAnalysisParams<<T as ComplexFloat>::Real>,
457
458 multidets: Vec<&'a MultiDeterminant<'a, T, B, SC>>,
460
461 symmetry_group: &'a SymmetryGroupDetectionResult,
464
465 sao: &'a Array2<T>,
470
471 #[builder(default = "None")]
478 sao_h: Option<&'a Array2<T>>,
479
480 angular_function_parameters: &'a AngularFunctionRepAnalysisParams,
482
483 #[builder(setter(skip), default = "None")]
485 result: Option<MultiDeterminantRepAnalysisResult<'a, G, T, B, SC>>,
486}
487
488impl<'a, G, T, B, SC> MultiDeterminantRepAnalysisDriverBuilder<'a, G, T, B, SC>
489where
490 G: SymmetryGroupProperties + Clone,
491 G::CharTab: SubspaceDecomposable<T>,
492 T: ComplexFloat + Lapack,
493 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug,
494 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
495 SC: StructureConstraint + Hash + Eq + fmt::Display,
496{
497 fn validate(&self) -> Result<(), String> {
498 let params = self.parameters.ok_or(
499 "No multi-determinantal wavefunction representation analysis parameters found."
500 .to_string(),
501 )?;
502
503 let sym_res = self
504 .symmetry_group
505 .ok_or("No symmetry group information found.".to_string())?;
506
507 let sao = self.sao.ok_or("No SAO matrix found.".to_string())?;
508
509 if let Some(sao_h) = self.sao_h.flatten() {
510 if sao_h.shape() != sao.shape() {
511 return Err("Mismatched shapes between `sao` and `sao_h`.".to_string());
512 }
513 }
514
515 let multidets = self
516 .multidets
517 .as_ref()
518 .ok_or("No multi-determinantal wavefunctions found.".to_string())?;
519 let mut n_spatial_set = multidets
520 .iter()
521 .flat_map(|multidet| {
522 multidet.basis().iter().map(|det_res| {
523 det_res.map(|det| {
524 (
525 det.bao().n_funcs(),
526 det.structure_constraint()
527 .n_explicit_comps_per_coefficient_matrix(),
528 )
529 })
530 })
531 })
532 .collect::<Result<HashSet<(usize, usize)>, _>>()
533 .map_err(|err| err.to_string())?;
534 let (n_spatial, n_comps) = if n_spatial_set.len() == 1 {
535 n_spatial_set
536 .drain()
537 .next()
538 .ok_or("Unable to retrieve the number of spatial AO basis functions and the number of explicit components.".to_string())
539 } else {
540 Err("Inconsistent numbers of spatial AO basis functions and/or explicit components across multi-determinantal wavefunctions.".to_string())
541 }?;
542
543 let sym = if params.use_magnetic_group.is_some() {
544 sym_res
545 .magnetic_symmetry
546 .as_ref()
547 .ok_or("Magnetic symmetry requested for representation analysis, but no magnetic symmetry found.")?
548 } else {
549 &sym_res.unitary_symmetry
550 };
551
552 if sym.is_infinite() && params.infinite_order_to_finite.is_none() {
553 Err(
554 format!(
555 "Representation analysis cannot be performed using the entirety of the infinite group `{}`. \
556 Consider setting the parameter `infinite_order_to_finite` to restrict to a finite subgroup instead.",
557 sym.group_name.as_ref().expect("No symmetry group name found.")
558 )
559 )
560 } else if (n_spatial != sao.nrows() || n_spatial != sao.ncols())
561 || (n_spatial * n_comps != sao.nrows() || n_spatial * n_comps != sao.ncols())
562 {
563 Err("The dimensions of the SAO matrix do not match either the number of spatial AO basis functions or the number of spatial AO basis functions multiplied by the number of explicit components per coefficient matrix.".to_string())
564 } else {
565 Ok(())
566 }
567 }
568}
569
570impl<'a, G, T, B, SC> MultiDeterminantRepAnalysisDriver<'a, G, T, B, SC>
578where
579 G: SymmetryGroupProperties + Clone,
580 G::CharTab: SubspaceDecomposable<T>,
581 T: ComplexFloat + Lapack,
582 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug,
583 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
584 SC: StructureConstraint + Clone + Hash + Eq + fmt::Display,
585{
586 pub fn builder() -> MultiDeterminantRepAnalysisDriverBuilder<'a, G, T, B, SC> {
588 MultiDeterminantRepAnalysisDriverBuilder::default()
589 }
590
591 fn construct_sao(&self) -> Result<(Array2<T>, Option<Array2<T>>), anyhow::Error> {
594 let mut structure_constraint_set = self
595 .multidets
596 .iter()
597 .map(|multidet| multidet.structure_constraint())
598 .collect::<HashSet<_>>();
599 let structure_constraint = if structure_constraint_set.len() == 1 {
600 structure_constraint_set.drain().next().ok_or(format_err!(
601 "Unable to retrieve the structure constraint of the multi-determinantal wavefunctions."
602 ))
603 } else {
604 Err(format_err!(
605 "Inconsistent structure constraints across multi-determinantal wavefunctions."
606 ))
607 }?;
608
609 let mut nbas_set = self
610 .multidets
611 .iter()
612 .map(|multidet| {
613 multidet
614 .basis()
615 .first()
616 .expect("Unable to obtain the first determinant in the basis.")
617 .bao()
618 .n_funcs()
619 })
620 .collect::<HashSet<_>>();
621 let nbas = if nbas_set.len() == 1 {
622 nbas_set.drain().next().ok_or(format_err!(
623 "Unable to retrieve the number of basis functions describing the multi-determinantal wavefunctions."
624 ))
625 } else {
626 Err(format_err!(
627 "Inconsistent numbers of basis functions across multi-determinantal wavefunctions."
628 ))
629 }?;
630 let ncomps = structure_constraint.n_explicit_comps_per_coefficient_matrix();
631 let provided_dim = self.sao.nrows();
632
633 if provided_dim == nbas {
634 let sao = {
635 let mut sao_mut = Array2::zeros((ncomps * nbas, ncomps * nbas));
636 (0..ncomps).for_each(|icomp| {
637 let start = icomp * nbas;
638 let end = (icomp + 1) * nbas;
639 sao_mut
640 .slice_mut(s![start..end, start..end])
641 .assign(self.sao);
642 });
643 sao_mut
644 };
645
646 let sao_h = self.sao_h.map(|sao_h| {
647 let mut sao_h_mut = Array2::zeros((ncomps * nbas, ncomps * nbas));
648 (0..ncomps).for_each(|icomp| {
649 let start = icomp * nbas;
650 let end = (icomp + 1) * nbas;
651 sao_h_mut
652 .slice_mut(s![start..end, start..end])
653 .assign(sao_h);
654 });
655 sao_h_mut
656 });
657
658 Ok((sao, sao_h))
659 } else {
660 assert_eq!(provided_dim, nbas * ncomps);
661 Ok((self.sao.clone(), self.sao_h.cloned()))
662 }
663 }
664}
665
666impl<'a, T, B, SC> MultiDeterminantRepAnalysisDriver<'a, UnitaryRepresentedSymmetryGroup, T, B, SC>
670where
671 T: ComplexFloat + Lapack + Sync + Send,
672 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug + Sync + Send,
673 for<'b> Complex<f64>: Mul<&'b T, Output = Complex<f64>>,
674 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
675 SC: StructureConstraint + Hash + Eq + fmt::Display,
676{
677 fn_construct_unitary_group!(
678 construct_unitary_group
681 );
682}
683
684impl<'a, T, B, SC> MultiDeterminantRepAnalysisDriver<'a, MagneticRepresentedSymmetryGroup, T, B, SC>
688where
689 T: ComplexFloat + Lapack + Sync + Send,
690 <T as ComplexFloat>::Real: From<f64> + Sync + Send + fmt::LowerExp + fmt::Debug,
691 for<'b> Complex<f64>: Mul<&'b T, Output = Complex<f64>>,
692 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
693 SC: StructureConstraint + Hash + Eq + fmt::Display,
694{
695 fn_construct_magnetic_group!(
696 construct_magnetic_group
699 );
700}
701
702#[duplicate_item(
706 duplicate!{
707 [ dtype_nested; [f64]; [Complex<f64>] ]
708 duplicate!{
709 [ sctype_nested; [SpinConstraint] ]
710 duplicate!{
711 [
712 [
713 btype_nested [OrbitBasis<'a, UnitaryRepresentedSymmetryGroup, SlaterDeterminant<'a, dtype_nested, sctype_nested>>]
714 calc_smat_nested [calc_smat_optimised]
715 ]
716 [
717 btype_nested [EagerBasis<SlaterDeterminant<'a, dtype_nested, sctype_nested>>]
718 calc_smat_nested [calc_smat]
719 ]
720 ]
721 [
722 gtype_ [ UnitaryRepresentedSymmetryGroup ]
723 dtype_ [ dtype_nested ]
724 btype_ [ btype_nested ]
725 sctype_ [ sctype_nested ]
726 doc_sub_ [ "Performs representation analysis using a unitary-represented group and stores the result." ]
727 analyse_fn_ [ analyse_representation ]
728 construct_group_ [ self.construct_unitary_group()? ]
729 calc_smat_ [ calc_smat_nested ]
730 calc_projections_ [
731 log_subtitle("Multi-determinantal wavefunction projection decompositions");
732 qsym2_output!("");
733 qsym2_output!(" Projections are defined w.r.t. the following inner product:");
734 qsym2_output!(" {}", multidet_orbit.origin().overlap_definition());
735 qsym2_output!("");
736 multidet_orbit
737 .projections_to_string(
738 &multidet_orbit.calc_projection_compositions()?,
739 params.integrality_threshold,
740 )
741 .log_output_display();
742 qsym2_output!("");
743 ]
744 ]
745 }
746 }
747 }
748 duplicate!{
749 [ dtype_nested; [f64]; [Complex<f64>] ]
750 duplicate!{
751 [ sctype_nested; [SpinConstraint] ]
752 duplicate!{
753 [
754 [
755 btype_nested [OrbitBasis<'a, MagneticRepresentedSymmetryGroup, SlaterDeterminant<'a, dtype_nested, sctype_nested>>]
756 calc_smat_nested [calc_smat_optimised]
757 ]
758 [
759 btype_nested [EagerBasis<SlaterDeterminant<'a, dtype_nested, sctype_nested>>]
760 calc_smat_nested [calc_smat]
761 ]
762 ]
763 [
764 gtype_ [ MagneticRepresentedSymmetryGroup ]
765 dtype_ [ dtype_nested ]
766 btype_ [ btype_nested ]
767 sctype_ [ sctype_nested ]
768 doc_sub_ [ "Performs corepresentation analysis using a magnetic-represented group and stores the result." ]
769 analyse_fn_ [ analyse_corepresentation ]
770 construct_group_ [ self.construct_magnetic_group()? ]
771 calc_smat_ [ calc_smat_nested ]
772 calc_projections_ [ ]
773 ]
774 }
775 }
776 }
777)]
778impl<'a> MultiDeterminantRepAnalysisDriver<'a, gtype_, dtype_, btype_, sctype_> {
779 #[doc = doc_sub_]
780 fn analyse_fn_(&mut self) -> Result<(), anyhow::Error> {
781 let params = self.parameters;
782 let (sao, sao_h) = self.construct_sao()?;
783 let group = construct_group_;
784 log_cc_transversal(&group);
785 let _ = find_angular_function_representation(&group, self.angular_function_parameters);
786 if group.is_double_group() {
787 let _ = find_spinor_function_representation(&group, self.angular_function_parameters);
788 }
789 if let Some(det) = self
790 .multidets
791 .get(0)
792 .and_then(|multidet| multidet.basis().first())
793 {
794 log_bao(det.bao());
795 }
796
797 let (multidet_symmetries, multidet_symmetries_thresholds): (Vec<_>, Vec<_>) = self.multidets
798 .iter()
799 .enumerate()
800 .map(|(i, multidet)| {
801 log_micsec_begin(&format!("Multi-determinantal wavefunction {i}"));
802 qsym2_output!("");
803 let res = MultiDeterminantSymmetryOrbit::builder()
804 .group(&group)
805 .origin(multidet)
806 .integrality_threshold(params.integrality_threshold)
807 .linear_independence_threshold(params.linear_independence_threshold)
808 .symmetry_transformation_kind(params.symmetry_transformation_kind.clone())
809 .eigenvalue_comparison_mode(params.eigenvalue_comparison_mode.clone())
810 .build()
811 .map_err(|err| format_err!(err))
812 .and_then(|mut multidet_orbit| {
813 multidet_orbit
814 .calc_smat_(Some(&sao), sao_h.as_ref(), params.use_cayley_table)?
815 .normalise_smat()?
816 .calc_xmat(false)?;
817 log_overlap_eigenvalues(
818 "Overlap eigenvalues",
819 multidet_orbit.smat_eigvals.as_ref().ok_or(format_err!("Orbit overlap eigenvalues not found."))?,
820 params.linear_independence_threshold,
821 ¶ms.eigenvalue_comparison_mode
822 );
823 qsym2_output!("");
824 let multidet_symmetry_thresholds = multidet_orbit
825 .smat_eigvals
826 .as_ref()
827 .map(|eigvals| {
828 let mut eigvals_vec = eigvals.iter().collect::<Vec<_>>();
829 match multidet_orbit.eigenvalue_comparison_mode() {
830 EigenvalueComparisonMode::Modulus => {
831 eigvals_vec.sort_by(|a, b| {
832 a.abs().partial_cmp(&b.abs()).expect("Unable to compare two eigenvalues based on their moduli.")
833 });
834 }
835 EigenvalueComparisonMode::Real => {
836 eigvals_vec.sort_by(|a, b| {
837 a.re().partial_cmp(&b.re()).expect("Unable to compare two eigenvalues based on their real parts.")
838 });
839 }
840 }
841 let eigval_above = match multidet_orbit.eigenvalue_comparison_mode() {
842 EigenvalueComparisonMode::Modulus => eigvals_vec
843 .iter()
844 .find(|val| {
845 val.abs() >= multidet_orbit.linear_independence_threshold
846 })
847 .copied()
848 .copied(),
849 EigenvalueComparisonMode::Real => eigvals_vec
850 .iter()
851 .find(|val| {
852 val.re() >= multidet_orbit.linear_independence_threshold
853 })
854 .copied()
855 .copied(),
856 };
857 eigvals_vec.reverse();
858 let eigval_below = match multidet_orbit.eigenvalue_comparison_mode() {
859 EigenvalueComparisonMode::Modulus => eigvals_vec
860 .iter()
861 .find(|val| {
862 val.abs() < multidet_orbit.linear_independence_threshold
863 })
864 .copied()
865 .copied(),
866 EigenvalueComparisonMode::Real => eigvals_vec
867 .iter()
868 .find(|val| {
869 val.re() < multidet_orbit.linear_independence_threshold
870 })
871 .copied()
872 .copied(),
873 };
874 (eigval_above, eigval_below)
875 })
876 .unwrap_or((None, None));
877 let multidet_sym = multidet_orbit.analyse_rep().map_err(|err| err.to_string());
878 { calc_projections_ }
879 Ok((multidet_sym, multidet_symmetry_thresholds))
880 })
881 .unwrap_or_else(|err| (Err(err.to_string()), (None, None)));
882 log_micsec_end(&format!("Multi-determinantal wavefunction {i}"));
883 qsym2_output!("");
884 res
885 }).unzip();
886
887 let result = MultiDeterminantRepAnalysisResult::builder()
888 .parameters(params)
889 .multidets(self.multidets.clone())
890 .group(group)
891 .multidet_symmetries(multidet_symmetries)
892 .multidet_symmetries_thresholds(multidet_symmetries_thresholds)
893 .build()?;
894 self.result = Some(result);
895
896 Ok(())
897 }
898}
899
900impl<'a, G, T, B, SC> fmt::Display for MultiDeterminantRepAnalysisDriver<'a, G, T, B, SC>
908where
909 G: SymmetryGroupProperties + Clone,
910 G::CharTab: SubspaceDecomposable<T>,
911 T: ComplexFloat + Lapack,
912 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug,
913 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
914 SC: StructureConstraint + Hash + Eq + fmt::Display,
915{
916 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
917 write_title(f, "Multi-determinantal Wavefunction Symmetry Analysis")?;
918 writeln!(f)?;
919 writeln!(f, "{}", self.parameters)?;
920 Ok(())
921 }
922}
923
924impl<'a, G, T, B, SC> fmt::Debug for MultiDeterminantRepAnalysisDriver<'a, G, T, B, SC>
925where
926 G: SymmetryGroupProperties + Clone,
927 G::CharTab: SubspaceDecomposable<T>,
928 T: ComplexFloat + Lapack,
929 <T as ComplexFloat>::Real: From<f64> + fmt::LowerExp + fmt::Debug,
930 B: Basis<SlaterDeterminant<'a, T, SC>> + Clone,
931 SC: StructureConstraint + Hash + Eq + fmt::Display,
932{
933 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
934 writeln!(f, "{self}")
935 }
936}
937
938#[duplicate_item(
942 duplicate!{
943 [ dtype_nested; [f64]; [Complex<f64>] ]
944 duplicate!{
945 [ sctype_nested; [SpinConstraint] ]
946 duplicate!{
947 [
948 btype_nested;
949 [OrbitBasis<'a, UnitaryRepresentedSymmetryGroup, SlaterDeterminant<'a, dtype_nested, sctype_nested>>];
950 [EagerBasis<SlaterDeterminant<'a, dtype_nested, sctype_nested>>]
951 ]
952 [
953 gtype_ [ UnitaryRepresentedSymmetryGroup ]
954 dtype_ [ dtype_nested ]
955 btype_ [ btype_nested ]
956 sctype_ [ sctype_nested ]
957 analyse_fn_ [ analyse_representation ]
958 ]
959 }
960 }
961 }
962 duplicate!{
963 [ dtype_nested; [f64]; [Complex<f64>] ]
964 duplicate!{
965 [ sctype_nested; [SpinConstraint] ]
966 duplicate!{
967 [
968 btype_nested;
969 [OrbitBasis<'a, MagneticRepresentedSymmetryGroup, SlaterDeterminant<'a, dtype_nested, sctype_nested>>];
970 [EagerBasis<SlaterDeterminant<'a, dtype_nested, sctype_nested>>]
971 ]
972 [
973 gtype_ [ MagneticRepresentedSymmetryGroup ]
974 dtype_ [ dtype_nested ]
975 btype_ [ btype_nested ]
976 sctype_ [ sctype_nested ]
977 analyse_fn_ [ analyse_corepresentation ]
978 ]
979 }
980 }
981 }
982)]
983impl<'a> QSym2Driver for MultiDeterminantRepAnalysisDriver<'a, gtype_, dtype_, btype_, sctype_> {
984 type Params = MultiDeterminantRepAnalysisParams<f64>;
985
986 type Outcome = MultiDeterminantRepAnalysisResult<'a, gtype_, dtype_, btype_, sctype_>;
987
988 fn result(&self) -> Result<&Self::Outcome, anyhow::Error> {
989 self.result.as_ref().ok_or_else(|| {
990 format_err!(
991 "No multi-determinantal wavefunction representation analysis results found."
992 )
993 })
994 }
995
996 fn run(&mut self) -> Result<(), anyhow::Error> {
997 self.log_output_display();
998 self.analyse_fn_()?;
999 self.result()?.log_output_display();
1000 Ok(())
1001 }
1002}