1use std::cmp::Ordering;
4use std::collections::HashSet;
5use std::fmt;
6use std::hash::{Hash, Hasher};
7use std::str::FromStr;
8
9use approx;
10use derive_builder::Builder;
11use indexmap::{IndexMap, IndexSet};
12use itertools::Itertools;
13use log;
14use nalgebra::Vector3;
15use ndarray::{Array2, ArrayView2, Axis};
16use num_traits::ToPrimitive;
17use phf::{phf_map, phf_set};
18use serde::{Deserialize, Serialize};
19
20use crate::chartab::character::Character;
21use crate::chartab::chartab_group::CharacterProperties;
22use crate::chartab::chartab_symbols::{
23 disambiguate_linspace_symbols, CollectionSymbol, DecomposedSymbol,
24 DecomposedSymbolBuilderError, GenericSymbol, GenericSymbolParsingError, LinearSpaceSymbol,
25 MathematicalSymbol, ReducibleLinearSpaceSymbol,
26};
27use crate::chartab::unityroot::UnityRoot;
28use crate::chartab::CharacterTable;
29use crate::group::FiniteOrder;
30use crate::symmetry::symmetry_element::symmetry_operation::SpecialSymmetryTransformation;
31use crate::symmetry::symmetry_element::SymmetryElement;
32use crate::symmetry::symmetry_element_order::ORDER_1;
33use crate::symmetry::symmetry_group::SymmetryGroupProperties;
34
35#[cfg(test)]
36#[path = "symmetry_symbols_tests.rs"]
37mod symmetry_symbols_tests;
38
39static MULLIKEN_IRREP_DEGENERACIES: phf::Map<&'static str, u64> = phf_map! {
45 "A" => 1u64,
46 "B" => 1u64,
47 "Σ" => 1u64,
48 "Γ" => 1u64,
49 "E" => 2u64,
50 "Π" => 2u64,
51 "Δ" => 2u64,
52 "Φ" => 2u64,
53 "T" => 3u64,
54 "F" => 4u64,
55 "H" => 5u64,
56 "I" => 6u64,
57 "J" => 7u64,
58 "K" => 8u64,
59 "L" => 9u64,
60 "M" => 10u64,
61 "A~" => 1u64,
62 "B~" => 1u64,
63 "Σ~" => 1u64,
64 "Γ~" => 1u64,
65 "E~" => 2u64,
66 "Π~" => 2u64,
67 "Δ~" => 2u64,
68 "Φ~" => 2u64,
69 "T~" => 3u64,
70 "F~" => 4u64,
71 "H~" => 5u64,
72 "I~" => 6u64,
73 "J~" => 7u64,
74 "K~" => 8u64,
75 "L~" => 9u64,
76 "M~" => 10u64,
77};
78
79static INV_MULLIKEN_IRREP_DEGENERACIES: phf::Map<u64, &'static str> = phf_map! {
81 2u64 => "E",
82 3u64 => "T",
83 4u64 => "F",
84 5u64 => "H",
85 6u64 => "I",
86 7u64 => "J",
87 8u64 => "K",
88 9u64 => "L",
89 10u64 => "M",
90};
91
92pub static FORCED_C3_PRINCIPAL_GROUPS: phf::Set<&'static str> = phf_set! {
94 "O",
95 "Oh",
96 "Td",
97 "O + θ·O",
98 "Oh + θ·Oh",
99 "Td + θ·Td",
100 "u[O + θ·O]",
101 "u[Oh + θ·Oh]",
102 "u[Td + θ·Td]",
103 "O*",
104 "Oh*",
105 "Td*",
106 "(O + θ·O)*",
107 "(Oh + θ·Oh)*",
108 "(Td + θ·Td)*",
109 "u[(O + θ·O)*]",
110 "u[(Oh + θ·Oh)*]",
111 "u[(Td + θ·Td)*]",
112};
113
114#[derive(Builder, Debug, Clone, PartialEq, Eq, Hash, PartialOrd, Serialize, Deserialize)]
124pub struct MullikenIrrepSymbol {
125 generic_symbol: GenericSymbol,
127}
128
129impl MullikenIrrepSymbol {
130 fn builder() -> MullikenIrrepSymbolBuilder {
131 MullikenIrrepSymbolBuilder::default()
132 }
133
134 pub fn new(symstr: &str) -> Result<Self, MullikenIrrepSymbolBuilderError> {
152 Self::from_str(symstr)
153 }
154}
155
156#[derive(Builder, Debug, Clone, PartialEq, Eq, Hash, PartialOrd, Serialize, Deserialize)]
162pub struct MullikenIrcorepSymbol {
163 inducing_irreps: DecomposedSymbol<MullikenIrrepSymbol>,
165}
166
167impl MullikenIrcorepSymbol {
168 fn builder() -> MullikenIrcorepSymbolBuilder {
170 MullikenIrcorepSymbolBuilder::default()
171 }
172
173 pub fn new(symstr: &str) -> Result<Self, MullikenIrcorepSymbolBuilderError> {
186 Self::from_str(symstr)
187 }
188}
189
190#[derive(Builder, Debug, Clone, Serialize, Deserialize)]
196pub struct SymmetryClassSymbol<R: Clone + Serialize> {
197 generic_symbol: GenericSymbol,
199
200 representatives: Option<Vec<R>>,
202}
203
204impl<R: Clone + Serialize> SymmetryClassSymbol<R> {
205 fn builder() -> SymmetryClassSymbolBuilder<R> {
206 SymmetryClassSymbolBuilder::default()
207 }
208
209 pub fn new(symstr: &str, rep: Option<Vec<R>>) -> Result<Self, GenericSymbolParsingError> {
240 let generic_symbol = GenericSymbol::from_str(symstr)?;
241 if generic_symbol.multiplicity().is_none() {
242 Err(GenericSymbolParsingError(format!(
243 "{symstr} contains no class multiplicity prefactor."
244 )))
245 } else {
246 Ok(Self::builder()
247 .generic_symbol(generic_symbol)
248 .representatives(rep)
249 .build()
250 .unwrap_or_else(|_| panic!("Unable to construct a class symbol from `{symstr}`.")))
251 }
252 }
253}
254
255impl MathematicalSymbol for MullikenIrrepSymbol {
264 fn main(&self) -> String {
266 self.generic_symbol.main()
267 }
268
269 fn presuper(&self) -> String {
272 self.generic_symbol.presuper()
273 }
274
275 fn presub(&self) -> String {
276 self.generic_symbol.presub()
277 }
278
279 fn postsuper(&self) -> String {
281 self.generic_symbol.postsuper()
282 }
283
284 fn postsub(&self) -> String {
287 self.generic_symbol.postsub()
288 }
289
290 fn prefactor(&self) -> String {
292 String::new()
293 }
294
295 fn postfactor(&self) -> String {
297 String::new()
298 }
299
300 fn multiplicity(&self) -> Option<usize> {
302 if let Some(&mult) = MULLIKEN_IRREP_DEGENERACIES.get(&self.main()) {
303 Some(
304 mult.try_into()
305 .unwrap_or_else(|_| panic!("Unable to convert {mult} to `usize`.")),
306 )
307 } else {
308 None
309 }
310 }
311}
312
313impl FromStr for MullikenIrrepSymbol {
314 type Err = MullikenIrrepSymbolBuilderError;
315
316 fn from_str(symstr: &str) -> Result<Self, Self::Err> {
341 let generic_symbol = GenericSymbol::from_str(symstr)
342 .unwrap_or_else(|_| panic!("Unable to parse {symstr} as a generic symbol."));
343 Self::builder().generic_symbol(generic_symbol).build()
344 }
345}
346
347impl LinearSpaceSymbol for MullikenIrrepSymbol {
348 fn dimensionality(&self) -> usize {
349 usize::try_from(
350 *MULLIKEN_IRREP_DEGENERACIES
351 .get(&self.main())
352 .unwrap_or_else(|| {
353 panic!(
354 "Unknown dimensionality for Mulliken symbol {}.",
355 self.main()
356 )
357 }),
358 )
359 .expect("Unable to convert the dimensionality of this irrep to `usize`.")
360 }
361
362 fn set_dimensionality(&mut self, dim: usize) -> bool {
363 let dim_u64 = u64::try_from(dim).unwrap_or_else(|err| {
364 log::error!("{err}");
365 panic!("Unable to convert `{dim}` to `u64`.")
366 });
367 let dim_from_current_main = *MULLIKEN_IRREP_DEGENERACIES.get(&self.main()).unwrap_or(&0);
368 if dim_from_current_main == dim_u64 {
369 log::debug!(
370 "The current main symbol `{}` already has the right dimension of `{dim}`. No new Mulliken symbols will be set.",
371 self.main()
372 );
373 false
374 } else {
375 let main_opt = INV_MULLIKEN_IRREP_DEGENERACIES.get(&dim_u64);
376 if let Some(main) = main_opt {
377 self.generic_symbol.set_main(main);
378 true
379 } else {
380 log::warn!("Unable to retrieve an unambiguous Mulliken symbol for dimensionality `{dim_u64}`. Main symbol of {self} will be kept unchanged.");
381 false
382 }
383 }
384 }
385}
386
387impl fmt::Display for MullikenIrrepSymbol {
388 fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
389 write!(f, "{}", self.generic_symbol)
390 }
391}
392
393impl From<DecomposedSymbolBuilderError> for MullikenIrcorepSymbolBuilderError {
398 fn from(value: DecomposedSymbolBuilderError) -> Self {
399 match value {
400 DecomposedSymbolBuilderError::UninitializedField(msg) => {
401 MullikenIrcorepSymbolBuilderError::UninitializedField(msg)
402 }
403 DecomposedSymbolBuilderError::ValidationError(msg) => {
404 MullikenIrcorepSymbolBuilderError::ValidationError(msg)
405 }
406 }
407 }
408}
409
410impl FromStr for MullikenIrcorepSymbol {
411 type Err = MullikenIrcorepSymbolBuilderError;
412
413 fn from_str(symstr: &str) -> Result<Self, Self::Err> {
433 MullikenIrcorepSymbol::builder()
434 .inducing_irreps(DecomposedSymbol::from_str(symstr)?)
435 .build()
436 }
437}
438
439impl ReducibleLinearSpaceSymbol for MullikenIrcorepSymbol {
440 type Subspace = MullikenIrrepSymbol;
441
442 fn from_subspaces(irreps: &[(Self::Subspace, usize)]) -> Self {
443 Self::builder()
444 .inducing_irreps(DecomposedSymbol::from_subspaces(irreps))
445 .build()
446 .expect("Unable to construct a Mulliken ircorep symbol from a slice of irrep symbols.")
447 }
448
449 fn subspaces(&self) -> Vec<(&Self::Subspace, &usize)> {
450 self.inducing_irreps.subspaces()
451 }
452}
453
454impl fmt::Display for MullikenIrcorepSymbol {
455 fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
456 write!(f, "D[{}]", self.main())
457 }
458}
459
460impl<R: Clone + Serialize> PartialEq for SymmetryClassSymbol<R> {
465 fn eq(&self, other: &Self) -> bool {
466 self.generic_symbol == other.generic_symbol
467 }
468}
469
470impl<R: Clone + Serialize> Eq for SymmetryClassSymbol<R> {}
471
472impl<R: Clone + Serialize> Hash for SymmetryClassSymbol<R> {
473 fn hash<H: Hasher>(&self, state: &mut H) {
474 self.generic_symbol.hash(state);
475 }
476}
477
478impl<R: Clone + Serialize> MathematicalSymbol for SymmetryClassSymbol<R> {
479 fn main(&self) -> String {
481 self.generic_symbol.main()
482 }
483
484 fn presuper(&self) -> String {
486 String::new()
487 }
488
489 fn presub(&self) -> String {
491 String::new()
492 }
493
494 fn postsuper(&self) -> String {
496 String::new()
497 }
498
499 fn postsub(&self) -> String {
501 String::new()
502 }
503
504 fn prefactor(&self) -> String {
506 self.generic_symbol.prefactor()
507 }
508
509 fn postfactor(&self) -> String {
511 String::new()
512 }
513
514 fn multiplicity(&self) -> Option<usize> {
517 self.generic_symbol.multiplicity()
518 }
519}
520
521impl<R: Clone + Serialize> CollectionSymbol for SymmetryClassSymbol<R> {
522 type CollectionElement = R;
523
524 fn from_reps(
525 symstr: &str,
526 reps: Option<Vec<Self::CollectionElement>>,
527 ) -> Result<Self, GenericSymbolParsingError> {
528 Self::new(symstr, reps)
529 }
530
531 fn representative(&self) -> Option<&Self::CollectionElement> {
532 self.representatives.as_ref().map(|reps| &reps[0])
533 }
534
535 fn representatives(&self) -> Option<&Vec<Self::CollectionElement>> {
536 self.representatives.as_ref()
537 }
538}
539
540impl<R: SpecialSymmetryTransformation + Clone + Serialize> SpecialSymmetryTransformation
541 for SymmetryClassSymbol<R>
542{
543 fn is_proper(&self) -> bool {
553 self.representative()
554 .as_ref()
555 .expect("No representative element found for this class.")
556 .is_proper()
557 }
558
559 fn is_spatial_identity(&self) -> bool {
560 self.representative()
561 .as_ref()
562 .expect("No representative element found for this class.")
563 .is_spatial_identity()
564 }
565
566 fn is_spatial_binary_rotation(&self) -> bool {
567 self.representative()
568 .as_ref()
569 .expect("No representative element found for this class.")
570 .is_spatial_binary_rotation()
571 }
572
573 fn is_spatial_inversion(&self) -> bool {
574 self.representative()
575 .as_ref()
576 .expect("No representative element found for this class.")
577 .is_spatial_inversion()
578 }
579
580 fn is_spatial_reflection(&self) -> bool {
581 self.representative()
582 .as_ref()
583 .expect("No representative element found for this class.")
584 .is_spatial_reflection()
585 }
586
587 fn contains_time_reversal(&self) -> bool {
597 self.representative()
598 .as_ref()
599 .expect("No representative element found for this class.")
600 .contains_time_reversal()
601 }
602
603 fn is_su2(&self) -> bool {
613 self.representative()
614 .as_ref()
615 .expect("No representative element found for this class.")
616 .is_su2()
617 }
618
619 fn is_su2_class_1(&self) -> bool {
625 self.representatives()
626 .expect("No representative element found for this class.")
627 .iter()
628 .all(|rep| rep.is_su2_class_1())
629 }
630}
631
632impl<R: FiniteOrder + Clone + Serialize> FiniteOrder for SymmetryClassSymbol<R> {
633 type Int = R::Int;
634
635 fn order(&self) -> Self::Int {
636 self.representative()
637 .as_ref()
638 .expect("No representative element found for this class.")
639 .order()
640 }
641}
642
643impl<R: Clone + Serialize> fmt::Display for SymmetryClassSymbol<R> {
644 fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
645 write!(f, "{}", self.generic_symbol)
646 }
647}
648
649pub(super) fn sort_irreps<R: Clone + Serialize + SpecialSymmetryTransformation>(
676 char_arr: &ArrayView2<Character>,
677 frobenius_schur_indicators: &[i8],
678 class_symbols: &IndexMap<SymmetryClassSymbol<R>, usize>,
679 principal_classes: &[SymmetryClassSymbol<R>],
680) -> (Array2<Character>, Vec<i8>) {
681 log::debug!("Sorting irreducible representations...");
682 let su2_0 = if class_symbols.keys().any(|cc_sym| cc_sym.is_su2()) {
683 "(Σ)"
684 } else {
685 ""
686 };
687 let class_e = class_symbols
688 .first()
689 .expect("No class symbols found.")
690 .0
691 .clone();
692 let class_e1: SymmetryClassSymbol<R> = SymmetryClassSymbol::new("1||E(QΣ)||", None)
693 .expect("Unable to construct class symbol `1||E(QΣ)||`.");
694 let class_i: SymmetryClassSymbol<R> = SymmetryClassSymbol::new(&format!("1||i{su2_0}||"), None)
695 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||i{su2_0}||`."));
696 let class_ti: SymmetryClassSymbol<R> =
697 SymmetryClassSymbol::new(&format!("1||θ·i{su2_0}||"), None)
698 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||θ·i{su2_0}||`."));
699 let class_s: SymmetryClassSymbol<R> =
700 SymmetryClassSymbol::new(&format!("1||σh{su2_0}||"), None)
701 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||σh{su2_0}||`."));
702 let class_s2: SymmetryClassSymbol<R> = SymmetryClassSymbol::new("1||σh(Σ), σh(QΣ)||", None)
703 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||σh(Σ), σh(QΣ)||`."));
704 let class_ts: SymmetryClassSymbol<R> =
705 SymmetryClassSymbol::new(&format!("1||θ·σh{su2_0}||"), None)
706 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||θ·σh{su2_0}||`."));
707 let class_t: SymmetryClassSymbol<R> = SymmetryClassSymbol::new(&format!("1||θ{su2_0}||"), None)
708 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||θ{su2_0}||`."));
709
710 let mut leading_classes: IndexSet<SymmetryClassSymbol<R>> = IndexSet::new();
711 let mut sign_only_classes: HashSet<SymmetryClassSymbol<R>> = HashSet::new();
712
713 if class_symbols.contains_key(&class_e1) {
715 leading_classes.insert(class_e1.clone());
716 sign_only_classes.insert(class_e1);
717 }
718
719 if class_symbols.contains_key(&class_t) {
721 leading_classes.insert(class_t.clone());
722 sign_only_classes.insert(class_t);
723 }
724
725 if class_symbols.contains_key(&class_i) {
728 leading_classes.insert(class_i.clone());
729 sign_only_classes.insert(class_i);
730 } else if class_symbols.contains_key(&class_ti) {
731 leading_classes.insert(class_ti.clone());
732 sign_only_classes.insert(class_ti);
733 } else if class_symbols.contains_key(&class_s) {
734 leading_classes.insert(class_s.clone());
735 sign_only_classes.insert(class_s);
736 } else if class_symbols.contains_key(&class_s2) {
737 leading_classes.insert(class_s2.clone());
738 sign_only_classes.insert(class_s2);
739 } else if class_symbols.contains_key(&class_ts) {
740 leading_classes.insert(class_ts.clone());
741 sign_only_classes.insert(class_ts);
742 };
743
744 leading_classes.insert(class_e);
746
747 leading_classes.extend(principal_classes.iter().cloned());
749
750 log::debug!("Irreducible representation sort order:");
751 for leading_cc in leading_classes.iter() {
752 log::debug!(
753 " {}{}",
754 leading_cc,
755 if sign_only_classes.contains(leading_cc) {
756 " (sign only)"
757 } else {
758 ""
759 }
760 );
761 }
762
763 let leading_idxs: IndexSet<usize> = leading_classes
764 .iter()
765 .map(|cc| {
766 *class_symbols
767 .get(cc)
768 .unwrap_or_else(|| panic!("Class `{cc}` cannot be found."))
769 })
770 .collect();
771
772 let n_rows = char_arr.nrows();
773 let mut col_idxs: Vec<usize> = Vec::with_capacity(n_rows);
774 col_idxs.extend(leading_idxs.iter());
775 col_idxs.extend((1..n_rows).filter(|i| !leading_idxs.contains(i)));
776 let mut sort_arr = char_arr.select(Axis(1), &col_idxs);
777 let one = Character::new(&[(UnityRoot::new(0, 2), 1)]);
778 let m_one = Character::new(&[(UnityRoot::new(1, 2), 1)]);
779 sign_only_classes.iter().for_each(|class| {
780 let col_idx = leading_classes
781 .get_index_of(class)
782 .unwrap_or_else(|| panic!("Unable to obtain the column index of class `{class}`."));
783 sort_arr.column_mut(col_idx).mapv_inplace(|character| {
784 let character_c = character.complex_value();
785 if approx::relative_eq!(
786 character_c.im,
787 0.0,
788 max_relative = character.threshold(),
789 epsilon = character.threshold()
790 ) {
791 if character_c.re > 0.0 {
792 one.clone()
793 } else {
794 m_one.clone()
795 }
796 } else if approx::relative_eq!(
797 character_c.re,
798 0.0,
799 max_relative = character.threshold(),
800 epsilon = character.threshold()
801 ) {
802 if character_c.im > 0.0 {
803 one.clone()
804 } else {
805 m_one.clone()
806 }
807 } else {
808 panic!("Character {character} is neither purely real nor purely imaginary for sign-only sorting.")
809 }
810 });
811 });
812 let sort_row_indices: Vec<_> = (0..n_rows)
832 .sorted_by(|&i, &j| {
833 let keys_i = sort_arr.row(i).iter().cloned().collect_vec();
834 let keys_j = sort_arr.row(j).iter().cloned().collect_vec();
835 keys_i
836 .partial_cmp(&keys_j)
837 .unwrap_or_else(|| panic!("`{keys_i:?}` and `{keys_j:?}` cannot be compared."))
838 })
839 .collect();
840 let char_arr = char_arr.select(Axis(0), &sort_row_indices);
841 let old_fs = frobenius_schur_indicators.iter().collect::<Vec<_>>();
842 let sorted_fs = sort_row_indices.iter().map(|&i| *old_fs[i]).collect_vec();
843 log::debug!("Sorting irreducible representations... Done.");
844 (char_arr, sorted_fs)
845}
846
847pub(super) fn deduce_principal_classes<R, P>(
871 class_symbols: &IndexMap<SymmetryClassSymbol<R>, usize>,
872 force_principal_predicate: Option<P>,
873 force_principal: Option<SymmetryClassSymbol<R>>,
874) -> Vec<SymmetryClassSymbol<R>>
875where
876 R: fmt::Debug + SpecialSymmetryTransformation + FiniteOrder + Clone + Serialize,
877 P: Copy + Fn(&SymmetryClassSymbol<R>) -> bool,
878{
879 log::debug!("Determining principal classes...");
880 let principal_classes = force_principal.map_or_else(|| {
881 let mut sorted_class_symbols = class_symbols
887 .clone()
888 .sorted_by(|cc1, _, cc2, _| {
889 PartialOrd::partial_cmp(
890 &(cc1.order(), !cc1.contains_time_reversal(), cc1.is_proper(), !cc1.is_su2_class_1()),
891 &(cc2.order(), !cc2.contains_time_reversal(), cc2.is_proper(), !cc2.is_su2_class_1()),
892 )
893 .expect("Unable to sort class symbols.")
894 })
895 .rev();
896 let (principal_rep, _) = if let Some(predicate) = force_principal_predicate {
897 sorted_class_symbols
900 .find(|(cc, _)| predicate(cc))
901 .expect("`Unable to find classes fulfilling the specified predicate.`")
902 } else {
903 sorted_class_symbols
906 .next()
907 .expect("Unexpected empty `sorted_class_symbols`.")
908 };
909
910 class_symbols
913 .keys()
914 .filter(|cc| {
915 cc.contains_time_reversal() == principal_rep.contains_time_reversal()
916 && cc.is_proper() == principal_rep.is_proper()
917 && cc.order() == principal_rep.order()
918 && cc.is_su2_class_1() == principal_rep.is_su2_class_1()
919 })
920 .cloned()
921 .collect::<Vec<_>>()
922 }, |principal_cc| {
923 assert!(
924 force_principal_predicate.is_none(),
925 "`force_principal_predicate` and `force_principal` cannot be both provided."
926 );
927 assert!(
928 class_symbols.contains_key(&principal_cc),
929 "Forcing principal-axis class to be `{principal_cc}`, but `{principal_cc}` is not a valid class in the group."
930 );
931
932 log::warn!(
933 "Principal-axis class forced to be `{principal_cc}`. Auto-detection of principal-axis classes will be skipped.",
934 );
935 vec![principal_cc]
936 });
937
938 assert!(!principal_classes.is_empty());
939 if principal_classes.len() == 1 {
940 log::debug!("Principal-axis class found: {}", principal_classes[0]);
941 } else {
942 log::debug!("Principal-axis classes found:");
943 for princc in &principal_classes {
944 log::debug!(" {}", princc);
945 }
946 }
947 log::debug!("Determining principal classes... Done.");
948 principal_classes
949}
950
951#[allow(clippy::too_many_lines)]
969pub(super) fn deduce_mulliken_irrep_symbols<R>(
970 char_arr: &ArrayView2<Character>,
971 class_symbols: &IndexMap<SymmetryClassSymbol<R>, usize>,
972 principal_classes: &[SymmetryClassSymbol<R>],
973) -> Vec<MullikenIrrepSymbol>
974where
975 R: fmt::Debug + SpecialSymmetryTransformation + FiniteOrder + Clone + Serialize,
976{
977 log::debug!("Generating Mulliken irreducible representation symbols...");
978
979 let su2_0 = if class_symbols.keys().any(|cc_sym| cc_sym.is_su2()) {
980 "(Σ)"
981 } else {
982 ""
983 };
984
985 let e_cc = class_symbols
986 .first()
987 .expect("No class symbols found.")
988 .0
989 .clone();
990 let e1_cc: SymmetryClassSymbol<R> = SymmetryClassSymbol::new("1||E(QΣ)||", None)
991 .expect("Unable to construct class symbol `1||E(QΣ)||`.");
992 let i_cc: SymmetryClassSymbol<R> = SymmetryClassSymbol::new(&format!("1||i{su2_0}||"), None)
993 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||i{su2_0}||`."));
994 let ti_cc: SymmetryClassSymbol<R> = SymmetryClassSymbol::new(&format!("1||θ·i{su2_0}||"), None)
995 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||θ·i{su2_0}||`."));
996 let s_cc: SymmetryClassSymbol<R> = SymmetryClassSymbol::new(&format!("1||σh{su2_0}||"), None)
997 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||σh{su2_0}||`."));
998 let s2_cc: SymmetryClassSymbol<R> = SymmetryClassSymbol::new("1||σh(Σ), σh(QΣ)||", None)
999 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||σh(Σ), σh(QΣ)||`."));
1000 let ts_cc: SymmetryClassSymbol<R> =
1001 SymmetryClassSymbol::new(&format!("1||θ·σh{su2_0}||"), None)
1002 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||θ·σh{su2_0}||`."));
1003 let t_cc: SymmetryClassSymbol<R> = SymmetryClassSymbol::new(&format!("1||θ{su2_0}||"), None)
1004 .unwrap_or_else(|_| panic!("Unable to construct class symbol `1||θ{su2_0}||`."));
1005
1006 let su2_1_parity = class_symbols.contains_key(&e1_cc);
1008
1009 let i_parity = class_symbols.contains_key(&i_cc);
1011
1012 let ti_parity = class_symbols.contains_key(&ti_cc);
1014
1015 let s_parity = class_symbols.contains_key(&s_cc) || class_symbols.contains_key(&s2_cc);
1017
1018 let ts_parity = class_symbols.contains_key(&ts_cc);
1020
1021 if su2_1_parity {
1023 log::debug!("E(QΣ) found. This will be used to determine projective representations.");
1024 }
1025
1026 if i_parity {
1028 log::debug!("Inversion centre found. This will be used for g/u ordering.");
1029 } else if ti_parity {
1030 log::debug!(
1031 "Time-reversed inversion centre found (but no inversion centre). This will be used for g/u ordering."
1032 );
1033 } else if s_parity {
1034 log::debug!(
1035 "Horizontal mirror plane found (but no inversion centre). This will be used for '/'' ordering."
1036 );
1037 } else if ts_parity {
1038 log::debug!(
1039 "Time-reversed horizontal mirror plane found (but no inversion centre). This will be used for '/'' ordering."
1040 );
1041 }
1042
1043 let t_parity = class_symbols.contains_key(&t_cc);
1045 if t_parity {
1046 log::debug!("Time reversal found. This will be used for magnetic ordering.");
1047 }
1048
1049 let e2p1 = UnityRoot::new(1u32, 2u32);
1050 let e2p2 = UnityRoot::new(2u32, 2u32);
1051 let char_p1 = Character::new(&[(e2p2, 1usize)]);
1052 let char_m1 = Character::new(&[(e2p1, 1usize)]);
1053
1054 log::debug!("First pass: assign symbols from rules");
1056
1057 let mut raw_irrep_symbols = char_arr.rows().into_iter().map(|irrep| {
1058 let dim = irrep[
1060 *class_symbols.get(&e_cc).unwrap_or_else(|| panic!("Class `{}` not found.", &e_cc))
1061 ].complex_value();
1062 assert!(
1063 approx::relative_eq!(dim.im, 0.0)
1064 && approx::relative_eq!(dim.re.round(), dim.re)
1065 && dim.re.round() > 0.0
1066 );
1067
1068 assert!(dim.re.round() > 0.0);
1069 #[allow(clippy::cast_sign_loss, clippy::cast_possible_truncation)]
1070 let dim = dim.re.round() as u64;
1071 assert!(dim > 0);
1072 let main = if dim >= 2 {
1073 INV_MULLIKEN_IRREP_DEGENERACIES.get(&dim).map_or_else(|| {
1075 log::warn!("{} cannot be assigned a standard dimensionality symbol. A generic 'Λ' will be used instead.", dim);
1076 "Λ"
1077 }, |sym| sym)
1078 } else {
1079 let complex = irrep.map(|character| character.complex_conjugate()) != irrep;
1081 if complex {
1082 "Γ"
1083 } else {
1084 let char_rots: HashSet<_> = principal_classes
1085 .iter()
1086 .map(|cc| {
1087 irrep[
1088 *class_symbols.get(cc).unwrap_or_else(|| panic!("Class `{cc}` not found."))
1089 ].clone()
1090 })
1091 .collect();
1092 if char_rots.len() == 1 && *char_rots.iter().next().expect("No rotation classes found.") == char_p1 {
1093 "A"
1094 } else if char_rots
1095 .iter()
1096 .all(|char_rot| char_rot.clone() == char_p1 || char_rot.clone() == char_m1)
1097 {
1098 "B"
1099 } else {
1100 panic!("");
1101 }
1105 }
1106 };
1107
1108 let projective = if su2_1_parity {
1109 let char_e1 = irrep[
1110 *class_symbols
1111 .get(&e1_cc)
1112 .unwrap_or_else(|| {
1113 panic!("The character for `{e1_cc}` could not be found.")
1114 })
1115 ].clone();
1116 let char_e1_c = char_e1.complex_value();
1117 assert!(
1118 approx::relative_eq!(
1119 char_e1_c.im,
1120 0.0,
1121 epsilon = char_e1.threshold(),
1122 max_relative = char_e1.threshold()
1123 ) && approx::relative_eq!(
1124 char_e1_c.re.round(),
1125 char_e1_c.re,
1126 epsilon = char_e1.threshold(),
1127 max_relative = char_e1.threshold()
1128 ),
1129 );
1130
1131 #[allow(clippy::cast_possible_truncation)]
1132 let char_e1_c = char_e1_c
1133 .re
1134 .round()
1135 .to_i32()
1136 .unwrap_or_else(|| panic!("Unable to convert the real part of the character for `{e1_cc}` to `i32`."));
1137 let dim_i32 = i32::try_from(dim).expect("Unable to convert the irrep dimensionality to `i32`.");
1138 if char_e1_c == dim_i32 {
1139 false
1140 } else {
1141 assert_eq!(char_e1_c, -dim_i32);
1142 true
1143 }
1144 } else {
1145 false
1146 };
1147 let projective_str = if projective { "~" } else { "" };
1148
1149 let (inv, mir) = if i_parity || ti_parity {
1150 let char_inv = irrep[
1153 *class_symbols
1154 .get(&i_cc)
1155 .unwrap_or_else(|| {
1156 class_symbols
1157 .get(&ti_cc)
1158 .unwrap_or_else(|| {
1159 panic!("Neither `{}` nor `{}` found.", &i_cc, &ti_cc)
1160 })
1161 })
1162 ].clone();
1163 let char_inv_c = char_inv.complex_value();
1164 assert!(
1165 approx::relative_eq!(
1166 char_inv_c.im,
1167 0.0,
1168 epsilon = char_inv.threshold(),
1169 max_relative = char_inv.threshold()
1170 ) && approx::relative_eq!(
1171 char_inv_c.re.round(),
1172 char_inv_c.re,
1173 epsilon = char_inv.threshold(),
1174 max_relative = char_inv.threshold()
1175 ),
1176 );
1177
1178 #[allow(clippy::cast_possible_truncation)]
1179 let char_inv_c = char_inv_c.re.round() as i32;
1180 match char_inv_c.cmp(&0) {
1181 Ordering::Greater => ("g", ""),
1182 Ordering::Less => ("u", ""),
1183 Ordering::Equal => panic!("Inversion character must not be zero."),
1184 }
1185 } else if !projective && (s_parity || ts_parity) {
1186 let char_ref = irrep[
1190 *class_symbols
1191 .get(&s_cc)
1192 .unwrap_or_else(|| {
1193 class_symbols
1194 .get(&s2_cc)
1195 .unwrap_or_else(|| {
1196 class_symbols
1197 .get(&ts_cc)
1198 .unwrap_or_else(|| {
1199 panic!("Neither `{}`, `{}`, nor `{}` found.", &s_cc, &s2_cc, &ts_cc)
1200 })
1201 })
1202 })
1203 ].clone();
1204 let char_ref_c = char_ref.complex_value();
1205 assert!(
1206 approx::relative_eq!(
1207 char_ref_c.im,
1208 0.0,
1209 epsilon = char_ref.threshold(),
1210 max_relative = char_ref.threshold()
1211 ) && approx::relative_eq!(
1212 char_ref_c.re.round(),
1213 char_ref_c.re,
1214 epsilon = char_ref.threshold(),
1215 max_relative = char_ref.threshold()
1216 ),
1217 );
1218
1219 #[allow(clippy::cast_possible_truncation)]
1220 let char_ref_c = char_ref_c.re.round() as i32;
1221 match char_ref_c.cmp(&0) {
1222 Ordering::Greater => ("", "'"),
1223 Ordering::Less => ("", "''"),
1224 Ordering::Equal => panic!("Reflection character must not be zero for linear irreducible representations."),
1225 }
1226 } else {
1227 ("", "")
1228 };
1229
1230
1231 let trev = if t_parity {
1232 let char_trev = irrep[
1234 *class_symbols.get(&t_cc).unwrap_or_else(|| {
1235 panic!("Class `{}` not found.", &t_cc)
1236 })
1237 ].clone();
1238 let char_trev_c = char_trev.complex_value();
1239 if approx::relative_eq!(
1240 char_trev_c.im,
1241 0.0,
1242 epsilon = char_trev.threshold(),
1243 max_relative = char_trev.threshold()
1244 ) && approx::relative_eq!(
1245 char_trev_c.re.round(),
1246 char_trev_c.re,
1247 epsilon = char_trev.threshold(),
1248 max_relative = char_trev.threshold()
1249 ) {
1250 #[allow(clippy::cast_possible_truncation)]
1252 let char_trev_c = char_trev_c.re.round() as i32;
1253 match char_trev_c.cmp(&0) {
1254 Ordering::Greater => "+",
1255 Ordering::Less => "-",
1256 Ordering::Equal => panic!("Real time-reversal character must not be zero."),
1257 }
1258 } else {
1259 ""
1261 }
1262 } else {
1263 ""
1264 };
1265
1266 MullikenIrrepSymbol::new(format!("|^({trev})|{main}{projective_str}|^({mir})_({inv})|").as_str())
1267 .unwrap_or_else(|_| {
1268 panic!(
1269 "Unable to construct symmetry symbol `|^({trev})|{main}{projective_str}|^({mir})_({inv})|`."
1270 )
1271 })
1272 }).collect_vec();
1273
1274 let mut complex_irrep_indices = char_arr
1275 .rows()
1276 .into_iter()
1277 .enumerate()
1278 .filter_map(|(i, irrep)| {
1279 let complex = irrep.map(|character| character.complex_conjugate()) != irrep;
1280 if complex {
1281 Some(i)
1282 } else {
1283 None
1284 }
1285 })
1286 .collect::<IndexSet<_>>();
1287 complex_irrep_indices.reverse();
1288
1289 let cc_pairs = if !complex_irrep_indices.is_empty() {
1290 log::debug!("Grouping pairs of complex-conjugate irreps...");
1291 let mut cc_pairs: Vec<(usize, usize)> = vec![];
1292 while !complex_irrep_indices.is_empty() {
1293 let complex_irrep_index = complex_irrep_indices.pop().unwrap();
1294 let complex_irrep = char_arr.row(complex_irrep_index);
1295 let complex_conj_irrep = complex_irrep.map(|character| character.complex_conjugate());
1296 let complex_conj_irrep_index = char_arr
1297 .rows()
1298 .into_iter()
1299 .enumerate()
1300 .find_map(|(i, irrep)| {
1301 if irrep == complex_conj_irrep {
1302 Some(i)
1303 } else {
1304 None
1305 }
1306 })
1307 .expect("Unable to find the complex-conjugate irrep.");
1308 complex_irrep_indices.shift_remove(&complex_conj_irrep_index);
1309 cc_pairs.push((complex_irrep_index, complex_conj_irrep_index));
1310 raw_irrep_symbols[complex_irrep_index]
1311 .generic_symbol
1312 .set_presub("a");
1313 raw_irrep_symbols[complex_conj_irrep_index]
1314 .generic_symbol
1315 .set_presub("b");
1316 }
1317 log::debug!(
1318 "There {} {} {} of complex-conjugate irreps.",
1319 if cc_pairs.len() == 1 { "is" } else { "are" },
1320 cc_pairs.len(),
1321 if cc_pairs.len() == 1 { "pair" } else { "pairs" },
1322 );
1323 log::debug!("Grouping pairs of complex-conjugate irreps... Done.");
1324 Some(cc_pairs)
1325 } else {
1326 None
1327 };
1328
1329 log::debug!("Second pass: disambiguate identical cases not distinguishable by rules");
1330 let mut irrep_symbols = disambiguate_linspace_symbols(raw_irrep_symbols.into_iter());
1331
1332 if let Some(cc_pairs_vec) = cc_pairs {
1333 log::debug!("Equalising post-subscripts in complex-conjugate pairs...");
1334 for (complex_irrep_index, complex_conj_irrep_index) in &cc_pairs_vec {
1335 let complex_irrep_postsub = irrep_symbols[*complex_irrep_index].postsub();
1336 irrep_symbols[*complex_conj_irrep_index]
1337 .generic_symbol
1338 .set_postsub(&complex_irrep_postsub);
1339 }
1340 log::debug!("Equalising post-subscripts in complex-conjugate pairs... Done.");
1341 }
1342
1343 log::debug!("Generating Mulliken irreducible representation symbols... Done.");
1344
1345 irrep_symbols
1346}
1347
1348#[must_use]
1362pub(super) fn deduce_sigma_symbol(
1363 sigma_axis: &Vector3<f64>,
1364 principal_element: &SymmetryElement,
1365 thresh: f64,
1366 force_d: bool,
1367) -> Option<String> {
1368 if approx::relative_eq!(
1369 principal_element.raw_axis().dot(sigma_axis).abs(),
1370 0.0,
1371 epsilon = thresh,
1372 max_relative = thresh
1373 ) && *principal_element.raw_proper_order() != ORDER_1
1374 {
1375 if force_d {
1377 Some("d".to_owned())
1378 } else {
1379 Some("v".to_owned())
1380 }
1381 } else if approx::relative_eq!(
1382 principal_element.raw_axis().cross(sigma_axis).norm(),
1383 0.0,
1384 epsilon = thresh,
1385 max_relative = thresh
1386 ) && *principal_element.raw_proper_order() != ORDER_1
1387 {
1388 Some("h".to_owned())
1390 } else {
1391 None
1392 }
1393}
1394
1395#[derive(Clone, PartialEq, Eq, Hash, Debug)]
1397pub enum MirrorParity {
1398 Even,
1400
1401 Odd,
1403
1404 Neither,
1406}
1407
1408impl fmt::Display for MirrorParity {
1409 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1410 match self {
1411 MirrorParity::Even => write!(f, "(+)"),
1412 MirrorParity::Odd => write!(f, "(-)"),
1413 MirrorParity::Neither => write!(f, "( )"),
1414 }
1415 }
1416}
1417
1418pub(crate) fn deduce_mirror_parities<G, R>(
1436 group: &G,
1437 rep: &R,
1438) -> IndexMap<<<G as CharacterProperties>::CharTab as CharacterTable>::ColSymbol, MirrorParity>
1439where
1440 G: SymmetryGroupProperties,
1441 R: ReducibleLinearSpaceSymbol<
1442 Subspace = <<G as CharacterProperties>::CharTab as CharacterTable>::RowSymbol,
1443 >,
1444{
1445 let id_cc = group
1446 .get_cc_symbol_of_index(0)
1447 .expect("Unable to obtain the conjugacy class symbol of the identity.");
1448 let mirrors = group.filter_cc_symbols(|cc| cc.is_spatial_reflection());
1449 mirrors
1450 .iter()
1451 .map(|sigma_cc| {
1452 let sigma_cc_chars = rep
1453 .subspaces()
1454 .iter()
1455 .map(|(irrep, _)| {
1456 let id_char = group.character_table().get_character(irrep, &id_cc);
1457 let char = group.character_table().get_character(irrep, sigma_cc);
1458 if *char == *id_char {
1459 MirrorParity::Even
1460 } else if *char == -id_char {
1461 MirrorParity::Odd
1462 } else {
1463 MirrorParity::Neither
1464 }
1465 })
1466 .collect::<HashSet<_>>();
1467 let sigma_cc_parity = if sigma_cc_chars.len() == 1 {
1468 sigma_cc_chars
1469 .into_iter()
1470 .next()
1471 .expect("Unable to extract the mirror parity.")
1472 } else {
1473 MirrorParity::Neither
1474 };
1475 (sigma_cc.clone(), sigma_cc_parity)
1476 })
1477 .collect::<IndexMap<_, _>>()
1478}