# Interpolatory rational model order reduction of parametric problems lacking
uniform inf-sup stability††thanks: This work was funded by the Swiss National
Science Foundation (SNF) through project no. 182236.

Davide Pradovera CSQI, EPFL, Lausanne, Switzerland (davide.pradovera@epfl.ch).

###### Abstract

We present a technique for the approximation of a class of Hilbert space-
valued maps which arise within the framework of Model Order Reduction for
parametric partial differential equations, whose solution map has a
meromorphic structure. Our MOR stategy consists in constructing an _explicit
rational_ approximation based on few snapshots of the solution, in an
interpolatory fashion. Under some restrictions on the structure of the
original problem, we describe _a priori_ convergence results for our
technique, hereafter called _minimal rational interpolation_ , which show its
ability to identify the main features (e.g. resonance locations) of the target
solution map. We also investigate some procedures to obtain _a posteriori_
error indicators, which may be employed to adapt the degree and the sampling
points of the minimal rational interpolant. Finally, some numerical
experiments are carried out to confirm the theoretical results and the
effectiveness of our technique.

Keywords: Model order reduction, rational approximation, Hilbert space-valued
meromorphic maps, non-coercive parametric PDEs, Helmholtz equation, frequency
response.

AMS Subject Classification: 30D30, 41A20, 41A25, 35J05, 35P15, 65D15.

##  1 Introduction

Mathematical models based on PDEs are employed to analyze numerically a wide
array of physical, financial, and engineering-related phenomena. In many
situations, such models depend on one or more parameters, either because of
uncertainties or for control/design purposes, and multiple solves of the model
have to be performed for different parameter values.

Often the “naive” approach of solving the original problem (which we will
refer to as _full order model_ , FOM) as many times as required is not
feasible, for instance because the number of evaluations is too large, or
because a real-time solution is required. In such cases, model order reduction
(MOR) is employed to reduce the computational effort needed for the solution
of the FOM. This is achieved by constructing a surrogate model, whose solution
is managed to be close to that of the FOM. The construction of the reduced
model often requires a considerable computational cost, which, however, can be
done _offline_ in a preliminary phase, whereas its evaluation at any parameter
value is quite cheap, and can be carried out _online_ in real-time.

Many MOR techniques have been proposed for general FOMs, with the most notable
and widely applied being the Reduced Basis (RB) method [15, 16, 17, 21, 26,
35, 36]. The RB method, in its simplest form, assumes that few samples of the
solution of the FOM (snapshots) are enough to capture the main features of the
solution manifold, i.e. the set of all solutions of the FOM for all parameter
values. Exploiting this idea, a surrogate model is built by projecting the FOM
onto a subspace generated by few pre-computed snapshots. The effectiveness of
the RB method has been verified in many cases. However, it is possible to find
some parametric problems for which projection-based MOR does not perform
optimally. Such examples are often characterized by some irregularity of the
FOM (e.g. nonsmooth dependence on the parameters, non-linearity, lack of
stability, etc. of the differential operator therein).

One simple case falling into this category, which is of primary importance in
our discussion, is the time-harmonic wave (Helmholtz) equation with parametric
wavenumber:

$\text{for }k\in K\subset\mathbb{C},\;\text{find }u(k)\in
V\;\text{s.t.}\;-\Delta u(k)-k^{2}u(k)\stackrel{{\scriptstyle(W)}}{{=}}f\in
W\text{,}$ (1)

where $V$ and $W\supset V$ are some suitable Banach function spaces. Depending
on the choice of $K$ and of the boundary conditions which complement Eq. 1,
the problem above may lack inf-sup stability over $K$ [32], due to the
presence of resonant wavenumbers. If one applies carelessly the RB method (in
both its main versions, POD and greedy [35, 36]) to Eq. 1, one may observe the
appearance of spurious (“non-physical”) resonances in the surrogate model.

Discussions concerning this effect can be found in works related to MOR
methods for dynamical systems (e.g. Krylov subspace methods [20, 22]), where
problems of similar form are quite common, due to the need for frequency-
domain computations. Usually, if the number of parameters is small (often,
just 1, the frequency), explicit rather than implicit approximants are
considered, with the surrogate model being built by enforcing interpolation or
moment matching conditions at the sample points in the parameter space,
possibly in a Least Squares (LS) way. The main representatives of this class
of methods are the Loewner framework [1, 2, 6, 7, 29] and the vector fitting
(VF) algorithm [18, 22, 23, 25, 30, 39], which are tailored to the specific
structure of the FOM: rational functions are employed, with the objective of
representing each resonance of the FOM by a pole of the approximant. This
turns out to be particularly useful when the resonances of the FOM have a
physical meaning, since their approximation is implicitly enclosed in the
surrogate model, from which it can be obtained with small or no computational
effort.

More recently, in the wake of these methods for dynamical systems, and somehow
trying to profit from their main advantages, univariate LS Padé approximants
have been introduced and studied, both in a standard [9, 11] and a fast [10]
version, in the context of a single parameter. Such techniques are based on
multiple solves of the FOM at a single parameter value, in the same spirit as
Krylov subspace methods, but yield an explicit rational approximant like VF.
Actually, LS Padé approximants are actually quite comparable to VF, the main
differences between the methods being:

  1. (i)

in VF, one usually must achieve a reasonably high “sampling density” [18],
namely use at least $2N$ samples to approximate $N$ resonances; conversely,
fast LS Padé approximants can get away with half as many samples;

  2. (ii)

the fast version of LS Padé approximants relies on the high-dimensionality of
the samples, namely on their linear independence, whereas VF can be applied
even for scalar outputs (this is not actually a limitation of LS Padé
approximants, since it suffices to build a rational approximation of the
solution map, and then to apply the desired linear functional to such
surrogate solution);

  3. (iii)

the error convergence of LS Padé approximants with respect to the number of
samples and to the degree of the denominator is well understood for a
reasonably broad class of problems [9, 10], while the convergence behavior of
the VF error has, to our knowledge, yet to be analyzed thoroughly [18];

  4. (iv)

LS Padé approximants are computed from moments of the solution at a single
parameter, the explicit computation of which may feature numerically
instabilities, see [12]; instead, VF uses distributed samples in a LS-Lagrange
interpolation fashion, for which ill-conditioning issues are much milder.

In this paper we wish to develop and analyze an extension of fast LS Padé
approximants which overcomes iv by allowing a distributed parameter sampling,
while keeping i and iii intact. In particular, in section 2 we introduce our
method, which we name _minimal rational interpolation_ , and we state more
precisely to which kinds of FOMs our theory applies to. Then, in section 3 we
describe the main convergence results concerning minimal rational
interpolation. Afterwards, in section 4 we discuss briefly the topic of _a
posteriori_ error/residual estimation for our technique. Possible extensions
to multi-parameter problems are described in section 5. In section 6 we show
two numerical examples to verify the effectiveness of our method and of the
presented error indicator, as well as the theoretical convergence rates. Some
concluding remarks are in section 7.

##  2 Description of the method

Let $\left(V,\langle\cdot,\cdot\rangle_{V}\right)$ be a Hilbert space over
$\mathbb{C}$, with induced norm $\left\|\cdot\right\|_{V}$, and
$K\subset\mathbb{C}$ be compact. Our task is to find a surrogate for a given
map $u:K\to V$. Notably, we seek an approach which builds a surrogate starting
only from few evaluations (snapshots) of $u$. Assuming $u$ to be the outcome
of a complex computational model, such technique qualifies as _non-intrusive_
, in the sense that it can employ any available solver as a “black box”.

We assume that $S>0$ samples of $u$ at the points
$\Xi_{S}=\\{\mu_{j}\\}_{j=1}^{S}\subset K$ are available. In particular, under
some additional regularity assumptions on $u$, we allow confluence of the
sample points: if $\mu_{j}$ appears $K_{j}+1$ times in $\Xi_{S}$, we assume to
have evaluated $u$, as well as its derivatives with respect to $\mu$ up to
order $K_{j}$, at the point $\mu_{j}$.

In order to set up our approximation, it is necessary to specify how many
“resonances” of $u$ we wish to approximate. To this aim, we consider the
integer $N\in\\{0,\ldots,S-1\\}$, as well as the set of normalized polynomials
of degree at most $N$

$\mathbb{P}_{N}^{\star}\left(\mathbb{C}\right)=\left\\{Q\in\mathbb{P}_{N}\left(\mathbb{C}\right),\;{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|Q\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}=1\right\\}\text{,}$
(2)

with $\mathbb{P}_{N}\left(\mathbb{C}\right)$ the space of $\mathbb{C}$-valued
polynomials of degree up to $N$, and
${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}$
an arbitrary (but fixed) Hilbertian norm on
$\mathbb{P}_{N}\left(\mathbb{C}\right)$. Our aim will be to approximate $N$
resonances of $u$ by the roots of some element of
$\mathbb{P}_{N}^{\star}\left(\mathbb{C}\right)$.

Let $I^{\Xi_{S}}$ be the polynomial interpolation operator111For stability
reasons it may be wise to replace $I^{\Xi_{S}}$ with a Least Squares-type
interpolator. With some effort, our results can be generalized to such case.
for $V$-valued functions based on samples at the points $\Xi_{S}$. Namely,
given $\phi:\Xi_{S}\to V$, the interpolant $I^{\Xi_{S}}(\phi)$ is the element
of the space of $V$-valued polynomials

$\mathbb{P}_{S-1}\left(\mathbb{C};V\right)=\left\\{\mu\mapsto\sum_{i=0}^{S-1}p_{i}\mu^{i},\;p_{i}\in
V\,\text{ for }i=0,\ldots,S-1\right\\}$ (3)

such that $\phi(\mu^{\prime})=I^{\Xi_{S}}(\phi)(\mu^{\prime})$ for
$\mu^{\prime}\in\Xi_{S}$. In case of confluent sample points, such
interpolation condition has to be interpreted in a Hermite sense, with
derivatives of the interpolant matching those of the original function.

######  Definition 1 (Minimal rational interpolants).

A minimal rational interpolant of $u$ of type $[S-1/N]$ with samples at
$\Xi_{S}$ is a rational function

$u_{N}^{\Xi_{S}}(\mu)=I^{\Xi_{S}}\left(uQ_{N}^{\Xi_{S}}\right)(\mu)\,{\Big{/}}\,Q_{N}^{\Xi_{S}}(\mu),\quad\mu\in\mathbb{C},Q_{N}^{\Xi_{S}}(\mu)\neq
0\text{,}$

where the denominator $Q_{N}^{\Xi_{S}}$ minimizes the functional

${\mathbb{P}_{N}\left(\mathbb{C}\right)\ni
Q\;\mapsto\;}j^{\Xi_{S}}(Q)^{2}=\left\|\frac{1}{(S-1)!}\frac{\operatorname{d}^{S-1}}{\operatorname{d}{\mu}^{S-1}}I^{\Xi_{S}}(uQ)\right\|_{V}^{2}\text{,}$
(4)

over $\mathbb{P}_{N}^{\star}\left(\mathbb{C}\right)$.

The specific expression of $j^{\Xi_{S}}$ is the natural generalization of the
corresponding functional in fast LS Padé approximation [10]. We will see that
this form is essential to prove all our convergence results, starting from the
central results described in Lemma 6.

We remark that, since $j^{\Xi_{S}}$ is convex, see Eq. 5, the existence of a
minimizer $Q_{N}^{\Xi_{S}}$ is guaranteed by the compactness of
$\mathbb{P}_{N}^{\star}\left(\mathbb{C}\right)$. From a more practical
perspective, the minimization of $j^{\Xi_{S}}$ can be carried out as follows:

  1. (i)

build a $V$-orthonormal basis $\left\\{\varphi_{i}\right\\}_{i=1}^{S}$ of the
span of the $S$ samples $\left\\{u(\mu_{j})\right\\}_{j=1}^{S}$, as well as
the corresponding “component map”

$\mathbf{w}(\mu_{j})=\left(\langle
u(\mu_{j}),\varphi_{i}\rangle_{V}\right)_{i=1}^{S},\quad\text{so that}\quad
u(\mu_{j})=\sum_{i=1}^{S}w_{i}(\mu_{j})\varphi_{i},\quad\text{for
}j=1,\ldots,S\text{;}$

in practice, the basis and the component map can be computed via the Gram-
Schmidt procedure or by a generalized QR decomposition of the snapshots;

  2. (ii)

fix a basis $\left\\{\psi_{l}\right\\}_{l=0}^{N}$ of
$\mathbb{P}_{N}\left(\mathbb{C}\right)$, orthonormal with respect to the
scalar product associated to
${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}$,
and employ it to represent any arbitrary
$Q\in\mathbb{P}_{N}^{\star}\left(\mathbb{C}\right)$ by its component vector
$\mathbf{q}$ (by construction, $\left\|\mathbf{q}\right\|_{2}=1$); for a more
detailed discussion on practical choices of the polynomial basis, we refer to
subsection 3.1;

  3. (iii)

build the Gramian factor

$\Psi=\left(\frac{1}{(S-1)!}\frac{\operatorname{d}^{S-1}}{\operatorname{d}{\mu}^{S-1}}I^{\Xi_{S}}(w_{i}\psi_{l})\right)_{i=1,\ldots,S,\,l=0,\ldots,N}\in\mathbb{C}^{S\times(N+1)}\text{,}$

which, by linearity, allows to express the target functional as a quadratic
form:

$\displaystyle j^{\Xi_{S}}(Q)^{2}=$
$\displaystyle\frac{1}{(S-1)!^{2}}\langle\frac{\operatorname{d}^{S-1}}{\operatorname{d}{\mu}^{S-1}}I^{\Xi_{S}}(uQ),\frac{\operatorname{d}^{S-1}}{\operatorname{d}{\mu}^{S-1}}I^{\Xi_{S}}(uQ)\rangle_{V}$
$\displaystyle=$
$\displaystyle\sum_{i,i^{\prime}=1}^{S}\sum_{l,l^{\prime}=0}^{N}\frac{1}{(S-1)!^{2}}q_{l}\frac{\operatorname{d}^{S-1}}{\operatorname{d}{\mu}^{S-1}}I^{\Xi_{S}}(w_{i}\psi_{l})\overline{q_{l^{\prime}}\frac{\operatorname{d}^{S-1}}{\operatorname{d}{\mu}^{S-1}}I^{\Xi_{S}}(w_{i^{\prime}}\psi_{l^{\prime}})}\langle\varphi_{i},\varphi_{i^{\prime}}\rangle_{V}$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{S}\sum_{l,l^{\prime}=0}^{N}q_{l}\overline{q_{l^{\prime}}}\,\Psi_{il}\overline{\Psi_{il^{\prime}}}=\mathbf{q}^{H}\Psi^{H}\Psi\mathbf{q}$
(5)

($\overline{\cdot}$ and $\cdot^{H}$ denote conjugation and conjugate
transposition, respectively);

  4. (iv)

by SVD, obtain the required minimizer, as the minimal right singular vector of
the Gramian factor $\Psi$.

We remark that introducing the orthogonal bases
$\left\\{\varphi_{i}\right\\}_{i=1}^{S}$ and
$\left\\{\psi_{l}\right\\}_{l=0}^{N}$ is not strictly necessary for the
minimization of $j^{\Xi_{S}}$. However, they allow to state the SVD problem in
the usual Euclidean norm, and to reduce it to size $S$, independently of the
dimension of $V$.

###  2.1 Parametric problem framework

While the map $u$ introduced above may be a very general Hilbert-space valued
function, we are often interested in approximating maps of a quite specific
type. Indeed, assume that we are given a parametric problem:

$\text{for }\mu\in K\subset\mathbb{C},\;\text{find
}u(\mu)\in\operatorname{Dom}(\mathcal{F}_{\mu})\;\text{s.t.}\;\mathcal{F}_{\mu}(u(\mu))=0\text{,}$
(6)

with $\mathcal{F}_{\mu}:V\supset\operatorname{Dom}(\mathcal{F}_{\mu})\to V$ a
family of densely defined operators, depending continuously on $\mu$. Provided
the set $K^{\star}\subset K$ of values of $\mu$, for which a unique $u(\mu)$
solving Eq. 6 exists, is non-empty, we focus on the _solution map_ associated
to Eq. 6:

$\displaystyle u:K^{\star}$ $\displaystyle\to V$ (7) $\displaystyle\mu$
$\displaystyle\mapsto u(\mu)\text{ solution of
}\lx@cref{creftype~refnum}{eq:parametric_problem}\text{.}$

In the following, we restrict our analysis to a special class of solution
maps, namely meromorphic $V$-valued functions over $\mathbb{C}$ with simple
poles:

$u(\mu)=\sum_{\lambda\in\Lambda}\frac{v_{\lambda}}{\lambda-\mu},\quad\mu\in\mathbb{C}\text{.}$
(8)

We denote by $\\{v_{\lambda}\\}_{\lambda\in\Lambda}\subset V$ the set of
_residues_ , which we assume to be $V$-square-summable:
$\sum_{\lambda\in\Lambda}\left\|v_{\lambda}\right\|_{V}^{2}<\infty$. Moreover,
we assume that the set $\Lambda$ of _poles_ , or _resonances_ , is countable,
with no finite limit point. If $\Lambda$ is not finite, Eq. 8 has to be
understood in the sense of $V$-convergence, for an arbitrary ordering of
$\Lambda$. Without loss of generality, we will assume that $\Lambda$ is
unbounded: if this is not the case, one can always extend $\Lambda$ by adding
(countably) infinitely many fictitious poles at $\infty$ with arbitrary
corresponding residues.

The bulk of our results will be obtained under the additional assumption that
the residues form an orthogonal family, i.e. that $\langle
v_{\lambda},v_{\lambda^{\prime}}\rangle_{V}=0$ whenever
$\lambda\neq\lambda^{\prime}$. This allows considerable simplifications in our
derivation: in particular, orthogonal residues have the great advantage of
introducing no “correlation” (in the $V$-topology) between components of the
solution corresponding to different resonances. Still, many of our conclusions
can be generalized to the non-orthogonal case: we defer a discussion of such
extensions to subsection 3.4.

To motivate our interest in this class of functions, we give here some
examples of operator families $\\{\mathcal{F}_{\mu}\\}_{\mu\in K}$ which lead
to meromorphic solution maps; for their approximation, it is reasonable to
employ rational functions, as minimal rational interpolants do. However, we
remark that the theory that we develop might not necessarily apply to all the
cases which we list here.

  * •

Frequency-domain dynamical systems: $\mathcal{F}_{\mu}(u)=(A-\mu E)u+B$, with
$A,E\in\mathbb{C}^{n\times n}$, and $B\in V=\mathbb{C}^{n\times n_{I}}$ (one
column per input). By employing a determinantal representation [5], one can
easily show that the solution map is a rational function; in particular, if
the pencil $(A,E)$ admits $n$ (not necessarily finite) eigenvalues counting
multiplicity, then $u$ is of the form Eq. 8 [4]. Depending on the spectral
properties of $A$, one may be able to ensure residue orthogonality: for
instance, it suffices for $E^{\dagger}A$ to be normal, where † denotes the
Moore-Penrose pseudo-inverse.

  * •

Polynomial eigenvalue-like problems: $\mathcal{F}_{\mu}(u)=A(\mu)u+B$, where
$A(\mu)\in\mathbb{C}^{n\times n}$ is a polynomial in $\mu$, and $B\in
V=\mathbb{C}^{n\times n_{I}}$. Using augmentation-based approaches [4, 24],
one can cast the problem as a frequency-domain dynamical system of size $nM$,
where $M$ is the degree of $A$. This implies that the solution map of the
original problem (of size $n$) is a rational function; in order to check
whether it is actually of the form Eq. 8, one has to investigate the spectral
properties of the augmented $(nM)\times(nM)$ problem.

  * •

Helmholtz-type problems:
$\mathcal{F}_{\mu}(u)=(\mathcal{L}-\mu\mathcal{I})u-f$, where
$\mathcal{L}:V\supset\operatorname{Dom}(\mathcal{L})\to V$ is densely defined
and bijective, with compact and normal inverse, and $f\in V$ (e.g.
$V=L^{2}(\Omega)$, $\mathcal{L}=-\Delta$, and
$\operatorname{Dom}(\mathcal{L})=\Delta^{-1}(V)$). The solution map is of the
form Eq. 8, and the normality of the resolvent $\mathcal{L}^{-1}$ ensures
orthogonal residues [10].

  * •

Analytic compact perturbations of bounded bijective operators:
$\mathcal{F}_{\mu}(u)=(\mathcal{J}-\mathcal{L}(\mu))u-f$, where
$\mathcal{J}:V\to V$ is bounded and bijective, $\mathcal{L}(\mu):V\to V$ is
compact for all $\mu\in\mathbb{C}$, and analytic with respect to $\mu$, and
$f\in V$; a proof of meromorphy is in [37]. Unfortunately, the form of the
solution map is, in general, not (8): for instance, it cannot be guaranteed
that the poles are simple, i.e. a $\lambda$-dependent exponent may appear in
the denominator of (8).

##  3 Convergence theory

In this section, we investigate the convergence properties of minimal rational
interpolants as the number of samples $S$ increases. To this aim, it becomes
necessary to introduce some assumptions on the asymptotic properties of the
sample set $\Xi_{S}$. In order to streamline the derivation, we choose to
follow a potential theory-inspired approach.

###  3.1 Preliminaries for convergence theory

Fix $K\subset\mathbb{C}$ either as the closure of the finite region bounded by
a Jordan curve or a line segment with positive length. There exist [38] a
(unique) positive real number $\operatorname{Cap}(K)$ and a continuous
bijective conformal mapping

$\phi_{K}:\overline{\mathbb{C}\setminus
K}\to\mathbb{C}\setminus\mathcal{B}(0,\operatorname{Cap}(K))=\\{z\in\mathbb{C},\left|z\right|\geq\operatorname{Cap}(K)\\}\text{,}$

holomorphic outside $K$, such that
$\left|\phi_{K}\right|=\operatorname{Cap}(K)$ over $\partial K$ and
$\lim_{\mu\to\infty}\left|\phi_{K}^{\prime}(\mu)\right|=1$.

The value $\operatorname{Cap}(K)$ is the _logarithmic capacity_ of $K$. It is
a powerful concept in complex analysis, algebraic geometry, and approximation
theory, and we refer to the monograph [27] for a thorough introduction to the
subject. In its more general framework, logarithmic capacity is well defined
for any compact set in the complex plane, and coincides with the definition
above under our assumptions.

Among its properties we remark some bounds for other metrics over complex
sets: for any compact set $A\subset\mathbb{C}$ we have

$\sqrt{\frac{\ell(A)}{\pi}}\leq\operatorname{Cap}(A)\leq\frac{1}{2}\max_{z,w\in
A}\left|z-w\right|\text{,}$ (9)

with $\ell$ denoting the Lebesgue measure over $\mathbb{C}$. We remark that
zero-measure sets are not guaranteed to have capacity zero: for instance,
$\operatorname{Cap}([0,4])=1$.

A crucial quantity for stating convergence results in approximation theory is
the Green’s potential, which is closely related to the complex magnitude of
$\phi_{K}$:

$\Phi_{K}(\mu)=\begin{cases}\left|\phi_{K}(\mu)\right|\quad&\text{if
}\mu\notin K,\\\ \operatorname{Cap}(K)\quad&\text{if }\mu\in
K\end{cases}\text{.}$ (10)

In particular, we use the Green’s potential to induce an ordering in the set
of poles $\Lambda=\\{\lambda_{j}\\}_{j=1,2,\ldots}$:

$\Phi_{K}(\lambda_{1})\leq\Phi_{K}(\lambda_{2})\leq\ldots\text{.}$ (11)

Moreover, we require that $\Phi_{K}(\lambda_{N+1})>\operatorname{Cap}(K)$,
i.e. that the number of poles within $K$ is at most $N$, see Eq. 10 and Eq.
11. Indeed, we will show in the following that minimal rational interpolants
fail to differentiate between poles of $u$ within $K$; hence, without the
assumption $\Phi_{K}(\lambda_{N+1})>\operatorname{Cap}(K)$, the
“unidentifiable” poles would distort the approximant and, in general, prevent
convergence to the target solution map.

Now, for a given finite point set $\Xi\subset\mathbb{C}$, we define the
corresponding nodal polynomial as
$\omega^{\Xi}(\mu)=\prod_{\mu^{\prime}\in\Xi}\left(\mu-\mu^{\prime}\right)$.
This allows us to state the result which motivates our interest in Green’s
potentials.

######  Theorem 2 (Fejér points [38]).

Let $K$ be either the closure of the finite region bounded by a Jordan curve
or a line segment with positive length. There exists a sequence of sample sets
$\big{(}\Xi_{S}^{K}\big{)}_{S=1,2,\ldots}$, whose elements are referred to as
the Fejér points for $K$, such that

$\displaystyle\limsup_{S\to\infty}\left|\omega^{\Xi_{S}^{K}}\right|^{1/S}\leq\Phi_{K}\quad\text{uniformly
over compact subsets of }\mathbb{C}\text{,}$ (12)
$\displaystyle\lim_{S\to\infty}\left|\omega^{\Xi_{S}^{K}}\right|^{1/S}=\Phi_{K}\quad\text{uniformly
over compact subsets of }\mathbb{C}\setminus\partial K\text{.}$ (13)

This means that the Green’s potential provides an asymptotic upperbound for
the magnitude of the nodal polynomial, as long as the sample points are chosen
in an optimal way. From here onward we assume that, given $S\in\mathbb{N}$,
the sample set $\Xi_{S}$ coincides with some set of Fejér points
$\Xi_{S}^{K}$, allowing us to drop the superscript.

######  Remark 3.

Fejér points are explicitly known only for very simple sets $K$; this is the
case, for instance, in two situations which are very common in applications:
the Fejér points of the unit disk are the roots of unity, and for a line
segment they are the Chebyshev nodes [38]. However, our theory applies even
for more general sampling strategies. The technical requirement is the
following: there exists $\Phi:\mathbb{C}\to\mathbb{R}$ such that

  * •

$\sup_{K}\Phi\leq\inf_{\mathbb{C}\setminus K}\Phi$,

  * •

$\limsup_{S\to\infty}\left|\omega^{\Xi_{S}}\right|^{1/S}\leq\Phi$ uniformly
over $K$,

  * •

$\lim_{S\to\infty}\left|\omega^{\Xi_{S}}\right|^{1/S}=\Phi$ uniformly over
compact subsets of $\mathbb{C}\setminus K$.

If this is the case, then one can replace the Green’s potential by $\Phi$ in
the entirety of our discussion; in particular, $\Phi$ should induce the pole
ordering, cf. Eq. 11.

For instance, one could take any sequence which approximates asymptotically
the Fejér points [13]. Another interesting case is that of repeated points,
employing Lagrange-Hermite interpolation: given
$\\{\mu_{i}^{0}\\}_{i=1}^{L}\subset\mathbb{C}$, one defines the sample points
as $\mu_{nL+i}=\mu_{i}^{0}$ for all $n=0,1,\ldots$ and $i=1,\ldots,L$. Then,
as long as $K$ is defined as the discrete set $K=\\{\mu_{i}^{0}\\}_{i=1}^{L}$,
the three conditions above hold with
$\Phi(\mu)=\prod_{i=1}^{L}\left|\mu-\mu_{i}^{0}\right|^{1/L}$. In particular,
this allows to analyze fast LS Padé approximants [10] as a special case of
minimal rational interpolants, with $L=1$ and
$\Phi(\mu)=\left|\mu-\mu_{1}^{0}\right|$.

A second class of considerations involves the denominator space
$\mathbb{P}_{N}^{\star}\left(\mathbb{C}\right)$, whose definition hinges upon
the choice of the norm
${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}$.
Depending on the type of result we wish to prove, we may need to introduce the
following assumption on the scaling of
${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}$
with respect to $N$.

######  Assumption 4.

There exist $\mu_{0}\in\mathbb{C}$ and three positive constants $R_{\mu_{0}}$,
$c_{\mu_{0}}$, and $C_{\mu_{0}}$, all independent of $N$, such that

$c_{\mu_{0}}\prod_{j=1}^{D}\left|\mu_{0}-z_{j}\right|\leq{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\prod_{j=1}^{D}\left(\,\cdot-
z_{j}\right)\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}\leq
C_{\mu_{0}}\prod_{j=1}^{D}\left(R_{\mu_{0}}+\left|\mu_{0}-z_{j}\right|\right)\text{,}$

for all $N\in\mathbb{N}$ and $\\{z_{j}\\}_{j=1}^{D}\in\mathbb{C}^{D}$, with
$0\leq D\leq N$.

In [10] it was shown that Assumption 4 holds (with
$R_{\mu_{0}}=c_{\mu_{0}}=C_{\mu_{0}}=1$) for the “monomial” norm
${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}={\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{\mu_{0},N}$
defined as

${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|Q\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{\mu_{0},N}=\sqrt{\sum_{j=0}^{N}\left|q_{j}\right|^{2}},\quad\text{where}\quad
Q=\sum_{j=0}^{N}q_{j}\left(\,\cdot-\mu_{0}\right)^{j}\text{.}$ (14)

We describe in Appendix A a broader framework where Assumption 4 can still be
shown to hold, with the (often ill-conditioned) monomials in Eq. 14 being
replaced by a more general polynomial basis.

Many of our results do not rely on Assumption 4, but exploit the following
weaker result.

######  Lemma 5.

For all $\mu_{0}\in\mathbb{C}$, there exist positive constants $R_{\mu_{0}}$
(independent of $N$), $c_{\mu_{0},N}$, and $C_{\mu_{0},N}$, such that

$c_{\mu_{0},N}\prod_{j=1}^{D}\left|\mu_{0}-z_{j}\right|\leq{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\prod_{j=1}^{D}\left(\,\cdot-
z_{j}\right)\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}\leq
C_{\mu_{0},N}\prod_{j=1}^{D}\left(R_{\mu_{0}}+\left|\mu_{0}-z_{j}\right|\right)\text{,}$

for all $N\in\mathbb{N}$ and $\\{z_{j}\\}_{j=1}^{D}\in\mathbb{C}^{D}$, with
$0\leq D\leq N$.

######  Proof.

As described above,
${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}={\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{\mu_{0},N}$
satisfies Lemma 5. The result follows by the equivalence of all norms on a
finite-dimensional space.

###  3.2 Convergence of minimal rational interpolant poles

In the present section we prove that the denominator $Q_{N}^{\Xi_{S}}$
described in Definition 1 performs well in identifying the approximate
position of some of the poles $\Lambda$ of the solution map, among which those
within the parameter domain $K$. As a preliminary step, it is useful to
exploit the definition of $Q_{N}^{\Xi_{S}}$ to derive a bound for the minimal
value of $j^{\Xi_{S}}$ over $\mathbb{P}_{N}^{\star}\left(\mathbb{C}\right)$.

######  Lemma 6.

Let $Q\in\mathbb{P}_{N}\left(\mathbb{C}\right)$. Then

$j^{\Xi_{S}}(Q)^{2}=\sum_{\lambda\in\Lambda}\frac{\left\|v_{\lambda}\right\|_{V}^{2}}{\left|\omega^{\Xi_{S}}(\lambda)\right|^{2}}\left|Q(\lambda)\right|^{2}\text{.}$

Moreover, for any subset $\Lambda_{0}$ of $\Lambda$ with cardinality at most
$N$, we have

$j^{\Xi_{S}}(Q_{N}^{\Xi_{S}})\leq\frac{\left(\sum_{\lambda\in\Lambda\setminus\Lambda_{0}}\left\|v_{\lambda}\right\|_{V}^{2}\right)^{1/2}}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\omega^{\Lambda_{0}}\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}}\sup_{\Lambda\setminus\Lambda_{0}}\left|\frac{\omega^{\Lambda_{0}}}{\omega^{\Xi_{S}}}\right|\text{.}$
(15)

In particular, for any $\varepsilon>0$, there exists $S_{\varepsilon,N}$ such
that

$j^{\Xi_{S}}(Q_{N}^{\Xi_{S}})\leq(1+\varepsilon)^{S}\frac{\left(\sum_{j>N}\left\|v_{\lambda_{j}}\right\|_{V}^{2}\right)^{1/2}}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\prod_{j=1}^{N}\left(\
\cdot-\lambda_{j}\right)\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}}\frac{\prod_{j=1}^{N}\left|\lambda_{N+1}-\lambda_{j}\right|}{\Phi_{K}(\lambda_{N+1})^{S}}\quad\text{for
}S>S_{\varepsilon,N}\text{,}$ (16)

provided $\Phi_{K}(\lambda_{N+1})<\Phi_{K}(\lambda_{N+2})$. (If this is not
true, there exists a reordering of the poles, still satisfying Eq. 11, for
which Eq. 16 holds.)

######  Proof.

The interpolation operator can be expressed in barycentric coordinates [8] as

${I^{\Xi_{S}}(\phi)(\mu)=\sum_{\mu^{\prime}\in\Xi_{S}}\frac{\phi(\mu^{\prime})\,\omega^{\Xi_{S}}(\mu)}{(\mu-\mu^{\prime})\left.\omega^{\Xi_{S}}\right.^{\prime}(\mu^{\prime})}\text{.}}$
(17)

Accordingly, since $Q$ is interpolated exactly,

$\displaystyle\frac{1}{(S-1)!}\frac{\operatorname{d}^{S-1}}{\operatorname{d}{\mu}^{S-1}}I^{\Xi_{S}}(uQ)=$
$\displaystyle\sum_{\mu^{\prime}\in\Xi_{S}}\frac{u(\mu^{\prime})Q(\mu^{\prime})}{\left.\omega^{\Xi_{S}}\right.^{\prime}(\mu^{\prime})}=\sum_{\lambda\in\Lambda}\sum_{\mu^{\prime}\in\Xi_{S}}\frac{v_{\lambda}Q(\mu^{\prime})}{(\lambda-\mu^{\prime})\left.\omega^{\Xi_{S}}\right.^{\prime}(\mu^{\prime})}$
$\displaystyle=$
$\displaystyle\sum_{\lambda\in\Lambda}\frac{v_{\lambda}}{\omega^{\Xi_{S}}(\lambda)}I^{\Xi_{S}}(Q)(\lambda)=\sum_{\lambda\in\Lambda}\frac{v_{\lambda}}{\omega^{\Xi_{S}}(\lambda)}Q(\lambda)\text{.}$

The first claim follows by orthonormality of
$\\{v_{\lambda}\\}_{\lambda\in\Lambda}$. From here, we remark that, by
construction, $j^{\Xi_{S}}(Q_{N}^{\Xi_{S}})\leq j^{\Xi_{S}}\left(Q\right)$ for
all $Q\in\mathbb{P}_{N}^{\star}\left(\mathbb{C}\right)$, so that, by choosing
$Q=\omega^{\Lambda_{0}}/{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\omega^{\Lambda_{0}}\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}$,
we obtain

$j^{\Xi_{S}}(Q_{N}^{\Xi_{S}})^{2}\leq
j^{\Xi_{S}}\left(\frac{\omega^{\Lambda_{0}}}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\omega^{\Lambda_{0}}\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}}\right)^{2}=\sum_{\lambda\in\Lambda\setminus\Lambda_{0}}\frac{\left\|v_{\lambda}\right\|_{V}^{2}}{\left|\omega^{\Xi_{S}}(\lambda)\right|^{2}}\frac{\left|\omega^{\Lambda_{0}}(\lambda)\right|^{2}}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\omega^{\Lambda_{0}}\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}^{2}}\text{,}$
(18)

from which Eq. 15 follows.

Now, let $\varepsilon>0$ and $N$ be fixed, and define
$\Lambda_{0}=\\{\lambda_{j}\\}_{j=1}^{N}$. It remains to find a bound for the
supremum of
$r^{\Lambda_{0},\Xi_{S}}(\lambda)=|\omega^{\Lambda_{0}}(\lambda)/\omega^{\Xi_{S}}(\lambda)|$
over $\lambda\in\Lambda\setminus\Lambda_{0}$. To this aim, let
$\lambda^{\prime}_{S}$ maximize $r^{\Lambda_{0},\Xi_{S}}$ over
$\Lambda\setminus\Lambda_{0}$ (for all $S$, the supremum must be attained
since $r^{\Lambda_{0},\Xi_{S}}(\lambda^{\prime})$ converges to 0 as
$\left|\lambda^{\prime}\right|$ diverges to infinity). The set of maximizers
$\\{\lambda^{\prime}_{S}\\}_{S>N}$ must be bounded, due to the polynomial
degree of $\omega^{\Lambda_{0}}$ being fixed, and smaller than that of
$\omega^{\Xi_{S}}$. Thus, by Eq. 13 there exists
$S^{\prime}=S^{\prime}(\varepsilon,N,\Lambda,K)$ such that, for
$S>S^{\prime}$,

$\sup_{\Lambda\setminus\Lambda_{0}}r^{\Lambda_{0},\Xi_{S}}=\sup_{\\{\lambda^{\prime}_{S}\\}_{S>N}}r^{\Lambda_{0},\Xi_{S}}\leq(1+\varepsilon)^{S}\sup_{\\{\lambda^{\prime}_{S}\\}_{S>N}}\frac{|\omega^{\Lambda_{0}}|}{\Phi_{K}^{S}}\leq(1+\varepsilon)^{S}\sup_{\Lambda\setminus\Lambda_{0}}\frac{|\omega^{\Lambda_{0}}|}{\Phi_{K}^{S}}\text{.}$

Let
$\Phi_{K}(\lambda_{N+1})=\Phi_{K}(\lambda_{N+\nu})<\Phi_{K}(\lambda_{N+\nu+1})$
for some $\nu\geq 1$, and assume that $\lambda_{N+1}$ maximizes
$|\omega^{\Lambda_{0}}(\lambda)|$ over
$\lambda\in\\{\lambda_{N+k}\\}_{k=1}^{\nu}$. Then, since
$|\omega^{\Lambda_{0}}|$ is independent of $S$, the last supremum is attained
at $\lambda_{N+1}$, provided $S$ is large enough. This yields Eq. 16.

Before proving convergence of the roots of the approximate denominator to
poles of the solution map, we show as a preliminary result that
$Q_{N}^{\Xi_{S}}$ takes small values close to some of the elements of
$\Lambda$.

######  Lemma 7.

Let $N$ be fixed. For any $\varepsilon>0$, there exists $S_{\varepsilon,N}$
such that

$\left|Q_{N}^{\Xi_{S}}\left(\lambda_{j}\right)\right|\leq
C_{j}(1+\varepsilon)^{4S}\left(\frac{\Phi_{K}(\lambda_{j})}{\Phi_{K}(\lambda_{N+1})}\right)^{2S}\quad\text{for
}S>S_{\varepsilon,N}\text{,}$ (19)

for $j=1,\ldots,N$, with $C_{j}$ independent of $S$ and $\varepsilon$.

######  Proof.

Through some isometry, we identify $\mathbb{P}_{N}\left(\mathbb{C}\right)$
with $\mathbb{C}^{N+1}$ (endowed with the Euclidean norm), so that, in
particular, $\mathbb{P}_{N}^{\star}\left(\mathbb{C}\right)$ corresponds to the
boundary of the unit sphere
$\partial\mathcal{B}_{N+1}(\bm{0},1)\subset\mathbb{C}^{N+1}$.

By duality, there exists a family
$\\{\bm{\omega}_{\lambda}\\}_{\lambda\in\Lambda}\subset\mathbb{C}^{N+1}$ such
that

$\mathbb{P}_{N}\left(\mathbb{C}\right)\ni
Q\leftrightsquigarrow\mathbf{q}\in\mathbb{C}^{N+1}\quad\text{implies}\quad\bm{\omega}_{\lambda}^{H}\mathbf{q}=Q(\lambda)\;\text{for
}\lambda\in\Lambda\text{.}$

Let $\mathbf{q}_{N}^{\Xi_{S}}\in\partial\mathcal{B}_{N+1}(\bm{0},1)$ be the
identification of $Q_{N}^{\Xi_{S}}$, and consider the Hermitian matrices
$G^{\Xi_{S}},G_{N}^{\Xi_{S}}\in\mathbb{C}^{(N+1)\times(N+1)}$ defined as

$G^{\Xi_{S}}=\sum_{\lambda\in\Lambda}\frac{\left\|v_{\lambda}\right\|_{V}^{2}}{\left|\omega^{\Xi_{S}}(\lambda)\right|^{2}}\bm{\omega}_{\lambda}\bm{\omega}_{\lambda}^{H}\quad\text{and}\quad
G_{N}^{\Xi_{S}}=\sum_{j=1}^{N}\frac{\left\|v_{\lambda_{j}}\right\|_{V}^{2}}{\left|\omega^{\Xi_{S}}(\lambda_{j})\right|^{2}}\bm{\omega}_{\lambda_{j}}\bm{\omega}_{\lambda_{j}}^{H}\text{,}$

From here onward, the proof resembles quite closely that of Lemma 5.4 in [10],
being essentially based on deriving a bound for the smallest singular value of
$G^{\Xi_{S}}$ and on finding a rank-$N$ decomposition of $G_{N}^{\Xi_{S}}$.
The final result that can be obtained is of the form: for $j=1,\ldots,N$, for
large $S$,

$\left|Q_{N}^{\Xi_{S}}\left(\lambda_{j}\right)\right|=\left|\bm{\omega}_{\lambda_{j}}^{H}\mathbf{q}_{N}^{\Xi_{S}}\right|\leq
C_{j}(1+\varepsilon)^{2S}\frac{\left|\omega^{\Xi_{S}}(\lambda_{j})\right|^{2}}{\Phi_{K}(\lambda_{N+1})^{2S}}\text{,}$

with $C_{j}$ independent of $S$ and $\varepsilon$. From here it suffices to
observe that $\left|\omega^{\Xi_{S}}(\lambda_{j})\right|$ is bounded by
$\left((1+\varepsilon)\Phi_{K}(\lambda_{j})\right)^{S}$, thanks to Eq. 12.

Finally, we are able to show that the poles of minimal rational approximants
converge, as $S$ increases, to some of the elements of $\Lambda$, namely the
$N$ poles which are the “closest” according to Eq. 11.

######  Theorem 8.

Let $N$ be fixed, and denote by $\\{\lambda_{j}^{\Xi_{S}}\\}_{j=1}^{N}$ the
roots of $Q_{N}^{\Xi_{S}}$ (in case of deficient degree, we assume the missing
roots to be $\infty$). For any $\varepsilon>0$, there exists
$S_{\varepsilon,N}$ such that

$\min_{i=1,\ldots,N}\left|\lambda_{i}^{\Xi_{S}}-\lambda_{j}\right|\leq
C_{j}(1+\varepsilon)^{4S}\left(\frac{\Phi_{K}(\lambda_{j})}{\Phi_{K}(\lambda_{N+1})}\right)^{2S}\quad\text{for
}S>S_{\varepsilon,N}\text{,}$ (20)

for $j=1,\ldots,N$, with $C_{j}$ independent of $S$ and $\varepsilon$.

######  Proof.

We start by observing that, by the normalization of $Q_{N}^{\Xi_{S}}$, it must
hold

$\left|Q_{N}^{\Xi_{S}}(\lambda)\right|=\frac{\prod_{j=1}^{S}\left|\lambda-\lambda_{j}^{\Xi_{S}}\right|}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\prod_{j=1}^{S}\left(\,\cdot-\lambda_{j}^{\Xi_{S}}\right)\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}}\geq\frac{1}{C_{\mu_{0},N}}\prod_{j=1}^{N}\frac{\left|\lambda-\lambda_{j}^{\Xi_{S}}\right|}{R_{\mu_{0}}+\left|\lambda-\mu_{0}\right|+\left|\lambda-\lambda_{j}^{\Xi_{S}}\right|}\\\
$ (21)

for any arbitrary $\mu_{0}\in\mathbb{C}$, see Lemma 5.

From here it suffices to apply the same proof as for Theorem 5.5 in [10], with
Eq. 19 replacing the corresponding bound for fast LS-Padé denominators. The
main steps are: (i) use Eq. 19 and Eq. 21 to show that, for each
$j=1,\ldots,N$, there is a sequence of approximate poles
$\\{\lambda_{i(S)}^{\Xi_{S}}\\}_{S>S_{\varepsilon,N}}$ which converges to
$\lambda_{j}$; (ii) exploiting the fact that the poles of $u$ are simple,
conclude that $Q_{N}^{\Xi_{S}}$ cannot have, for large $S$, two roots too
close to each other; (iii) apply Eq. 19 locally to derive Eq. 20.

Going back to our original problem formulation, let $0\leq\nu\leq N$ be such
that
$\operatorname{Cap}(K)=\Phi_{K}(\lambda_{\nu})<\Phi_{K}(\lambda_{\nu+1})$.
Theorem 8 shows that, for fixed $N$ and large enough $S$, minimal rational
interpolants approximate well the $\nu$ poles $\\{\lambda_{j}\\}_{j=1}^{\nu}$
within $K$, with (exactly) one root of $Q_{N}^{\Xi_{S}}$ converging to each of
them at rate

$\mathcal{O}\left(\left(\operatorname{Cap}(K)\,\big{/}\,\Phi_{K}(\lambda_{N+1})\right)^{2S}\right)\text{.}$
(22)

####  3.2.1 Pole convergence with respect to denominator degree

The statements of Lemma 7 and Theorem 8 may seem somewhat limited due to the
fact that the denominator degree $N$ enters the convergence bounds in a
complex way. This makes it impossible to determine how the pole approximation
error behaves if $N$ is increased together with $S$, namely if it diverges to
infinity. In this section we try to solve this shortcoming. The main result we
prove is the following.

######  Theorem 9.

Take two diverging sequences of integers $(N_{k})_{k\in\mathbb{N}}$ and
$(S_{k})_{k\in\mathbb{N}}$, increasing and strictly increasing respectively,
such that $N_{k}<S_{k}$ for all $k$. Let $\Xi^{(k)}$ be the shorthand for
$\Xi_{S_{k}}$, so that $\smash{Q_{N_{k}}^{\Xi^{(k)}}}$ is the denominator of
the $\smash{[S_{k}-1/N_{k}]}$ minimal rational interpolant computed from
samples of $u$ at $\Xi_{S_{k}}$. We denote by
$\smash{\\{\lambda_{j}^{(k)}\\}_{j=1}^{N_{k}}}$ the roots of
$\smash{Q_{N_{k}}^{\Xi^{(k)}}}$. Also, let Assumption 4 be satisfied. For all
$j=1,2,\ldots$, we have

$\lim_{k\to\infty}\min_{i=1,\ldots,N_{k}}\left|\lambda_{i}^{(k)}-\lambda_{j}\right|=0\text{.}$

######  Proof.

Since we are interested in asymptotic properties of the poles, in the present
proof we employ a slightly different ordering of the resonances
$\Lambda=\\{\sigma_{j}\\}_{j=1,2,\ldots}$, such that

$\left|\sigma_{1}-\mu_{0}\right|\leq\left|\sigma_{2}-\mu_{0}\right|\leq\ldots\text{.}$
(23)

Let $j\in\mathbb{N}$ be arbitrary. Thanks to Assumption 4, we can obtain the
lower bound

$\displaystyle\left|Q_{N_{k}}^{\Xi^{(k)}}(\sigma_{j})\right|=\frac{\prod_{i=1}^{N_{k}}\left|\sigma_{j}-\lambda_{i}^{(k)}\right|}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\prod_{i=1}^{N_{k}}\left(\,\cdot-\lambda_{i}^{(k)}\right)\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N_{k}}}\geq$
$\displaystyle\frac{1}{C_{\mu_{0}}}\prod_{i=1}^{N_{k}}\frac{\left|\sigma_{j}-\lambda_{i}^{(k)}\right|}{R_{\mu_{0}}+\left|\sigma_{j}-\mu_{0}\right|+\left|\sigma_{j}-\lambda_{i}^{(k)}\right|}$
$\displaystyle=$
$\displaystyle\frac{1}{C_{\mu_{0}}}\prod_{i=1}^{N_{k}}\varphi_{j}\left(\left|\sigma_{j}-\lambda_{i}^{(k)}\right|\right)\text{,}$

with $\varphi_{j}=\cdot/(R_{\mu_{0}}+\left|\sigma_{j}-\mu_{0}\right|+\cdot)$
non-negative and strictly increasing over $\mathbb{R}^{+}$. On the other hand,
due to Assumption 4 and Lemma 6 with
$\Lambda_{0}=\\{\sigma_{i}\\}_{i=1}^{N_{k}}$, one has

$\displaystyle\left|Q_{N_{k}}^{\Xi^{(k)}}(\sigma_{j})\right|\leq$
$\displaystyle\frac{\left|\omega^{\Xi^{(k)}}(\sigma_{j})\right|}{\left\|v_{\sigma_{j}}\right\|_{V}}j^{\Xi^{(k)}}(Q_{N_{k}}^{\Xi^{(k)}})$
$\displaystyle\leq$
$\displaystyle\frac{\left|\omega^{\Xi^{(k)}}(\sigma_{j})\right|}{\left\|v_{\sigma_{j}}\right\|_{V}}\frac{\left(\sum_{\lambda\in\Lambda\setminus\Lambda_{0}}\left\|v_{\lambda}\right\|_{V}^{2}\right)^{1/2}}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\prod_{i=1}^{N_{k}}\left(\
\cdot-\sigma_{i}\right)\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}}\sup_{l>N_{k}}\frac{\prod_{i=1}^{N_{k}}\left|\sigma_{l}-\sigma_{i}\right|}{\left|\omega^{\Xi^{(k)}}(\sigma_{l})\right|}$
$\displaystyle\leq$
$\displaystyle\frac{\left(\sum_{\lambda\in\Lambda}\left\|v_{\lambda}\right\|_{V}^{2}\right)^{1/2}}{c_{\mu_{0}}\left\|v_{\sigma_{j}}\right\|_{V}}\sup_{l>N_{k}}\left(\left|\frac{\omega^{\Xi^{(k)}}(\sigma_{j})}{\omega^{\Xi^{(k)}}(\sigma_{l})}\right|\prod_{i=1}^{N_{k}}\left|\frac{\sigma_{l}-\sigma_{i}}{\mu_{0}-\sigma_{i}}\right|\right)\text{.}$

The two bounds above can be combined with the monotonicity of $\varphi_{j}$ to
obtain

$\varphi_{j}\left(\min_{i=1,\ldots,N_{k}}\left|\sigma_{j}-\lambda_{i}^{(k)}\right|\right)\leq\left(C_{j}^{\prime}\sup_{l>N_{k}}\left(\left|\frac{\omega^{\Xi^{(k)}}(\sigma_{j})}{\omega^{\Xi^{(k)}}(\sigma_{l})}\right|\prod_{i=1}^{N_{k}}\left|\frac{\sigma_{l}-\sigma_{i}}{\mu_{0}-\sigma_{i}}\right|\right)\right)^{1/N_{k}}\text{,}$
(24)

with
$C_{j}^{\prime}=C_{\mu_{0}}(\sum_{\lambda\in\Lambda}\left\|v_{\lambda}\right\|_{V}^{2})^{1/2}/(c_{\mu_{0}}\left\|v_{\sigma_{j}}\right\|_{V})$.

Thanks to the bijectivity of $\varphi_{j}$, which maps 0 to itself, the claim
follows if we can show that the right hand side of Eq. 24 converges to 0 as
$k$ increases. To this aim, a first step is to remove the supremum, whose
argument can be simplified by introducing the bounds

$\prod_{i=1}^{N_{k}}\left|\frac{\sigma_{l}-\sigma_{i}}{\mu_{0}-\sigma_{i}}\right|\leq\prod_{i=1}^{N_{k}}\frac{\left|\mu_{0}-\sigma_{l}\right|+\left|\mu_{0}-\sigma_{i}\right|}{\left|\mu_{0}-\sigma_{i}\right|}\leq\prod_{i=1}^{N_{k}}\frac{2\left|\mu_{0}-\sigma_{l}\right|}{\left|\mu_{0}-\sigma_{i}\right|}=\frac{2^{N_{k}}\left|\mu_{0}-\sigma_{l}\right|^{N_{k}}}{\prod_{i=1}^{N_{k}}\left|\mu_{0}-\sigma_{i}\right|}$

and

$\displaystyle\left|\frac{\omega^{\Xi^{(k)}}(\sigma_{j})}{\omega^{\Xi^{(k)}}(\sigma_{l})}\right|=\prod_{\mu^{\prime}\in\Xi^{(k)}}\left|\frac{\mu^{\prime}-\sigma_{j}}{\mu^{\prime}-\sigma_{l}}\right|\leq\prod_{\mu^{\prime}\in\Xi^{(k)}}\frac{\left|\mu_{0}-\sigma_{j}\right|+\left|\mu_{0}-\mu^{\prime}\right|}{\left|\left|\mu_{0}-\sigma_{l}\right|-\left|\mu_{0}-\mu^{\prime}\right|\right|}\text{,}$

where the term $\left|\mu_{0}-\mu^{\prime}\right|$ is bounded from above by
$R_{\mu_{0},K}=\max_{K}\left|\,\cdot-\mu_{0}\right|$.

Now, let $k$ be large enough, so that
$\left|\mu_{0}-\sigma_{N_{k}+1}\right|>R_{\mu_{0},K}$. Then the absolute value
in the last bound can be omitted, and it must hold

$\left|\frac{\omega^{\Xi^{(k)}}(\sigma_{j})}{\omega^{\Xi^{(k)}}(\sigma_{l})}\right|\prod_{i=1}^{N_{k}}\left|\frac{\sigma_{l}-\sigma_{i}}{\mu_{0}-\sigma_{i}}\right|\leq(C_{l,k,j,\mu_{0}})^{N_{k}}\prod_{i=1}^{N_{k}}\frac{\left|\mu_{0}-\sigma_{j}\right|+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{i}\right|}\text{,}$
(25)

with

$C_{l,k,j,\mu_{0}}=\frac{2\left|\mu_{0}-\sigma_{l}\right|}{\left|\mu_{0}-\sigma_{l}\right|-R_{\mu_{0},K}}\left(\frac{\left|\mu_{0}-\sigma_{j}\right|+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{l}\right|-R_{\mu_{0},K}}\right)^{S_{k}/N_{k}-1}\text{.}$

Going back to the supremum, it is easy to see that all the terms depending on
$l$ in the right hand side of Eq. 25 are maximized (over $l>N_{k}$) by the
choice $l=N_{k}+1$, thanks to Eq. 23. This, after a suitable rearrangement of
the terms, leads to

$\varphi_{j}\left(\min_{i=1,\ldots,N_{k}}\left|\sigma_{j}-\lambda_{i}^{(k)}\right|\right)\leq\left.C_{j}^{\prime}\right.^{1/N_{k}}C_{N_{k}+1,k,j,\mu_{0}}\prod_{i=1}^{N_{k}}\left(\frac{\left|\mu_{0}-\sigma_{j}\right|+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{i}\right|}\right)^{1/N_{k}}\text{.}$
(26)

All of the factors appearing in Eq. 26, except for the last one, can be easily
shown to be bounded as $k$ (and consequently $N_{k}$) goes to $\infty$. Hence,
it suffices to show that

$\lim_{k\to\infty}\prod_{i=1}^{N_{k}}\left(\frac{\left|\mu_{0}-\sigma_{j}\right|+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{i}\right|}\right)^{1/N_{k}}=0\text{.}$

To this aim, the Stolz-Cesàro theorem [3] can be applied to prove that

$\lim_{N\to\infty}\frac{1}{N}\sum_{i=1}^{N}\log\left(\frac{\left|\mu_{0}-\sigma_{j}\right|+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{i}\right|}\right)=\lim_{N\to\infty}\log\left(\frac{\left|\mu_{0}-\sigma_{j}\right|+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{N+1}\right|}\right)=-\infty\text{,}$

yielding the claim.

We have thus shown that, by minimal rational interpolation, we can approximate
an arbitrary number of poles of $u$, as long as we allow the denominator
degree $N$ to become sufficiently large. However, this more general result
does not provide any information on the rate of convergence, at least not a
meaningful one for typical applications. Indeed, most of the steps in the
proof above introduce gross simplifications, which make it impossible for the
intermediate bounds Eq. 24 and Eq. 26 to be reasonably sharp.

###  3.3 Convergence of minimal rational interpolants

As a complement to the previous results, in this section we investigate the
approximation properties of the whole minimal rational interpolant
$u_{N}^{\Xi_{S}}$, the reference being the map $u$. Similarly to the analysis
of the denominator $Q_{N}^{\Xi_{S}}$, we start by fixing the denominator
degree $N$, and by proving convergence with respect to the number of samples
$S$.

######  Theorem 10.

Fix $N$ and $\varepsilon>0$, and let $E$ be an arbitrary compact subset of
$\\{\mu\in\mathbb{C},\Phi_{K}(\mu)<\Phi_{K}(\lambda_{N+1})\\}\setminus\Lambda$.
There exists $S_{\varepsilon,N,E}$ such that

$\left\|u(\mu)-u_{N}^{\Xi_{S}}(\mu)\right\|_{V}\leq
C_{E}(1+\varepsilon)^{2S}\left(\frac{\Phi_{K}(\mu)}{\Phi_{K}(\lambda_{N+1})}\right)^{S}\quad\text{for
all }\mu\in E,\,S>S_{\varepsilon,N,E}\text{,}$

with $C_{E}$ independent of $S$ and $\varepsilon$.

######  Proof.

We start by deriving a useful identity for the error $u-u_{N}^{\Xi_{S}}$:
since $Q_{N}^{\Xi_{S}}$ is interpolated exactly by $I^{\Xi_{S}}$, we can
exploit identity Eq. 17 to obtain

$\displaystyle Q_{N}^{\Xi_{S}}(\mu)\left(u(\mu)-u_{N}^{\Xi_{S}}(\mu)\right)=$
$\displaystyle
u(\mu)I^{\Xi_{S}}(Q_{N}^{\Xi_{S}})(\mu)-I^{\Xi_{S}}(uQ_{N}^{\Xi_{S}})(\mu)$
(27) $\displaystyle=$
$\displaystyle\sum_{\mu^{\prime}\in\Xi_{S}}\sum_{\lambda\in\Lambda}\frac{Q_{N}^{\Xi_{S}}(\mu^{\prime})\omega^{\Xi_{S}}(\mu)}{(\mu-\mu^{\prime})\left.\omega^{\Xi_{S}}\right.^{\prime}(\mu^{\prime})}\left(\frac{v_{\lambda}}{\lambda-\mu}-\frac{v_{\lambda}}{\lambda-\mu^{\prime}}\right)$
$\displaystyle=$
$\displaystyle{\omega^{\Xi_{S}}(\mu)\sum_{\lambda\in\Lambda}\frac{v_{\lambda}}{\left(\lambda-\mu\right)\omega^{\Xi_{S}}(\lambda)}\underbrace{\sum_{\mu^{\prime}\in\Xi_{S}}\frac{Q_{N}^{\Xi_{S}}(\mu^{\prime})\omega^{\Xi_{S}}(\lambda)}{(\lambda-\mu^{\prime})\left.\omega^{\Xi_{S}}\right.^{\prime}(\mu^{\prime})}}_{=I^{\Xi_{S}}(Q_{N}^{\Xi_{S}})(\lambda)=Q_{N}^{\Xi_{S}}(\lambda)}}\text{.}$
(28)

After dividing by $Q_{N}^{\Xi_{S}}(\mu)$, taking the norm, and applying Lemma
6, we obtain

$\displaystyle\left\|u(\mu)-u_{N}^{\Xi_{S}}(\mu)\right\|_{V}=$
$\displaystyle\left|\frac{\omega^{\Xi_{S}}(\mu)}{Q_{N}^{\Xi_{S}}(\mu)}\right|\left(\sum_{\lambda\in\Lambda}\frac{\left\|v_{\lambda}\right\|_{V}^{2}}{\left|\lambda-\mu\right|^{2}\left|\omega^{\Xi_{S}}(\lambda)\right|^{2}}\left|Q_{N}^{\Xi_{S}}(\lambda)\right|^{2}\right)^{1/2}$
$\displaystyle\leq$
$\displaystyle\left|\frac{\omega^{\Xi_{S}}(\mu)}{Q_{N}^{\Xi_{S}}(\mu)}\right|\frac{1}{\min_{\Lambda}\left|\,\cdot-\mu\right|}\,j^{\Xi_{S}}(Q_{N}^{\Xi_{S}})$
(29) $\displaystyle\leq$
$\displaystyle\left|\frac{\omega^{\Xi_{S}}(\mu)}{Q_{N}^{\Xi_{S}}(\mu)}\right|\frac{1}{\min_{\Lambda}\left|\,\cdot-\mu\right|}(1+\varepsilon)^{S}\frac{C^{\prime}}{\Phi_{K}(\lambda_{N+1})^{S}}\text{,}$

for $S$ large enough. In particular, the constant $C^{\prime}$ is independent
of $S$ and $\varepsilon$.

Now, due to Eq. 12, the term $\left|\omega^{\Xi_{S}}\right|$ is bounded by
$(1+\varepsilon)^{S}\Phi_{K}^{S}$ for $S$ large enough. Thus, it just remains
to show that $|Q_{N}^{\Xi_{S}}(\mu)|\min_{\Lambda}\left|\,\cdot-\mu\right|$ is
bounded away from 0 for $\mu\in E$. Since $E$ is a compact subset of
$\mathbb{C}\setminus\Lambda$, the second factor is not troublesome. The first
term is more problematic, since $Q_{N}^{\Xi_{S}}$ may have roots inside $E$.
However, Theorem 8 and the triangular inequality ensure that, for $S$ large
enough, each root of $Q_{N}^{\Xi_{S}}$ has distance from $E$ bounded away from
0. This yields the claim.

If we restrict our interest to compact subsets $E$ of $K\setminus\Lambda$, we
can obtain the simplified result

$\max_{\mu\in E}\left\|u(\mu)-u_{N}^{\Xi_{S}}(\mu)\right\|_{V}\leq
C_{E}\left(\frac{(1+\varepsilon)^{2}\operatorname{Cap}(K)}{\Phi_{K}(\lambda_{N+1})}\right)^{S}\text{.}$
(30)

Thus, the rate of convergence of the approximation error is the same one that
we could expect from polynomial approximants for a function holomorphic over
$\\{\mu\in\mathbb{C},\Phi_{K}(\mu)<\Phi_{K}(\lambda_{N+1})\\}$, having a pole
or singularity at $\lambda_{N+1}$.

####  3.3.1 Error convergence with respect to denominator degree

Just like subsection 3.2.1 describes an extension for variable $N$ of the
convergence theory for approximate poles, here we consider a generalization of
the error convergence theory in the case of diverging denominator degree. The
result we are able to show is of a flavor similar to Theorem 9, and is
summarized in the following.

######  Theorem 11.

Let $(N_{k})_{k\in\mathbb{N}}$, $(S_{k})_{k\in\mathbb{N}}$, and $\Xi^{(k)}$ be
as in Theorem 9. Fix $E\subset\mathbb{C}$ compact, and let Assumption 4 be
satisfied for some fixed $\mu_{0}\in\mathbb{C}$. Then, for all
$\varepsilon>0$, we have

$\lim_{k\to\infty}\operatorname{Cap}\left(\left\\{\mu\in
E,\,\left\|u(\mu)-u_{N_{k}}^{\Xi^{(k)}}(\mu)\right\|_{V}\geq\varepsilon\right\\}\right)=0\text{,}$

with $\operatorname{Cap}$ denoting the logarithmic capacity, see subsection
3.1.

######  Proof.

We rely on some concepts introduced in the proofs of Theorem 9 and Theorem 10,
namely the pole ordering Eq. 23 and the intermediate bound Eq. 29. To start,
we follow the same steps as in the proof of Theorem 9, exploiting Eq. 29 and
Eq. 15 with $\Lambda_{0}=\\{\sigma_{j}\\}_{j=1}^{N_{k}}$. The resulting bound
is

$\left\|u(\mu)-u_{N_{k}}^{\Xi^{(k)}}(\mu)\right\|_{V}\leq\frac{C_{k,\mu_{0},R}\left(\sum_{\lambda\in\Lambda}\left\|v_{\lambda}\right\|_{V}^{2}\right)^{1/2}}{c_{\mu_{0}}\left|Q_{N_{k}}^{\Xi^{(k)}}(\mu)\right|\min_{\Lambda}\left|\,\cdot-\mu\right|}\prod_{i=1}^{N_{k}}\frac{R+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{i}\right|}\text{,}$
(31)

with $R=\max_{E}\left|\,\cdot-\mu_{0}\right|$ and

$C_{k,\mu_{0},R}=\left(\frac{R+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{N_{k}+1}\right|-R_{\mu_{0},K}}\right)^{S_{k}-N_{k}}\left(\frac{2\left|\mu_{0}-\sigma_{N_{k}+1}\right|}{\left|\mu_{0}-\sigma_{N_{k}+1}\right|-R_{\mu_{0},K}}\right)^{N_{k}}\text{.}$

In order to obtain the desired result, we need to manage carefully the two
troublesome $\mu$-dependent terms in the denominator of Eq. 31. To this aim,
let $\Lambda^{(k)}=\smash{\\{\lambda_{j}^{(k)}\\}_{j=1}^{N_{k}}}$ be the set
of roots of $Q_{N_{k}}^{\Xi^{(k)}}$. We partition $\Lambda^{(k)}$ into two
(potentially empty) sets $\Lambda_{in}^{(k)}$ and $\Lambda_{out}^{(k)}$
according to the characterization

$\lambda_{j}^{(k)}\in\Lambda_{in}^{(k)}\quad\text{if and only
if}\quad\left|\lambda_{j}^{(k)}-\mu_{0}\right|\leq 2R\text{.}$

Given $\\#\Lambda_{in}^{(k)}$ and $\overline{N}$ the cardinalities of
$\Lambda_{in}^{(k)}$ and $\Lambda\cap E$ respectively, we define the family of
lemniscates

$\mathcal{E}_{k,\delta}=\left\\{\mu\in\mathbb{C},\,\left(\prod_{\lambda^{\prime}\in\Lambda_{in}^{(k)}}\left|\mu-\lambda^{\prime}\right|\right)\left(\prod_{\lambda\in\Lambda\cap
E}\left|\mu-\lambda\right|\right)\leq\delta^{\\#\Lambda_{in}^{(k)}+\overline{N}}\right\\}\text{,}$
(32)

which depends on the index $k$ and on the (small) value $\delta>0$.

We remark that the exponent of $\delta$ in Eq. 32 is equal to the degree of
the monic polynomial (in $\mu$) whose magnitude is the left hand side of the
inequality. Thus [5, Theorem 6.6.3], the logarithmic capacity of
$\mathcal{E}_{k,\delta}$ is equal to $\delta$. (We are assuming without loss
of generality that $\smash{\\#\Lambda_{in}^{(k)}+\overline{N}}>0$. If this is
not the case, $\mathcal{E}_{k,\delta}$ is empty for $\delta<1$, and the claim
holds quite trivially.)

The main structure of the remainder of the proof is the following:

  1. (i)

we define explicitly a sequence $(\delta_{k})_{k=1,2,\ldots}$ converging to 0;

  2. (ii)

we show that
$\left\|u(\mu)-u_{N_{k}}^{\Xi^{(k)}}(\mu)\right\|_{V}<\varepsilon$ for all
$\mu\in E\setminus\mathcal{E}_{k,\delta_{k}}$, for large $k$.

Then the claim follows, since, by inclusion,

$0\leq\lim_{k\to\infty}\operatorname{Cap}\left(\left\\{\mu\in
E,\,\left\|u(\mu)-u_{N_{k}}^{\Xi^{(k)}}(\mu)\right\|_{V}\geq\varepsilon\right\\}\right)\leq\lim_{k\to\infty}\operatorname{Cap}\left(\mathcal{E}_{k,\delta_{k}}\right)=\lim_{k\to\infty}\delta_{k}=0\text{.}$

However, before we can proceed with either step, it is necessary to introduce
some additional bounds in Eq. 31. First, since $Q_{N_{k}}^{\Xi^{(k)}}$ is
normalized, by Assumption 4 it must hold

$\left|Q_{N_{k}}^{\Xi^{(k)}}(\mu)\right|\geq\frac{1}{C_{\mu_{0}}}\left(\prod_{\lambda^{\prime}\in\Lambda_{in}^{(k)}}\frac{\left|\mu-\lambda^{\prime}\right|}{R_{\mu_{0}}+\left|\mu_{0}-\lambda^{\prime}\right|}\right)\left(\prod_{\lambda^{\prime\prime}\in\Lambda_{out}^{(k)}}\frac{\left|\mu-\lambda^{\prime\prime}\right|}{R_{\mu_{0}}+\left|\mu_{0}-\lambda^{\prime\prime}\right|}\right)\text{.}$

In the first group of factors, $\left|\mu_{0}-\lambda^{\prime}\right|$ is
bounded from above by $2R$. In the second group, since
$\left|\mu-\mu_{0}\right|\leq R$ for all $\mu\in E$ by definition, we have

$\frac{\left|\mu-\lambda^{\prime\prime}\right|}{R_{\mu_{0}}+\left|\mu_{0}-\lambda^{\prime\prime}\right|}\geq\frac{\left|\mu_{0}-\lambda^{\prime\prime}\right|}{R_{\mu_{0}}+\left|\mu_{0}-\lambda^{\prime\prime}\right|}-\frac{\left|\mu-\mu_{0}\right|}{R_{\mu_{0}}+\left|\mu_{0}-\lambda^{\prime\prime}\right|}\geq\frac{R}{R_{\mu_{0}}+2R}\text{.}$

Hence,

$\left|Q_{N_{k}}^{\Xi^{(k)}}(\mu)\right|\geq\frac{1}{C_{\mu_{0}}}\frac{R^{N_{k}-\\#\Lambda_{in}^{(k)}}}{\left(R_{\mu_{0}}+2R\right)^{N_{k}}}\prod_{\lambda^{\prime}\in\Lambda_{in}^{(k)}}\left|\mu-\lambda^{\prime}\right|\text{,}$

which, provided $\delta\leq R$, for all $\mu\in
E\setminus\mathcal{E}_{k,\delta}$, implies

$\displaystyle\left|Q_{N_{k}}^{\Xi^{(k)}}(\mu)\right|\prod_{\lambda\in\Lambda\cap
E}\left|\mu-\lambda\right|\geq$
$\displaystyle\frac{R^{N_{k}-\\#\Lambda_{in}^{(k)}}}{C_{\mu_{0}}\left(R_{\mu_{0}}+2R\right)^{N_{k}}}\left(\prod_{\lambda^{\prime}\in\Lambda_{in}^{(k)}}\left|\mu-\lambda^{\prime}\right|\right)\left(\prod_{\lambda\in\Lambda\cap
E}\left|\mu-\lambda\right|\right)$ $\displaystyle>$
$\displaystyle\frac{R^{N_{k}-\\#\Lambda_{in}^{(k)}}\delta^{\\#\Lambda_{in}^{(k)}+\overline{N}}}{C_{\mu_{0}}\left(R_{\mu_{0}}+2R\right)^{N_{k}}}\geq\frac{\delta^{N_{k}+\overline{N}}}{C_{\mu_{0}}\left(R_{\mu_{0}}+2R\right)^{N_{k}}}\text{.}$

Accordingly, we can simplify Eq. 31 for all $\mu\in
E\setminus\mathcal{E}_{k,\delta}$:

$\left\|u(\mu)-u_{N_{k}}^{\Xi^{(k)}}(\mu)\right\|_{V}<\frac{C_{\mu_{0}}\left(\sum_{\lambda\in\Lambda}\left\|v_{\lambda}\right\|_{V}^{2}\right)^{1/2}\prod_{\lambda\in\Lambda\cap
E}\left|\mu-\lambda\right|}{c_{\mu_{0}}\min_{\Lambda}\left|\,\cdot-\mu\right|}\times\\\
\times
C_{k,\mu_{0},R}\frac{\left(R_{\mu_{0}}+2R\right)^{N_{k}}}{\delta^{N_{k}+\overline{N}}}\prod_{i=1}^{N_{k}}\frac{R+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{i}\right|}\text{.}$

Since $E$ is compact and $\Lambda$ has no finite limit point, the distance
between $E$ and $\Lambda\setminus E$ is strictly positive. Hence, the quantity

$C_{\mu_{0}}\left(\sum_{\lambda\in\Lambda}\left\|v_{\lambda}\right\|_{V}^{2}\right)^{1/2}\prod_{\lambda\in\Lambda\cap
E}\left|\mu-\lambda\right|\,\Big{/}\,\left(c_{\mu_{0}}\min_{\Lambda}\left|\,\cdot-\mu\right|\right)$

is bounded, for all $\mu\in E$, by a constant $C^{\prime}$ depending only on
$E$, $\mu_{0}$, and $\Lambda$.

Now it is trivial to achieve ii by setting

$\delta=\delta_{k}=\bigg{(}\frac{C^{\prime}C_{k,\mu_{0},R}\left(R_{\mu_{0}}+2R\right)^{N_{k}}}{\varepsilon}\prod_{i=1}^{N_{k}}\frac{R+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{i}\right|}\bigg{)}^{1/(N_{k}+\overline{N})}\text{.}$
(33)

(We remark that we still need to guarantee that the working assumption
$\delta\leq R$ is satisfied. However, since we are planning to prove i, such
condition will trivially hold, provided $k$ is large enough.)

At this point, it suffices to apply the same strategy employed in the final
part of the proof of Theorem 9. In particular, we have

$\lim_{k\to\infty}\left(\prod_{i=1}^{N_{k}}\frac{R+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{i}\right|}\right)^{1/(N_{k}+\overline{N})}=\lim_{N\to\infty}\frac{R+R_{\mu_{0},K}}{\left|\mu_{0}-\sigma_{N+1}\right|}=0\text{,}$

while the other factors in Eq. 33 stay bounded as $k$ increases.

Due to Eq. 9, Theorem 11 can be weakened to prove the convergence in
probability of $\smash{u_{N_{k}}^{\Xi^{(k)}}}$ to $u$, when both are
interpreted as functions from the probability space
$(E,\mathcal{B}(E),\ell/\ell(E))$ to the Banach space
$(V,\left\|\cdot\right\|_{V})$.

Still, convergence in capacity is somewhat stronger than that in probability:
for instance, it ensures that the set of parameter values for which the
approximant yields an error above $\varepsilon$ in the limit, namely

$\bigcap_{l\in\mathbb{N}}\bigcup_{k\geq l}\left\\{\mu\in
E\,:\,\left\|u(\mu)-u_{N_{k}}^{\Xi^{(k)}}(\mu)\right\|_{V}\geq\varepsilon\right\\}\text{,}$

cannot include any curve in $\mathbb{C}$.

###  3.4 Generalizations to non-orthogonal residues

In the previous sections, we have been working under the assumption that the
residues $\\{v_{\lambda}\\}_{\lambda\in\Lambda}$ form a $V$-orthogonal set.
However, most of our results for fixed $N$ can be extended to the general
case, in which we impose no condition on the angles between residues. As such,
in this section we remove the assumption of residue orthogonality.

We start by extending Lemma 6.

######  Lemma 12.

Let $Q$, $\Lambda_{0}$, and $\varepsilon$ be the same as in Lemma 6. Then

$j^{\Xi_{S}}(Q)^{2}=\sum_{\lambda,\lambda^{\prime}\in\Lambda}\frac{\langle
v_{\lambda},v_{\lambda^{\prime}}\rangle_{V}}{\omega^{\Xi_{S}}(\lambda)\overline{\omega^{\Xi_{S}}(\lambda^{\prime})}}Q(\lambda)\overline{Q(\lambda^{\prime})}\text{.}$

Moreover,

$j^{\Xi_{S}}(Q_{N}^{\Xi_{S}})\leq\frac{\sum_{\lambda\in\Lambda\setminus\Lambda_{0}}\left\|v_{\lambda}\right\|_{V}}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\omega^{\Lambda_{0}}\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}}\sup_{\Lambda\setminus\Lambda_{0}}\left|\frac{\omega^{\Lambda_{0}}}{\omega^{\Xi_{S}}}\right|\text{,}$

and there exists $S_{\varepsilon,N}$ such that

$j^{\Xi_{S}}(Q_{N}^{\Xi_{S}})\leq(1+\varepsilon)^{S}\frac{\sum_{j>N}\left\|v_{\lambda_{j}}\right\|_{V}}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\prod_{j=1}^{N}\left(\
\cdot-\lambda_{j}\right)\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}}\frac{\prod_{j=1}^{N}\left|\lambda_{N+1}-\lambda_{j}\right|}{\Phi_{K}(\lambda_{N+1})^{S}}\quad\text{for
}S>S_{\varepsilon,N}\text{,}$

provided $\Phi_{K}(\lambda_{N+1})<\Phi_{K}(\lambda_{N+2})$.

######  Proof (sketch).

One can follow most of the proof of Lemma 6. However, Eq. 18 does not hold due
to the lack of orthogonality. Still, we can replace it by

$j^{\Xi_{S}}(Q_{N}^{\Xi_{S}})^{2}\leq
j^{\Xi_{S}}\left(\frac{\omega^{\Lambda_{0}}}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\omega^{\Lambda_{0}}\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}}\right)^{2}\leq\frac{1}{{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\omega^{\Lambda_{0}}\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}^{2}}\left(\sum_{\lambda\in\Lambda\setminus\Lambda_{0}}\left\|v_{\lambda}\right\|_{V}\left|\frac{\omega^{\Lambda_{0}}(\lambda)}{\omega^{\Xi_{S}}(\lambda)}\right|\right)^{2}\text{,}$

which holds by triangular inequality.

We remark the the discrete 2-norm of the residue $V$-norms in Lemma 6 has been
replaced by the stronger 1-norm. In order to derive convergence properties, we
need to introduce the assumption that the $V$-norms of the residues are
summable, so that the inequalities in Lemma 12 are nontrivial.

######  Theorem 13.

Let $N$, $\\{\lambda_{j}^{\Xi_{S}}\\}_{j=1}^{N}$, and $\varepsilon$ be the
same as in Theorem 8. If
$\sum_{\lambda\in\Lambda}\left\|v_{\lambda}\right\|_{V}<\infty$, there exists
$S_{\varepsilon,N}$ such that

$\left|Q_{N}^{\Xi_{S}}\left(\lambda_{j}\right)\right|\leq
C_{j}(1+\varepsilon)^{4S}\left(\frac{\Phi_{K}(\lambda_{j})}{\Phi_{K}(\lambda_{N+1})}\right)^{S}\quad\text{for
}S>S_{\varepsilon,N}$

and

$\min_{i=1,\ldots,N}\left|\lambda_{i}^{\Xi_{S}}-\lambda_{j}\right|\leq
C_{j}(1+\varepsilon)^{4S}\left(\frac{\Phi_{K}(\lambda_{j})}{\Phi_{K}(\lambda_{N+1})}\right)^{S}\quad\text{for
}S>S_{\varepsilon,N}\text{,}$

for $j=1,\ldots,N$, with $C_{j}$ independent of $S$ and $\varepsilon$.

######  Proof (sketch).

The proof is similar to that of Lemma 7 and Theorem 8. However, we can observe
that the exponents in the convergence rates are halved. This is caused by
“interference” between the modes: the difference between the representative
matrices $G^{\Xi_{S}}$ and $G_{N}^{\Xi_{S}}$ is

$G^{\Xi_{S}}-G_{N}^{\Xi_{S}}=\sum_{\begin{subarray}{c}\lambda,\lambda^{\prime}\in\Lambda\\\
(\lambda,\lambda^{\prime})\notin\\{1,\ldots,N\\}^{2}\end{subarray}}\frac{\langle
v_{\lambda},v_{\lambda^{\prime}}\rangle_{V}}{\omega^{\Xi_{S}}(\lambda)\overline{\omega^{\Xi_{S}}(\lambda^{\prime})}}\bm{\omega}_{\lambda^{\prime}}\bm{\omega}_{\lambda}^{H}\text{,}$

which, in particular, contains “off-diagonal” terms, for instance the one
corresponding to $\lambda=1$ and $\lambda^{\prime}=N+1$. When computing an
upperbound for a term of this form, we cannot hope to obtain more than a
_single_ factor of $\Phi_{K}(\lambda_{N+1})^{-S}$, whereas two such factors
can be gained in the orthogonal case. Similarly, the correlation between the
first $N$ resonances causes one less factor of
$\left|\omega^{\Xi_{S}}(\lambda_{j})\right|$ to be obtainable from
$G_{N}^{\Xi_{S}}$.

######  Theorem 14.

Let $N$, $\varepsilon$, and $E$ be as in Theorem 10. If
$\sum_{\lambda\in\Lambda}\left\|v_{\lambda}\right\|_{V}<\infty$, there exists
$S_{\varepsilon,N,E}$ such that

$\left\|u(\mu)-u_{N}^{\Xi_{S}}(\mu)\right\|_{V}\leq
C_{E}(1+\varepsilon)^{4S}\left(\frac{\Phi_{K}(\mu)}{\Phi_{K}(\lambda_{N+1})}\right)^{S}\quad\text{for
all }\mu\in E,\,S>S_{\varepsilon,N,E}\text{,}$

with $C_{E}$ independent of $S$ and $\varepsilon$.

######  Proof (sketch).

Since Lemma 6 does not hold, it is impossible to apply the middle part of the
proof of Theorem 10. Instead, we can split the sum in Eq. 28 into two parts:
for the first $N$ poles $\lambda_{1},\ldots,\lambda_{N}$, Theorem 13 can be
applied to bound $|Q_{N}^{\Xi_{S}}/\omega^{\Xi_{S}}|$ from above; for the
remaining resonances, it suffices to use the argument from the last portion of
Lemma 6, with some modifications to account for the presence of
$Q_{N}^{\Xi_{S}}$.

Unfortunately, convergence results for variable $N$ do not appear to be
possible in the non-orthogonal case, without some strong hypotheses on the
“amount of correlation” between resonances: for instance, we believe that one
could work under the assumption that the worst-case angle between residues

$\alpha(\lambda)=\inf_{\lambda^{\prime}\in\Lambda\setminus\\{\lambda\\}}\frac{\langle
v_{\lambda},v_{\lambda^{\prime}}\rangle_{V}}{\left\|v_{\lambda}\right\|_{V}\left\|v_{\lambda^{\prime}}\right\|_{V}}$

is bounded away from 0 for all poles $\lambda$ in the region of interest.

##  4 A posteriori error indicators

All convergence results presented above are _a priori_ estimates. In
particular, they contain quantities which are often not available, namely the
number and locations of the poles inside (and of some of the ones outside) the
parameter domain $K$. Moreover, the results in Theorems 8 and 11 hold only
(with some extensions) if the solution map is of the form Eq. 8, whereas
minimal rational interpolants can, in principle, be applied to more general
parametric problems and sample points sets.

A typical MOR problem can be cast in the following form: let a certain
parameter set $K$ and a parametric problem as in Eq. 6 be given, along with
some fixed tolerance $\varepsilon>0$. The task is to obtain an approximate
solution map $\widetilde{u}$ such that

$\max_{\mu\in
K}\left\|u(\mu)-\widetilde{u}(\mu)\right\|_{V}\leq\varepsilon\quad\text{or}\quad\max_{\mu\in
K}\frac{\left\|u(\mu)-\widetilde{u}(\mu)\right\|_{V}}{\left\|u(\mu)\right\|_{V}}\leq\varepsilon\text{.}$

In our framework, the presence of singularities in the true and, possibly, the
approximate solution map within $K$ may make the conditions above meaningless.
In such cases, it may be preferable [17] to consider a residual-based accuracy
requirement

$\max_{\mu\in
K}\left\|\mathcal{F}_{\mu}\left(\widetilde{u}(\mu)\right)\right\|_{W}\leq\varepsilon\quad\text{or}\quad\max_{\mu\in
K}\frac{\left\|\mathcal{F}_{\mu}\left(\widetilde{u}(\mu)\right)\right\|_{W}}{\left\|\mathcal{F}_{\mu}\left(0\right)\right\|_{W}}\leq\varepsilon\text{,}$
(34)

with the $\left\|\cdot\right\|_{W}$-norm replacing $\left\|\cdot\right\|_{V}$
to account for regularity differences between residual and solution.

Now, let us take as approximate map $\widetilde{u}$ the minimal rational
interpolant $u_{N}^{Z_{S}}$ with a certain denominator degree $N$ and samples
at the $S$ points $Z_{S}$, which need not be Fejér points for $K$. We face the
task of understanding whether the required tolerance $\varepsilon$ is
achieved, in the sense of Eq. 34.

In particular, let us consider the problem of computing _a posteriori_ the
norm of the residual
$\left\|\mathcal{F}_{\mu}\left(\widetilde{u}(\mu)\right)\right\|_{W}$ at some
given point $\mu\in K$. In order to proceed, we will assume that the operator
$\mathcal{F}_{\mu}$ is linear and has a separable form in $\mu$, i.e. that
there exist two families of complex-valued functions
$\\{\theta_{i}^{F}\\}_{i=1}^{n_{F}}$ and $\\{\theta_{i}^{f}\\}_{i=1}^{n_{f}}$,
a family of operators $\\{F_{i}\\}_{i=1}^{n_{F}}$ over (suitable subsets of)
$V$, and a family $\\{f_{i}\\}_{i=1}^{n_{f}}$ of elements of $V$, such that

$\mathcal{F}_{\mu}(v)=F(\mu)v-f(\mu)=\sum_{i=1}^{n_{F}}\theta_{i}^{F}(\mu)F_{i}v-\sum_{i=1}^{n_{f}}\theta_{i}^{f}(\mu)f_{i}\text{.}$
(35)

This is the ideal situation to employ RB methods [16, 26, 35, 36], and many
strategies have been devised to approximate general parametric problems into
the form Eq. 35, e.g. the (D)EIM [15, 21] and hyper-reduction techniques [14].

Now, $F(\mu)$ can be applied to both sides of Eq. 27 to obtain

$Q_{N}^{Z_{S}}(\mu)\mathcal{F}_{\mu}\left(u_{N}^{Z_{S}}(\mu)\right)=\sum_{\mu^{\prime}\in
Z_{S}}\frac{Q_{N}^{Z_{S}}(\mu^{\prime})\,\omega^{Z_{S}}(\mu)}{(\mu-\mu^{\prime})\left.\omega^{Z_{S}}\right.^{\prime}(\mu^{\prime})}\mathcal{F}_{\mu}\left(u(\mu^{\prime})\right)\\\
=\sum_{i=1}^{n_{F}}F_{i}\Delta^{Z_{S}}\left(\theta_{i}^{F}I^{Z_{S}}\left(uQ_{N}^{Z_{S}}\right)\right)(\mu)-\sum_{i=1}^{n_{f}}\Delta^{Z_{S}}\left(\theta_{i}^{f}Q_{N}^{Z_{S}}\right)(\mu)\,f_{i}\text{.}$
(36)

with $\Delta^{Z_{S}}(\phi)=\phi-I^{Z_{S}}(\phi)$ the interpolation error for a
function $\phi$.

This, after dividing by $\smash{Q_{N}^{Z_{S}}}(\mu)$, provides an affine
decomposition of the residual with $n_{F}+n_{f}$ terms. Hence, its norm
(squared) admits an affine decomposition with $(n_{F}+n_{f})^{2}$ terms. In
particular, assuming an offline-online framework [26, 35, 36], the separable
terms of the residual can be precomputed offline, i.e. together with the $S$
expensive evaluations of the target solution map $u$, without increasing the
complexity of the overall MOR procedure.

Of course, the discussion above can intrinsically be applied only in an
intrusive fashion, since we require access to the operators defining the
parametric problem. In general, obtaining an _a posteriori_ non-intrusive
error indicator can be quite tricky. Here we propose a simplified approach,
which is rigorous only for very specific problem structures, namely

$\mathcal{F}_{\mu}(u)=F(\mu)u-f=(F_{0}+\mu F_{1})u-f\text{.}$ (37)

For such problems, thanks to Eq. 17 and Eq. 36, it must hold

$\displaystyle
Q_{N}^{Z_{S}}(\mu)\mathcal{F}_{\mu}\left(u_{N}^{Z_{S}}(\mu)\right)=$
$\displaystyle\omega^{Z_{S}}(\mu)F_{1}\sum_{\mu^{\prime}\in
Z_{S}}\frac{Q_{N}^{Z_{S}}(\mu^{\prime})u(\mu^{\prime})}{\left.\omega^{Z_{S}}\right.^{\prime}(\mu^{\prime})}$
$\displaystyle=$
$\displaystyle\frac{\omega^{Z_{S}}(\mu)}{(S-1)!}F_{1}\frac{\operatorname{d}^{S-1}}{\operatorname{d}{\mu}^{S-1}}I^{Z_{S}}(uQ_{N}^{Z_{S}})\text{.}$

Now it suffices to take the $\left\|\cdot\right\|_{W}$-norm and exploit Eq. 4
to obtain

$\displaystyle\left\|\mathcal{F}_{\mu}\left(u_{N}^{Z_{S}}(\mu)\right)\right\|_{W}=$
$\displaystyle\frac{1}{(S-1)!}\left\|F_{1}\frac{\operatorname{d}^{S-1}}{\operatorname{d}{\mu}^{S-1}}I^{Z_{S}}(uQ_{N}^{Z_{S}})\right\|_{W}\left|\frac{\omega^{Z_{S}}(\mu)}{Q_{N}^{Z_{S}}(\mu)}\right|$
(38) $\displaystyle\leq$ $\displaystyle\left(\sup_{v\in
V\setminus\\{0\\}}\frac{\left\|F_{1}v\right\|_{W}}{\left\|v\right\|_{V}}\right)j^{Z_{S}}\left(Q_{N}^{Z_{S}}\right)\left|\frac{\omega^{Z_{S}}(\mu)}{Q_{N}^{Z_{S}}(\mu)}\right|\text{.}$
(39)

One obvious limitation of Eq. 38 and Eq. 39 is the need to evaluate or bound
expensive norms (either the $W$ one or the $V\to W$ operator one) involving
$F_{1}$. We remark that this drawback can be interpreted as an issue in
identifying the exact scaling of the error estimator, a common problem even
for stable parametric problems, where the inf-sup (or coercivity) constant
must be estimated to link residual and error [28, 35, 36]. If one can overcome
such limitation, see subsection 6.2, Eq. 38 and Eq. 39 are certainly quite
appealing from a computational standpoint, since their evaluation consists
only of few scalar computations.

If Eq. 37 does not hold, we can replace $F(\mu)$ by a suitable linear
approximation (e.g. its first-order Taylor series around some parameter) in a
spirit similar to that of the Empirical Interpolation Method. Then Eq. 38 and
Eq. 39 can be applied heuristically, with an accuracy which may depend quite
sharply on the smoothness of $\mathcal{F}_{\mu}$, in particular on its second-
order variations over the parameter domain.

As a supplement to our discussion, we would like to note that _a posteriori_
error indicators are often used in RB procedures not only to determine whether
the required accuracy has been achieved, but also to drive the selection of
the sample points. For instance, in the weak-greedy RB approach [35, 36], one
keeps adding a new snapshot at the parameter value which maximizes the _a
posteriori_ error indicator of choice, thus updating the surrogate model (and
the error indicator as well), until convergence. We envision that a similar
adaptive procedure could be devised also for minimal rational interpolation,
with new parameter values being added to the sample set according to the same
logic.

##  5 Extensions to multiple parameters

Throughout our discussion, we have only dealt with the approximation of
univariate meromorphic solution maps. However, it is often of interest to
approximate the solution map of multi-parameter problems: for example, in the
analysis of parametric dynamical systems, one may need to approximate the
solution of

$(A(\mathbf{p})-\mu
E(\mathbf{p}))u(\mu,\mathbf{p})=B(\mathbf{p}),\quad\text{with
}\mathbf{p}\in\mathbb{C}^{d-1},$ (40)

where $\mathbf{p}$ is a collection of physical or geometric parameters. For
such problems, it is often quite complicated to prove properties of $u$,
jointly in $\mu$ and $\mathbf{p}$; however, the analysis is much simpler if
the extra parameters $\mathbf{p}$ are kept fixed (“marginalized out”) and one
looks at the univariate solution map $u_{\mathbf{p}}(\mu):\mu\mapsto
u(\mu,\mathbf{p})$.

This has lead to the development of _parametric MOR_ (pMOR) [6, 19, 33, 40,
41], whose core idea is the following: take a set of parameter values
$\\{\mathbf{p}_{1},\ldots,\mathbf{p}_{T}\\}\subset\mathbb{C}^{d-1}$ and, for
each $j=1,\ldots,T$, build a surrogate $\widetilde{u}_{\mathbf{p}_{j}}$ of
$u_{\mathbf{p}_{j}}$ (this involves computing the snapshots
$u(\mu_{i},\mathbf{p}_{j})$ for $i=1,\ldots,S_{j}$, and applying a single-
parameter MOR strategy, with respect to $\mu$ only); then, the joint reduced
model $\widetilde{u}$ of $u$ is constructed as a suitable combination of
$\widetilde{u}_{\mathbf{p}_{1}},\ldots,\widetilde{u}_{\mathbf{p}_{T}}$.

In a pMOR framework, one could apply minimal rational approximation to obtain
the surrogates in frequency only
$\\{\widetilde{u}_{\mathbf{p}_{j}}\\}_{j=1}^{T}$. In particular, the results
presented in the previous sections apply: for instance, the residual estimator
in section 4 is still valid, allowing for a greedy approach with respect to
$\mu$. This line of approximation for multi-parameter problems is currently
under investigation.

There are several reasons to prefer pMOR over global MOR (in $d$ parameters)
for problems like Eq. 40. Many of them are related to the curse of
dimensionality, which forces to take a large number of snapshots, and to build
a surrogate model of considerable size, undermining the cost-effectiveness of
MOR. We refer to the general pMOR literature cited above for a broader
discussion.

Here, we wish to remark that these disadvantages affect also the natural
extension of minimal rational interpolation to $d$ parameters, which we could
define as follows: assume that the sets of polynomials
$\mathbb{P}_{N}\left(\mathbb{C}^{d};\mathbb{C}\right)$ and
$\mathbb{P}_{N}^{\star}\left(\mathbb{C}^{d};\mathbb{C}\right)$, as well as the
interpolation operator $I^{\Xi_{S}}$, are suitably defined (for instance,
$\mathbb{P}_{N}\left(\mathbb{C}^{d};\mathbb{C}\right)$ could denote the space
of polynomials of _total_ degree at most $N$); the optimization problem which
yields the surrogate denominator $Q_{N}^{\Xi_{S}}$ is now: find
$Q_{N}^{\Xi_{S}}\in\mathbb{P}_{N}^{\star}\left(\mathbb{C}^{d};\mathbb{C}\right)$
which minimizes over
$\mathbb{P}_{N}^{\star}\left(\mathbb{C}^{d};\mathbb{C}\right)$ the convex
functional

$Q\;\mapsto\;j^{\Xi_{S}}(Q)^{2}=\stackrel{{\scriptstyle[}}{{\alpha}}_{1}+\ldots+\alpha_{d}=M]{}{\sum_{\alpha_{1}=0,1\ldots}\cdots\sum_{\alpha_{d}=0,1\ldots}}\left\|\frac{1}{\alpha_{1}!\cdots\alpha_{d}!}\frac{\operatorname{d}^{\alpha_{1}}}{\operatorname{d}{p}_{1}^{\alpha_{1}}}\cdots\frac{\operatorname{d}^{\alpha_{d-1}}}{\operatorname{d}{p}_{d-1}^{\alpha_{d-1}}}\frac{\operatorname{d}^{\alpha_{d}}}{\operatorname{d}{\mu}^{\alpha_{d}}}I^{\Xi_{S}}(uQ)\right\|_{V}^{2}\text{,}$

where $M$ is the maximal total degree of the interpolant $I^{\Xi_{S}}(uQ)$, so
that, in general, $S\sim M^{d}$. In practice, an algorithm for minimizing this
functional can be found as a generalization of the strategy proposed in
section 2; while we omit the details here, we remark that the Gramian matrix
$\Psi$ representing $j^{\Xi_{S}}$, see Eq. 5 for the univariate case, becomes
of size approximately $S^{2-1/d}\times N^{d}$.

In addition to the other limitations of multivariate minimal rational
interpolation, we also observe that it is impossible to extend directly our
single-parameter theory to the multivariate case, since many of the core ideas
(e.g. the barycentric expansion Eq. 17) do not generalize nicely to more than
one dimension. Overall, also considering the numerical issues intrinsic to
high-dimensional interpolation, we believe minimal rational approximants to be
an appropriate MOR technique only in low dimension: extensions to many
parameters should rely on pMOR or similar approaches.

##  6 Numerical examples

In this section, we perform some numerical experiments to show the usefulness
of minimal rational interpolation, and to verify the theoretical convergence
rates and error indicator. The first example satisfies the assumption under
which our results have been derived; in particular, the residues are
orthogonal. Instead, the solution map in the second example can be proven
meromorphic, but not necessarily of the form Eq. 8. For the sake of
reproducibility, the code used to run the simulations is available at [34].

###  6.1 Normal eigenvalue problem

Let $n=100$ and take $A\in\mathbb{C}^{n\times n}$ a normal matrix whose
eigenvalues are randomly generated according to a uniform distribution over
$\\{x+\iota y,(x,y)\in[-5,5]^{2}\\}$, and whose eigenvectors are prescribed by
orthonormalizing a matrix with random (complex) standard gaussian entries. In
particular, we order the eigenvalues $\\{\lambda_{j}\\}_{j=1}^{n}$ according
to their distance from 0.

Our task is to estimate the eigenvalues of $A$ which lie within the unit disk
$\smash{K=\overline{\mathcal{B}(0,1)}}$. To this aim, we fix
$\mathbf{v}\in\mathbb{C}^{n}$ a random gaussian vector, and we approximate by
minimal rational interpolants the solution map of the parametric problem

$\text{for }\mu\in K,\;\text{find
}\mathbf{u}(\mu)\in\mathbb{C}^{n}\;\text{s.t.}\;\left(A-\mu
I\right)\mathbf{u}(\mu)=\mathbf{v}\text{,}$ (41)

with $I\in\mathbb{C}^{n\times n}$ the identity matrix. We endow
$\mathbb{C}^{n}$ with the Euclidean scalar product. As sample points
$\Xi_{S}$, we take the $S$-th roots of unity, and we choose to employ the
polynomial norm
${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{0,N}$,
see Eq. 14. This norm is induced by monomials, which are $L^{2}(\partial
K)$-orthogonal. Moreover, the action of the interpolation operator
$\smash{I^{\Xi_{S}}}$ can be evaluated quite efficiently as a Fast Fourier
Transform.

Let us assume that our interest lies in approximating just the eigenvalues of
$A$ inside $K$. As their exact number is unknown, we can choose the
denominator degree $N$ according to different strategies:

  1. •

we fix $N$ (in our case, we set $N=10$) and increase progressively the number
of samples $S$; to verify if the starting guess on $N$ was large enough, we
can simply check if at least one of the approximate poles has converged (for
large $S$) to some point outside $K$, see Theorem 8.

  2. •

we can increase $N$ along with the number of samples $S$, for instance by
choosing $N=S-1$ for all $S$; in this case we are sure to approximate all
poles, see Theorem 9.

We compare visually the two approaches in Figure 2, where we show the error
norm obtained by approximants of type $[20/10]$ and $[20/20]$, both based on
$S=21$ samples at roots of unity. We can see that the error increases much
quicker near the boundary of the sampled square in the case $N=10$. This is to
be expected, as the region of convergence increases together with $N$, see
Theorem 10. A (possibly) more interesting observation is that the error
appears to be globally smaller when we choose $N=20$. Overall, we can conclude
that the choice $N=S-1$ leads to a better approximant.

Figure 1: Natural logarithm of the relative error norm
$\|\mathbf{u}-\mathbf{u}_{N}^{\Xi_{S}}\|_{2}/\left\|\mathbf{u}\right\|_{2}$
achieved by minimal rational interpolants based on samples at
$\Xi_{S}=\\{e^{2j\pi\iota/21}\\}_{j=1}^{21}$ (black crosses) with $N=10$
(left) and $N=20$ (right). The exact and approximated poles are represented in
blue and white respectively.

Figure 2: Pole approximation error for the poles inside $K$. On the left
$N=10$ is kept fixed, whereas $N=S-1$ on the right. In black the theoretical
rate (22). The dashed lines show the results obtained by projection on
subspaces of dimension $N$ computed through POD.

Now it just remains to check how well each approach manages to capture the
location of the eigenvalues of $A$. To this aim, we increase $S$ from $11$ to
$30$, and compute the pole approximation error, defined in the left hand side
of Eq. 20, for the 5 poles within $K$. The results are shown in Figure 2. For
fixed $N=10$, the rate predicted by Eq. 22 can be observed (Green’s potential
for a centered disk coincides with the complex magnitude outside the disk).
Remarkably, a similar, if not slightly better, rate of convergence is obtained
for increasing $N=S-1$, despite the absence of theoretical results in this
regard.

In Figure 2 we also show the pole approximation error obtained through a RB-
like projection technique: given the same sample set $\Xi_{S}$, we perform a
Proper Orthogonal Decomposition (POD) [35] of the samples, identifying the $N$
dominant “sample directions”; then we project orthogonally the original
problem on the subspace spanned by these $N$ directions, and we use the
solution of the resulting $N\times N$ eigenvalue problem to approximate the
original one. We can observe that the projection technique performs much worse
than minimal rational interpolants for fixed $N=10$, whereas for $N=S-1$ the
two methods achieve very similar results.

Figure 3: Logarithm of the relative error norm
$\|\mathbf{u}-\mathbf{u}_{N}^{\Xi_{S}}\|_{2}/\left\|\mathbf{u}\right\|_{2}$
achieved by minimal rational interpolants based on samples (white crosses) at
the origin (left) and at quasi-random points (right) with $N=20$. The exact
and approximated poles are represented in blue and white respectively.

Figure 4: Pole approximation error for the poles inside $K$. The samples are
all at the origin on the left, whereas on the right the samples are at quasi-
random points. In both cases we keep $N=S-1$.

Lastly, we wish to test two other approximation strategies: a fast LS Padé
approximant [10] centered at the origin, and a minimal rational interpolant
based on quasi-random (Halton) points in the unit disk. We show the error
norms and the approximated poles in Figure 4 for $N=S-1=20$. Qualitatively, it
appears that, for both choices, the accuracy in the approximation of the
resonances is similar to that achieved with samples at the roots of unity.
With respect to the $[20/20]$ minimal rational interpolant in Figure 2, the
errors for the two new approximants appear smaller within $K$, but larger
outside $K$; this is to be expected, since the new sample points are closer to
the origin.

More quantitatively, we show in Figure 4 the convergence of the pole
approximation error with respect to $S$ (here as well, we keep $N=S-1$). For
fast LS Padé approximants, when $S$ is small, the error converges quite
quickly (the error is uniformly smaller than the one obtained with samples at
the roots of unity), but then it stagnates, and even increases starting from
$S=20$; this is a symptom of the instability which affects the construction of
this centered approximant, namely the computation of high-order derivatives of
a meromorphic function.

The results for quasi-random sample points appear relatively similar, with the
error flattening out starting from $S=20$. In this case, the non-monotone
behavior of the error for large $S$ is due to the interaction between two
conflicting effects: the addition of sample points close to the resonances
(which would, on its own, improve the accuracy of the approximation) and the
ill-conditioning of the interpolation problem (for $S=30$, the Lebesgue
constant is $\sim 10^{7}$, as opposed to $\sim 3$ for the $30$-th roots of
unity).

###  6.2 Time-harmonic vibrations of a tuning fork

Let $\Omega\subset\mathbb{R}^{3}$ be the region occupied by a tuning fork.
Assume that the tuning fork is clamped on a portion $\Gamma_{D}$ of the
handle. Let $\smash{\widetilde{\mathbf{g}}}(t,\mathbf{x})=e^{\iota
2\pi\overline{\nu}t}\mathbf{g}(\mathbf{x})$ be a pressure pulse (sinusoidal in
time and gaussian in space) impinging on a portion $\Gamma_{\mathbf{g}}$ of
the exterior of the tuning fork, which we assume to be sound-hard. See Figure
5 for a graphical representation of the setup.

\begin{overpic}[height=51.21504pt,trim=35.56593pt 42.67912pt 25.6073pt
55.48277pt,clip]{RMOR-figDomB} \put(2.0,5.0){$\Gamma_{D}$}
\put(32.0,13.5){$\Gamma_{free}$} \put(80.0,10.5){$\Gamma_{\mathbf{g}}$}
\end{overpic}\begin{overpic}[height=51.21504pt,trim=35.56593pt 42.67912pt
25.6073pt 49.79231pt,clip]{RMOR-figDom} \put(100.5,-0.75){\scriptsize$0$\hskip
24.18501pt[m]} \put(100.5,5.125){\scriptsize$9.54\cdot 10^{-8}$}
\put(100.5,11.75){\scriptsize$1.91\cdot 10^{-7}$}
\put(100.5,17.875){\scriptsize$2.86\cdot 10^{-7}$}
\put(100.5,24.25){\scriptsize$3.82\cdot 10^{-7}$} \end{overpic}

Figure 5: Reference domain (left) and domain warped by the real part of the
displacement, magnified by a factor of $3\cdot 10^{3}$ (right). The
displacement is the solution of (42) for frequency $\nu=10$ kHz and damping
coefficient $\eta=0$ Hz. The total number of degrees of freedom of the
discretized problem is 36330.

We are interested in computing the time-harmonic deformation $\mathbf{u}(\nu)$
with frequency $\nu$ which the tuning fork undergoes. Assuming the tuning fork
to have density $\rho$ and linear stress constitutive relation
$\sigma(\cdot)$, the displacement $\mathbf{u}(\nu)$ can be found by solving a
(damped) time-harmonic linear elasticity problem

$\begin{cases}-\nabla\cdot\sigma(\mathbf{u}(\nu))-2\pi\nu(2\pi\nu-\eta\iota)\rho\mathbf{u}(\nu)=\bm{0}\quad&\text{in
}\Omega\\\ \mathbf{u}(\nu)=\bm{0}\quad&\text{on }\Gamma_{D}\\\
\sigma(\mathbf{u}(\nu))\mathbf{n}=\mathbf{g}\quad&\text{on
}\Gamma_{\mathbf{g}}\\\ \sigma(\mathbf{u}(\nu))\mathbf{n}=\bm{0}\quad&\text{on
}\Gamma_{free}=\partial\Omega\setminus\Gamma_{D}\setminus\Gamma_{\mathbf{g}}\end{cases}\text{.}$
(42)

As $V$-inner product, we choose the elastic energy inner product:

$\langle\mathbf{u},\mathbf{v}\rangle_{V}=\int_{\Omega}\sigma(\mathbf{u}):\varepsilon(\mathbf{v})\textrm{d}\mathbf{x}=\int_{\Omega}\sigma(\mathbf{u}):\frac{\nabla\mathbf{v}+\nabla\mathbf{v}^{\top}}{2}\textrm{d}\mathbf{x}.$

In particular, we choose to approximate problem Eq. 42 using $\mathbb{P}_{1}$
FE on a sufficiently fine tetrahedral discretization of $\Omega$. The FE
discretization of Eq. 42 defines a parametric problem of the form Eq. 6, whose
solution map $\mathbf{u}$ can be proven meromorphic using compactness
arguments [11, 31]. Thus, we decide to approximate it for $\nu\in K=[10,40]$
kHz with minimal rational interpolants relying on samples at the Chebyshev
nodes of $K$. The chosen polynomial norm
${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}$
is the one induced by the family of Chebyshev polynomials over $K$, see
subsection 3.1.

Figure 6: Convergence of the approximation error in the energy norm with
respect to the number of samples at $\nu=13$ kHz (left) and $\nu=31$ kHz
(right).

Figure 7: Energy norm of the solution map $\mathbf{u}$ of (42) in blue, for
$\eta\in\\{0,500\\}$ Hz. In green and red the norm of the $[20/20]$ and
$[40/40]$ minimal rational interpolant of $\mathbf{u}$.

Figure 8: Relative $L^{2}(\Omega)^{3}$-norm of the residual for problem (42)
and $\widetilde{\mathbf{u}}=\mathbf{u}_{20}^{\Xi_{21}}$ (blue) and
$\widetilde{\mathbf{u}}=\mathbf{u}_{40}^{\Xi_{41}}$ (black). In green and red
the heuristic approximations (38) of the residual. The dashed grey line is the
residual of the RB approximation of $\mathbf{u}$ with samples at $\Xi_{21}$.

First, we investigate the behavior of the relative approximation error
$\frac{\left\|\widetilde{\mathbf{u}}(\nu)-\mathbf{u}(\nu)\right\|_{V}}{\left\|\mathbf{u}(\nu)\right\|_{V}}$
at $\nu=13$ kHz and $\nu=31$ kHz, for $\eta=0$ Hz. We apply different
approximation strategies: for $S=1,\ldots,51$, we consider $[S-1/S-1]$ fast LS
Padé approximants [10] with $S$ samples at the center of $K$, $[S-1/N]$
minimal rational approximants for fixed $N=10,20,30$, and diagonal $[S-1/S-1]$
minimal rational approximants.

The results are shown in Figure 8. As expected from Theorem 10, minimal
rational approximants of type $[S-1/10]$ do not converge ($K$ contains 16
resonances), whereas the ones of type $[S-1/20]$ and $[S-1/30]$ do; we can
observe that the diagonal approximant $[S-1/S-1]$ yields the best results
overall.

The error achieved by fast LS Padé approximants converges extremely quickly at
$\nu=31$ kHz for small $S$, but no convergence can be observed at $\nu=13$
kHz; this is reasonable, since this type of approximant, relying on
information at a single point, can deliver an extremely accurate surrogate
close to the center, but the quality of the approximation degrades as we move
further away. Also, we observe that the convergence plots for fast LS Padé
approximants stop decreasing already at $S=11$, and terminate abruptly at
$S=35$: for $S\geq 36$, the construction of the surrogate model fails, due to
conditioning issues in the evaluation of the derivatives of the solution map,
and in the application of the (Taylor) interpolation operator [10, 12].

The norm of the solution map is shown in Figure 8 for $\eta\in\\{0,500\\}$ Hz,
along with the approximated norm profiles, obtained by minimal rational
interpolants of type $[20/20]$ and $[40/40]$. If $S=21$, for both values of
$\eta$, we can observe that the approximation of
$\left\|\mathbf{u}\right\|_{V}$ is quite accurate except for the low-frequency
region, where the approximant for $\eta=0$ Hz actually appears quite unstable.
This behavior can be improved by increasing the number of samples to $S=41$.

In a (more) realistic application, in order to determine whether there is
actual need for more than 21 samples, one could rely on residual-based
estimators such as Eq. 34. To this aim, we can apply Eq. 36 to obtain an
expression for the residual in the $L^{2}(\Omega)^{3}$-norm. The evaluation of
Eq. 36 on a fine grid is shown in Figure 8.

In the same Figure, we also show the values obtained through the simplified
residual estimator Eq. 39. In particular, in order to evaluate the estimator,
we need to compute or approximate the operator norm
$\left\|F_{1}\right\|_{\mathcal{L}(V;L^{2}(\Omega)^{3})}$, with $F_{1}$ the
derivative of the operator in Eq. 42 at some point $\mu^{\prime}$, see section
4. We manage to avoid this, thus keeping the estimator as non-intrusive as
possible, by the following argument: if we assume that Eq. 38 holds
approximately (the problem is quadratic and not linear in the parameter), we
can identify a constant $C$, independent of $\mu$, for which, employing the
notation of section 4,

$\left\|\mathcal{F}_{\mu}\left(u_{N}^{Z_{S}}(\mu)\right)\right\|_{W}\approx
C\left|\frac{\omega^{Z_{S}}(\mu)}{Q_{N}^{Z_{S}}(\mu)}\right|\text{.}$

In practice, the value of $C$ can be easily found by evaluating (offline and
at a relatively low computational cost) the exact value of the residual at
some point $\mu^{\prime}\in K\setminus Z_{S}$, so that the value of $C$ can be
set to

$C=\left\|\mathcal{F}_{\mu}\left(u_{N}^{Z_{S}}(\mu^{\prime})\right)\right\|_{W}\left|\frac{Q_{N}^{Z_{S}}(\mu^{\prime})}{\omega^{Z_{S}}(\mu^{\prime})}\right|\text{.}$

In our example we choose $\mu^{\prime}=25.5$ kHz, and we can observe the
corresponding residual estimator to be extremely accurate. Indeed, the values
obtained with the two estimators never differ by more than 1%, and the
corresponding curves cannot be distinguished in the plot.

Finally, in Figure 8 we also perform a comparison between minimal rational
interpolants and RBs, by showing the norm of the residual obtained with the RB
approximant of Eq. 42 computed from samples at $\Xi_{21}$. We can observe that
the two residuals are quite similar in terms of both order of magnitude and
behavior with respect to $\nu$. In particular, simplifying grossly, minimal
rational interpolants seem to have a very slight edge on RBs at higher
frequencies, while the opposite appears to be true at lower frequencies.

##  7 Conclusions

In this paper we have presented minimal rational interpolants as the natural
extension of fast LS-Padé approximants [10], and we have shown that most of
the theoretical properties of this latter technique generalize nicely to the
former. The new method can be considered superior to the original one for two
main reasons:

  1. •

the stability is greatly improved, since computing derivatives of the solution
map accumulates numerical errors; the Arnoldi-like reorthogonalization
technique described in [10] can somehow limit but not totally eliminate this
issue;

  2. •

the distributed sampling helps to achieve a globally (e.g. using the
$L^{\infty}(K)$-metric) small approximation error, which is the target in most
applications; of course, an interpolatory approach cannot be as accurate close
to the samples as a Taylor-type approximation, assuming the number of samples
to be the same.

Our numerical example has shown minimal rational interpolants to preform well
in the MOR of a normal eigenproblem and of a time-harmonic problem. In
particular, we have shown the effectiveness of a heuristic and cheap to
compute _a posteriori_ residual estimator. Accordingly, we deem of interest to
investigate a possible greedy-type algorithm, in the same spirit as the one
for the RB method, but with minimal rational interpolants replacing Galerkin
projection. In particular, it is still unclear how well the approximation
properties of the method are preserved if the samples are not computed at
optimal (Fejér) points.

### Acknowledgments

The author gratefully thanks Fabio Nobile for the many fruitful discussions
and for his helpful advice at several stages of this work.

##  Appendix A Polynomial norm bounds

Let $\Gamma\subset\mathbb{C}$ be a Jordan curve, and consider
$\\{\varphi_{j}\\}_{j=0,1,\ldots}$ a hierarchical (degree$(\varphi_{j})=j$ for
all $j$) family of polynomials, orthonormal over $L^{2}(\Gamma)$:

$\frac{\int_{\Gamma}\varphi_{j}\overline{\varphi_{i}}\left|\text{d}l\right|}{\int_{\Gamma}\left|\text{d}l\right|}=\delta_{ij}\quad\forall(i,j)\in\\{0,1,\ldots\\}^{2}\text{.}$
(43)

For instance, if $\Gamma=\partial\mathcal{B}(\mu_{0},1)$, we can choose
$\varphi_{j}(\mu)=\mu^{j}$ for all $j$.

This, for all $N\in\mathbb{N}$, induces a norm over
$\mathbb{P}_{N}\left(\mathbb{C}\right)$ in the following sense:

${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|Q\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}=\sqrt{\sum_{j=0}^{N}\left|q_{j}\right|^{2}},\quad\text{where}\quad
Q=\sum_{j=0}^{N}q_{j}\varphi_{j}.$ (44)

The following holds.

######  Lemma 15.

The norm
${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|\cdot\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}$
in Eq. 44 satisfies Assumption 4 for any $\mu_{0}$ chosen within the interior
of $\Gamma$.

######  Proof.

For a given $N\in\mathbb{N}$, let $0\leq D\leq N$,
$\\{z_{j}\\}_{j=1}^{D}\in\mathbb{C}^{D}$, and $Q=\prod_{j=1}^{D}\left(\,\cdot-
z_{j}\right)$. Due to orthonormality, it holds

${\left|\kern-1.07639pt\left|\kern-1.07639pt\left|Q\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}^{2}\int_{\Gamma}\left|\text{d}l\right|=\int_{\Gamma}\left|Q\right|^{2}\left|\text{d}l\right|=\int_{\Gamma}\prod_{j=1}^{D}\left|\,\cdot-
z_{j}\right|^{2}\left|\text{d}l\right|\text{.}$ (45)

On one side, let $R_{\mu_{0}}=\max_{\Gamma}\left|\,\cdot-\mu_{0}\right|$. We
can use the triangular inequality to obtain

$\displaystyle{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|Q\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}^{2}\int_{\Gamma}\left|\text{d}l\right|\leq$
$\displaystyle\int_{\Gamma}\prod_{j=1}^{D}\left(\left|\,\cdot-\mu_{0}\right|+\left|\mu_{0}-z_{j}\right|\right)^{2}\left|\text{d}l\right|$
$\displaystyle\leq$
$\displaystyle\prod_{j=1}^{D}\left(R_{\mu_{0}}+\left|\mu_{0}-z_{j}\right|\right)^{2}\int_{\Gamma}\left|\text{d}l\right|\text{.}$
(46)

On the other hand, let
$c_{\mu_{0}}=\left(2\pi\min_{\Gamma}\left|\,\cdot-\mu_{0}\right|/\int_{\Gamma}\left|\textup{d}l\right|\right)^{2}$.
We can exploit the Cauchy-Schwarz inequality to show that

$\left(\int_{\Gamma}\prod_{j=1}^{D}\left|\,\cdot-
z_{j}\right|^{2}\left|\text{d}l\right|\right)\left(\int_{\Gamma}\frac{\left|\text{d}l\right|}{\left|\,\cdot-\mu_{0}\right|^{2}}\right)\geq\left|\int_{\Gamma}\frac{\prod_{j=1}^{D}\left(\,\cdot-
z_{j}\right)}{\,\cdot-\mu_{0}}\text{d}l\right|^{2}\text{,}$

so that, due to the Cauchy integral formula,

$\displaystyle{\left|\kern-1.07639pt\left|\kern-1.07639pt\left|Q\right|\kern-1.07639pt\right|\kern-1.07639pt\right|}_{N}^{2}\int_{\Gamma}\left|\text{d}l\right|\geq$
$\displaystyle\left|\int_{\Gamma}\frac{\prod_{j=1}^{D}\left(\,\cdot-
z_{j}\right)}{\,\cdot-\mu_{0}}\text{d}l\right|^{2}\left(\int_{\Gamma}\frac{\left|\text{d}l\right|}{\left|\,\cdot-\mu_{0}\right|^{2}}\right)^{-1}$
$\displaystyle\geq$ $\displaystyle
4\pi^{2}\prod_{j=1}^{D}\left|\mu_{0}-z_{j}\right|^{2}\left(\min_{\Gamma}\left|\,\cdot-\mu_{0}\right|\right)^{2}\left(\int_{\Gamma}\left|\text{d}l\right|\right)^{-1}$
$\displaystyle\geq$ $\displaystyle
c_{\mu_{0}}\prod_{j=1}^{D}\left|\mu_{0}-z_{j}\right|^{2}\int_{\Gamma}\left|\text{d}l\right|\text{.}$
(47)

Dividing Eq. 46 and Eq. 47 by $\int_{\Gamma}\left|\text{d}l\right|$ yields the
proof.

## References

  * [1] A.C. Antoulas.  Approximation of Large-Scale Dynamical Systems: An Overview.  IFAC Proceedings Volumes, 37(11):19–28, 2017.
  * [2] A.C. Antoulas, C.A. Beattie, and S. Gugercin.  Interpolatory Model Reduction of Large-Scale Dynamical Systems, pages 3–58.  Springer US, Boston, MA, 2010.
  * [3] J.M. Ash, A. Berele, and S. Catoiu.  Plausible and genuine extensions of L’Hospital’s rule.  Mathematics Magazine, 85(1):52–60, 2012.
  * [4] Z. Bai, J. Demmel, J. Dongarra, A. Ruhe, and H. van der Vorst.  Templates for the Solution of Algebraic Eigenvalue Problems.  Society for Industrial and Applied Mathematics, 2000.
  * [5] G. A. Baker and P. R. Graves-Morris.  Padé approximants.  Cambridge University Press, 1996.
  * [6] U. Baur, C. Beattie, P. Benner, and S. Gugercin.  Interpolatory Projection Methods for Parameterized Model Reduction.  SIAM Journal on Scientific Computing, 33(5):2489–2518, 2011.
  * [7] C.A. Beattie and S. Gugercin.  Model Reduction by Rational Interpolation.  Comput. Sci. Engrg, 15:297–334, 2014.
  * [8] J.-P. Berrut and L.N. Trefethen.  Barycentric Lagrange Interpolation.  SIAM Review, 46(3):501–517, 2004.
  * [9] F. Bonizzoni, F. Nobile, and I. Perugia.  Convergence analysis of Padé approximations for Helmholtz frequency response problems.  ESAIM: M2AN., 2017.  DOI: 10.1051/m2an/2017050.
  * [10] F. Bonizzoni, F. Nobile, I. Perugia, and D. Pradovera.  Fast Least-Squares Padé approximation of self-adjoint problems with meromorphic structure.  ArXiv e-prints, 2018.  DOI: arXiv/1808.03112.
  * [11] F. Bonizzoni, F. Nobile, I. Perugia, and D. Pradovera.  Least-Squares Padé approximation of parametric and stochastic Helmholtz maps.  ArXiv e-prints, 2018.  DOI: arXiv/1805.05031.
  * [12] F. Bonizzoni and D. Pradovera.  Distributed sampling for rational approximation of the acoustic scattering of an airfoil.  PAMM, 19(1), 2019.
  * [13] L.P. Bos and N. Levenberg.  On the calculation of approximate fekete points: The univariate case.  Electronic Transactions on Numerical Analysis, 30:377–397, 2008\.
  * [14] K. Carlberg, C. Farhat, J. Cortial, and D. Amsallem.  The GNAT method for nonlinear model reduction: Effective implementation and application to computational fluid dynamics and turbulent flows.  Journal of Computational Physics, 242:623–647, 2013.
  * [15] S. Chaturantabut and D.C. Sorensen.  Nonlinear Model Reduction via Discrete Empirical Interpolation.  SIAM Journal on Scientific Computing, 32(5):2737–2764, 2010.
  * [16] Y. Chen, J.S. Hesthaven, Y. Maday, and J. Rodríguez.  Certified reduced basis methods and output bounds for the harmonic Maxwell’s equations.  SIAM J. Sci. Comput., 32(2):970–996, 2010.
  * [17] M. Chevreuil and A. Nouy.  Model order reduction based on proper generalized decomposition for the propagation of uncertainties in structural dynamics.  International Journal for Numerical Methods in Engineering, 89(2):241–268, 2012.
  * [18] Z. Drmač, S. Gugercin, and C. Beattie.  Quadrature-Based Vector Fitting for Discretized $\mathcal{H}_{2}$ Approximation.  SIAM Journal on Scientific Computing, 37(2):625–652, 2015.
  * [19] Francesco Ferranti, Luc Knockaert, and Tom Dhaene.  Passivity-preserving parametric macromodeling by means of scaled and shifted state-space systems.  IEEE Transactions on Microwave Theory and Techniques, 59(10 PART 1):2394–2403, 2011.
  * [20] Roland W. Freund.  Model reduction methods based on Krylov subspaces.  Acta Numerica, 12(2003):267–319, 2003.
  * [21] M.A. Grepl, Y. Maday, N.C. Nguyen, and A.T. Patera.  Efficient reduced-basis treatment of nonaffine and nonlinear partial differential equations.  ESAIM: Mathematical Modelling and Numerical Analysis, 41(3):575–605, 2007.
  * [22] E. Grimme.  Krylov projection methods for model reduction.  Electrical engineering, University of Illinois at Urbana-Champaign, 1997\.
  * [23] S. Grivet-Talocia and B. Gustavsen.  Passive Macromodeling: Theory and Applications.  Wiley series in microwave and optical engineering. John Wiley & Sons, Inc, Hoboken New Jersey, 2015.
  * [24] P. Guillaume.  Nonlinear eigenproblems.  SIAM Journal on Matrix Analysis and Applications, 20(3):575–595, 1999.
  * [25] B. Gustavsen and A. Semlyen.  Rational approximation of frequency domain responses by vector fitting.  IEEE Transactions on Power Delivery, 14(3):1052–1061, jul 1999\.
  * [26] J.S. Hesthaven, G. Rozza, and B. Stamm.  Certified reduced basis methods for parametrized partial differential equations.  SpringerBriefs in Mathematics. Springer Cham; BCAM Basque Center for Applied Mathematics, Bilbao, 2016.
  * [27] E Hille.  Analytic function theory.  Chelsea Publishing Company, New York, N.Y, second ed. edition, 1965.
  * [28] D.B.P. Huynh, G. Rozza, S. Sen, and A.T. Patera.  A successive constraint linear optimization method for lower bounds of parametric coercivity and inf-sup stability constants.  Comptes Rendus Mathematique, 345(8):473–478, oct 2007.
  * [29] A.C. Ionita and A.C. Antoulas.  Data-Driven Parametrized Model Reduction in the Loewner Framework.  SIAM Journal on Scientific Computing, 36(3):A984–A1007, 2014.
  * [30] C. Lei, Y. Wang, Q. Chen, and N. Wong.  On Vector Fitting Methods in Signal / Power Integrity Applications.  II, 2010.
  * [31] M. Lenoir, M. Vullierme-Ledard, and C. Hazard.  Variational Formulations for the Determination of Resonant States in Scattering Problems.  SIAM Journal on Mathematical Analysis, 23(3):579–608, 1992.
  * [32] P. Monk.  Finite Element Methods for Maxwell’s Equations.  Clarendon, Oxford, 2003.
  * [33] H. Panzer, J. Mohring, R. Eid, and B. Lohmann.  Parametric model order reduction by matrix interpolation.  At-Automatisierungstechnik, 58(8):475–484, 2010.
  * [34] D. Pradovera.  Interpolatory rational model order reduction of parametric problems lacking uniform inf-sup stability - numerical tests.  Feb 2020.
  * [35] A. Quarteroni, A. Manzoni, and F. Negri.  Reduced Basis Methods for Partial Differential Equations: An Introduction.  UNITEXT. Springer International Publishing, 2015.
  * [36] G. Rozza, D. B.P. Huynh, and A. T. Patera.  Reduced basis approximation and a posteriori error estimation for affinely parametrized elliptic coercive partial differential equations: Application to transport and continuum mechanics.  Archives of Computational Methods in Engineering, 15(3):229–275, 2008.
  * [37] S. Steinberg.  Meromorphic families of compact operators.  Archive for Rational Mechanics and Analysis, 31:372–379, dec 1968\.
  * [38] J L Walsh.  Interpolation and Approximation by Rational Functions in the Complex Domain.  American Mathematical Society, 1969.
  * [39] D.S. Weile, E. Michielssen, E. Grimme, and K. Gallivan.  A method for generating rational interpolant reduced order models of two-parameter linear systems.  Applied Mathematics Letters, 12(5):93–102, 1999.
  * [40] Yi Qing Xiao, Stefano Grivet-Talocia, Paolo Manfredi, and Roni Khazaka.  A Novel Framework for Parametric Loewner Matrix Interpolation.  IEEE Transactions on Components, Packaging and Manufacturing Technology, 9(12):1–1, 2019.
  * [41] Y. Yue, L. Feng, and P. Benner.  Reduced-order modelling of parametric systems via interpolation of heterogeneous surrogates.  Advanced Modeling and Simulation in Engineering Sciences, 0123456789, 2019.