Core

Core#

diffaaable.core.aaa(Z, F, tol=1e-13, mmax=100, return_errors=False)[source]#

This is a wrapped version of aaa as provided by baryrat, providing a custom jvp to enable differentiability.

For detailed information on the usage of aaa please refer to the original documentation:

Compute a rational approximation of `F` over the points `Z` using the
AAA algorithm.

Arguments:
    Z (array): the sampling points of the function. Unlike for interpolation
        algorithms, where a small number of nodes is preferred, since the
        AAA algorithm chooses its support points adaptively, it is better
        to provide a finer mesh over the support.
    F: the function to be approximated; can be given as a function or as an
        array of function values over ``Z``.
    tol: the approximation tolerance
    mmax: the maximum number of iterations/degree of the resulting approximant
    return_errors: if `True`, also return the history of the errors over
        all iterations

Returns:
    BarycentricRational: an object which can be called to evaluate the
    rational function, and can be queried for the poles, residues, and
    zeros of the function.

For more information, see the paper

  | The AAA Algorithm for Rational Approximation
  | Yuji Nakatsukasa, Olivier Sete, and Lloyd N. Trefethen
  | SIAM Journal on Scientific Computing 2018 40:3, A1494-A1522
  | https://doi.org/10.1137/16M1106122

as well as the Chebfun package <http://www.chebfun.org>. This code is an
almost direct port of the Chebfun implementation of aaa to Python.

Attention

Returns nodes, values, weights and poles, in contrast to the baryrat implementation that returns the BarycentricRational ready to be evaluated. This is done to facilitate differentiability.

Parameters:

z_k (array (M,)) – the sampling points of the function. Unlike for interpolation algorithms, where a small number of nodes is preferred, since the AAA algorithm chooses its support points adaptively, it is better to provide a finer mesh over the support.
f_k (array (M,)) – the function to be approximated; can be given as a callable function or as an array of function values over z_k.
tol (float) – the approximation tolerance
mmax (int) – the maximum number of iterations/degree of the resulting approximant

Returns:

z_j (array (m,)) – nodes of the barycentric approximant
f_j (array (m,)) – values of the barycentric approximant
w_j (array (m,)) – weights of the barycentric approximant
z_n (array (m-1,)) – poles of the barycentric approximant (for convenience)

diffaaable.core.aaa_jvp(primals, tangents)[source]#: Derivatives according to [1]. The implemented matrix expressions are motivated in the appendix of [2]:

AAA Derivatives

Here we will briefly elaborate how the derivatives introduced in [1] are implemented as JAX compatible Jacobian Vector products (JVPs) in diffaaable.

Given the tangents \(\frac{\partial f_k}{\partial p}\) we will use the chain rule on \(r(w_j, f_j, z)\) along its weights \(w_j\) and values \(f_j\) (the nodes \(z_j\) are treated as independent of \(p\)) (Equation 4 of [1]): \[ \frac{\partial f_k}{\partial p} \approx \frac{\partial r_k}{\partial p}= \sum_{j=1}^m\frac{\partial r_k}{\partial f_j}\frac{\partial f_j}{\partial p}+\sum_{j=1}^m\frac{\partial r_k}{\partial w_j}\frac{\partial w_j}{\partial p} \]

To solve this system of equations for \(\frac{\partial w_j}{\partial p}\) we express it in matrix form:

\[\mathbf{A}\mathbf{w}^\prime = \mathbf{b}\]

where \(\mathbf{w}^\prime\) is the column vector containing \(\frac{\partial w_j}{\partial p}\) and \(\mathbf{b}\) and \(\mathbf{A}\) are defined element wise:

\[\begin{split}\begin{aligned} b_k &= \frac{\partial f_k}{\partial p} - \sum_{j=1}^m\frac{\partial r_k}{\partial f_j}\frac{\partial f_j}{\partial p}\\ A_{kj} &= \frac{\partial r_k}{\partial w_j} \end{aligned}\end{split}\]

These are augmented with Equation 5 of \cite{betzEfficientRationalApproximation2024} which removes the ambiguity in \(\mathbf{w}^\prime\) associated with a shared phase of all weights.

The expressions for the derivatives of \(r_k\) are found in the definition of \(r(z)\) (Equation 1 in [2]): \[ \frac{\partial r_k}{\partial f_j}= \frac{1}{d(z_k)} \frac{w_j}{z_k-z_j}\] and \[ \frac{\partial r_k}{\partial w_j}= \frac{1}{d(z_k)} \frac{f_j-r_k}{z_k-z_j} \approx \frac{1}{d(z_k)} \frac{f_j-f_k}{z_k-z_j} \]

diffaaable.core.poles(z_j, w_j)[source]#

The poles of a barycentric rational with given nodes and weights. Poles lifted by zeros of the nominator are included. Thus the values \(f_j\) do not contribute and don’t need to be provided The implementation was modified from baryrat to support JAX AD.

Parameters:

z_j (array (m,)) – nodes of the barycentric rational
w_j (array (m,)) – weights of the barycentric rational

Returns:

z_n – poles of the barycentric rational (more strictly zeros of the denominator)

Return type:

array (m-1,)

diffaaable.core.residues(z_j, f_j, w_j, z_n)[source]#

Residues for given poles via formula for simple poles of quotients of analytic functions. The implementation was modified from baryrat to support JAX AD.

Parameters:

z_j (array (m,)) – nodes of the barycentric rational
w_j (array (m,)) – weights of the barycentric rational
z_n (array (n,)) – poles of interest of the barycentric rational (n<=m-1)

Returns:

r_n – residues of poles z_n

Return type:

array (n,)

Core

Contents

Core#