Docstrings

This page lists all docstrings in Polynomials4ML.jl including for functions that are not part of the public API. Please check with the Public API which functionality for which we aim to guarantee semver-stability.

Polynomials4ML.AbstractP4MLBasisType

abstract type AbstractP4MLBasis end

Annotates types that map a low-dimensional input, scalar or SVector, to a vector of scalars (feature vector, embedding, basis...). Can be used as a Lux layer.

source
Polynomials4ML.CTrigBasisType

Complex trigonometric polynomials up to degree N (inclusive). The basis is constructed in the order

[1, exp(im*θ), exp(-im*θ), exp(2im*θ), exp(-2im*θ), ..., 
                                exp(N*im*θ), exp(-N*im*θ) ]

where θ is input variable.

source
Polynomials4ML.ChebBasisType

ChebBasis(N):

Chebyshev polynomials up to degree N-1 (inclusive). i.e basis with length N. The basis is ordered as

\[[1, x, 2x^2-1, 4x^3-3x, ..., 2xT_{N-1}(x)-T_{N-2}(x)]\]

where x is input variable.

The differences between ChebBasis and chebyshev_basis is that ChebBasis computes the basis on the fly when it is compiled and it does not store the recursion coefficients as in chebyshev_basis. There might be a small performance benefit from this.

Secondly, ChebBasis and chebyshev_basis use different normalization.

source
Polynomials4ML.MonoBasisType

Standard Monomials basis. This should very rarely be used. Possibly useful in combination with a transformation of the inputs, e.g. exponential.

source
Polynomials4ML.OrthPolyBasis1D3TType

OrthPolyBasis1D3T: defines a basis of polynomials in terms of a 3-term recursion,

\[\begin{aligned} P_1(x) &= A_1 \\ P_2 &= A_2 x + B_2 \\ P_{n} &= (A_n x + B_n) P_{n-1}(x) + C_n P_{n-2}(x) \end{aligned}\]

Typically (but not necessarily) such bases are obtained by orthogonalizing the monomials with respect to a user-specified distribution, which can be either continuous or discrete but must have a density function. See also

  • legendre_basis
  • chebyshev_basis
  • jacobi_basis
source
Polynomials4ML.RTrigBasisType

RTrigBasis(N):

Real trigonometric polynomials up to degree N (inclusive). The basis is ordered as

[1, cos(θ), sin(θ), cos(2θ), sin(2θ), ..., cos(Nθ), sin(Nθ) ]

where θ is input variable.

source
Polynomials4ML.StaticBatchType

StaticBatch{N,T} : an auxiliary StaticArray type that is distinct from SVector{N,T}. It can be used to create a batch of inputs of static size N. It is in used to convert function calls with single inputs to function calls with a batch of inputs.

source
Polynomials4ML.TransformedBasisType

struct TransformedBasis

Basically a two-stage chain, consisting of an input transformation, and basis evaluation. Constructor:

TransformedBasis(trans, basis)

The point of this structure is to provide a transformed basis that behaves exactly as all other P4ML bases.

Comments

  • the "natural indices" will be the same as for basis
  • _generate_input is not implemented for general input transforms; to

implement it for an input transform of type TIN one should monkey-patch

Polynomials4ML._generate_input(::TransformedBasis{TIN}, T::Type) where {TIN} = ...
  • _valtype is implemented but unclear how well it behaves, might be necessary

to monkey-patch it as well

source
Polynomials4ML._gradtypeFunction

_gradtype(basis, x)

If the intention is that P, dP = evaluate_ed(basis, x) then then _gradtype(basis, x) should return etype(dP).

source
Polynomials4ML._valtypeFunction

_valtype(basis, x)

If the intention is that P = basis(x) where P is a Vector{T} then _valtype(basis, x) should return T.

Here, x can be a single input, a batch or a type. A new basis type TB only needs to implement _valtype(::TB, x::Type).

source
Polynomials4ML.chebyshev_basisMethod

chebyshev_basis(N::Integer):

Constructs an OrthPolyBasis1D3T object representing a possibly rescaled version of the basis of Chebyshev polynomials of the first kind. N is the length of the basis, not the degree.

Careful: the normalisation may be non-standard.

source
Polynomials4ML.construct_basisMethod
construct_basis(RadialDecay, ζ_raw, D_raw, decay)

Construct a RadialDecay object from raw matrix data ζ_raw, D_raw and a DecayFunction(f, df) representing the decay form and its derivative. All input is converted to statically-sized SMatrix for efficiency.

source
Polynomials4ML.indexFunction

index(basis, k) -> Integer

Given a "natural description" of a basis element return the index of that basis element in the computed vector of basis function values. For example, for Chebyshev polynomials, index(basis, n) returns n+1.

source
Polynomials4ML.jacobi_basisMethod

jacobi_basis(N::Integer, α::Real, β::Real):

Constructs an OrthPolyBasis1D3T object representing a possibly rescaled version of the basis of Jacobi polynomials Jαβ. N is the length of the basis, not the degree.

Careful: the normalisation may be non-standard.

source
Polynomials4ML.legendre_basisMethod

legendre_basis(N::Integer):

Constructs an OrthPolyBasis1D3T object representing a possibly rescaled version of the basis of Legendre polynomials (L2 orthonormal on [-1, 1]). N is the length of the basis, not the degree.

Careful: the normalisation may be non-standard.

source
Polynomials4ML.natural_indicesFunction

natural_indices(basis) -> AbstractVector

Returns an abstract vector of "natural" descriptions of the basis functions in the order that they are stored in the computed vector of basis function values. For example, for Chebyshev polynomials, natural_indices(basis) returns 0:N, where N+1 is the length of the basis. For Spherical Harmmonics, a natural description requires two indices (l, m), so the output will be a vector of tuples.

At the moment, this function is used only for inspection and testing so no strict format is enforced.

source
Polynomials4ML.orthpolybasisMethod

function orthpolybasis(...) : construct a univariate orthogonal polynomial basis with respect to some specified inner product. For the standard 3-term recursion polynomials, use legendre_basis, jacobi_basis or chebyshev_basis.

The orthpolybasis currently implements orthogonal polynomials for discrete weights with the following constructors:

orthpolybasis(N::Integer, W::DiscreteWeights{TW}; TX = Float64)
orthpolybasis(N::Integer, X::AbstractVector{<: Real}, 
              W::AbstractVector{<: Real}, normalizeW=false; kwargs...)
source
Polynomials4ML.Utils.gensparseMethod

gensparse(...) : utility function to generate high-dimensional sparse grids which are downsets. All arguments are keyword arguments (with defaults):

  • NU : maximum correlation order
  • minvv = 0 : minvv[i] gives the minimum value forvv[i]`
  • maxvv = Inf : maxvv[i] gives the minimum value forvv[i]`
  • tup2b = vv -> vv :
  • admissible = _ -> false : determines whether a tuple belongs to the downset
  • filter = _ -> true : a callable object that returns true of tuple is to be kept and

false otherwise (whether or not it is part of the downset!) This is used, e.g. to enfore conditions such as ∑ lₐ = even or |∑ mₐ| ≦ M

  • INT = Int : integer type to be used
  • ordered = false : whether only ordered tuples are produced; ordered tuples

correspond to permutation-invariant basis functions

source