FastTransforms.jl Documentation

Introduction

FastTransforms.jl allows the user to conveniently work with orthogonal polynomials with degrees well into the millions.

This package provides a Julia wrapper for the C library of the same name. Additionally, all three types of nonuniform fast Fourier transforms available, as well as the Padua transform.

Fast orthogonal polynomial transforms

For this documentation, please see the documentation for FastTransforms. Most transforms have separate forward and inverse plans. In some instances, however, the inverse is in the sense of least-squares, and therefore only the forward transform is planned.

Fast Cholesky factorization of the Gram matrix

FastTransforms.GramMatrix — Type

GramMatrix(W::AbstractMatrix, X::AbstractMatrix)

Construct a symmetric positive-definite Gram matrix with data stored in $W$. Given a family of orthogonal polynomials ${\bf P}(x) = \{p_0(x), p_1(x),\ldots\}$ and a continuous inner product $\langle f, g\rangle$, the Gram matrix is defined by:

\[W_{i,j} = \langle p_{i-1}, p_{j-1}\rangle.\]

Moreover, given $X$, the transposed Jacobi matrix that satisfies $x {\bf P}(x) = {\bf P}(x) X$, the Gram matrix satisfies the skew-symmetric rank-2 displacement equation ($X = X_{1:n, 1:n}$):

\[X^\top W - WX = GJG^\top,\]

where $J = \begin{pmatrix} 0 & 1\\ -1 & 0\end{pmatrix}$ and where:

\[G_{:, 1} = e_n,\quad{\rm and}\quad G_{:, 2} = W_{n-1, :}X_{n-1, n} - X^\top W_{:, n}.\]

Fast ($O(n^2)$) Cholesky factorization of the Gram matrix returns the connection coefficients between ${\bf P}(x)$ and the polynomials ${\bf Q}(x)$ orthogonal in the modified inner product, ${\bf P}(x) = {\bf Q}(x) R$.

FastTransforms.ChebyshevGramMatrix — Type

ChebyshevGramMatrix(μ::AbstractVector)

Construct a Chebyshev–Gram matrix of size (length(μ)+1)÷2 with entries:

\[W_{i,j} = \frac{\mu_{|i-j|+1} +\mu_{i+j-1}}{2}.\]

Due to the linearization of a product of two first-kind Chebyshev polynomials, the Chebyshev–Gram matrix can be constructed from modified Chebyshev moments:

\[\mu_{n} = \langle T_{n-1}, 1\rangle.\]

Specialized construction and Cholesky factorization is given for this type.

See also GramMatrix for the general case.

Nonuniform fast Fourier transforms

FastTransforms.nufft1 — Function

Computes a nonuniform fast Fourier transform of type I:

\[f_j = \sum_{k=0}^{N-1} c_k e^{-2\pi{\rm i} \frac{j}{N} \omega_k},\quad{\rm for}\quad 0 \le j \le N-1.\]

Computes a 2D nonuniform fast Fourier transform of type I-I:

\[F_{i,j} = \sum_{k=0}^{M-1}\sum_{\ell=0}^{N-1} C_{k,\ell} e^{-2\pi{\rm i} (\frac{i}{M} \omega_k + \frac{j}{N} \pi_{\ell})},\quad{\rm for}\quad 0 \le i \le M-1,\quad 0 \le j \le N-1.\]

FastTransforms.nufft2 — Function

Computes a nonuniform fast Fourier transform of type II:

\[f_j = \sum_{k=0}^{N-1} c_k e^{-2\pi{\rm i} x_j k},\quad{\rm for}\quad 0 \le j \le N-1.\]

Computes a 2D nonuniform fast Fourier transform of type II-II:

\[F_{i,j} = \sum_{k=0}^{M-1}\sum_{\ell=0}^{N-1} C_{k,\ell} e^{-2\pi{\rm i} (x_i k + y_j \ell)},\quad{\rm for}\quad 0 \le i \le M-1,\quad 0 \le j \le N-1.\]

FastTransforms.nufft3 — Function

Computes a nonuniform fast Fourier transform of type III:

\[f_j = \sum_{k=0}^{N-1} c_k e^{-2\pi{\rm i} x_j \omega_k},\quad{\rm for}\quad 0 \le j \le N-1.\]

FastTransforms.inufft1 — Function

Computes an inverse nonuniform fast Fourier transform of type I.

FastTransforms.inufft2 — Function

Computes an inverse nonuniform fast Fourier transform of type II.

FastTransforms.paduatransform — Function

Padua Transform maps from interpolant values at the Padua points to the 2D Chebyshev coefficients.

FastTransforms.ipaduatransform — Function

Inverse Padua Transform maps the 2D Chebyshev coefficients to the values of the interpolation polynomial at the Padua points.

Other Exported Methods

FastTransforms.gaunt — Function

Calculates the Gaunt coefficients, defined by:

\[a(m,n,\mu,\nu,q) = \frac{2(n+\nu-2q)+1}{2} \frac{(n+\nu-2q-m-\mu)!}{(n+\nu-2q+m+\mu)!} \int_{-1}^{+1} P_n^m(x) P_\nu^\mu(x) P_{n+\nu-2q}^{m+\mu}(x) {\rm\,d}x.\]

or defined by:

\[P_n^m(x) P_\nu^\mu(x) = \sum_{q=0}^{q_{\rm max}} a(m,n,\mu,\nu,q) P_{n+\nu-2q}^{m+\mu}(x)\]

This is a Julia implementation of the stable recurrence described in:

Y.-l. Xu, Fast evaluation of Gaunt coefficients: recursive approach, J. Comp. Appl. Math., 85:53–65, 1997.

Calculates the Gaunt coefficients in 64-bit floating-point arithmetic.

FastTransforms.paduapoints — Function

Returns coordinates of the $(n+1)(n+2)/2$ Padua points.

FastTransforms.sphevaluate — Function

Pointwise evaluation of real orthonormal spherical harmonic:

\[Y_\ell^m(\theta,\varphi) = (-1)^{|m|}\sqrt{(\ell+\frac{1}{2})\frac{(\ell-|m|)!}{(\ell+|m|)!}} P_\ell^{|m|}(\cos\theta) \sqrt{\frac{2-\delta_{m,0}}{2\pi}} \left\{\begin{array}{ccc} \cos m\varphi & {\rm for} & m \ge 0,\\ \sin(-m\varphi) & {\rm for} & m < 0.\end{array}\right.\]

Internal Methods

Miscellaneous Special Functions

FastTransforms.half — Function

Compute a typed 0.5.

FastTransforms.two — Function

Compute a typed 2.

FastTransforms.δ — Function

The Kronecker $\delta$ function:

\[\delta_{k,j} = \left\{\begin{array}{ccc} 1 & {\rm for} & k = j,\\ 0 & {\rm for} & k \ne j.\end{array}\right.\]

FastTransforms.Λ — Function

The Lambda function $\Lambda(z) = \frac{\Gamma(z+\frac{1}{2})}{\Gamma(z+1)}$ for the ratio of gamma functions.

For 64-bit floating-point arithmetic, the Lambda function uses the asymptotic series for $\tau$ in Appendix B of

I. Bogaert and B. Michiels and J. Fostier, 𝒪(1) computation of Legendre polynomials and Gauss–Legendre nodes and weights for parallel computing, SIAM J. Sci. Comput., 34:C83–C101, 2012.

The Lambda function $\Lambda(z,λ₁,λ₂) = \frac{\Gamma(z+\lambda_1)}{Γ(z+\lambda_2)}$ for the ratio of gamma functions.

FastTransforms.lambertw — Function

The principal branch of the Lambert-W function, defined by $x = W_0(x) e^{W_0(x)}$, computed using Halley's method for $x \in [-e^{-1},\infty)$.

FastTransforms.pochhammer — Function

Pochhammer symbol $(x)_n = \frac{\Gamma(x+n)}{\Gamma(x)}$ for the rising factorial.

FastTransforms.stirlingseries — Function

Stirling's asymptotic series for $\Gamma(z)$.

Modified Chebyshev Moment-Based Quadrature

FastTransforms.clenshawcurtisnodes — Function

Compute nodes of the Clenshaw—Curtis quadrature rule.

FastTransforms.clenshawcurtisweights — Function

Compute weights of the Clenshaw—Curtis quadrature rule with modified Chebyshev moments of the first kind $\mu$.

FastTransforms.fejernodes1 — Function

Compute nodes of Fejer's first quadrature rule.

FastTransforms.fejerweights1 — Function

Compute weights of Fejer's first quadrature rule with modified Chebyshev moments of the first kind $\mu$.

FastTransforms.fejernodes2 — Function

Compute nodes of Fejer's second quadrature rule.

FastTransforms.fejerweights2 — Function

Compute weights of Fejer's second quadrature rule with modified Chebyshev moments of the second kind $\mu$.

FastTransforms.chebyshevmoments1 — Function

Modified Chebyshev moments of the first kind:

\[ \int_{-1}^{+1} T_n(x) {\rm\,d}x.\]

FastTransforms.chebyshevjacobimoments1 — Function

Modified Chebyshev moments of the first kind with respect to the Jacobi weight:

\[ \int_{-1}^{+1} T_n(x) (1-x)^\alpha(1+x)^\beta{\rm\,d}x.\]

FastTransforms.chebyshevlogmoments1 — Function

Modified Chebyshev moments of the first kind with respect to the logarithmic weight:

\[ \int_{-1}^{+1} T_n(x) \log\left(\frac{1-x}{2}\right){\rm\,d}x.\]

FastTransforms.chebyshevmoments2 — Function

Modified Chebyshev moments of the second kind:

\[ \int_{-1}^{+1} U_n(x) {\rm\,d}x.\]

FastTransforms.chebyshevjacobimoments2 — Function

Modified Chebyshev moments of the second kind with respect to the Jacobi weight:

\[ \int_{-1}^{+1} U_n(x) (1-x)^\alpha(1+x)^\beta{\rm\,d}x.\]

FastTransforms.chebyshevlogmoments2 — Function

Modified Chebyshev moments of the second kind with respect to the logarithmic weight:

\[ \int_{-1}^{+1} U_n(x) \log\left(\frac{1-x}{2}\right){\rm\,d}x.\]

Elliptic

FastTransforms.Elliptic — Module

FastTransforms submodule for the computation of some elliptic integrals and functions.

Complete elliptic integrals of the first and second kinds:

\[K(k) = \int_0^{\frac{\pi}{2}} \frac{{\rm d}\theta}{\sqrt{1-k^2\sin^2\theta}},\quad{\rm and},\]

\[E(k) = \int_0^{\frac{\pi}{2}} \sqrt{1-k^2\sin^2\theta} {\rm\,d}\theta.\]

Jacobian elliptic functions:

\[x = \int_0^{\operatorname{sn}(x,k)} \frac{{\rm d}t}{\sqrt{(1-t^2)(1-k^2t^2)}},\]

\[x = \int_{\operatorname{cn}(x,k)}^1 \frac{{\rm d}t}{\sqrt{(1-t^2)[1-k^2(1-t^2)]}},\]

\[x = \int_{\operatorname{dn}(x,k)}^1 \frac{{\rm d}t}{\sqrt{(1-t^2)(t^2-1+k^2)}},\]

and the remaining nine are defined by:

\[\operatorname{pq}(x,k) = \frac{\operatorname{pr}(x,k)}{\operatorname{qr}(x,k)} = \frac{1}{\operatorname{qp}(x,k)}.\]