API Reference

Here we present the API reference for all functions and types within the module. The end user must only use the exported objects but private objects are also documented for completeness

Public Modules

QuasinormalModes.QuasinormalModesModule

This package contains routines for computing eigenvalues of second order ordinary differential equations and in particular the quasinormal modes (QNMs) of black holes in General Relativity using the "Asymptotic Iteration Method" [1] using the implementation based on the "improved" version of the AIM, described in [2].

References:

1 2

source

Public types

QuasinormalModes.AIMCacheType

Cache of coefficient arrays for the AIM. To each AIM problem corresponds a cache. As long as the problem doesn't change, the cache can be reused.

Members

  • icda::Array{T,1}: Hold the initial c data, i.e., c^i_0.
  • ccda::Array{T,1}: Hold the coefficients for the current aim step, c^i_n.
  • pcda::Array{T,1}: Hold the coefficients for the previous aim step, c^i_{n-1}.
  • bcda::Array{T,1}: The work buffer used to actually compute the c coefficients in parallel.
  • idda::Array{T,1}: Hold the initial d data, i.e., c^i_0.
  • cdda::Array{T,1}: Hold the coefficients for the current aim step, d^i_n.
  • pdda::Array{T,1}: Hold the coefficients for the previous aim step, d^i_{n-1}.
  • bdda::Array{T,1}: The work buffer used to actually compute the d coefficients in parallel.
  • size::N: The size of the arrays in the cache.
source
QuasinormalModes.AIMCacheMethod
AIMCache(p::NumericAIMProblem{N,T}) where {N <: Unsigned, T <: Number}

Create an AIMCache object suitable for Numeric Eigenvalue Problems.

Input

  • p::NumericAIMProblem: The problem data.

Output

An AIMCache{N,T} object.

source
QuasinormalModes.AIMCacheMethod
AIMCache(p::QuadraticEigenvalueProblem{N,T}) where {N <: Unsigned, T <: Number}

Create an AIMCache object suitable for Quadratic Eigenvalue Problems.

Input

  • p::QuadraticEigenvalueProblem: The problem data.

Output

An AIMCache{N,Polynomial{T}} object.

source

Public functions

QuasinormalModes.S0Method

All problem types must implement a S0 function. This behavior is enforced by the default implementations.

source
QuasinormalModes.computeDelta!Method
computeDelta!(m::AIMSteppingMethod, p::NumericAIMProblem{N,T}, c::AIMCache{N,T}, ω::T) where {N <: Unsigned, T <: Number}

Compute and return the AIM "quantization condition".

Input

  • m::AIMSteppingMethod: The stepping method to use.
  • p::QuadraticEigenvalueProblem: A quadratic frequency problem.
  • c::AIMCache: An AIM cache object created from p.
  • ω::T: Point to evaluate the quantization condition.

Output

An object of type T which represents the AIM quantization condition at point ω.

source
QuasinormalModes.computeDelta!Method
computeDelta!(p::QuadraticEigenvalueProblem{N,T}, c::AIMCache{N,Polynomial{T}}) where {N <: Unsigned, T <: Number}

Compute and return the AIM "quantization condition".

Input

  • m::AIMSteppingMethod: The stepping method to use.
  • p::QuadraticEigenvalueProblem: A quadratic frequency problem.
  • c::AIMCache: An AIM cache object created from p.

Output

An object of type Polynomial{T} whose roots are the problem's eigenvalues.

source
QuasinormalModes.computeEigenvaluesMethod
computeEigenvalues(
    m::AIMSteppingMethod,
    p::NumericAIMProblem{N,T},
    c::AIMCache{N,T},
    guess::T;
    nls_xtol::Real = convert(T, 1.0e-10),
    nls_ftol::Real = convert(T, 1.0e-10),
    nls_iterations::Int = 1000
    ) where {N <: Unsigned, T <: Complex}

Compute a single eigenvalue for the problem p with corresponding cache c.

Input

  • m::AIMSteppingMethod: The stepping method to use.
  • p::NumericAIMProblem: The previously defined problem data.
  • c::AIMCache: The cache constructed from p.
  • nls_xtol::Real: Norm difference in x between two successive iterates under which convergence is declared.
  • nls_ftol::Real: Infinite norm of residuals under which convergence is declared.
  • nls_iterations::Int: Maximum number of iterations performed by NLsolve.

Output

An object of type SolverResults returned by nlsolve. See NLsolve.jl for further details.

source
QuasinormalModes.computeEigenvaluesMethod
computeEigenvalues(
    m::AIMSteppingMethod,
    p::NumericAIMProblem{N,T},
    c::AIMCache{N,T},
    guess::T;
    roots_atol::Real = 1.0e-10,
    roots_rtol::Real = 1.0e-10,
    roots_xatol::Real = 1.0e-10,
    roots_xrtol::Real = 1.0e-10,
    roots_maxevals::Int = 100
    ) where {N <: Unsigned, T <: Real}

Compute a single eigenvalue for the problem p with corresponding cache c. For details on convergence settings see Roots.jl.

Input

  • m::AIMSteppingMethod: The stepping method to use.
  • p::NumericAIMProblem{N,T}: The previously defined problem data.
  • c::AIMCache{N,T}: The cache constructed from p.
  • grid::Tuple{T,T}: A tuple consisting of (start point, end point).
  • roots_atol::Real: Absolute tolerance.
  • roots_rtol::Real: Relative tolerance.
  • roots_xatol::Real: Floating point comparison absolute tolerance.
  • roots_xrtol::Real: Floating point comparison relative tolerance.
  • roots_maxevals::Int: Number of algorithm iterations performed.

Output

An object of type T containing the found eigenvalue.

source
QuasinormalModes.computeEigenvaluesMethod
computeEigenvalues(
    m::AIMSteppingMethod,
    p::QuadraticEigenvalueProblem{N,T},
    c::AIMCache{N,Polynomial{T}};
    plr_polish::Bool = true, 
    plr_epsilon::Real = convert(T, 1.0e-10)
    ) where {N <: Unsigned, T <: Number}

Compute the eigenvalues for the problem p with corresponding cache c.

Input

  • m::AIMSteppingMethod: The stepping method to use.
  • p::QuadraticEigenvalueProblem: The previously defined problem data.
  • c::AIMCache: The cache constructed from p.
  • plr_polish::Bool: Tell PolynomialRoots to divide the original polynomial by each root found and polish the results using the full polynomial.
  • plr_epsilon::Real: The stopping criterion described in Skowron & Gould paper. This is not the precision with which the roots will be calculated.

Output

An object of type Array{T,1} containing the computed eigenvalues.

source
QuasinormalModes.eigenvaluesInGridMethod
eigenvaluesInGrid(
    m::AIMSteppingMethod,
    p::NumericAIMProblem{N,T},
    c::AIMCache{N,T},
    grid::Tuple{T,T,Int64,Int64};
    xtol::Real = 1.0e-10,
    ftol::Real = 1.0e-10,
    iterations::Int = 1000
    ) where {N <: Unsigned, T <: Complex}

Attempts to find eigenvalues using a grid of complex plane data points as initial guesses passed to nlsolve.

Input

  • m::AIMSteppingMethod: The stepping method to use.
  • p::NumericAIMProblem{N,T}: The previously defined problem data.
  • c::AIMCache{N,T}: The cache constructed from p.
  • grid::Tuple{T,T,Int64,Int64}: A tuple consisting of (start point, end point, num. of real pts., num. of imag. pots.).
  • xtol::Real: Norm difference in x between two successive iterates under which convergence is declared.
  • ftol::Real: Infinite norm of residuals under which convergence is declared.
  • iterations::Int: Maximum number of iterations performed by NLsolve.

Output

An object of type Array{T,1} containing the modes found within the grid.

source
QuasinormalModes.eigenvaluesInGridMethod
eigenvaluesInGrid(
    m::AIMSteppingMethod,
    p::NumericAIMProblem{N,T},
    c::AIMCache{N,T},
    grid::Tuple{T,T};
    roots_atol::Real = 1.0e-10,
    roots_rtol::Real = 1.0e-10,
    roots_xatol::Real = 1.0e-10,
    roots_xrtol::Real = 1.0e-10,
    roots_maxevals::Int = 100
    ) where {N <: Unsigned, T <: Real}

Attempts to find eigenvalues using a range of real data points as a search region to find_zeros. For details on convergence settings see Roots.jl.

Input

  • m::AIMSteppingMethod: The stepping method to use.
  • p::NumericAIMProblem{N,T}: The previously defined problem data.
  • c::AIMCache{N,T}: The cache constructed from p.
  • grid::Tuple{T,T}: A tuple consisting of (start point, end point).
  • roots_atol::Real: Absolute tolerance.
  • roots_rtol::Real: Relative tolerance.
  • roots_xatol::Real: Floating point comparison absolute tolerance.
  • roots_xrtol::Real: Floating point comparison relative tolerance.
  • roots_maxevals::Int: Number of algorithm iterations performed.

Output

An object of type Array{T,1} containing the eigenvalues found within the grid.

source
QuasinormalModes.λ0Method

All problem types must implement a λ0 function. This behavior is enforced by the default implementations.

source

Private types

Private functions

QuasinormalModes.AIMStep!Method
AIMStep!(::Serial, p::AIMProblem{N,T}, c::AIMCache{N,U}) where {N <: Unsigned, T <: Number, U <: Any}

Performs a single step of the AIM algorithm sequentially:

  1. The initial data arrays are not altered.
  2. The previous arrays receive the values of the current arrays.
  3. The results of the next step computed using the initial and current arrays. Results are stored in the buffer arrays.
  4. The current arrays receive the contents of the buffer arrays.

Input

  • ::Serial: Instance of Serial object.
  • p::AIMProblem: The problem data to use in the computation.
  • c::AIMCache: The cache of arrays that corresponds to the problem p.

Output

nothing
source
QuasinormalModes.AIMStep!Method
AIMStep!(::Threaded, p::AIMProblem{N,T}, c::AIMCache{N,U}) where {N <: Unsigned, T <: Number, U <: Any}

Performs a single step of the AIM algorithm in parallel using threads:

  1. The initial data arrays are not altered.
  2. The previous arrays receive the values of the current arrays.
  3. The results of the next step computed using the initial and current arrays. Results are stored in the buffer arrays.
  4. The current arrays receive the contents of the buffer arrays.

Input

  • ::Threaded: Instance of the Threaded object.
  • p::AIMProblem: The problem data to use in the computation.
  • c::AIMCache: The cache of arrays that corresponds to the problem p.

Output

nothing
source
QuasinormalModes.computePolynomialFactorsMethod
computePolynomialFactors(p::QuadraticEigenvalueProblem{N,T}, n::N, f::Function) where {N <: Unsigned, T <: Number}

Create a second order Polynomial object in the variable ω by computing derivatives of λ0 or S0.

Input

  • p::QuadraticEigenvalueProblem: The problem data with the expressions to compute the derivative.
  • n::Unsigned: The order of the derivative.
  • f::Function: the function to extract the polynomial from. Either λ0 or S0.

Output

An object of type Polynomial{T} containing the polynomial resulting from the derivation of the expression.

source
QuasinormalModes.createPolyMethod
createPoly(p::QuadraticEigenvalueProblem{N,T}, n::N, f::Function)

Compute the n-th coefficient of the Taylor expansion around x0 for the functions λ0 or S0

Input

  • p::QuadraticEigenvalueProblem: The problem data with the expressions to compute the derivative.
  • n::Unsigned: The order of the derivative.
  • f::Function: the function to extract the polynomial from. Either λ0 or S0.

Output

An object of type Polynomial{T} containing the polynomial Taylor coefficient.

source
QuasinormalModes.dnxMethod
dnx(p::AnalyticAIMProblem{N,T}, n::N, f::Function, v::Function) where {N <: Unsigned, T <: Number}

Computes the n-th derivative of the AIM expressions with respect to ODE's variable. This function is only a thin wrapper around SymEngine's own diff function. It works as a barrier function that produces a type stable Basic result.

Input

  • p::AnalyticAIMProblem: The problem data with the expressions to compute the derivative.
  • n::Unsigned: The order of the derivative.
  • f::Function: The actual expression to compute the derivative. Either λ0 or S0.
  • v::Function: The variable with respec to which the derivative will be computed. Either getODEvar or getODEeigen.

Output

A SymEngine.Basic object with the derived expression.

source
QuasinormalModes.dnxMethod
dnx(n::T, f::Basic, v::Basic) where {T <: Unsigned}

Computes the n-th derivative of the the function f with respect to the variable v. This function re implements SymEngine's own diff function using an early quitting strategy.

Input

  • n::T: The order of the derivative.
  • f::Basic: The expression to compute the derivative.
  • v::Basic: The variable with respect to which the derivative will be computed.

Output

A SymEngine.Basic object with the derived expression.

source
QuasinormalModes.recomputeInitials!Method
recomputeInitials(p::NumericAIMProblem{N,T}, c::AIMCache{N,T}, ω::T) where {N <: Unsigned, T <: Number}

Reevaluate (in-place) the initial data arrays. The initial data array elements are the Taylor expansion coefficients of λ0 and S0 in the ODE variable x around x0 of order get_niter(p) at a point ω.

Input

  • p::NumericAIMProblem: The problem data.
  • c::AIMCache: The problem data associated cache.
  • ω::T: The value of the eigenvalue to evaluate the arrays in.

Output

nothing
source