qsym2/symmetry/symmetry_element/mod.rs
1//! Geometrical symmetry elements.
2
3use std::convert::TryInto;
4use std::fmt;
5use std::hash::{Hash, Hasher};
6
7use approx;
8use derive_builder::Builder;
9use fraction;
10use log;
11use nalgebra::Vector3;
12use num::integer::gcd;
13use num_traits::{ToPrimitive, Zero};
14use serde::{Deserialize, Serialize};
15
16use crate::auxiliary::geometry;
17use crate::auxiliary::misc::{self, HashableFloat};
18use crate::symmetry::symmetry_element_order::ElementOrder;
19
20type F = fraction::GenericFraction<u32>;
21
22pub mod symmetry_operation;
23pub use symmetry_operation::*;
24
25#[cfg(test)]
26mod symmetry_element_tests;
27
28// ====================================
29// Enum definitions and implementations
30// ====================================
31
32/// Enumerated type to classify the type of the antiunitary term that contributes to a symmetry
33/// element.
34#[derive(Copy, Clone, Debug, PartialEq, Eq, Hash, Serialize, Deserialize)]
35pub enum AntiunitaryKind {
36 /// Variant for the antiunitary term being a complex-conjugation operation.
37 ComplexConjugation,
38
39 /// Variant for the antiunitary term being a time-reversal operation.
40 TimeReversal,
41}
42
43/// Enumerated type to classify the types of symmetry element.
44#[derive(Copy, Clone, Debug, PartialEq, Eq, Hash, Serialize, Deserialize)]
45pub enum SymmetryElementKind {
46 /// Proper symmetry element which consists of just a proper rotation axis.
47 ///
48 /// The associated option indicates the type of the associated antiunitary operation, if any.
49 Proper(Option<AntiunitaryKind>),
50
51 /// Improper symmetry element in the mirror-plane convention, which consists of a proper
52 /// rotation axis and an orthogonal mirror plane.
53 ///
54 /// The associated option indicates the type of the associated antiunitary operation, if any.
55 ImproperMirrorPlane(Option<AntiunitaryKind>),
56
57 /// Improper symmetry element in the inversion-centre convention, which consists of a proper
58 /// rotation axis and an inversion centre.
59 ///
60 /// The associated option indicates the type of the associated antiunitary operation, if any.
61 ImproperInversionCentre(Option<AntiunitaryKind>),
62}
63
64impl SymmetryElementKind {
65 /// Indicates if the antiunitary part of this element contains a time-reversal operation.
66 #[must_use]
67 pub fn contains_time_reversal(&self) -> bool {
68 match self {
69 Self::Proper(tr)
70 | Self::ImproperMirrorPlane(tr)
71 | Self::ImproperInversionCentre(tr) => *tr == Some(AntiunitaryKind::TimeReversal),
72 }
73 }
74
75 /// Indicates if an antiunitary operation is associated with this element.
76 #[must_use]
77 pub fn contains_antiunitary(&self) -> bool {
78 match self {
79 Self::Proper(au)
80 | Self::ImproperMirrorPlane(au)
81 | Self::ImproperInversionCentre(au) => au.is_some(),
82 }
83 }
84
85 /// Converts the current kind to the desired time-reversal form.
86 ///
87 /// # Arguments
88 ///
89 /// * `tr` - A boolean indicating whether time reversal is included (`true`) or removed
90 /// (`false`).
91 ///
92 /// # Returns
93 ///
94 /// A copy of the current kind with or without the time-reversal antiunitary operation.
95 #[must_use]
96 pub fn to_tr(&self, tr: bool) -> Self {
97 if tr {
98 self.to_antiunitary(Some(AntiunitaryKind::TimeReversal))
99 } else {
100 self.to_antiunitary(None)
101 }
102 }
103
104 /// Converts the associated antiunitary operation to the desired kind.
105 ///
106 /// # Arguments
107 ///
108 /// * `au` - An option containing the desired antiunitary kind.
109 ///
110 /// # Returns
111 ///
112 /// A new symmetry element kind with the desired antiunitary kind.
113 pub fn to_antiunitary(&self, au: Option<AntiunitaryKind>) -> Self {
114 match self {
115 Self::Proper(_) => Self::Proper(au),
116 Self::ImproperMirrorPlane(_) => Self::ImproperMirrorPlane(au),
117 Self::ImproperInversionCentre(_) => Self::ImproperInversionCentre(au),
118 }
119 }
120}
121
122/// Enumerated type to indicate the rotation group to which the unitary proper rotation part of a
123/// symmetry element belongs.
124#[derive(Clone, Hash, PartialEq, Eq, Debug, Serialize, Deserialize)]
125pub enum RotationGroup {
126 /// Variant indicating that the proper part of the symmetry element generates rotations in
127 /// $`\mathsf{SO}(3)`$.
128 SO3,
129
130 /// Variant indicating that the proper part of the symmetry element generates rotations in
131 /// $`\mathsf{SU}(2)`$. The associated boolean indicates whether the proper part of the element
132 /// itself is connected to the identity via a homotopy path of class 0 (`true`) or class 1
133 /// (`false`).
134 SU2(bool),
135}
136
137impl RotationGroup {
138 /// Indicates if the rotation is in $`\mathsf{SU}(2)`$.
139 pub fn is_su2(&self) -> bool {
140 matches!(self, RotationGroup::SU2(_))
141 }
142
143 /// Indicates if the rotation is in $`\mathsf{SU}(2)`$ and connected to the
144 /// identity via a homotopy path of class 1.
145 pub fn is_su2_class_1(&self) -> bool {
146 matches!(self, RotationGroup::SU2(false))
147 }
148}
149
150// ======================================
151// Struct definitions and implementations
152// ======================================
153
154pub struct SymmetryElementKindConversionError(String);
155
156impl fmt::Debug for SymmetryElementKindConversionError {
157 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
158 f.debug_struct("SymmetryElementKindConversionError")
159 .field("Message", &self.0)
160 .finish()
161 }
162}
163
164impl fmt::Display for SymmetryElementKindConversionError {
165 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
166 write!(
167 f,
168 "SymmetryElementKindConversionError with message: {}",
169 &self.0
170 )
171 }
172}
173
174impl std::error::Error for SymmetryElementKindConversionError {}
175
176impl TryInto<geometry::ImproperRotationKind> for SymmetryElementKind {
177 type Error = SymmetryElementKindConversionError;
178
179 fn try_into(self) -> Result<geometry::ImproperRotationKind, Self::Error> {
180 match self {
181 Self::Proper(_) => Err(SymmetryElementKindConversionError(
182 "Unable to convert a proper element to an `ImproperRotationKind`.".to_string(),
183 )),
184 Self::ImproperMirrorPlane(_) => Ok(geometry::ImproperRotationKind::MirrorPlane),
185 Self::ImproperInversionCentre(_) => Ok(geometry::ImproperRotationKind::InversionCentre),
186 }
187 }
188}
189
190impl fmt::Display for SymmetryElementKind {
191 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
192 match self {
193 Self::Proper(au) => match au {
194 None => write!(f, "Proper"),
195 Some(AntiunitaryKind::TimeReversal) => write!(f, "Time-reversed proper"),
196 Some(AntiunitaryKind::ComplexConjugation) => write!(f, "Complex-conjugated proper"),
197 },
198 Self::ImproperMirrorPlane(au) => match au {
199 None => write!(f, "Improper (mirror-plane convention)"),
200 Some(AntiunitaryKind::TimeReversal) => {
201 write!(f, "Time-reversed improper (mirror-plane convention)")
202 }
203 Some(AntiunitaryKind::ComplexConjugation) => {
204 write!(f, "Complex-conjugated improper (mirror-plane convention)")
205 }
206 },
207 Self::ImproperInversionCentre(au) => match au {
208 None => write!(f, "Improper (inversion-centre convention)"),
209 Some(AntiunitaryKind::TimeReversal) => {
210 write!(f, "Time-reversed improper (inversion-centre convention)")
211 }
212 Some(AntiunitaryKind::ComplexConjugation) => {
213 write!(
214 f,
215 "Complex-conjugated improper (inversion-centre convention)"
216 )
217 }
218 },
219 }
220 }
221}
222
223/// Structure for storing and managing symmetry elements.
224///
225/// Each symmetry element is a geometrical object in $`\mathbb{R}^3`$ that encodes the following
226/// pieces of information:
227///
228/// * the axis of rotation $`\hat{\mathbf{n}}`$,
229/// * the angle of rotation $`\phi`$,
230/// * the associated improper operation, if any, and
231/// * the associated antiunitary operation, if any.
232///
233/// These pieces of information can be stored in the following representation of a symmetry element
234/// $`\hat{g}`$:
235///
236/// ```math
237/// \hat{g} = \hat{\alpha} \hat{\gamma} \hat{C}_n^k,
238/// ```
239///
240/// where
241/// * $`n \in \mathbb{N}_{+}`$, $`k \in \mathbb{Z}/n\mathbb{Z}`$ such that
242/// $`\lfloor -n/2 \rfloor < k \le \lfloor n/2 \rfloor`$,
243/// * $`\hat{\gamma}`$ is either the identity $`\hat{e}`$, the inversion operation $`\hat{i}`$, or
244/// a reflection operation $`\hat{\sigma}`$ perpendicular to the axis of rotation,
245/// * $`\hat{\alpha}`$ is either the identity $`\hat{e}`$, the complex conjugation $`\hat{K}`$, or
246/// the time reversal $`\hat{\theta}`$.
247///
248/// We shall refer to $`\hat{C}_n^k`$ as the *unitary proper rotation part*, $`\hat{\gamma}`$ the
249/// *improper rotation part*, and $`\hat{\alpha}`$ the *antiunitary part* of the symmetry element.
250///
251/// With this definition, the above pieces of information required to specify a geometrical symmetry
252/// element are given as follows:
253///
254/// * the axis of rotation $`\hat{\mathbf{n}}`$ is given by the axis of $`\hat{C}_n^k`$,
255/// * the angle of rotation $`\phi = 2\pi k/n \in (-\pi, \pi]`$,
256/// * the improper contribution $`\hat{\gamma}`$,
257/// * the antiunitary contribution $`\hat{\alpha}`$.
258///
259/// This definition also allows the unitary part of $`\hat{g}`$ (*i.e.* the composition
260/// $`\hat{\gamma} \hat{C}_n^k`$) to be interpreted as an element of either $`\mathsf{O}(3)`$ or
261/// $`\mathsf{SU}'(2)`$, which means that the unitary part of $`\hat{g}`$ is also a symmetry
262/// operation in the corresponding group, and a rather special one that can be used to generate
263/// other symmetry operations of the group. $`\hat{g}`$ thus serves as a bridge between molecular
264/// symmetry and abstract group theory.
265///
266/// There is one small caveat: for infinite-order elements, $`n`$ and $`k`$ can no longer be used
267/// to give the angle of rotation. There must thus be a mechanism to allow for infinite-order
268/// elements to be interpreted as an arbitrary finite-order one. An explicit specification of the
269/// angle of rotation $`\phi`$ seems to be the best way to do this. In other words, the angle of
270/// rotation of each element is specified by either a tuple of integers $`(k, n)`$ or a
271/// floating-point number $`\phi`$.
272#[derive(Builder, Clone, Serialize, Deserialize)]
273pub struct SymmetryElement {
274 /// The rotational order $`n`$ of the unitary proper rotation part of the symmetry element. This
275 /// can be finite or infinite, and will determine whether the proper power is `None` or
276 /// contains an integer value.
277 #[builder(setter(name = "proper_order"))]
278 raw_proper_order: ElementOrder,
279
280 /// The power $`k \in \mathbb{Z}/n\mathbb{Z}`$ of the unitary proper rotation part of the
281 /// symmetry element such that $`\lfloor -n/2 \rfloor < k <= \lfloor n/2 \rfloor`$. This is only
282 /// defined if [`Self::raw_proper_order`] is finite.
283 ///
284 /// The unitary proper rotation does not include any additional unitary rotations introduced by
285 /// the antiunitary part of the symmetry element (*e.g.* time reversal).
286 #[builder(setter(custom, name = "proper_power"), default = "None")]
287 raw_proper_power: Option<i32>,
288
289 /// The normalised axis of the unitary proper rotation part of the symmetry element whose
290 /// direction is as specified at construction.
291 #[builder(setter(custom))]
292 raw_axis: Vector3<f64>,
293
294 /// The spatial and antiunitary kind of the symmetry element.
295 #[builder(default = "SymmetryElementKind::Proper(None)")]
296 kind: SymmetryElementKind,
297
298 /// The rotation group in which the unitary proper rotation part of the symmetry element shall
299 /// be interpreted.
300 rotation_group: RotationGroup,
301
302 /// A boolean indicating whether the symmetry element is a generator of the group to which it
303 /// belongs.
304 #[builder(default = "false")]
305 generator: bool,
306
307 /// A threshold for approximate equality comparisons.
308 #[builder(setter(custom))]
309 threshold: f64,
310
311 /// An additional superscript for distinguishing symmetry elements.
312 #[builder(default = "String::new()")]
313 pub(crate) additional_superscript: String,
314
315 /// An additional subscript for distinguishing symmetry elements.
316 #[builder(default = "String::new()")]
317 pub(crate) additional_subscript: String,
318
319 /// The fraction $`k/n \in (-1/2, 1/2]`$ of the unitary proper rotation, represented exactly
320 /// for hashing and comparison purposes.
321 ///
322 /// The unitary proper rotation does not include any additional unitary rotations introduced by
323 /// the antiunitary part of the symmetry element (*e.g.* time reversal).
324 ///
325 /// This is not defined for infinite-order elements and cannot be set arbitrarily.
326 #[builder(setter(skip), default = "self.calc_proper_fraction()")]
327 proper_fraction: Option<F>,
328
329 /// The normalised proper angle $`\phi \in (-\pi, \pi]`$ corresponding to the unitary proper
330 /// rotation.
331 ///
332 /// The unitary proper rotation does not include any additional unitary rotations introduced by
333 /// the antiunitary part of the symmetry element (*e.g.* time reversal).
334 ///
335 /// This can be set arbitrarily only for infinite-order elements.
336 #[builder(setter(custom), default = "self.calc_proper_angle()")]
337 proper_angle: Option<f64>,
338}
339
340impl SymmetryElementBuilder {
341 /// Sets the power of the unitary proper rotation part of the element.
342 ///
343 /// # Arguments
344 ///
345 /// * `prop_pow` - A proper power to be set. This will be folded into the interval
346 /// $`(\lfloor -n/2 \rfloor, \lfloor n/2 \rfloor]`$.
347 pub fn proper_power(&mut self, prop_pow: i32) -> &mut Self {
348 let raw_proper_order = self
349 .raw_proper_order
350 .as_ref()
351 .expect("Proper order has not been set.");
352 self.raw_proper_power = match raw_proper_order {
353 ElementOrder::Int(io) => {
354 let io_i32 =
355 i32::try_from(*io).expect("Unable to convert the integer order to `i32`.");
356 let residual = prop_pow.rem_euclid(io_i32);
357 if residual > io_i32.div_euclid(2) {
358 Some(Some(residual - io_i32))
359 } else {
360 Some(Some(residual))
361 }
362 }
363 ElementOrder::Inf => None,
364 };
365 self
366 }
367
368 /// Sets the rotation angle of the unitary proper rotation part of the infinite-order element.
369 ///
370 /// # Arguments
371 ///
372 /// * `ang` - A proper rotation angle to be set. This will be folded into the interval
373 /// $`(-\pi, \pi]`$.
374 ///
375 /// # Panics
376 ///
377 /// Panics when `self` is of finite order.
378 pub fn proper_angle(&mut self, ang: f64) -> &mut Self {
379 let proper_order = self
380 .raw_proper_order
381 .as_ref()
382 .expect("Proper order has not been set.");
383 self.proper_angle = match proper_order {
384 ElementOrder::Int(_) => panic!(
385 "Arbitrary proper rotation angles can only be set for infinite-order elements."
386 ),
387 ElementOrder::Inf => {
388 let (normalised_rotation_angle, _) = geometry::normalise_rotation_angle(
389 ang,
390 self.threshold.expect("Threshold value has not been set."),
391 );
392 Some(Some(normalised_rotation_angle))
393 }
394 };
395 self
396 }
397
398 /// Sets the raw axis of the element.
399 ///
400 /// # Arguments
401 ///
402 /// * `axs` - The raw axis which will be normalised.
403 pub fn raw_axis(&mut self, axs: Vector3<f64>) -> &mut Self {
404 let thresh = self.threshold.expect("Threshold value has not been set.");
405 if approx::abs_diff_eq!(axs.norm(), 1.0, epsilon = thresh) {
406 self.raw_axis = Some(axs);
407 } else {
408 log::warn!("Axis not normalised. Normalising...");
409 self.raw_axis = Some(axs.normalize());
410 }
411 self
412 }
413
414 /// Sets the comparison threshold of the element.
415 ///
416 /// # Arguments
417 ///
418 /// * `thresh` - The comparison threshold..
419 pub fn threshold(&mut self, thresh: f64) -> &mut Self {
420 if thresh >= 0.0 {
421 self.threshold = Some(thresh);
422 } else {
423 log::error!(
424 "Threshold value `{}` is invalid. Threshold must be non-negative.",
425 thresh
426 );
427 self.threshold = None;
428 }
429 self
430 }
431
432 /// Calculates the fraction $`k/n \in (-1/2, 1/2]`$ of the unitary proper rotation, represented
433 /// exactly for hashing and comparison purposes.
434 ///
435 /// The unitary proper rotation does not include any additional unitary rotations introduced by
436 /// the antiunitary part of the symmetry element (*e.g.* time reversal).
437 ///
438 /// This is not defined for infinite-order elements and cannot be set arbitrarily.
439 fn calc_proper_fraction(&self) -> Option<F> {
440 let raw_proper_order = self
441 .raw_proper_order
442 .as_ref()
443 .expect("Proper order has not been set.");
444 match raw_proper_order {
445 ElementOrder::Int(io) => {
446 // The generating element has a proper fraction, pp/n.
447 //
448 // If pp/n > 1/2, we seek a positive integer x such that
449 // -1/2 < pp/n - x <= 1/2.
450 // It turns out that x ∈ [pp/n - 1/2, pp/n + 1/2).
451 //
452 // If pp/n <= -1/2, we seek a positive integer x such that
453 // -1/2 < pp/n + x <= 1/2.
454 // It turns out that x ∈ (-pp/n - 1/2, -pp/n + 1/2].
455 //
456 // x is then used to bring pp/n back into the (-1/2, 1/2] interval.
457 //
458 // See S.L. Altmann, Rotations, Quaternions, and Double Groups (Dover
459 // Publications, Inc., New York, 2005) for further information.
460 let pp = self
461 .raw_proper_power
462 .expect("Proper power has not been set.")
463 .expect("No proper powers found.");
464 let total_proper_fraction = if pp >= 0 {
465 F::new(pp.unsigned_abs(), *io)
466 } else {
467 F::new_neg(pp.unsigned_abs(), *io)
468 };
469 Some(total_proper_fraction)
470 }
471 ElementOrder::Inf => None,
472 }
473 }
474
475 /// Calculates the normalised proper angle $`\phi \in (-\pi, \pi]`$ corresponding to the unitary
476 /// proper rotation.
477 ///
478 /// The unitary proper rotation does not include any additional unitary rotations introduced by
479 /// the antiunitary part of the symmetry element (*e.g.* time reversal).
480 ///
481 /// This can be set arbitrarily only for infinite-order elements.
482 fn calc_proper_angle(&self) -> Option<f64> {
483 let raw_proper_order = self
484 .raw_proper_order
485 .as_ref()
486 .expect("Proper order has not been set.");
487 match raw_proper_order {
488 ElementOrder::Int(io) => {
489 let pp = self
490 .raw_proper_power
491 .expect("Proper power has not been set.")
492 .expect("No proper powers found.");
493 let total_proper_fraction = if pp >= 0 {
494 F::new(pp.unsigned_abs(), *io)
495 } else {
496 F::new_neg(pp.unsigned_abs(), *io)
497 };
498 Some(
499 total_proper_fraction
500 .to_f64()
501 .expect("Unable to convert the proper fraction to `f64`.")
502 * 2.0
503 * std::f64::consts::PI,
504 )
505 }
506 ElementOrder::Inf => self.proper_angle.unwrap_or(None),
507 }
508 }
509}
510
511impl SymmetryElement {
512 /// Returns a builder to construct a new symmetry element.
513 ///
514 /// # Returns
515 ///
516 /// A builder to construct a new symmetry element.
517 #[must_use]
518 pub fn builder() -> SymmetryElementBuilder {
519 SymmetryElementBuilder::default()
520 }
521
522 /// Returns the raw order of the unitary proper rotation. This might not be equal to the value of
523 /// $`n`$ if the fraction $`k/n`$ has been reduced.
524 pub fn raw_proper_order(&self) -> &ElementOrder {
525 &self.raw_proper_order
526 }
527
528 /// Returns the raw power of the unitary proper rotation. This might not be equal to the value of
529 /// $`k`$ if the fraction $`k/n`$ has been reduced.
530 pub fn raw_proper_power(&self) -> Option<&i32> {
531 self.raw_proper_power.as_ref()
532 }
533
534 /// Returns the raw axis of the unitary proper rotation.
535 pub fn raw_axis(&self) -> &Vector3<f64> {
536 &self.raw_axis
537 }
538
539 /// Returns the axis of the unitary proper rotation in the standard positive hemisphere.
540 pub fn standard_positive_axis(&self) -> Vector3<f64> {
541 geometry::get_standard_positive_pole(&self.raw_axis, self.threshold)
542 }
543
544 /// Returns the axis of the unitary proper rotation multiplied by the sign of the rotation angle.
545 /// If the proper rotation is a binary rotation, then the positive axis is always returned.
546 pub fn signed_axis(&self) -> Vector3<f64> {
547 let au = self.contains_antiunitary();
548 if self.is_o3_binary_rotation_axis(au) || self.is_o3_mirror_plane(au) {
549 self.standard_positive_axis()
550 } else {
551 self.proper_fraction
552 .map(|frac| {
553 frac.signum()
554 .to_f64()
555 .expect("Unable to obtain the sign of the proper fraction.")
556 })
557 .or_else(|| self.proper_angle.map(|proper_angle| proper_angle.signum())).map(|signum| signum * self.raw_axis)
558 .unwrap_or_else(|| {
559 log::warn!("No rotation signs could be obtained. The positive axis will be used for the signed axis.");
560 self.standard_positive_axis()
561 })
562 }
563 }
564
565 /// Returns the pole of the unitary proper rotation part of the element while leaving any improper
566 /// and antiunitary contributions intact.
567 ///
568 /// If the unitary proper rotation part if a binary rotation, the pole is always in the standard
569 /// positive hemisphere.
570 ///
571 /// # Returns
572 ///
573 /// The position vector of the proper rotation pole.
574 pub fn proper_rotation_pole(&self) -> Vector3<f64> {
575 match *self.raw_proper_order() {
576 ElementOrder::Int(_) => {
577 let frac_1_2 = F::new(1u32, 2u32);
578 let proper_fraction = self.proper_fraction.expect("No proper fractions found.");
579 if proper_fraction == frac_1_2 {
580 // Binary proper rotations
581 geometry::get_standard_positive_pole(&self.raw_axis, self.threshold)
582 } else if proper_fraction > F::zero() {
583 // Positive proper rotation angles
584 self.raw_axis
585 } else if proper_fraction < F::zero() {
586 // Negative proper rotation angles
587 -self.raw_axis
588 } else {
589 // Identity or inversion
590 assert!(proper_fraction.is_zero());
591 Vector3::zeros()
592 }
593 }
594 ElementOrder::Inf => {
595 if approx::abs_diff_eq!(
596 self.proper_angle.expect("No proper angles found."),
597 std::f64::consts::PI,
598 epsilon = self.threshold
599 ) {
600 // Binary proper rotations
601 geometry::get_standard_positive_pole(&self.raw_axis, self.threshold)
602 } else if approx::abs_diff_ne!(
603 self.proper_angle.expect("No proper angles found."),
604 0.0,
605 epsilon = self.threshold
606 ) {
607 self.proper_angle.expect("No proper angles found.").signum() * self.raw_axis
608 } else {
609 approx::assert_abs_diff_eq!(
610 self.proper_angle.expect("No proper angles found."),
611 0.0,
612 epsilon = self.threshold
613 );
614 Vector3::zeros()
615 }
616 }
617 }
618 }
619
620 /// Returns the unitary proper fraction for this element, if any.
621 ///
622 /// The element lacks a unitary proper fraction if it is infinite-order.
623 pub fn proper_fraction(&self) -> Option<&F> {
624 self.proper_fraction.as_ref()
625 }
626
627 /// Returns the unitary proper angle for this element, if any.
628 ///
629 /// The element lacks a unitary proper angle if it is infinite-order and the rotation angle has
630 /// not been set.
631 pub fn proper_angle(&self) -> Option<f64> {
632 self.proper_angle
633 }
634
635 /// Returns the spatial and antiunitary kind of this element.
636 pub fn kind(&self) -> &SymmetryElementKind {
637 &self.kind
638 }
639
640 /// Returns the rotation group and possibly the identity-connected homotopy class in which the
641 /// proper rotation part of this element is to be interpreted.
642 pub fn rotation_group(&self) -> &RotationGroup {
643 &self.rotation_group
644 }
645
646 /// Returns a boolean indicating if the element is a generator of a group.
647 pub fn is_generator(&self) -> bool {
648 self.generator
649 }
650
651 /// Returns the threshold for approximate comparisons.
652 pub fn threshold(&self) -> f64 {
653 self.threshold
654 }
655
656 /// Checks if the symmetry element contains a time-reversal operator as the antiunitary part.
657 ///
658 /// # Returns
659 ///
660 /// A boolean indicating if the symmetry element contains a time-reversal operator as the
661 /// antiunitary part.
662 #[must_use]
663 pub fn contains_time_reversal(&self) -> bool {
664 self.kind.contains_time_reversal()
665 }
666
667 /// Checks if the symmetry element contains an antiunitary part.
668 ///
669 /// # Returns
670 ///
671 /// Returns `None` if the symmetry element has no antiunitary parts, or `Some` wrapping around
672 /// the antiunitary kind if the symmetry element contains an antiunitary part.
673 pub fn contains_antiunitary(&self) -> Option<AntiunitaryKind> {
674 match self.kind {
675 SymmetryElementKind::Proper(au)
676 | SymmetryElementKind::ImproperMirrorPlane(au)
677 | SymmetryElementKind::ImproperInversionCentre(au) => au,
678 }
679 }
680
681 /// Checks if the unitary proper rotation part of the element is in $`\mathsf{SU}(2)`$.
682 #[must_use]
683 pub fn is_su2(&self) -> bool {
684 self.rotation_group.is_su2()
685 }
686
687 /// Checks if the unitary proper rotation part of the element is in $`\mathsf{SU}(2)`$ and
688 /// connected to the identity via a homotopy path of class 1.
689 ///
690 /// See S.L. Altmann, Rotations, Quaternions, and Double Groups (Dover Publications, Inc., New
691 /// York, 2005) for further information.
692 #[must_use]
693 pub fn is_su2_class_1(&self) -> bool {
694 self.rotation_group.is_su2_class_1()
695 }
696
697 /// Checks if the spatial part of the symmetry element is proper and has the specified
698 /// antiunitary attribute.
699 ///
700 /// More specifically, this checks if the combination $`\hat{\gamma} \hat{C}_n^k`$ is proper,
701 /// and if $`\hat{\alpha}`$ is as specified by `au`.
702 ///
703 /// # Arguments
704 ///
705 /// * `au` - An `Option` for the desired antiunitary kind.
706 ///
707 /// # Returns
708 ///
709 /// A boolean indicating if the symmetry element is proper and has the specified antiunitary
710 /// attribute.
711 #[must_use]
712 pub fn is_o3_proper(&self, au: Option<AntiunitaryKind>) -> bool {
713 self.kind == SymmetryElementKind::Proper(au)
714 }
715
716 /// Checks if the symmetry element is spatially an identity element and has the specified
717 /// antiunitary attribute.
718 ///
719 /// More specifically, this checks if the combination $`\hat{\gamma} \hat{C}_n^k`$ is the
720 /// identity, and if $`\hat{\alpha}`$ is as specified by `au`.
721 ///
722 /// # Arguments
723 ///
724 /// * `au` - An `Option` for the desired antiunitary kind.
725 ///
726 /// # Returns
727 ///
728 /// A boolean indicating if this symmetry element is spatially an identity element and has the
729 /// specified antiunitary attribute.
730 #[must_use]
731 pub fn is_o3_identity(&self, au: Option<AntiunitaryKind>) -> bool {
732 self.kind == SymmetryElementKind::Proper(au)
733 && self
734 .proper_fraction
735 .map(|frac| frac.is_zero())
736 .or_else(|| {
737 self.proper_angle.map(|proper_angle| {
738 approx::abs_diff_eq!(proper_angle, 0.0, epsilon = self.threshold,)
739 })
740 })
741 .unwrap_or(false)
742 }
743
744 /// Checks if the symmetry element is spatially an inversion centre and has the specified
745 /// antiunitary attribute.
746 ///
747 /// More specifically, this checks if the combination $`\hat{\gamma} \hat{C}_n^k`$ is the
748 /// spatial inversion, and if $`\hat{\alpha}`$ is as specified by `au`.
749 ///
750 /// # Arguments
751 ///
752 /// * `au` - An `Option` for the desired antiunitary kind.
753 ///
754 /// # Returns
755 ///
756 /// A boolean indicating if this symmetry element is spatially an inversion centre and has the
757 /// specified antiunitary attribute.
758 #[must_use]
759 pub fn is_o3_inversion_centre(&self, au: Option<AntiunitaryKind>) -> bool {
760 match self.kind {
761 SymmetryElementKind::ImproperMirrorPlane(sig_au) => {
762 sig_au == au
763 && self
764 .proper_fraction
765 .map(|frac| frac == F::new(1u32, 2u32))
766 .or_else(|| {
767 self.proper_angle.map(|proper_angle| {
768 approx::abs_diff_eq!(
769 proper_angle,
770 std::f64::consts::PI,
771 epsilon = self.threshold,
772 )
773 })
774 })
775 .unwrap_or(false)
776 }
777 SymmetryElementKind::ImproperInversionCentre(inv_au) => {
778 inv_au == au
779 && self
780 .proper_fraction
781 .map(|frac| frac.is_zero())
782 .or_else(|| {
783 self.proper_angle.map(|proper_angle| {
784 approx::abs_diff_eq!(proper_angle, 0.0, epsilon = self.threshold,)
785 })
786 })
787 .unwrap_or(false)
788 }
789 _ => false,
790 }
791 }
792
793 /// Checks if the symmetry element is spatially a binary rotation axis and has the specified
794 /// antiunitary attribute.
795 ///
796 /// More specifically, this checks if the combination $`\hat{\gamma} \hat{C}_n^k`$ is a binary
797 /// rotation, and if $`\hat{\alpha}`$ is as specified by `au`.
798 ///
799 /// # Arguments
800 ///
801 /// * `au` - An `Option` for the desired antiunitary kind.
802 ///
803 /// # Returns
804 ///
805 /// A boolean indicating if this symmetry element is spatially a binary rotation axis and has
806 /// the specified antiunitary attribute.
807 #[must_use]
808 pub fn is_o3_binary_rotation_axis(&self, au: Option<AntiunitaryKind>) -> bool {
809 self.kind == SymmetryElementKind::Proper(au)
810 && self
811 .proper_fraction
812 .map(|frac| frac == F::new(1u32, 2u32))
813 .or_else(|| {
814 self.proper_angle.map(|proper_angle| {
815 approx::abs_diff_eq!(
816 proper_angle,
817 std::f64::consts::PI,
818 epsilon = self.threshold,
819 )
820 })
821 })
822 .unwrap_or(false)
823 }
824
825 /// Checks if the symmetry element is spatially a mirror plane and has the specified
826 /// antiunitary attribute.
827 ///
828 /// More specifically, this checks if the combination $`\hat{\gamma} \hat{C}_n^k`$ is a
829 /// reflection, and if $`\hat{\alpha}`$ is as specified by `au`.
830 ///
831 /// # Arguments
832 ///
833 /// * `au` - An `Option` for the desired antiunitary kind.
834 ///
835 /// # Returns
836 ///
837 /// A boolean indicating if this symmetry element is spatially a mirror plane and has the
838 /// specified antiunitary attribute.
839 #[must_use]
840 pub fn is_o3_mirror_plane(&self, au: Option<AntiunitaryKind>) -> bool {
841 match self.kind {
842 SymmetryElementKind::ImproperMirrorPlane(sig_au) => {
843 sig_au == au
844 && self
845 .proper_fraction
846 .map(|frac| frac.is_zero())
847 .or_else(|| {
848 self.proper_angle.map(|proper_angle| {
849 approx::abs_diff_eq!(proper_angle, 0.0, epsilon = self.threshold,)
850 })
851 })
852 .unwrap_or(false)
853 }
854 SymmetryElementKind::ImproperInversionCentre(inv_au) => {
855 inv_au == au
856 && self
857 .proper_fraction
858 .map(|frac| frac == F::new(1u32, 2u32))
859 .or_else(|| {
860 self.proper_angle.map(|proper_angle| {
861 approx::abs_diff_eq!(
862 proper_angle,
863 std::f64::consts::PI,
864 epsilon = self.threshold,
865 )
866 })
867 })
868 .unwrap_or(false)
869 }
870 _ => false,
871 }
872 }
873
874 /// Returns the full symbol for this symmetry element, which does not classify certain
875 /// improper rotation axes into inversion centres or mirror planes, but which does simplify
876 /// the power/order ratio, and which displays only the absolute value of the power since
877 /// symmetry elements do not distinguish senses of rotations since rotations of oposite
878 /// directions are inverses of each other, both of which must exist in the group.
879 ///
880 /// Some additional symbols that can be unconventional include:
881 ///
882 /// * `θ`: time reversal,
883 /// * `(Σ)`: the spatial part is in homotopy class 0 of $`\mathsf{SU}'(2)`$,
884 /// * `(QΣ)`: the spatial part is in homotopy class 1 of $`\mathsf{SU}'(2)`$.
885 ///
886 /// See [`RotationGroup`] for further information.
887 ///
888 /// # Returns
889 ///
890 /// The full symbol for this symmetry element.
891 #[must_use]
892 pub fn get_full_symbol(&self) -> String {
893 let tr_sym = if self.kind.contains_time_reversal() {
894 "θ·"
895 } else {
896 ""
897 };
898 let main_symbol: String = match self.kind {
899 SymmetryElementKind::Proper(_) => {
900 format!("{tr_sym}C")
901 }
902 SymmetryElementKind::ImproperMirrorPlane(_) => {
903 if *self.raw_proper_order() != ElementOrder::Inf
904 && (*self
905 .proper_fraction
906 .expect("No proper fractions found for a finite-order element.")
907 .numer()
908 .expect("Unable to extract the numerator of the proper fraction.")
909 == 1
910 || self
911 .proper_fraction
912 .expect("No proper fractions found for a finite-order element.")
913 .is_zero())
914 {
915 format!("{tr_sym}S")
916 } else {
917 format!("{tr_sym}σC")
918 }
919 }
920 SymmetryElementKind::ImproperInversionCentre(_) => {
921 if *self.raw_proper_order() != ElementOrder::Inf
922 && (*self
923 .proper_fraction
924 .expect("No proper fractions found for a finite-order element.")
925 .numer()
926 .expect("Unable to extract the numerator of the proper fraction.")
927 == 1
928 || self
929 .proper_fraction
930 .expect("No proper fractions found for a finite-order element.")
931 .is_zero())
932 {
933 format!("{tr_sym}Ṡ")
934 } else {
935 format!("{tr_sym}iC")
936 }
937 }
938 };
939
940 let su2_sym = if self.is_su2_class_1() {
941 "(QΣ)"
942 } else if self.is_su2() {
943 "(Σ)"
944 } else {
945 ""
946 };
947
948 if let Some(proper_fraction) = self.proper_fraction {
949 if proper_fraction.is_zero() {
950 format!("{main_symbol}1{su2_sym}")
951 } else {
952 let proper_order = proper_fraction
953 .denom()
954 .expect("Unable to extract the denominator of the proper fraction.")
955 .to_string();
956 let proper_power = {
957 let pow = *proper_fraction
958 .numer()
959 .expect("Unable to extract the numerator of the proper fraction.");
960 if pow > 1 {
961 format!("^{pow}")
962 } else {
963 String::new()
964 }
965 };
966 format!("{main_symbol}{proper_order}{proper_power}{su2_sym}")
967 }
968 } else {
969 assert_eq!(*self.raw_proper_order(), ElementOrder::Inf);
970 let proper_angle = if let Some(proper_angle) = self.proper_angle {
971 format!("({:+.3})", proper_angle.abs())
972 } else {
973 String::new()
974 };
975 format!(
976 "{main_symbol}{}{proper_angle}{su2_sym}",
977 self.raw_proper_order()
978 )
979 }
980 }
981
982 /// Returns the simplified symbol for this symmetry element, which classifies special symmetry
983 /// elements (identity, inversion centre, time reversal, mirror planes), and which simplifies
984 /// the power/order ratio and displays only the absolute value of the power since symmetry
985 /// elements do not distinguish senses of rotations, as rotations of oposite directions are
986 /// inverses of each other, both of which must exist in the group.
987 ///
988 /// # Returns
989 ///
990 /// The simplified symbol for this symmetry element.
991 #[must_use]
992 pub fn get_simplified_symbol(&self) -> String {
993 let (main_symbol, needs_power) = match self.kind {
994 SymmetryElementKind::Proper(au) => {
995 if self.is_o3_identity(au) {
996 match au {
997 None => ("E".to_owned(), false),
998 Some(AntiunitaryKind::TimeReversal) => ("θ".to_owned(), false),
999 Some(AntiunitaryKind::ComplexConjugation) => ("K".to_owned(), false),
1000 }
1001 } else {
1002 let au_sym = match au {
1003 None => "",
1004 Some(AntiunitaryKind::TimeReversal) => "θ·",
1005 Some(AntiunitaryKind::ComplexConjugation) => "K·",
1006 };
1007 (format!("{au_sym}C"), true)
1008 }
1009 }
1010 SymmetryElementKind::ImproperMirrorPlane(au) => {
1011 let au_sym = match au {
1012 None => "",
1013 Some(AntiunitaryKind::TimeReversal) => "θ·",
1014 Some(AntiunitaryKind::ComplexConjugation) => "K·",
1015 };
1016 if self.is_o3_mirror_plane(au) {
1017 (format!("{au_sym}σ"), false)
1018 } else if self.is_o3_inversion_centre(au) {
1019 (format!("{au_sym}i"), false)
1020 } else if *self.raw_proper_order() == ElementOrder::Inf
1021 || *self
1022 .proper_fraction
1023 .expect("No proper fractions found for a finite-order element.")
1024 .numer()
1025 .expect("Unable to extract the numerator of the proper fraction.")
1026 == 1
1027 {
1028 (format!("{au_sym}S"), false)
1029 } else {
1030 (format!("{au_sym}σC"), true)
1031 }
1032 }
1033 SymmetryElementKind::ImproperInversionCentre(au) => {
1034 let au_sym = match au {
1035 None => "",
1036 Some(AntiunitaryKind::TimeReversal) => "θ·",
1037 Some(AntiunitaryKind::ComplexConjugation) => "K·",
1038 };
1039 if self.is_o3_mirror_plane(au) {
1040 (format!("{au_sym}σ"), false)
1041 } else if self.is_o3_inversion_centre(au) {
1042 (format!("{au_sym}i"), false)
1043 } else if *self.raw_proper_order() == ElementOrder::Inf
1044 || *self
1045 .proper_fraction
1046 .expect("No proper fractions found for a finite-order element.")
1047 .numer()
1048 .expect("Unable to extract the numerator of the proper fraction.")
1049 == 1
1050 {
1051 (format!("{au_sym}Ṡ"), false)
1052 } else {
1053 (format!("{au_sym}iC"), true)
1054 }
1055 }
1056 };
1057
1058 let su2_sym = if self.is_su2_class_1() {
1059 "(QΣ)"
1060 } else if self.is_su2() {
1061 "(Σ)"
1062 } else {
1063 ""
1064 };
1065
1066 if let Some(proper_fraction) = self.proper_fraction {
1067 let au = self.contains_antiunitary();
1068 let proper_order = if self.is_o3_identity(au)
1069 || self.is_o3_inversion_centre(au)
1070 || self.is_o3_mirror_plane(au)
1071 {
1072 String::new()
1073 } else {
1074 proper_fraction
1075 .denom()
1076 .expect("Unable to extract the denominator of the proper fraction.")
1077 .to_string()
1078 };
1079
1080 let proper_power = if needs_power {
1081 let pow = *proper_fraction
1082 .numer()
1083 .expect("Unable to extract the numerator of the proper fraction.");
1084 if pow > 1 {
1085 format!("^{pow}")
1086 } else {
1087 String::new()
1088 }
1089 } else {
1090 String::new()
1091 };
1092 format!(
1093 "{main_symbol}{}{proper_order}{proper_power}{}{su2_sym}",
1094 self.additional_superscript, self.additional_subscript
1095 )
1096 } else {
1097 assert_eq!(*self.raw_proper_order(), ElementOrder::Inf);
1098 let proper_angle = if let Some(proper_angle) = self.proper_angle {
1099 format!("({:+.3})", proper_angle.abs())
1100 } else {
1101 String::new()
1102 };
1103 format!(
1104 "{main_symbol}{}{}{proper_angle}{}{su2_sym}",
1105 self.additional_superscript,
1106 *self.raw_proper_order(),
1107 self.additional_subscript
1108 )
1109 }
1110 }
1111
1112 /// Returns the simplified symbol for this symmetry element, which classifies special symmetry
1113 /// elements (identity, inversion centre, mirror planes), and which simplifies the power/order
1114 /// ratio and displays only the absolute value of the power since symmetry elements do not
1115 /// distinguish senses of rotations. Rotations of oposite directions are inverses of each
1116 /// other, both of which must exist in the group.
1117 ///
1118 /// # Returns
1119 ///
1120 /// The simplified symbol for this symmetry element.
1121 #[must_use]
1122 pub fn get_simplified_symbol_signed_power(&self) -> String {
1123 let (main_symbol, needs_power) = match self.kind {
1124 SymmetryElementKind::Proper(au) => {
1125 if self.is_o3_identity(au) {
1126 match au {
1127 None => ("E".to_owned(), false),
1128 Some(AntiunitaryKind::TimeReversal) => ("θ".to_owned(), false),
1129 Some(AntiunitaryKind::ComplexConjugation) => ("K".to_owned(), false),
1130 }
1131 } else {
1132 let au_sym = match au {
1133 None => "",
1134 Some(AntiunitaryKind::TimeReversal) => "θ·",
1135 Some(AntiunitaryKind::ComplexConjugation) => "K·",
1136 };
1137 (format!("{au_sym}C"), true)
1138 }
1139 }
1140 SymmetryElementKind::ImproperMirrorPlane(au) => {
1141 let au_sym = match au {
1142 None => "",
1143 Some(AntiunitaryKind::TimeReversal) => "θ·",
1144 Some(AntiunitaryKind::ComplexConjugation) => "K·",
1145 };
1146 if self.is_o3_mirror_plane(au) {
1147 (format!("{au_sym}σ"), false)
1148 } else if self.is_o3_inversion_centre(au) {
1149 (format!("{au_sym}i"), false)
1150 } else if *self.raw_proper_order() == ElementOrder::Inf
1151 || *self
1152 .proper_fraction
1153 .expect("No proper fractions found for a finite-order element.")
1154 .numer()
1155 .expect("Unable to extract the numerator of the proper fraction.")
1156 == 1
1157 {
1158 (format!("{au_sym}S"), true)
1159 } else {
1160 (format!("{au_sym}σC"), true)
1161 }
1162 }
1163 SymmetryElementKind::ImproperInversionCentre(au) => {
1164 let au_sym = match au {
1165 None => "",
1166 Some(AntiunitaryKind::TimeReversal) => "θ·",
1167 Some(AntiunitaryKind::ComplexConjugation) => "K·",
1168 };
1169 if self.is_o3_mirror_plane(au) {
1170 (format!("{au_sym}σ"), false)
1171 } else if self.is_o3_inversion_centre(au) {
1172 (format!("{au_sym}i"), false)
1173 } else if *self.raw_proper_order() == ElementOrder::Inf
1174 || *self
1175 .proper_fraction
1176 .expect("No proper fractions found for a finite-order element.")
1177 .numer()
1178 .expect("Unable to extract the numerator of the proper fraction.")
1179 == 1
1180 {
1181 (format!("{au_sym}Ṡ"), true)
1182 } else {
1183 (format!("{au_sym}iC"), true)
1184 }
1185 }
1186 };
1187
1188 let su2_sym = if self.is_su2_class_1() {
1189 "(QΣ)"
1190 } else if self.is_su2() {
1191 "(Σ)"
1192 } else {
1193 ""
1194 };
1195
1196 if let Some(proper_fraction) = self.proper_fraction {
1197 let au = self.contains_antiunitary();
1198 let proper_order = if self.is_o3_identity(au)
1199 || self.is_o3_inversion_centre(au)
1200 || self.is_o3_mirror_plane(au)
1201 {
1202 String::new()
1203 } else {
1204 proper_fraction
1205 .denom()
1206 .expect("Unable to extract the denominator of the proper fraction.")
1207 .to_string()
1208 };
1209
1210 let proper_power = if needs_power {
1211 let pow = *proper_fraction
1212 .numer()
1213 .expect("Unable to extract the numerator of the proper fraction.");
1214 if !geometry::check_standard_positive_pole(
1215 &self.proper_rotation_pole(),
1216 self.threshold,
1217 ) {
1218 format!("^(-{pow})")
1219 } else if pow > 1 {
1220 format!("^{pow}")
1221 } else {
1222 String::new()
1223 }
1224 } else {
1225 String::new()
1226 };
1227 format!(
1228 "{main_symbol}{}{proper_order}{proper_power}{}{su2_sym}",
1229 self.additional_superscript, self.additional_subscript
1230 )
1231 } else {
1232 assert_eq!(*self.raw_proper_order(), ElementOrder::Inf);
1233 let proper_angle = if let Some(proper_angle) = self.proper_angle {
1234 format!("({:+.3})", proper_angle.abs())
1235 } else {
1236 String::new()
1237 };
1238 format!(
1239 "{main_symbol}{}{}{proper_angle}{}{su2_sym}",
1240 self.additional_superscript,
1241 *self.raw_proper_order(),
1242 self.additional_subscript
1243 )
1244 }
1245 }
1246
1247 /// Returns a copy of the current improper symmetry element that has been converted to the
1248 /// required improper kind. For $`\mathsf{SU}'(2)`$ elements, the conversion will be carried
1249 /// out in the same homotopy class. The antiunitary part will be kept unchanged.
1250 ///
1251 /// To convert between the two improper kinds, we essentially seek integers
1252 /// $`n, n' \in \mathbb{N}_{+}`$ and $`k \in \mathbb{Z}/n\mathbb{Z}`$,
1253 /// $`k' \in \mathbb{Z}/n'\mathbb{Z}`$, such that
1254 ///
1255 /// ```math
1256 /// \sigma C_n^k = i C_{n'}^{k'},
1257 /// ```
1258 ///
1259 /// where the axes of all involved elements are parallel. By noting that
1260 /// $`\sigma = i C_2`$ and that $`k`$ and $`k'`$ must have opposite signs, we can easily show
1261 /// that, for $`k \ge 0, k' < 0`$,
1262 ///
1263 /// ```math
1264 /// \begin{aligned}
1265 /// n' &= \frac{2n}{\operatorname{gcd}(2n, n - 2k)},\\
1266 /// k' &= -\frac{n - 2k}{\operatorname{gcd}(2n, n - 2k)},
1267 /// \end{aligned}
1268 /// ```
1269 ///
1270 /// whereas for $`k < 0, k' \ge 0`$,
1271 ///
1272 /// ```math
1273 /// \begin{aligned}
1274 /// n' &= \frac{2n}{\operatorname{gcd}(2n, n + 2k)},\\
1275 /// k' &= \frac{n + 2k}{\operatorname{gcd}(2n, n + 2k)}.
1276 /// \end{aligned}
1277 /// ```
1278 ///
1279 /// The above relations are self-inverse. It can be further shown that
1280 /// $`\operatorname{gcd}(n', k') = 1`$. Hence, for symmetry *element* conversions, we can simply
1281 /// take $`k' = 1`$. This is because a symmetry element plays the role of a generator, and the
1282 /// coprimality of $`n'`$ and $`k'`$ means that $`i C_{n'}^{1}`$ is as valid a generator as
1283 /// $`i C_{n'}^{k'}`$.
1284 ///
1285 /// # Arguments
1286 ///
1287 /// * `improper_kind` - The improper kind to which `self` is to be converted. There is no need
1288 /// to make sure the time reversal specification in `improper_kind` matches that of `self` as
1289 /// the conversion will take care of this.
1290 /// * `preserves_power` - Boolean indicating if the proper rotation power $`k'`$ should be
1291 /// preserved or should be set to $`1`$.
1292 ///
1293 /// # Returns
1294 ///
1295 /// A copy of the current improper symmetry element that has been converted.
1296 ///
1297 /// # Panics
1298 ///
1299 /// Panics when `self` is not an improper element, or when `improper_kind` is not one of the
1300 /// improper variants.
1301 #[must_use]
1302 pub fn convert_to_improper_kind(
1303 &self,
1304 improper_kind: &SymmetryElementKind,
1305 preserves_power: bool,
1306 ) -> Self {
1307 let au = self.contains_antiunitary();
1308 assert!(
1309 !self.is_o3_proper(au),
1310 "Only improper elements can be converted."
1311 );
1312 let improper_kind = improper_kind.to_antiunitary(au);
1313 assert!(
1314 !matches!(improper_kind, SymmetryElementKind::Proper(_)),
1315 "`improper_kind` must be one of the improper variants."
1316 );
1317
1318 // self.kind and improper_kind must now have the same antiunitary part.
1319
1320 if self.kind == improper_kind {
1321 return self.clone();
1322 }
1323
1324 let (dest_order, dest_proper_power) = match *self.raw_proper_order() {
1325 ElementOrder::Int(_) => {
1326 let proper_fraction = self
1327 .proper_fraction
1328 .expect("No proper fractions found for an element with integer order.");
1329 let n = *proper_fraction.denom().unwrap();
1330 let k = if proper_fraction.is_sign_negative() {
1331 -i32::try_from(*proper_fraction.numer().unwrap_or_else(|| {
1332 panic!("Unable to retrieve the numerator of {proper_fraction:?}.")
1333 }))
1334 .expect("Unable to convert the numerator of the proper fraction to `i32`.")
1335 } else {
1336 i32::try_from(*proper_fraction.numer().unwrap_or_else(|| {
1337 panic!("Unable to retrieve the numerator of {proper_fraction:?}.")
1338 }))
1339 .expect("Unable to convert the numerator of the proper fraction to `i32`.")
1340 };
1341
1342 if k >= 0 {
1343 // k >= 0, k2 < 0
1344 let n_m_2k = n
1345 .checked_sub(2 * k.unsigned_abs())
1346 .expect("The value of `n - 2k` is negative.");
1347 let n2 = ElementOrder::Int(2 * n / (gcd(2 * n, n_m_2k)));
1348 let k2: i32 = if preserves_power {
1349 -i32::try_from(n_m_2k / gcd(2 * n, n_m_2k))
1350 .expect("Unable to convert `k'` to `i32`.")
1351 } else {
1352 1
1353 };
1354 (n2, k2)
1355 } else {
1356 // k < 0, k2 >= 0
1357 let n_p_2k = n
1358 .checked_sub(2 * k.unsigned_abs())
1359 .expect("The value of `n + 2k` is negative.");
1360 let n2 = ElementOrder::Int(2 * n / (gcd(2 * n, n_p_2k)));
1361 let k2: i32 = if preserves_power {
1362 i32::try_from(n_p_2k / (gcd(2 * n, n_p_2k)))
1363 .expect("Unable to convert `k'` to `i32`.")
1364 } else {
1365 1
1366 };
1367 (n2, k2)
1368 }
1369 }
1370 ElementOrder::Inf => (ElementOrder::Inf, 1),
1371 };
1372
1373 match dest_order {
1374 ElementOrder::Int(_) => Self::builder()
1375 .threshold(self.threshold)
1376 .proper_order(dest_order)
1377 .proper_power(dest_proper_power)
1378 .raw_axis(self.raw_axis)
1379 .kind(improper_kind)
1380 .rotation_group(self.rotation_group.clone())
1381 .generator(self.generator)
1382 .additional_superscript(self.additional_superscript.clone())
1383 .additional_subscript(self.additional_subscript.clone())
1384 .build()
1385 .expect("Unable to construct a symmetry element."),
1386
1387 ElementOrder::Inf => {
1388 if let Some(ang) = self.proper_angle {
1389 Self::builder()
1390 .threshold(self.threshold)
1391 .proper_order(dest_order)
1392 .proper_power(dest_proper_power)
1393 .proper_angle(-std::f64::consts::PI + ang)
1394 .raw_axis(self.raw_axis)
1395 .kind(improper_kind)
1396 .rotation_group(self.rotation_group.clone())
1397 .generator(self.generator)
1398 .additional_superscript(self.additional_superscript.clone())
1399 .additional_subscript(self.additional_subscript.clone())
1400 .build()
1401 .expect("Unable to construct a symmetry element.")
1402 } else {
1403 Self::builder()
1404 .threshold(self.threshold)
1405 .proper_order(dest_order)
1406 .proper_power(dest_proper_power)
1407 .raw_axis(self.raw_axis)
1408 .kind(improper_kind)
1409 .rotation_group(self.rotation_group.clone())
1410 .generator(self.generator)
1411 .additional_superscript(self.additional_superscript.clone())
1412 .additional_subscript(self.additional_subscript.clone())
1413 .build()
1414 .expect("Unable to construct a symmetry element.")
1415 }
1416 }
1417 }
1418 }
1419
1420 /// Converts the symmetry element to one with the desired time-reversal property.
1421 ///
1422 /// # Arguments
1423 ///
1424 /// * `tr` - A boolean indicating if time reversal is to be included.
1425 ///
1426 /// # Returns
1427 ///
1428 /// A new symmetry element with or without time reversal as indicated by `tr`.
1429 pub fn to_tr(&self, tr: bool) -> Self {
1430 let mut c_self = self.clone();
1431 c_self.kind = c_self.kind.to_tr(tr);
1432 c_self
1433 }
1434
1435 /// Convert the proper rotation of the current element to one in $`\mathsf{SU}(2)`$.
1436 ///
1437 /// # Arguments
1438 ///
1439 /// * `normal` - A boolean indicating whether the resultant $`\mathsf{SU}(2)`$ proper rotation
1440 /// is of homotopy class 0 (`true`) or 1 (`false`) when connected to the identity.
1441 ///
1442 /// # Returns
1443 ///
1444 /// A symmetry element in $`\mathsf{SU}(2)`$, or `None` if the current symmetry element
1445 /// is already in $`\mathsf{SU}(2)`$.
1446 pub fn to_su2(&self, normal: bool) -> Option<Self> {
1447 if self.is_su2() {
1448 None
1449 } else {
1450 let mut element = self.clone();
1451 element.rotation_group = RotationGroup::SU2(normal);
1452 Some(element)
1453 }
1454 }
1455
1456 /// The closeness of the symmetry element's axis to one of the three space-fixed Cartesian axes.
1457 ///
1458 /// # Returns
1459 ///
1460 /// A tuple of two values:
1461 /// - A value $`\gamma \in [0, 1-1/\sqrt{3}]`$ indicating how close the axis is to one of the
1462 /// three Cartesian axes. The closer $`\gamma`$ is to $`0`$, the closer the alignment.
1463 /// - An index for the closest axis: `0` for $`z`$, `1` for $`y`$, `2` for $`x`$.
1464 ///
1465 /// # Panics
1466 ///
1467 /// Panics when $`\gamma`$ is outside the required closed interval $`[0, 1-1/\sqrt{3}]`$ by
1468 /// more than the threshold value in `self`.
1469 #[must_use]
1470 pub fn closeness_to_cartesian_axes(&self) -> (f64, usize) {
1471 let pos_axis = self.standard_positive_axis();
1472 let rev_pos_axis = Vector3::new(pos_axis[2], pos_axis[1], pos_axis[0]);
1473 let (amax_arg, amax_val) = rev_pos_axis.abs().argmax();
1474 let axis_closeness = 1.0 - amax_val;
1475 let thresh = self.threshold;
1476 assert!(
1477 -thresh <= axis_closeness && axis_closeness <= (1.0 - 1.0 / 3.0f64.sqrt() + thresh)
1478 );
1479
1480 // closest axis: 0 for z, 1 for y, 2 for x
1481 // This is so that z axes are preferred.
1482 let closest_axis = amax_arg;
1483 (axis_closeness, closest_axis)
1484 }
1485}
1486
1487impl fmt::Debug for SymmetryElement {
1488 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1489 let signed_axis = self.signed_axis();
1490 write!(
1491 f,
1492 "{}({:+.3}, {:+.3}, {:+.3})",
1493 self.get_full_symbol(),
1494 signed_axis[0] + 0.0,
1495 signed_axis[1] + 0.0,
1496 signed_axis[2] + 0.0
1497 )
1498 }
1499}
1500
1501impl fmt::Display for SymmetryElement {
1502 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1503 let au = self.contains_antiunitary();
1504 if self.is_o3_identity(au) || self.is_o3_inversion_centre(au) {
1505 write!(f, "{}", self.get_simplified_symbol())
1506 } else {
1507 let signed_axis = self.signed_axis();
1508 write!(
1509 f,
1510 "{}({:+.3}, {:+.3}, {:+.3})",
1511 self.get_simplified_symbol(),
1512 signed_axis[0] + 0.0,
1513 signed_axis[1] + 0.0,
1514 signed_axis[2] + 0.0
1515 )
1516 }
1517 }
1518}
1519
1520impl PartialEq for SymmetryElement {
1521 /// Two symmetry elements are equal if and only if the following conditions are all satisfied:
1522 ///
1523 /// * they are both in the same rotation group and belong to the same homotopy class;
1524 /// * they are both proper or improper;
1525 /// * they both have the same antiunitary properties;
1526 /// * their axes are either parallel or anti-parallel;
1527 /// * their proper rotation angles have equal absolute values.
1528 ///
1529 /// For improper elements, proper rotation angles are taken in the inversion centre convention.
1530 ///
1531 /// Thus, symmetry element equality is less strict than symmetry operation equality. This is so
1532 /// that parallel or anti-parallel symmetry elements with the same spatial and time-reversal
1533 /// parities and angle of rotation are deemed identical, thus facilitating symmetry detection
1534 /// where one does not yet care much about directions of rotations.
1535 #[allow(clippy::too_many_lines)]
1536 fn eq(&self, other: &Self) -> bool {
1537 if self.rotation_group != other.rotation_group {
1538 // Different rotation groups or homotopy classes.
1539 return false;
1540 }
1541
1542 if self.contains_antiunitary() != other.contains_antiunitary() {
1543 // Different anti-unitary kinds.
1544 return false;
1545 }
1546
1547 let au = self.contains_antiunitary();
1548
1549 if self.is_o3_proper(au) != other.is_o3_proper(au) {
1550 // Different spatial parities.
1551 return false;
1552 }
1553
1554 if self.is_o3_identity(au) && other.is_o3_identity(au) {
1555 // Both are spatial identity.
1556 // assert_eq!(
1557 // misc::calculate_hash(self),
1558 // misc::calculate_hash(other),
1559 // "{self} and {other} have unequal hashes."
1560 // );
1561 // return true;
1562 return misc::calculate_hash(self) == misc::calculate_hash(other);
1563 }
1564
1565 if self.is_o3_inversion_centre(au) && other.is_o3_inversion_centre(au) {
1566 // Both are spatial inversion centre.
1567 // assert_eq!(
1568 // misc::calculate_hash(self),
1569 // misc::calculate_hash(other),
1570 // "{self} and {other} have unequal hashes."
1571 // );
1572 // return true;
1573 return misc::calculate_hash(self) == misc::calculate_hash(other);
1574 }
1575
1576 let thresh = (self.threshold * other.threshold).sqrt();
1577
1578 let result = if self.is_o3_proper(au) {
1579 // Proper.
1580
1581 // Parallel or anti-parallel axes.
1582 let similar_poles = approx::relative_eq!(
1583 geometry::get_standard_positive_pole(&self.raw_axis, thresh),
1584 geometry::get_standard_positive_pole(&other.raw_axis, thresh),
1585 epsilon = thresh,
1586 max_relative = thresh
1587 );
1588
1589 // Same angle of rotation (irrespective of signs).
1590 let similar_angles = match (*self.raw_proper_order(), *other.raw_proper_order()) {
1591 (ElementOrder::Inf, ElementOrder::Inf) => {
1592 match (self.proper_angle, other.proper_angle) {
1593 (Some(s_angle), Some(o_angle)) => {
1594 approx::relative_eq!(
1595 s_angle.abs(),
1596 o_angle.abs(),
1597 epsilon = thresh,
1598 max_relative = thresh
1599 )
1600 }
1601 (None, None) => similar_poles,
1602 _ => false,
1603 }
1604 }
1605 (ElementOrder::Int(_), ElementOrder::Int(_)) => {
1606 let c_proper_fraction = self
1607 .proper_fraction
1608 .expect("Proper fraction for `self` not found.");
1609 let o_proper_fraction = other
1610 .proper_fraction
1611 .expect("Proper fraction for `other` not found.");
1612 c_proper_fraction.abs() == o_proper_fraction.abs()
1613 }
1614 _ => false,
1615 };
1616
1617 similar_poles && similar_angles
1618 } else {
1619 // Improper => convert to inversion-centre convention.
1620 let inv_au = SymmetryElementKind::ImproperInversionCentre(au);
1621 let c_self = self.convert_to_improper_kind(&inv_au, false);
1622 let c_other = other.convert_to_improper_kind(&inv_au, false);
1623
1624 // Parallel or anti-parallel axes.
1625 let similar_poles = approx::relative_eq!(
1626 geometry::get_standard_positive_pole(&c_self.raw_axis, thresh),
1627 geometry::get_standard_positive_pole(&c_other.raw_axis, thresh),
1628 epsilon = thresh,
1629 max_relative = thresh
1630 );
1631
1632 // Same angle of rotation (irrespective of signs).
1633 let similar_angles = match (*c_self.raw_proper_order(), *c_other.raw_proper_order()) {
1634 (ElementOrder::Inf, ElementOrder::Inf) => {
1635 match (c_self.proper_angle, c_other.proper_angle) {
1636 (Some(s_angle), Some(o_angle)) => {
1637 approx::relative_eq!(
1638 s_angle.abs(),
1639 o_angle.abs(),
1640 epsilon = thresh,
1641 max_relative = thresh
1642 )
1643 }
1644 (None, None) => similar_poles,
1645 _ => false,
1646 }
1647 }
1648 (ElementOrder::Int(_), ElementOrder::Int(_)) => {
1649 let c_proper_fraction = c_self
1650 .proper_fraction
1651 .expect("Proper fraction for `c_self` not found.");
1652 let o_proper_fraction = c_other
1653 .proper_fraction
1654 .expect("Proper fraction for `c_other` not found.");
1655 c_proper_fraction.abs() == o_proper_fraction.abs()
1656 }
1657 _ => false,
1658 };
1659
1660 similar_poles && similar_angles
1661 };
1662
1663 // if result {
1664 // assert_eq!(
1665 // misc::calculate_hash(self),
1666 // misc::calculate_hash(other),
1667 // "`{self}` and `{other}` have unequal hashes."
1668 // );
1669 // }
1670 result && (misc::calculate_hash(self) == misc::calculate_hash(other))
1671 }
1672}
1673
1674impl Eq for SymmetryElement {}
1675
1676impl Hash for SymmetryElement {
1677 fn hash<H: Hasher>(&self, state: &mut H) {
1678 self.rotation_group.hash(state);
1679
1680 let au = self.contains_antiunitary();
1681 au.hash(state);
1682
1683 self.is_o3_proper(au).hash(state);
1684
1685 if self.is_o3_identity(au) || self.is_o3_inversion_centre(au) {
1686 true.hash(state);
1687 } else if self.kind == SymmetryElementKind::ImproperMirrorPlane(au) {
1688 let c_self = self
1689 .convert_to_improper_kind(&SymmetryElementKind::ImproperInversionCentre(au), false);
1690 let pole = geometry::get_standard_positive_pole(&c_self.raw_axis, c_self.threshold);
1691 pole[0]
1692 .round_factor(self.threshold)
1693 .integer_decode()
1694 .hash(state);
1695 pole[1]
1696 .round_factor(self.threshold)
1697 .integer_decode()
1698 .hash(state);
1699 pole[2]
1700 .round_factor(self.threshold)
1701 .integer_decode()
1702 .hash(state);
1703 if let ElementOrder::Inf = *c_self.raw_proper_order() {
1704 if let Some(angle) = c_self.proper_angle {
1705 angle
1706 .abs()
1707 .round_factor(self.threshold)
1708 .integer_decode()
1709 .hash(state);
1710 } else {
1711 0.hash(state);
1712 }
1713 } else {
1714 c_self
1715 .proper_fraction
1716 .expect("No proper fractions for `c_self` found.")
1717 .abs()
1718 .hash(state);
1719 };
1720 } else {
1721 let pole = geometry::get_standard_positive_pole(&self.raw_axis, self.threshold);
1722 pole[0]
1723 .round_factor(self.threshold)
1724 .integer_decode()
1725 .hash(state);
1726 pole[1]
1727 .round_factor(self.threshold)
1728 .integer_decode()
1729 .hash(state);
1730 pole[2]
1731 .round_factor(self.threshold)
1732 .integer_decode()
1733 .hash(state);
1734 if let ElementOrder::Inf = *self.raw_proper_order() {
1735 if let Some(angle) = self.proper_angle {
1736 angle
1737 .abs()
1738 .round_factor(self.threshold)
1739 .integer_decode()
1740 .hash(state);
1741 } else {
1742 0.hash(state);
1743 }
1744 } else {
1745 self.proper_fraction
1746 .expect("No proper fractions for `self` found.")
1747 .abs()
1748 .hash(state);
1749 };
1750 };
1751 }
1752}
1753
1754/// Time-reversal antiunitary kind.
1755pub const TR: AntiunitaryKind = AntiunitaryKind::TimeReversal;
1756
1757/// Proper rotation symmetry element kind.
1758pub const ROT: SymmetryElementKind = SymmetryElementKind::Proper(None);
1759
1760/// Improper symmetry element kind in the mirror-plane convention.
1761pub const SIG: SymmetryElementKind = SymmetryElementKind::ImproperMirrorPlane(None);
1762
1763/// Improper symmetry element kind in the inversion-centre convention.
1764pub const INV: SymmetryElementKind = SymmetryElementKind::ImproperInversionCentre(None);
1765
1766/// Time-reversed proper rotation symmetry element kind.
1767pub const TRROT: SymmetryElementKind = SymmetryElementKind::Proper(Some(TR));
1768
1769/// Time-reversed improper symmetry element kind in the mirror-plane convention.
1770pub const TRSIG: SymmetryElementKind = SymmetryElementKind::ImproperMirrorPlane(Some(TR));
1771
1772/// Time-reversed improper symmetry element kind in the inversion-centre convention.
1773pub const TRINV: SymmetryElementKind = SymmetryElementKind::ImproperInversionCentre(Some(TR));
1774
1775/// Rotation group $`\mathsf{SO}(3)`$.
1776pub const SO3: RotationGroup = RotationGroup::SO3;
1777
1778/// Rotation group $`\mathsf{SU}(2)`$, homotopy path of class 0.
1779pub const SU2_0: RotationGroup = RotationGroup::SU2(true);
1780
1781/// Rotation group $`\mathsf{SU}(2)`$, homotopy path of class 1.
1782pub const SU2_1: RotationGroup = RotationGroup::SU2(false);