GitHub

General Utilities

Einstein.Utils.BarycentricInterpolation — Type
BarycentricInterpolation{TF<:AbstractFloat}(points::AbstractVector{TF}, weights::AbstractVector{TF})

A structure representing barycentric interpolation with precomputed weights.

Fields

  • points::AbstractVector{TF}: Vector of interpolation points (typically Chebyshev points)
  • weights::AbstractVector{TF}: Vector of barycentric weights

Methods

(itp::BarycentricInterpolation{TF})(values::AbstractVector{TFC}, x::TF) where {TF<:AbstractFloat,TFC<:Union{TF,Complex{TF}}}

Evaluate the interpolant at point x for function values.

Reference

source
Einstein.Utils.SWSFun — Type
SWSFun{TR<:AbstractFloat}()
(swsf::SWSFun{TR})(θ::TR)
coefficients(swsf::SWSFun)

Spherical-weighted spheroidal harmonics. The eigenvectors contain the $C$ coefficients in the equation:

\[{}_s S_{\ell m}(x; c)=\sum_{\ell^{\prime}=\ell_{\min }}^{\infty} C_{\ell^{\prime} \ell m}(c) {}_s S_{\ell^{\prime} m}(x; 0)\]

where $C$ is normalized by

\[\sum_{\ell^{\prime}=\ell_{\text {min }}}^{\ell_{\max }}\left|C_{\ell^{\prime} \ell m}(c)\right|^2=1\]

and the phase is chosen such that $C_{\ell^{\prime} \ell m}(c)$ is real for $\ell^{\prime}=\ell$ [CZ14]. The spin-weighted spheroidal harmonics are normalized such that

\[\int_0^\pi {}_s S_{\ell m}(x; c) {}_s S^*_{\ell m}(x; c) \sin(\theta) d\theta = 1\]

Arguments

  • s::Integer: spin
  • c::Complex{TR}: oblateness parameter
  • m::Integer: azimuthal number
  • l::Integer: angular number
  • l_max::Integer: maximum angular number
source
Einstein.Utils.barycentric_differentiation_matrix — Method
barycentric_differentiation_matrix(x::AbstractVector{TF}, w::AbstractVector{TF}, k::Integer=1, t::Union{AbstractVector{TF},Nothing}=nothing) where {TF<:AbstractFloat}

Compute the barycentric differentiation matrix.

Arguments

  • x::AbstractVector{TF} : Vector of interpolation points
  • w::AbstractVector{TF} : Barycentric weights of the interpolation points
  • k::Integer : Order of the derivative (default: 1)
  • t::AbstractVector{TF} : Vector of angles (default: nothing)

References:

source
Einstein.Utils.barycentric_weights — Method
barycentric_weights(x::AbstractVector{TF}) where {TF<:AbstractFloat}
barycentric_weights(x::AbstractRange{TF}) where {TF<:AbstractFloat}
barycentric_weights([TF=Float64], order::Integer) where {TF<:AbstractFloat}

Compute normalized barycentric weights for interpolation nodes. These weights are used in barycentric Lagrange interpolation. For Chebyshev-Gauss nodes or Chebyshev-Lobatto nodes, chebgaussbarycentric_weights and cheblobattobarycentric_weights are more efficient implementations.

Arguments

  • x::AbstractVector{TF}: Vector of distinct interpolation nodes, for equidistant nodes, use range is recommended x = x_min:dx:x_max
  • order::Integer: order of the interpolation (order = length(x) - 1)

Returns

  • Vector{TF}: Barycentric weights scaled such that their infinity norm equals 1

Details

The barycentric weights $w_j$ for nodes $x_j$ are computed as:

\[w_j = \frac{1}{\prod_{k \neq j} (x_j - x_k)}\]

For equidistant nodes, a more efficient formula is used based on binomial coefficients.

Examples

# For arbitrary nodes
x = [0.0, 1.0, 2.0]
w = barycentric_weights(x)

# For equidistant nodes
w = barycentric_weights(2)  # 3 nodes: 0, 1, 2

# or
x = 0.0:0.1:1.0
w = barycentric_weights(x)

Notes

  • The weights are scaled to have unit infinity norm for numerical stability
  • Returns NaN weights if input points are not distinct
  • log-sum-exp trick is used to prevent underflow/overflow
  • For equidistant nodes, a more efficient algorithm based on binomial coefficients is used

References

source
Einstein.Utils.dot_kahan — Method
dot_kahan(v1::StridedVector{T}, v2::StridedVector{T}) where {T<:Number}

Neumaier's variant of Kahan summation algorithm to reduce numerical errors.

Note

  • Slower than dot_xsum for large vectors, but faster for small vectors.
  • Similar performance to dot_kahan
  • Uses loop unrolling for better performance while maintaining Kahan summation's

numerical stability. Processes two elements per iteration when possible.

References

source
Einstein.Utils.dot_xsum — Method
dot_xsum(x::StridedVector{T}, y::StridedVector{T}) where {T<:Real}

Compute the dot product of two vectors using extended precision accumulation. Uses the xsum package for improved numerical accuracy.

Note

  • Very fast for large vectors, but a bit slower than Kahan summation for small vectors.
  • Both input vectors must have the same length, which is not checked for performance reasons.

References

source
Einstein.Utils.sum_kahan — Method
sum_kahan(v::AbstractVector{T}) where {T<:Number}

Neumaier's variant of Kahan summation algorithm to reduce numerical errors.

Note

  • Slower than sum_xsum for large vectors, but faster for small vectors.
  • Uses loop unrolling for better performance while maintaining Kahan summation's

numerical stability. Processes two elements per iteration when possible.

References

source
Einstein.Utils.sws_A0 — Method
sws_A0(s::Integer, l::Integer) where {TR<:AbstractFloat}

Calculate angular separation constant at $a = 0$.

\[{}_s A_{\ell m}(0)=l(l+1)-s(s+1)\]

Arguments

  • s::Integer: spin
  • l::Integer: angular number
source
Einstein.Utils.sws_eigM — Method
sws_eigM(::Type{TR}, s::Integer, c::Complex{TR}, m::Integer, l_max::Integer) where {TR<:AbstractFloat}
sws_eigM!(M::AbstractMatrix{TR}, s::Integer, c::Complex, m::Integer, l_max::Integer) where {TR<:AbstractFloat}

Construct the spherical-spheroidal decomposition matrix truncated at l_max.

Arguments

  • TR: Type for floating point conversion
  • s::Integer: spin
  • c::Complex: oblateness parameter
  • m::Integer: azimuthal number
  • l_max::Integer: maximum angular number

References

source
Einstein.Utils.sws_eigen — Method
sws_eigen(::Type{TR}, s::Integer, c::Complex{TR}, m::Integer, l_max::Integer) where {TR<:AbstractFloat}

Calculate eigenvalues and eigenvectors of the spherical-spheroidal decomposition matrix.

Arguments

  • TR: Type for floating point conversion
  • s::Integer: spin
  • c::Complex: oblateness parameter
  • m::Integer: azimuthal number
  • l_max::Integer: maximum angular number

References

source
Einstein.Utils.sws_l_min — Method
sws_l_min(s::Integer, m::Integer)

Minimum allowed value of l for given s, m. Returns max(|s|, |m|).

Arguments

  • s::Integer: spin
  • m::Integer: azimuthal number
source