treams_core/lattice/
wave.rs

1//! The family of a sum, the scalar wave it adds up, and its angular factors.
2//!
3//! Upstream: the `lsumsw*` (spherical) and `lsumcw*` (cylindrical) families of
4//! `treams.lattice`.
5
6use crate::{Complex, Error, MAX_DEGREE, Result, numerics::Jet, special};
7
8/// The family of a lattice sum: spherical or cylindrical.
9///
10/// Each variant is the scalar outgoing wave the sum adds up. `vectorwaves::Family`
11/// selects vector waves instead; code outside each module names both with their module.
12///
13/// Upstream: the `lsumsw*` and `lsumcw*` functions of `treams.lattice`.
14#[derive(Clone, Copy, Debug)]
15pub enum Family {
16    /// `h_l^(1)(kr) Y_lm(theta, phi)`.
17    Spherical {
18        /// Degree.
19        l: i32,
20        /// Order.
21        m: i32,
22    },
23    /// `H_m^(1)(kr) exp(i m phi)`.
24    Cylindrical {
25        /// Signed order.
26        m: i32,
27    },
28}
29
30impl Family {
31    /// Checks the labels, and for cylindrical sums the lattice dimension: spherical sums
32    /// take `0 <= l <= MAX_DEGREE` and `|m| <= l`, cylindrical sums take
33    /// `|m| <= MAX_DEGREE` and `dim <= 2`.
34    ///
35    /// The bound keeps the Kambe integrals of the sums within the orders that
36    /// `special::MAX_LABEL = 2 MAX_DEGREE + 4` admits: spherical real-space terms and
37    /// their derivatives take orders `2 l` and `2 l + 2`, and the reciprocal terms of 1D
38    /// spherical sums orders down to `-(2 l + 3)`.
39    pub(super) fn validate(self, dim: usize) -> Result<()> {
40        let valid = match self {
41            Self::Spherical { l, m } => {
42                (0..=MAX_DEGREE).contains(&l) && m.unsigned_abs() <= l.unsigned_abs()
43            }
44            Self::Cylindrical { m } => (-MAX_DEGREE..=MAX_DEGREE).contains(&m) && dim <= 2,
45        };
46        if !valid {
47            return Err(Error::InvalidInput(
48                "invalid lattice wave degree or dimension".into(),
49            ));
50        }
51        Ok(())
52    }
53
54    /// The Cartesian axis of each lattice coordinate: z for 1D spherical sums, and x, y,
55    /// z in order otherwise (x for 1D cylindrical sums, xy for 2D sums).
56    pub(super) fn axes(self, dim: usize) -> [usize; 3] {
57        if matches!(self, Self::Spherical { .. }) && dim == 1 {
58            [2, 0, 1]
59        } else {
60            [0, 1, 2]
61        }
62    }
63}
64
65/// `x + i y`, or `x - i y` for negative `m`: the azimuthal factor whose `|m|`-th power
66/// carries the order of cylindrical waves and of the lattice polynomials.
67pub(super) fn azimuthal<const N: usize>(x: Jet<N>, y: Jet<N>, m: i32) -> Jet<N> {
68    x + Complex::i() * if m < 0 { -y } else { y }
69}
70
71/// `(-1)^m` for negative orders, relating `H_m` to `H_|m|`, and one otherwise.
72pub(super) const fn negative_order_sign(m: i32) -> f64 {
73    if m < 0 && m % 2 != 0 { -1.0 } else { 1.0 }
74}
75
76/// The solid harmonic `r^l P_l^m(cos theta) e^(i m phi)` of the jet vector `r`, with its
77/// derivatives by the chain rule through the Cartesian gradient. Times
78/// `harmonic_normalization(l, m)` it is `|r|^l Y_lm(r / |r|)`, a polynomial in `r` that
79/// stays regular at `r = 0`.
80pub(super) fn solid_jet<const N: usize>(l: i32, m: i32, r: [Jet<N>; 3]) -> Jet<N> {
81    let solid = special::solid::<false>(l, m, r.map(|r| r.value.re));
82    Jet {
83        value: solid.value,
84        derivative: std::array::from_fn(|i| {
85            (0..3).map(|j| solid.gradient[j] * r[j].derivative[i]).sum()
86        }),
87    }
88}
treams_core/lattice/wave.rs

treams_core/lattice/
wave.rs
