Introduction {#Sec1}
============

Many complete search methods for problems in global optimization and constraint satisfaction hinge on the ability to compute enclosures for the range of nonconvex functions as well as sets defined by multiple equalities and/or inequalities \[[@CR44]\]. A variety of parameterizations can be used to describe or approximate compact sets, including convex sets such as intervals, ellipsoids and polytopes as well as nonconvex sets for instance in the form of the image of a multivariate polynomial. Chachuat et al. \[[@CR13]\] have argued that the foregoing parameterizations and many others can be cast as affine set-parameterizations with either convex or nonconvex basis sets. In particular, a number of arithmetics have been developed over the years that enable propagation of certain affine-set parameterizations through the atom operations of factorable functions in order to determine rigorous enclosures of their image sets. This includes interval arithmetic \[[@CR40], [@CR41]\], ellipsoidal calculus \[[@CR28], [@CR61]\], polyhedral relaxations \[[@CR54], [@CR55]\], and Taylor model arithmetic \[[@CR11], [@CR36]\], to name but a few. Other established techniques to relax sets or programs defined by nonconvex factorable functions involve constructing a pair of convex/concave bounds for each function, such as McCormick relaxations \[[@CR38], [@CR39]\] and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\alpha $$\end{document}$BB/$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\gamma $$\end{document}$BB relaxations \[[@CR1], [@CR3], [@CR23]\]. In the case of semialgebraic problems, hierarchies of linear matrix inequality (LMI) relaxations can also be constructed based on a semidefinite programming formulation of the sum of squares decomposition for multivariate polynomials \[[@CR29], [@CR46]\].

The focus in this paper is on *polynomial models* of multivariable functions, namely enclosures in the form of a polynomial approximant of given order and an interval remainder bound on the approximation error. Taylor models (or Taylor forms) are a special case of polynomial models, whereby the polynomial approximation matches the multivariate Taylor expansion of the function at a given point in the variable domain \[[@CR43]\]. The idea of such interval polynomial enclosures dates back to Moore \[[@CR40]\], and methods for computing the remainder bounds concurrently with the polynomial were later developed in the early 1980s \[[@CR17], [@CR51]\] and popularized from the mid 1990s by Berz et al. \[[@CR7], [@CR35], [@CR36]\]. Similar in essence to interval arithmetic, a well-developed Taylor model arithmetic is available \[[@CR8], [@CR34]\], encompassing rules for binary sums, binary products and outer-compositions with a library of univariate functions such as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\exp (\cdot )$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\log (\cdot )$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\sqrt{\cdot }$$\end{document}$. For sufficiently smooth functions, the diameter of the remainder interval constructed with (*q*th-order) Taylor model arithmetic can be proved to be a high-order power ($\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$) of the diameter of the variables domain \[[@CR11]\]. This gives Taylor models a clear advantage over traditional interval extensions or centered forms for sufficiently narrow domains, but conversely it may result in a large overestimation or may even be poorer than naive interval evaluation over wider domains. Nevertheless, this approach has proved successful in computing tight enclosures for the solutions of differential equations and implicit algebraic equations \[[@CR24], [@CR33], [@CR42], [@CR49], [@CR50], [@CR53], [@CR60]\], and it has enabled complete search for a range of global optimization or constraint satisfaction problems that could not be tackled using interval techniques alone (see, e.g., \[[@CR4], [@CR9], [@CR31], [@CR32], [@CR47], [@CR52]\]). Such higher-order inclusion techniques are indeed appealing in complete search applications based on branching or subdivision, where they can mitigate the clustering effect \[[@CR15], [@CR62]\].

A rather natural idea involves replacing the Taylor series approximation with a better one, ideally the minimax polynomial approximation, or near-minimax approximations such as Chebyshev series. The development of affine arithmetics \[[@CR19]\] yielding first-order approximations can be seen as a precursor of this approach. Regarding univariate functions, improved variants based on the Chebyshev polynomials, Bernstein polynomials and others were widely used during the early 1980s under the name of ultra-arithmetic \[[@CR26]\]. A comparison of univariate Taylor forms and Chebyshev forms in \[[@CR27]\] showed that expansions in Chebyshev series may be orders of magnitude more accurate than expansions in Taylor series. The use of Chebyshev series expansion as an arithmetic is also the philosophy behind the package Chebfun by Trefethen et al. \[[@CR6], [@CR57]\], which has recently been extended to functions in two variables as well \[[@CR56]\]. It should be noted that Chebfun relies on computations with functions to 15-digit accuracy, thus removing the need for propagating a remainder term, but the results are not validated per se. In contrast, Brisebarre and Joldeş \[[@CR12], [@CR25]\] developed an arithmetic for univariate functions based on Chebyshev interpolation polynomials and a remainder term. More recently, Dzetkulič \[[@CR16]\] extended the previous work to enable enclosure of multivariate functions using truncated Chebyshev series and an interval remainder term, namely Chebyshev models. This latter work makes use of Chebyshev models as a means of propagating bounds for differential equations, but does not study the convergence properties of this arithmetic.

Building upon these existing contributions, this article presents a Chebyshev model arithmetic in order to compute enclosures for an inclusive class of factorable functions (Sect. [3](#Sec3){ref-type="sec"}). Besides formalizing the rules of addition, multiplication and composition as well as the polynomial range bounders, we describe an improvement of the remainder term in composition operations with certain univariate functions. We also present a convergence analysis of Chebyshev models, which investigates the propagation of local convergence properties through the addition, multiplication and composition operations (Sect. [4](#Sec11){ref-type="sec"}), and we discuss the global convergence of this arithmetic as the polynomial expansion order increases. Finally, we present numerical case studies illustrating the convergence properties of Chebyshev models and comparing their tightness and computational performance with Taylor models (Sect. [5](#Sec15){ref-type="sec"}). The underlying implementation of Chebyshev model arithmetic is part of the library MC++, which is made freely available at: <http://omega-icl.bitbucket.org/mcpp/>.

Notation and background {#Sec2}
=======================

The set of *n*-dimensional interval vectors is denoted by $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathbb {IR}^{n}$$\end{document}$. The midpoint and radius of an interval vector $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$X:=[x^{\mathrm{L}},x^{\mathrm{U}}]\in \mathbb {IR}^n$$\end{document}$ are defined as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$${{\mathrm{mid\,}}}X := \frac{1}{2} (x^{\mathrm {L}}+ x^{\mathrm {U}})\in \mathbb {R}^n$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$${{\mathrm{rad\,}}}X := \frac{1}{2} (x^{\mathrm {U}}- x^{\mathrm {L}})\in \mathbb {R}^n$$\end{document}$, respectively. The image of the set $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$X\in \mathbb {IR}^n$$\end{document}$ under a continuous function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f:\mathbb {R}^n\rightarrow \mathbb {R}$$\end{document}$ is given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \bar{f}(X) := \{f(x)|x\in X\}\in \mathbb {IR}. \end{aligned}$$\end{document}$$The total variation of a (possibly discontinuous) function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$g:\mathbb {R}\rightarrow \mathbb {R}$$\end{document}$ on the interval $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Y:=[y^{\mathrm {L}},y^{\mathrm {U}}]\in \mathbb {IR}$$\end{document}$ is the quantity$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left. V[g]_Y := \sup \left\{ \sum _{i=1}^{N} \left| g(\xi _{i})-g(\xi _{i-1}) \right| \right| N>0,\quad y^{\mathrm {L}}=\xi _0<\xi _1<\cdots <\xi _N=y^{\mathrm {U}}\right\} , \end{aligned}$$\end{document}$$where $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\xi _0,\ldots ,\xi _N$$\end{document}$ is any finite partition of *Y*, and the supremum norm of *g* on *Y* is given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \Vert g\Vert _Y := \sup \left\{ |g(\xi )||\xi \in Y \right\} \end{aligned}$$\end{document}$$Moreover, by a slight abuse of the notation, we use *V*\[*g*\] and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\Vert g\Vert $$\end{document}$ when $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Y=[-1,1]$$\end{document}$.

Chebyshev polynomials of the first kind, denoted by $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$T_k$$\end{document}$, are defined by the three-term recurrence relation$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} T_0(x) := 1, \quad T_1(x) := x, \quad T_{k+1}(x) := 2xT_k(x)-T_{k-1}(x). \end{aligned}$$\end{document}$$An explicit expression of these polynomials on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$ is$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} T_k (x) = \cos (k \arccos (x)). \end{aligned}$$\end{document}$$The Chebyshev polynomials are extremal polynomials with respect to the property that, for each $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$k>0$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$T_k(x)$$\end{document}$ is the polynomial with the largest possible leading coefficient subject to the condition that its range is in $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$ for $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x\in [-1,1]$$\end{document}$.

Any polynomial function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}^q$$\end{document}$ of degree *q* can be written in terms of the Chebyshev polynomials as \[[@CR37]\]$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}^q(x) = \sum _{k=0}^{q} a_k T_k(x). \end{aligned}$$\end{document}$$An efficient way of evaluating such a linear combination of Chebyshev polynomials is using the Clenshaw algorithm \[[@CR14]\],$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}^q(x) = a_0 + xb_1(x) - b_2(x), \end{aligned}$$\end{document}$$where the $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$b_k$$\end{document}$'s are obtained using the (reverse) recurrence formula$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} b_k(x) = a_k + 2xb_{k+1}(x) - b_{k+2}(x)\,\quad \text {with}\quad b_{q+1}(x)=b_{q+2}(x)=0. \end{aligned}$$\end{document}$$Moreover, the product between two Chebyshev polynomials expands to$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} T_j (x)T_k (x)&= \frac{1}{2}\left( T_{j+k}(x)+T_{|j-k|}(x)\right) . \end{aligned}$$\end{document}$$The Chebyshev polynomials form a sequence of orthogonal polynomials with respect to the weight $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\frac{1}{\sqrt{1-x^2}}$$\end{document}$ on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$,$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \int _{-1}^1 \frac{T_j(x)T_k(x)}{\sqrt{1-x^2}}dx = \left\{ \begin{array}{ll} \pi &{}\quad \text {if}\,j=k=0\\ \pi /2 &{}\quad \text {if }\, j=k>0\\ 0 &{}\quad \text {otherwise.}\end{array}\right. \end{aligned}$$\end{document}$$Any Lipschitz continuous function *f* on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$ has a unique representation as an absolutely and uniformly convergent series \[[@CR58]\],where the notation $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\sum ^{\prime }$$\end{document}$ indicates that the first term is halved. In general, the Chebyshev coefficients $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a_k$$\end{document}$ can be approximated using a numerical quadrature, possibly in combination with the fast Fourier cosine transform to improve the computation speed \[[@CR21]\].

For practical purposes, a *q*th-order polynomial approximant $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_f^q$$\end{document}$ for *f* can be obtained from the partial sums of ([3](#Equ3){ref-type=""}) asThe smoother the function *f* on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$, the faster the convergence of its approximants $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_f^q$$\end{document}$ as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\rightarrow \infty $$\end{document}$. The following bounds are available with regards to the accuracy of the partial sum $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_f^q$$\end{document}$, depending on the regularity of *f*:If *f* and its derivatives through $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f^{(s-1)}$$\end{document}$ for a given $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$s\ge 0$$\end{document}$ are absolutely continuous on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$ and if $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f^{(s)}$$\end{document}$ has finite total variation, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$V[f^{(s)}]<\infty $$\end{document}$, then \[[@CR58]\] $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall q\ge s+1, \quad \left\| f-\mathcal {P}_f^q\right\| \le \frac{2}{\pi }\frac{V[f^{(s)}]}{s(q-s)^s}. \end{aligned}$$\end{document}$$If *f* and its derivatives through $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f^{(q)}$$\end{document}$ are continuously differentiable on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$, then \[[@CR18]\] $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left\| f-\mathcal {P}_f^q\right\| \le \frac{\Vert f^{(q+1)}\Vert }{2^{q}(q+1)!}. \end{aligned}$$\end{document}$$An alternative approach to obtaining a *q*th-order polynomial approximant is by constructing the polynomial interpolant $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\widehat{\mathcal {P}}_f^q$$\end{document}$ at Chebyshev points,The polynomial $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\widehat{\mathcal {P}}_f^q$$\end{document}$ yields a near-minimax approximation, and counterparts to the error estimates in ([5](#Equ5){ref-type=""}) and ([6](#Equ6){ref-type=""}) are, respectively, \[[@CR18], [@CR58]\]$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall q\ge s+1,\quad \left\| f-\widehat{\mathcal {P}}_f^q\right\|&\le \frac{4}{\pi }\frac{V[f^{(s)}]}{s(q-s)^s},\end{aligned}$$\end{document}$$ $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left\| f-\widehat{\mathcal {P}}_f^q\right\|&\le \frac{\Vert f^{(q+1)}\Vert }{2^{q}(q+1)!}. \end{aligned}$$\end{document}$$While the truncated Chebyshev series expansion $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_f^q$$\end{document}$ is generally a better approximation than the Chebyshev interpolation polynomial $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\widehat{\mathcal {P}}_f^q$$\end{document}$, the latter comes with the inherent advantage that the coefficients $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\widehat{a}_k$$\end{document}$ in ([7](#Equ7){ref-type=""}) can be obtained from a few function evaluations only.

The (over-)approximation of multivariate functions using Chebyshev polynomials is at the heart of the present paper. Like in the univariate case, any *q*th-order polynomial function in *n* variables can be written in terms of the Chebyshev polynomials as$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}^q(x) = \sum _{|\kappa |\le q} a_\kappa T_\kappa (x), \end{aligned}$$\end{document}$$where $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa \in \mathbb {N}^n$$\end{document}$ is a multi-index; the order of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa $$\end{document}$ is given by $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|\kappa | := \sum \nolimits _{i=1}^n \kappa _i$$\end{document}$; and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$T_\kappa (x)$$\end{document}$ is a shorthand notation for the expression $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\prod \nolimits _{i=1}^n T_{\kappa _i}(x_i)$$\end{document}$ for any $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x\in \mathbb {R}^n$$\end{document}$.

Chebyshev model arithmetic {#Sec3}
==========================

Similar to Taylor models \[[@CR11], [@CR43]\], Chebyshev models are estimators made by the sum of a multivariate polynomial in Chebyshev form and a remainder interval providing a bound on the approximation error. A formal definition follows:

Definition 1 {#FPar1}
------------

Let the function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f:Z\rightarrow \mathbb {R}$$\end{document}$ be defined on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Z\subseteq \mathbb {R}^n$$\end{document}$. For every interval vector $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Y\subset Z$$\end{document}$, let the multivariate polynomial $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f,Y}^{q}:[-1,1]^n\rightarrow \mathbb {R}$$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f,Y}^{q}(\xi ) := \sum \nolimits _{|\kappa |\le q} a_\kappa T_\kappa (\xi )$$\end{document}$ and the symmetric interval $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {R}_{f,Y}^{q}\in \mathbb {IR}$$\end{document}$ be such that$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \xi \in [-1,1]^n, \quad f({{\mathrm{mid\,}}}Y +\xi \circ {{\mathrm{rad\,}}}Y) - \mathcal {P}_{f,Y}^{q}({\xi }) \in \mathcal {R}_{f,Y}^{q}, \end{aligned}$$\end{document}$$where $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\circ $$\end{document}$ denotes the entrywise (or Hadamard) product. We call the pair $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})$$\end{document}$ a *qth-order Chebyshev model* of *f* on *Y*, and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})_{Y\subset Z}$$\end{document}$ a *scheme of qth-order Chebyshev models* for *f* in *Z*.

Notice that, in the sense of Definition [1](#FPar1){ref-type="sec"}, the polynomial part $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f,Y}^{q}$$\end{document}$ of a Chebyshev model may be different from the actual truncated Chebyshev expansion of *f*. Notice also the need for scaling the variables in $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$ for consistency with the definition and properties of Chebyshev polynomials of the first kind. The definition of schemes of Chebyshev models is in analogy with the scheme of estimators introduced by Bompadre et al. \[[@CR10], [@CR11]\] to support convergence analysis.

Similar in essence to interval (or Taylor model) arithmetic, our focus in the following subsections is on the development of an arithmetic that enables the construction of Chebyshev models for factorable functions, namely functions which can be represented by a finite number of binary sums, binary products and outer compositions with a univariate function. As well as the specification of these three "elementary" operations, developing such an arithmetic also calls for an initialization procedure for variables and constants, and the ability to bound the range of Chebyshev models.

Initialization {#Sec4}
--------------

Given a set $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Z\subseteq \mathbb {R}^n$$\end{document}$, each variable $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$z_i\in Z_i$$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$i\in \{1,\ldots ,n\}$$\end{document}$ can be represented by a scheme of *q*th-order Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{z_i,Y}^{q},\mathcal {R}_{z_i,Y}^{q})_{Y\subset Z}$$\end{document}$ with$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \xi \in [-1,1]^n,\quad \mathcal {P}_{z_i,Y}^{q}(\xi )&:= \left\{ \begin{array}{ll} {{\mathrm{mid\,}}}Y_i \, T_0(\xi _i) &{}\quad \text {if}\,q=0\\ {{\mathrm{mid\,}}}Y_i \, T_0(\xi _i) + {{\mathrm{rad\,}}}Y_i \, T_1(\xi _i)&{}\quad \text {if}\,q>0 \end{array}\right. \end{aligned}$$\end{document}$$ $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {R}_{z_i,Y}^{q}&:= \left\{ \begin{array}{ll} Y_i-{{\mathrm{mid\,}}}Y_i &{}\quad \text {if}\,q=0\\ {[}0,0] &{}\quad \text {if}\,q>0. \end{array}\right. \end{aligned}$$\end{document}$$Likewise, any real constant $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r\in \mathbb {R}$$\end{document}$ is trivially represented by the scheme of Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(r \, T_0,[0,0])_{Y\subset Z}$$\end{document}$.

Binary sum operations {#Sec5}
---------------------

Consider two functions $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1,f_2:Z\subseteq \mathbb {R}^n\rightarrow \mathbb {R}$$\end{document}$, and suppose that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1,Y}^{q},\mathcal {R}_{f_1,Y}^{q})_{Y\subset Z}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_2,Y}^{q},\mathcal {R}_{f_2,Y}^{q})_{Y\subset Z}$$\end{document}$ are corresponding schemes of Chebyshev models. Then, a scheme of Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1\pm f_2,Y}^{q},\mathcal {R}_{f_1\pm f_2,Y}^{q})_{Y\subset Z}$$\end{document}$ for $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1\pm f_2$$\end{document}$ in *Z* is simply obtained by addition/subtraction of the polynomial parts and of the remainder parts,$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned}&\displaystyle \forall \xi \in [-1,1]^n,\quad \mathcal {P}_{f_1\pm f_2,Y}^{q}(\xi ) := \mathcal {P}_{f_1,Y}^{q}(\xi ) \pm \mathcal {P}_{f_2,Y}^{q}(\xi ) \end{aligned}$$\end{document}$$ $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned}&\displaystyle \mathcal {R}_{f_1\pm f_2,Y}^{q} := \mathcal {R}_{f_1,Y}^{q} \pm \mathcal {R}_{f_2,Y}^{q}. \end{aligned}$$\end{document}$$

Binary product operations {#Sec6}
-------------------------

Consider again two functions $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1,f_2:Z\subseteq \mathbb {R}^n\rightarrow \mathbb {R}$$\end{document}$, and suppose that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1,Y}^{q},\mathcal {R}_{f_1,Y}^{q})_{Y\subset Z}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_2,Y}^{q},\mathcal {R}_{f_2,Y}^{q})_{Y\subset Z}$$\end{document}$ are corresponding schemes of Chebyshev models, with$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}_{f_1,Y}^{q}(\xi ) := \sum _{|\kappa |\le q} a^q_{Y,\kappa } T_\kappa (\xi ) \quad \text {and} \quad \mathcal {P}_{f_2,Y}^{q}(\xi ) := \sum _{|\kappa |\le q} b^q_{Y,\kappa } T_\kappa (\xi ). \end{aligned}$$\end{document}$$In order to arrive at a valid scheme of *q*th-order Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1 f_2,Y}^{q},\mathcal {R}_{f_1 f_2,Y}^{q})_{Y\subset Z}$$\end{document}$ for the product function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1 f_2$$\end{document}$ in *Z*, we start by expressing the product between the respective polynomial parts,$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}_{f_1,Y}^{q}(\xi ) \mathcal {P}_{f_2,Y}^{q}(\xi ) = \sum _{|\lambda |\le q} \sum _{|\mu |\le q} a^q_{Y,\lambda } b^q_{Y,\mu }\, T_{\lambda }(\xi )\, T_{\mu }(\xi ). \end{aligned}$$\end{document}$$A complication with polynomial multiplication in Chebyshev basis is that, according to the property ([2](#Equ2){ref-type=""}), each product term $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$T_{\lambda }(\xi ) T_{\mu }(\xi )$$\end{document}$ generates $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$2^{N(\lambda ,\mu )}$$\end{document}$ terms with$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} N(\lambda ,\mu ) := \mathrm{card}\{ i\in \{1,\ldots ,n\} \mid \lambda _i\mu _i>0 \} \le n. \end{aligned}$$\end{document}$$See also Sect. [5](#Sec15){ref-type="sec"} for further discussions and related implementation considerations.

Introducing the map $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathbb {P}_q:\mathbb {N}^n\rightarrow \mathbb {N}^n\times \mathbb {N}^n$$\end{document}$ given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathbb {P}_q(\kappa ) := \left\{ (\lambda ,\mu )\in \mathbb {N}^{2n} \Big | \begin{array}{l} \forall i\in \{1,\ldots ,n\},\quad \lambda _i+\mu _i=\kappa _i\ \vee \ |\lambda _i-\mu _i|=\kappa _i\\ |\lambda |\le q,\quad |\mu |\le q \end{array}\right\} , \end{aligned}$$\end{document}$$it follows from ([2](#Equ2){ref-type=""}) that$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}_{f_1,Y}^{q}(\xi ) \mathcal {P}_{f_2,Y}^{q}(\xi ) = \sum _{|\kappa |\le 2q} c^q_{Y,\kappa } T_{\kappa }(\xi ) \quad \text {with}\quad c^q_{Y,\kappa } := \sum _{(\lambda ,\mu )\in \mathbb {P}_q(\kappa )} \frac{a^q_{Y,\lambda } b^q_{Y,\mu }}{2^{N(\lambda ,\mu )}}. \end{aligned}$$\end{document}$$Observe that the resulting product polynomial is of order 2*q* and that the index sets $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathbb {P}_q(\kappa )$$\end{document}$ are non-empty for $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|\kappa |\le 2q$$\end{document}$. A simple way of defining the polynomial part $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f_1 f_2,Y}^{q}$$\end{document}$ in the resulting Chebyshev model is therefore by retaining all the terms of order no larger than *q* as$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \xi \in [-1,1]^n,\quad \mathcal {P}_{f_1 f_2,Y}^{q}(\xi ) := \sum _{|\kappa |\le q} c^q_{Y,\kappa } T_{\kappa }(\xi ). \end{aligned}$$\end{document}$$Then, the corresponding remainder term $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {R}_{f_1 f_2,Y}^{q}$$\end{document}$ can be obtained by assembling four terms, namely the multiplication of the two remainders, the multiplication of each remainder with a bound on the other polynomial range, and a bound on all the terms of order larger than *q* from the product of the two polynomials,$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {R}_{f_1 f_2,Y}^{q}&:= \mathcal {R}_{f_1,Y}^{q} \mathcal {R}_{f_2,Y}^{q} + \mathcal {R}_{f_1,Y}^{q} \left[ \mathcal {P}_{f_2,Y}^{q}\right] + \left[ \mathcal {P}_{f_1,Y}^{q}\right] \mathcal {R}_{f_2,Y}^{q} \nonumber \\&\qquad + \sum _{q<|\kappa |\le 2q} \left| c^q_{Y,\kappa }\right| [-1,1]. \end{aligned}$$\end{document}$$The latter hinges on the ability to compute bounds on the polynomial parts of the Chebyshev model operands, denoted by $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[\mathcal {P}_{f_1,Y}^{q}]\supseteq \{\mathcal {P}_{f_1,Y}^{q}(\xi ) \mid \xi \in [-1,1]^n\}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[\mathcal {P}_{f_2,Y}^{q}]\supseteq \{\mathcal {P}_{f_2,Y}^{q}(\xi ) \mid \xi \in [-1,1]^n\}$$\end{document}$ here. Such range bounders for Chebyshev models are discussed later on in Sect. [3.5](#Sec8){ref-type="sec"}.

An important aspect about polynomial multiplication in Chebyshev basis is that the computed *q*th-order polynomial $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f_1 f_2,Y}^{q}$$\end{document}$ in ([17](#Equ17){ref-type=""}) will not correspond to the *q*th-order Chebyshev expansion of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1f_2$$\end{document}$ in general, even when $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f_1,Y}^{q}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f_2,Y}^{q}$$\end{document}$ are themselves *q*th-order Chebyshev expansion of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_2$$\end{document}$, respectively. The reason behind this discrepancy is that the Chebyshev coefficients in the *q*th-order expansion of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1f_2$$\end{document}$ depend on terms of order greater than *q* in the Chebyshev expansions of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_2$$\end{document}$, which are missing here because of the truncation. Dzetkulič \[[@CR16]\] argue that the resulting Chebyshev models nevertheless provide a tighter approximation than by multiplying $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f_1,Y}^{q}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f_2,Y}^{q}$$\end{document}$ in monomial basis.

Univariate outer-composition operations {#Sec7}
---------------------------------------

Consider the functions $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f:Z\subseteq \mathbb {R}^n\rightarrow \mathbb {R}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F:X\subseteq \mathbb {R}\rightarrow \mathbb {R}$$\end{document}$. Suppose that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})_{Y\subset Z}$$\end{document}$ is a scheme of *q*th-order Chebyshev models for *f*, such that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$B^q_{Y}:=[\mathcal {P}_{f,Y}^{q}]+\mathcal {R}_{f,Y}^{q} \subseteq X$$\end{document}$ for all $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Y\subseteq Z$$\end{document}$. Let the scaling function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\varPhi ^q_Y:Y\rightarrow [-1,1]$$\end{document}$ be given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \varPhi ^q_Y(y) := \frac{f(y)-{{\mathrm{mid\,}}}B^q_{Y}}{{{\mathrm{rad\,}}}B^q_{Y}} \quad \left( :=0\quad \text { if}\,{{\mathrm{rad\,}}}B^q_Y=0\right) , \end{aligned}$$\end{document}$$and denote by $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\varPhi ,Y}^{q},\mathcal {R}_{\varPhi ,Y}^{q})_{Y\subset Z}$$\end{document}$ the corresponding scheme of Chebyshev models for $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\varPhi ^q_Y$$\end{document}$,$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left( \mathcal {P}_{\varPhi ,Y}^{q},\mathcal {R}_{\varPhi ,Y}^{q}\right)&= \frac{(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})-{{\mathrm{mid\,}}}B^q_{Y}}{{{\mathrm{rad\,}}}B^q_{Y}} \quad \left( =(0,[0,0])\quad \text {if}\,{{\mathrm{rad\,}}}B^q_Y=0\right) , \end{aligned}$$\end{document}$$with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{\varPhi ,Y}^{q}(0)=0$$\end{document}$. Suppose also that a scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{F,W}^{q},\mathcal {R}_{F,W}^{q})_{W\subset X}$$\end{document}$ of Chebyshev models for *F* is available, with$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \zeta \in [-1,1], \quad \mathcal {P}_{F,W}^{q}(\zeta ) := \sum _{k=0}^q \varphi ^q_{W,k}\, T_k(\zeta ). \end{aligned}$$\end{document}$$We will discuss the computation of such a scheme at the end of this subsection.

In order to construct a valid scheme of *q*th-order Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{F\circ f,Y}^{q},\mathcal {R}_{F\circ f,Y}^{q})_{Y\subset Z}$$\end{document}$ for the composite function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F\circ f$$\end{document}$ in *Z*, we start by applying the Clenshaw recurrence formula ([1](#Equ1){ref-type=""}) in order to determine a scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\mathcal {P}_F\circ \varPhi , Y}^{q},\mathcal {R}_{\mathcal {P}_F\circ \varPhi , Y}^{q})_{Y\subset Z}$$\end{document}$ of Chebyshev models for the (parameterized) composite function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{F,B^q_{Y}}^{q}\circ \varPhi ^q_Y$$\end{document}$,$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} (\mathcal {P}_{\mathcal {P}_F\circ \varPhi , Y}^{q},\mathcal {R}_{\mathcal {P}_F\circ \varPhi ,Y}^{q})&:= \varphi ^q_{B^q_{Y},0} + \left( \mathcal {P}_{\varPhi ,Y}^{q},\mathcal {R}_{\varPhi ,Y}^{q}\right) \left( \mathcal {P}_{\beta _1,Y}^{q},\mathcal {R}_{\beta _1,Y}^{q}\right) \nonumber \\&\quad \quad -\left( \mathcal {P}_{\beta _2,Y}^{q},\mathcal {R}_{\beta _2,Y}^{q}\right) \end{aligned}$$\end{document}$$ $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \text {with:}\qquad \left( \mathcal {P}_{\beta _k,Y}^{q},\mathcal {R}_{\beta _k,Y}^{q}\right)&= \varphi ^q_{B^q_{Y},k} + 2 \left( \mathcal {P}_{\varPhi ,Y}^{q},\mathcal {R}_{\varPhi ,Y}^{q}\right) \left( \mathcal {P}_{\beta _{k+1},Y}^{q},\mathcal {R}_{\beta _{k+1},Y}^{q}\right) \nonumber \\&\qquad -\left( \mathcal {P}_{\beta _{k+2},Y}^{q},\mathcal {R}_{\beta _{k+2},Y}^{q}\right) , \quad k=0,\ldots ,q\end{aligned}$$\end{document}$$ $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left( \mathcal {P}_{\beta _{q+1},Y}^{q},\mathcal {R}_{\beta _{q+1},Y}^{q}\right)&= \left( \mathcal {P}_{\beta _{q+2},Y}^{q},\mathcal {R}_{\beta _{q+2},Y}^{q}\right) = (0,[0,0]). \end{aligned}$$\end{document}$$Observe that this recurrence consists of binary sum and product operations only, for which rules have already been established in Sects. [3.2](#Sec5){ref-type="sec"} and [3.3](#Sec6){ref-type="sec"}. Observe also that when $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$${{\mathrm{rad\,}}}B^q_Y=0$$\end{document}$, which occurs if either *f* is constant on *Y* or $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$${{\mathrm{rad\,}}}Y=0$$\end{document}$, we have$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left( \mathcal {P}_{\mathcal {P}_F\circ \varPhi , Y}^{q},\mathcal {R}_{\mathcal {P}_F\circ \varPhi ,Y}^{q}\right)&:= \left( \varphi ^q_{B^q_{Y},0}, [0,0]\right) . \end{aligned}$$\end{document}$$Finally, the polynomial and remainder parts of the desired Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{F\circ f,Y}^{q},\mathcal {R}_{F\circ f,Y}^{q})_{Y\subset Z}$$\end{document}$ are given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned}&\displaystyle \forall \xi \in [-1,1]^n,\quad \mathcal {P}_{F\circ f,Y}^{q}(\xi ) := \mathcal {P}_{\mathcal {P}_F\circ \varPhi ,Y}^{q}(\xi ) \end{aligned}$$\end{document}$$ $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned}&\displaystyle \mathcal {R}_{F\circ f,Y}^{q} := \mathcal {R}_{\mathcal {P}_F\circ \varPhi ,Y}^{q} + \mathcal {R}_{F,B^q_{Y}}^{q}. \end{aligned}$$\end{document}$$In order to use this approach, one needs the ability to construct Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{F,W}^{q},\mathcal {R}_{F,W}^{q})$$\end{document}$ corresponding to a given library of univariate functions, such as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F(x)=\exp (x)$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\log (x)$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\sqrt{x}$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\frac{1}{x}$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\sin (x)$$\end{document}$, \|*x*\|, etc.A systematic way of doing this is by computing a truncated Chebyshev expansion of *F* of order *q* as Then, if *F* is $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$-times continuously differentiable on *W*, the truncation error can be bounded as $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {R}_{F,W}^{q} = \frac{\left\| F^{(q+1)}\right\| _W ({{\mathrm{rad\,}}}W)^{q+1}}{2^{q}(q+1)!} [-1,1]. \end{aligned}$$\end{document}$$ Otherwise, if *F* and its derivatives through $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F^{(s-1)}$$\end{document}$ are absolutely continuous on *W* and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F^{(s)}$$\end{document}$ has finite total variation $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$V[F^{(s)}]_W<\infty $$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$s<q$$\end{document}$, a bound on the truncation error is given by $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {R}_{F,W}^{q} = \frac{2}{\pi }\frac{V\left[ F^{(s)}\right] _W ({{\mathrm{rad\,}}}W)^s}{s(q-s)^s} [-1,1]. \end{aligned}$$\end{document}$$Alternatively, a Chebyshev interpolating polynomial of order *q* may be computed as The error bound $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {R}_{F,W}^{q}$$\end{document}$ in ([26](#Equ26){ref-type=""}) remains valid for any $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$-times continuously differentiable function *F* on *W*. Otherwise, the following bound may be used instead of ([27](#Equ27){ref-type=""}) if *F* and its derivatives through $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F^{(s-1)}$$\end{document}$ are absolutely continuous on *W* and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F^{(s)}$$\end{document}$ has finite total variation $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$V[F^{(s)}]_W<\infty $$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$s<q$$\end{document}$, $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {R}_{F,W}^{q} = \frac{4}{\pi }\frac{V\left[ F^{(s)}\right] _W ({{\mathrm{rad\,}}}W)^s}{s(q-s)^s} [-1,1]. \end{aligned}$$\end{document}$$In comparison with their Taylor expansion counterparts for sufficiently smooth functions, we note that the width of the remainder is smaller by a factor of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$2^q$$\end{document}$. Moreover, the Chebyshev approach remains applicable for continuous, yet nonsmooth, functions. For instance, since the total variation of the derivative of \|*x*\| on the interval $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-\rho ,\rho ]$$\end{document}$ is finite, equal to 2, a bound for the interpolation error based on ([29](#Equ29){ref-type=""}) is $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {R}_{|\cdot |,[-\rho ,\rho ]}^{q} = [-\frac{8\rho }{\pi (q-1)},\frac{8\rho }{\pi (q-1)}]$$\end{document}$, for any $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q>1$$\end{document}$.

In practice, the bounds $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {R}_{F,W}^{q}$$\end{document}$ on the truncation error may be tightened further if *F* meets certain regularity and monotonicity conditions. Under these extra conditions, it can be established that the maximum absolute error between *F* and its polynomial approximant $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{F,W}^{q}$$\end{document}$---as given by ([25](#Equ25){ref-type=""}) or ([28](#Equ28){ref-type=""})---always occurs at one of the endpoints of *W*. This result is formalized in the following lemma:

### Lemma 1 {#FPar2}

Let the function *F* have an analytic extension on the closed unit disk in the complex plane and have all its successive derivatives $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F^{(k)}$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$k\ge 1$$\end{document}$ either of the same or alternating sign in $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$. Then, we have$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \Vert F-\mathcal {P}_{F}^q\Vert&= \max \left\{ \left| F(-1)-\mathcal {P}_{F}^q(-1)\right| , \left| F(1)-\mathcal {P}_{F}^q(1)\right| \right\} ,\\ \text {and} \quad \Vert F-\widehat{\mathcal {P}}_{F}^q\Vert&= \max \left\{ \left| F(-1)-\widehat{\mathcal {P}}_{F}^q(-1)\right| , \left| F(1)-\widehat{\mathcal {P}}_{F}^q(1)\right| \right\} \end{aligned}$$\end{document}$$for any $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\ge 0$$\end{document}$, with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{F}^q$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\widehat{\mathcal {P}}_{F}^q$$\end{document}$ given by ([4](#Equ4){ref-type=""}) and ([7](#Equ7){ref-type=""}), respectively.

### Proof {#FPar3}

See "Appendix 1". $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\square $$\end{document}$

This method of calculating the remainder is particularly useful as functions satisfying the conditions of Lemma [1](#FPar2){ref-type="sec"} include the exponential, logarithm, inverse and square-root functions on any finite interval in their domains of definition.

Range bounding {#Sec8}
--------------

A range bounder of the Chebyshev model $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})$$\end{document}$ of a function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f:Z\subseteq \mathbb {R}^n\rightarrow \mathbb {R}$$\end{document}$ on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Y\subset Z$$\end{document}$ is the set $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[\mathcal {P}_{f,Y}^{q}]+\mathcal {R}_{f,Y}^{q}\in \mathbb {IR}$$\end{document}$, with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[\mathcal {P}_{f,Y}^{q}] \supseteq \{\mathcal {P}^q_f(\xi ) \mid \xi \in [-1,1]^n\}$$\end{document}$. Such bounders are needed for computing binary product and univariate operations in Chebyshev model arithmetic---see Sects. [3.3](#Sec6){ref-type="sec"} and [3.4](#Sec7){ref-type="sec"} above, thus directly affecting the tightness of the computed Chebyshev models. Unfortunately, exact range bounding for multivariate polynomials is NP-hard, so one has to resort to some sort of over-approximation in practice.

A number of methods exist for bounding the range of multivariate polynomials in monomial form \[[@CR20], [@CR30], [@CR36], [@CR43]\], which can be readily adapted to their Chebyshev counterparts. Apart from the naive bounding of each term separately,$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left[ \mathcal {P}_{f,Y}^{q}\right] = a^q_{Y,0} + \sum _{0<|\kappa |\le q} \left| a^q_{Y,\kappa }\right| [-1,1], \end{aligned}$$\end{document}$$which may produce weak bounds, this subsection presents adaptations of: (i) the methods in \[[@CR33]\] that involves exact bounding of the first- and diagonal second-order terms; and (ii) the method of Bernstein.

### Exact bounding of first- and diagonal second-order terms {#Sec9}

We start by rewriting the polynomial part in the form$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}^q_{f,Y}({\xi })&= \sum _{j=1}^{n} \left[ a^q_{Y,\kappa ^{2,j}} T_2(\xi _j) + a^q_{Y,\kappa ^{1,j}} T_1(\xi _j)\right] + \mathcal {H}^q_{f,Y}(\xi ), \end{aligned}$$\end{document}$$where $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa ^{1,j}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa ^{2,j}$$\end{document}$ are the multi-indices given by $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa ^{1,j}_i := \delta _{i,j}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa ^{2,j}_i := 2\delta _{i,j}$$\end{document}$ for each $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$j\in \{1,\ldots ,n\}$$\end{document}$, using the Kronecker $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\delta $$\end{document}$ notation; and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {H}^q_{f,Y}$$\end{document}$ is the multivariate polynomial containing the same terms as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}^q_{f,Y}$$\end{document}$ apart from those indexed by $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa ^1_j$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa ^2_j$$\end{document}$. Then, a similar rearrangement to the one proposed in \[[@CR33]\] can be obtained by using the property that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$T_2(\xi )=2T_1(\xi )^2-1$$\end{document}$,$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}^q_{f,Y}({\xi })&= \sum _{j=1}^{n} \left[ 2a^q_{Y,\kappa ^{2,j}} \left( T_1(\xi _j)+\frac{a^q_{Y,\kappa ^{1,j}}}{4a^q_{Y,\kappa ^{2,j}}} \right) ^2 - a^q_{Y,\kappa ^{2,j}} - \frac{\left( a^q_{Y,\kappa ^{1,j}}\right) ^2}{8a^q_{Y,\kappa ^{2,j}}} \right] + \mathcal {H}^q_{f,Y}(\xi ). \end{aligned}$$\end{document}$$This way, an interval enclosure is obtained as$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left[ \mathcal {P}^q_{f,Y}\right] = \sum _{j=1}^{n} \left[ 2a^q_{Y,\kappa ^{2,j}} \left( [-1,1]+\frac{a^q_{Y,\kappa ^{1,j}}}{4a^q_{Y,\kappa ^{2,j}}} \right) ^2 - a^q_{Y,\kappa ^{2,j}} - \frac{\left( a^q_{Y,\kappa ^{1,j}}\right) ^2}{8a^q_{Y,\kappa ^{2,j}}} \right] + \left[ \mathcal {H}^q_{f,Y}\right] , \end{aligned}$$\end{document}$$with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[\mathcal {H}^q_{f,Y}]$$\end{document}$ computed directly from ([30](#Equ30){ref-type=""}). This bounder may only be considered when $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|a^q_{Y,\kappa ^{2,j}}|\ge \epsilon $$\end{document}$, where $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\epsilon $$\end{document}$ is a small positive number; otherwise, separate bounding of the linear and quadratic terms for the corresponding component *j* as in ([30](#Equ30){ref-type=""}) may be used.

### Bounding method of Bernstein {#Sec10}

The Bernstein algorithm is a well-established tool for computing bounds on the range of multivariate polynomials over a box; see, e.g., \[[@CR20], [@CR30]\] and references therein. The procedure involves rewriting the multivariate polynomial from Chebyshev form into Bernstein form:$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \xi \in [-1,1]^n, \quad \mathcal {P}^q_{f,Y}(\xi ) = \sum _{\kappa \in \{0,\ldots ,r\}^{n}} b^q_{Y,\kappa } \prod _{i=1}^{n} B^r_{\kappa _i}\left( \frac{\xi _i+1}{2}\right) , \end{aligned}$$\end{document}$$with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r\ge q$$\end{document}$, and where $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$B^r_i$$\end{document}$, with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$i\le r$$\end{document}$, denotes the *i*th Bernstein polynomial of order *r* on \[0, 1\], given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} B_i^r(z) = \left( {\begin{array}{c}r\\ i\end{array}}\right) z^i(1-z)^{k-i}. \end{aligned}$$\end{document}$$In particular, the Bernstein coefficient $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$b^q_{Y,\kappa }$$\end{document}$ in ([31](#Equ31){ref-type=""}) can be expressed in terms of the Chebyshev coefficients $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a^q_{Y,\kappa }$$\end{document}$ as follows \[[@CR48]\]$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \kappa \in \{0,\ldots ,r\}^{n}, \quad b^q_{Y,\kappa }&= \sum _{|\nu |\le q} a^q_{Y,\nu } \prod _{i=1}^{n} M_{\kappa _i,\gamma _i}\\ \text {with:} \quad M_{k,j}&:= \left( {\begin{array}{c}r\\ k\end{array}}\right) ^{-1} \sum _{i=\max \{0,j+k-r\}}^{\min \{j,k\}} (-1)^{j-i} \left( {\begin{array}{c}2j\\ 2i\end{array}}\right) \left( {\begin{array}{c}r-j\\ k-i\end{array}}\right) . \end{aligned}$$\end{document}$$This transformation turns out to be particularly well-conditioned from a numerical stability standpoint \[[@CR48]\]. At this point, an interval enclosure of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}^q_{f,Y}$$\end{document}$ on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]^n$$\end{document}$ is simply obtained as$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left[ \mathcal {P}^q_{f,Y}\right] = \left[ \min _{\nu \in \{0,\ldots ,r\}^{n}} b^q_{Y,\nu }, \max _{\nu \in \{0,\ldots ,r\}^{n}} b^q_{Y,\nu } \right] . \end{aligned}$$\end{document}$$The main advantage of this approach is that the foregoing enclosure converges to the actual polynomial range in the Hausdorff sense as the order *r* increases \[[@CR20]\]. Moreover, it is possible to bound the actual under- or over-estimation for a given order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r\ge q$$\end{document}$.

Convergence analysis {#Sec11}
====================

This section presents an analysis of the local convergence rate of Chebyshev models constructed from the application of the arithmetic operations described in Sect. [3](#Sec3){ref-type="sec"}. The global convergence of this arithmetic as the polynomial expansion order *q* increases is also discussed at the end of the section.

Definition 2 {#FPar4}
------------

Let the function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f:Z\rightarrow \mathbb {R}$$\end{document}$ be defined on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Z\subseteq \mathbb {R}^n$$\end{document}$. The scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})_{Y\subset Z}$$\end{document}$ of *q*th-order Chebyshev models for *f*, with$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}_{f,Y}^{q}(\xi ) := \sum _{|\kappa |\le q} a^{q}_{Y,\kappa } T_\kappa (\xi ), \end{aligned}$$\end{document}$$is said to have *local convergence order r on Z* if$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \kappa \in \mathbb {N}^n\quad \text {with}\,|\kappa |\le q, \quad a^{q}_{Y,\kappa }&= \mathcal O\left( ({{\mathrm{rad\,}}}Y)^{\min \{|\kappa |,r\}}\right) \end{aligned}$$\end{document}$$ $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \text {and}\quad {{\mathrm{rad\,}}}\mathcal {R}_{f,Y}^{q}&= \mathcal O\left( ({{\mathrm{rad\,}}}Y)^{r}\right) , \end{aligned}$$\end{document}$$for all interval vectors $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Y\subset Z$$\end{document}$ with sufficiently small $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$${{\mathrm{rad\,}}}Y$$\end{document}$.

Notice the extra conditions ([32](#Equ32){ref-type=""}) imposed on the coefficients of the polynomial approximant in the previous definition. While these conditions are not strictly necessary to establish local convergence, the propagation of (local) convergence order for the polynomial coefficients through binary sum, binary product and univariate outer-composition operations provides insight into the propagation of the convergence order for the remainder term in the following analysis.

As far as initialization is concerned, it is clear that a scheme of *q*th-order Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{z_i,Y}^{q},\mathcal {R}_{z_i,Y}^{q})_{Y\subset Z}$$\end{document}$ for the variable $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$z_i\in Z_i$$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$i\in \{1,\ldots ,n\}$$\end{document}$, as given by ([11](#Equ11){ref-type=""}) and ([12](#Equ12){ref-type=""}), has local convergence order 1 if $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=0$$\end{document}$, and infinite order if $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\ge 1$$\end{document}$. In both cases, the schemes of *q*th-order Chebyshev models for the variables thus have convergence order no less than $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$. Regarding the range bounders introduced in Sect. [3.5](#Sec8){ref-type="sec"}, the convergence condition ([32](#Equ32){ref-type=""}) imposes that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$${{\mathrm{rad\,}}}([\mathcal {P}_{f,Y}^{q}]+\mathcal {R}_{f,Y}^{q})=\mathcal O({{\mathrm{rad\,}}}Y)$$\end{document}$ whenever the scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})_{Y\subset Z}$$\end{document}$ has convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r\ge 1$$\end{document}$.

Next, we investigate the local convergence order of schemes of Chebyshev models as propagated through factorable functions, with the corresponding binary sum, binary product and univariate outer-composition operations derived in Sects. [3.2](#Sec5){ref-type="sec"}--[3.4](#Sec7){ref-type="sec"}. With regards to binary sum operations given by ([13](#Equ13){ref-type=""})--([14](#Equ14){ref-type=""}), the sum of two schemes of *q*th-order Chebyshev models trivially preserves the least convergence order of the operands, be they smaller than, equal to, or larger than $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$. In particular, the addition/subtraction of two schemes with local convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$ is itself a scheme of order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$. The propagation of convergence orders through binary product and univariate outer-composition operations is somewhat more involved, and detailed in Sects. [4.1](#Sec12){ref-type="sec"} and [4.2](#Sec13){ref-type="sec"}.

Local convergence rate of binary product operations {#Sec12}
---------------------------------------------------

Adopting the notation introduced in Sect. [3.3](#Sec6){ref-type="sec"}, we assume here that the schemes of *q*th-order Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1,Y}^{q},\mathcal {R}_{f_1,Y}^{q})_{Y\subset Z}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_2,Y}^{q},\mathcal {R}_{f_2,Y}^{q})_{Y\subset Z}$$\end{document}$ have local convergence orders $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r_1\ge 1$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r_2\ge 1$$\end{document}$ on *Z*, respectively.

The following lemma is instrumental to prove convergence of the Chebyshev coefficients in the product scheme:

### Lemma 2 {#FPar5}

For any $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\ge 0$$\end{document}$, any $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa \in \mathbb {N}^q$$\end{document}$, and any pair $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\lambda ,\mu )\in \mathbb {P}_q(\kappa )$$\end{document}$, we have $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|\lambda +\mu |\ge |\kappa |$$\end{document}$.

### Proof {#FPar6}

By construction of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathbb {P}_q(\kappa )$$\end{document}$ in ([15](#Equ15){ref-type=""}), we either have $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\lambda _i+\mu _i=\kappa _i$$\end{document}$ or $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|\lambda _i-\mu _i|=\kappa _i$$\end{document}$, for each $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$i\in \{1,\ldots ,n\}$$\end{document}$. Therefore, the result follows by noting that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\lambda _i+\mu _i\ge |\lambda _i-\mu _i|$$\end{document}$. $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\square $$\end{document}$

By the construction of the product polynomial in ([16](#Equ16){ref-type=""}) and by Lemma [2](#FPar5){ref-type="sec"}, we have$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left| c^q_{Y,\kappa }\right| \le \sum _{(\lambda ,\mu )\in \mathbb {P}_q(\kappa )} \frac{\left| a^q_{Y,\lambda }\right| \left| b^q_{Y,\mu }\right| }{2^{N(\lambda ,\mu )}} = \mathcal O\left( ({{\mathrm{rad\,}}}Y)^{\min \{r_1,r_2,|\kappa |\}}\right) , \end{aligned}$$\end{document}$$for every $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa \in \mathbb {N}^n$$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|\kappa |\le 2q$$\end{document}$. Therefore, the product rule ([17](#Equ17){ref-type=""}) propagates the least convergence order of the operands through the polynomial coefficients, up to order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|\kappa |$$\end{document}$. Moreover, we have $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[\mathcal {P}_{f_1,Y}^{q}] = \mathcal O(1)$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[\mathcal {P}_{f_2,Y}^{q}] = \mathcal O(1)$$\end{document}$, and$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned}&\sum _{q<|\kappa |\le 2q} \left| c^q_{Y,\kappa }\right| [-1,1] = \mathcal O\left( ({{\mathrm{rad\,}}}Y)^{\min \{r_1,r_2,q+1\}}\right) , \end{aligned}$$\end{document}$$so it follows from the remainder propagation rule ([18](#Equ18){ref-type=""}) that$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} {{\mathrm{rad\,}}}\mathcal {R}_{f_1 f_2,Y}^{q}&= \mathcal O\left( ({{\mathrm{rad\,}}}Y)^{\min \left\{ r_1,r_2,q+1\right\} }\right) . \end{aligned}$$\end{document}$$Overall, the least convergence order of the operands thus propagates through the remainder, up to order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$. In particular, the product scheme has local convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$ whenever both operands have convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$. These results are summarized in the following theorem:

### Theorem 1 {#FPar7}

Let $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1,Y}^{q},\mathcal {R}_{f_1,Y}^{q})_{Y\subset Z}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_2,Y}^{q},\mathcal {R}_{f_2,Y}^{q})_{Y\subset Z}$$\end{document}$ be schemes of *q*th-order Chebyshev models, with local convergence orders $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r_1\ge 1$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r_2\ge 1$$\end{document}$, respectively, for the functions $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1,f_2:Z\subseteq \mathbb {R}^n\rightarrow \mathbb {R}$$\end{document}$. Then, the product scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1 f_2,Y}^{q},\mathcal {R}_{f_1 f_2,Y}^{q})_{Y\subset Z}$$\end{document}$ given by ([17](#Equ17){ref-type=""}) and ([18](#Equ18){ref-type=""}) has local convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\min \{r_1,r_2,q+1\}$$\end{document}$.

Local convergence rate of univariate outer-composition operations {#Sec13}
-----------------------------------------------------------------

Using the notation introduced in Sect. [3.4](#Sec7){ref-type="sec"}, we assume that the outer functions $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F:X\subseteq \mathbb {R}\rightarrow \mathbb {R}$$\end{document}$ is *s*-times continuously differentiable, and we consider the scheme of *q*th-order Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{F,W}^q,\mathcal R_{F,W}^q)_{W\subset X}$$\end{document}$, as constructed from ([25](#Equ25){ref-type=""}), ([26](#Equ26){ref-type=""}) \[resp. from ([26](#Equ26){ref-type=""}), ([28](#Equ28){ref-type=""})\] if $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$s\ge q+1$$\end{document}$; or otherwise constructed from ([25](#Equ25){ref-type=""}), ([27](#Equ27){ref-type=""}) \[resp. from ([28](#Equ28){ref-type=""}), ([29](#Equ29){ref-type=""})\] if $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F^{(s)}$$\end{document}$ has a finite total variation on *X*.

### Lemma 3 {#FPar8}

The scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{F,W}^q,\mathcal R_{F,W}^q)_{W\subset X}$$\end{document}$ has local convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$${\min \{s,q+1\}}$$\end{document}$ on *X*.

### Proof {#FPar9}

See "Appendix 2". $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\square $$\end{document}$

To carry out the analysis, we also assume that the scheme of *q*th-order Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})_{Y\subset Z}$$\end{document}$ for the inner functions $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f:Z\subseteq \mathbb {R}^n\rightarrow \mathbb {R}$$\end{document}$ has local convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r\ge 1$$\end{document}$. Notice that the Clenshaw summations ([20](#Equ20){ref-type=""})--([22](#Equ22){ref-type=""}) consists of binary sum and product operations only and is initialized with constants. However, by construction of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\varPhi ,Y}^{q},\mathcal {R}_{\varPhi ,Y}^{q})_{Y\subset Z}$$\end{document}$ in ([19](#Equ19){ref-type=""}), we have$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned}&a^q_{Y,0}=0,\qquad \forall \kappa \in \mathbb {N}^n\quad \text {with}\, 1\le |\kappa |\le q,\quad a^q_{Y,\kappa }=\mathcal O\left( ({{\mathrm{rad\,}}}Y)^{\min \{|\kappa |,r\}-1}\right) ,\\&\quad \quad \text {and}\quad {{\mathrm{rad\,}}}\mathcal {R}_{\varPhi ,Y}^{q}=\mathcal O\left( ({{\mathrm{rad\,}}}Y)^{r-1}\right) , \end{aligned}$$\end{document}$$where $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a^q_{Y,\kappa }$$\end{document}$ denote the coefficients of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{\varPhi ,Y}^{q}$$\end{document}$. Even though $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\varPhi ,Y}^{q},\mathcal {R}_{\varPhi ,Y}^{q})_{Y\subset Z}$$\end{document}$ has local convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r-1$$\end{document}$ only, the following lemma shows that the original order *r* of the inner operand nonetheless propagates through the scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\mathcal {P}_F\circ \varPhi ,Y}^{q},\mathcal {R}_{\mathcal {P}_F\circ \varPhi ,Y}^{q})_{Y\subset Z}$$\end{document}$, up to order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$:

### Lemma 4 {#FPar10}

The scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\mathcal {P}_F\circ \varPhi ,Y}^{q},\mathcal {R}_{\mathcal {P}_F\circ \varPhi ,Y}^{q})_{Y\subset Z}$$\end{document}$ generated through the Clenshaw recursion ([20](#Equ20){ref-type=""})--([22](#Equ22){ref-type=""}) satisfies$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \kappa \in \mathbb {N}^n\quad \text {with}\, |\kappa |\le q, \quad c^q_{Y,\kappa }&= \mathcal O\left( ({{\mathrm{rad\,}}}Y)^{\min \{|\kappa |,r,s\}}\right) ,\\ {{\mathrm{rad\,}}}\mathcal {R}_{\mathcal {P}_F\circ \varPhi ,Y}^{q}&= \mathcal O\left( ({{\mathrm{rad\,}}}Y)^{\min \{r,q+1\}}\right) , \end{aligned}$$\end{document}$$where $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$c^q_{Y,\kappa }$$\end{document}$ denote the coefficients of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{\mathcal {P}_F\circ \varPhi ,Y}^{q}$$\end{document}$.

### Proof {#FPar11}

See "Appendix 3". $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\square $$\end{document}$

A direct consequence of Lemmata [3](#FPar8){ref-type="sec"} and [4](#FPar10){ref-type="sec"} is that the remainder term $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {R}_{F\circ f,Y}^{q}$$\end{document}$, as defined in ([24](#Equ24){ref-type=""}), has convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\min \{r,s,q+1\}$$\end{document}$. Therefore, the convergence order of the inner operand may propagate through the remainder, up to order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$, in a *q*th-order scheme of Chebyshev models. Moreover, if the outer operand is *s*-times continuously differentiable, the convergence order of the composition scheme may not be larger than *s*. In particular, the composition scheme has convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$ whenever the inner scheme has convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$ and the outer function is $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$-times continuously differentiable. These results are summarized in the following theorem:

### Theorem 2 {#FPar12}

Let $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{F,W}^{q},\mathcal {R}_{F,W}^{q})_{W\subset X}$$\end{document}$ be a scheme of *q*th-order Chebyshev models for the *s*-times continuously-differentiable function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F:X\subseteq \mathbb {R}\rightarrow \mathbb {R}$$\end{document}$, as constructed from ([25](#Equ25){ref-type=""}), ([26](#Equ26){ref-type=""}) \[resp. from ([26](#Equ26){ref-type=""}), ([28](#Equ28){ref-type=""})\] if $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$s\ge q+1$$\end{document}$; or otherwise constructed from ([25](#Equ25){ref-type=""}), ([27](#Equ27){ref-type=""}) \[resp. from ([28](#Equ28){ref-type=""}), ([29](#Equ29){ref-type=""})\] if $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F^{(s)}$$\end{document}$ has a finite total variation on *X*. Let also $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})_{Y\subset Z}$$\end{document}$ be a scheme of *q*th-order Chebyshev models for the function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f:Z\subseteq \mathbb {R}^n\rightarrow \mathbb {R}$$\end{document}$ with local convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$r\ge 1$$\end{document}$. Then, the composition scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{F\circ f,Y}^{q},\mathcal {R}_{F\circ f,Y}^{q})_{Y\subset Z}$$\end{document}$ given by ([23](#Equ23){ref-type=""}) and ([24](#Equ24){ref-type=""}) has convergence order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\min \{r,s,q+1\}$$\end{document}$ on *Z*.

Global convergence of Chebyshev model arithmetic {#Sec14}
------------------------------------------------

The previous convergence analysis investigates the rate at which the remainder bound shrinks as the domain *Y* of the variables is progressively reduced to a point, for a given expansion order *q*. This section discusses a different type of convergence, namely the property of the remainder bound to shrink to zero upon increasing the expansion order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\rightarrow \infty $$\end{document}$, for a fixed variable domain *Y*.

### Definition 3 {#FPar13}

Let the function $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f:Z\rightarrow \mathbb {R}$$\end{document}$ be defined on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Z\subseteq \mathbb {R}^n$$\end{document}$, and let $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathbb {IR}^n\ni Y\subset Z$$\end{document}$. The scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f,Y}^{q},\mathcal {R}_{f,Y}^{q})_{q\ge 0}$$\end{document}$ of *q*th-order Chebyshev models for *f*, with$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \mathcal {P}_{f,Y}^{q}(\xi ) := \sum _{|\kappa |\le q} a^{q}_{Y,\kappa } T_\kappa (\xi ), \end{aligned}$$\end{document}$$is said to be *globally convergent on Y* if$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \lim _{q\rightarrow \infty } \sum _{|\kappa |\le q} \left| a^{q}_{Y,\kappa }\right|&= \bar{A}_Y < \infty \quad \text {and}\quad \lim _{q\rightarrow \infty } \mathcal {R}_{f,Y}^{q} = \{0\}. \end{aligned}$$\end{document}$$

As far as initialization is concerned, it is clear that a scheme of *q*th-order Chebyshev models $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{y_i,Y}^{q},\mathcal {R}_{y_i,Y}^{q})_{q\ge 0}$$\end{document}$ for the variable $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$y_i\in Y_i$$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$i\in \{1,\ldots ,n\}$$\end{document}$, as given by ([11](#Equ11){ref-type=""}) and ([12](#Equ12){ref-type=""}), is globally convergent on *Y*. Regarding the range bounders introduced in Sect. [3.5](#Sec8){ref-type="sec"}, it is also worth noting that the sequence $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\{{{\mathrm{rad\,}}}[\mathcal {P}_{f,Y}^{q}]\}_{q\ge 0}$$\end{document}$ is bounded when ([36](#Equ36){ref-type=""}) holds.

Next, we adopt the notation introduced in Sects. [3.2](#Sec5){ref-type="sec"} and [3.3](#Sec6){ref-type="sec"}, and we assume that the schemes $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1,Y}^{q},\mathcal {R}_{f_1,Y}^{q})_{q\ge 0}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_2,Y}^{q},\mathcal {R}_{f_2,Y}^{q})_{q\ge 0}$$\end{document}$ are globally convergent on *Y*. The summation/subtraction scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1\pm f_2,Y}^{q},\mathcal {R}_{f_1\pm f_2,Y}^{q})_{q\ge 0}$$\end{document}$, as given in ([14](#Equ14){ref-type=""}), is trivially globally convergent by the triangle inequality. Concerning the product operation in Sect. [3.3](#Sec6){ref-type="sec"}, we have from ([16](#Equ16){ref-type=""}) that$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \sum _{|\kappa |\le q} \left| c^{q}_{Y,\kappa }\right|&= \sum _{|\kappa |\le q} \left| \sum _{(\lambda ,\mu )\in \mathbb {P}_q(\kappa )} \frac{a^q_{Y,\lambda } b^q_{Y,\mu }}{2^{N(\lambda ,\mu )}}\right| \\&\le \sum _{|\kappa |\le q} \sum _{(\lambda ,\mu )\in \mathbb {P}_q(\kappa )} \frac{\left| a^q_{Y,\lambda }\right| \left| b^q_{Y,\mu }\right| }{2^{N(\lambda ,\mu )}}\\&= \sum _{|\lambda |\le q} \sum _{|\mu |\le q} \left| a^q_{Y,\lambda }\right| \left| b^q_{Y,\mu }\right| , \end{aligned}$$\end{document}$$so$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \lim _{q\rightarrow \infty } \sum _{|\kappa |\le q} \left| c^{q}_{Y,\kappa }\right| \le \lim _{q\rightarrow \infty } \sum _{|\lambda |\le q} \left| a^{q}_{Y,\lambda }\right| \sum _{|\mu |\le q}\left| b^{q}_{Y,\mu }\right| < \infty , \end{aligned}$$\end{document}$$and, in particular,$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \lim _{q\rightarrow \infty } \sum _{q<|\kappa |\le 2q} \left| c^{q}_{Y,\kappa }\right| \le \lim _{q\rightarrow \infty }\sum _{|\kappa |> q} \left| c^{q}_{Y,\kappa }\right| = 0. \end{aligned}$$\end{document}$$It then follows from ([18](#Equ18){ref-type=""}) that the product scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1 f_2,Y}^{q},\mathcal {R}_{f_1 f_2,Y}^{q})_{q\ge 0}$$\end{document}$ is globally convergent on *Y*, as long as the sequences $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\{{{\mathrm{rad\,}}}[\mathcal {P}_{f_1,Y}^{q}]\}_{q\ge 0}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\{{{\mathrm{rad\,}}}[\mathcal {P}_{f_2,Y}^{q}]\}_{q\ge 0}$$\end{document}$ are bounded. These considerations are summarized in the following theorem:

### Theorem 3 {#FPar14}

Let $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1,Y}^{q},\mathcal {R}_{f_1,Y}^{q})_{Y\subset Z}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_2,Y}^{q},\mathcal {R}_{f_2,Y}^{q})_{Y\subset Z}$$\end{document}$ be globally convergent schemes of Chebyshev models for the functions $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$f_1,f_2:Z\subseteq \mathbb {R}^n\rightarrow \mathbb {R}$$\end{document}$, respectively. Then, the summation/subtraction scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1\pm f_2,Y}^{q},\mathcal {R}_{f_1\pm f_2,Y}^{q})_{q\ge 0}$$\end{document}$ given by ([13](#Equ13){ref-type=""}) and ([14](#Equ14){ref-type=""}) is globally convergent on *Y*; and so is the product scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{f_1 f_2,Y}^{q},\mathcal {R}_{f_1 f_2,Y}^{q})_{q\ge 0}$$\end{document}$ given by ([17](#Equ17){ref-type=""}) and ([18](#Equ18){ref-type=""}) for the range bounding methods in Sect. [3.5](#Sec8){ref-type="sec"}.

In order to formally establish global convergence of the Chebyshev model arithmetic presented in Sect. [3](#Sec3){ref-type="sec"}, one also needs to establish sufficient conditions on both the inner and outer operands of a composition operation in order for the resulting composition scheme to be globally convergent. The development of such conditions calls for an analysis of the propagation of the global convergence property through the Clenshaw recurrence, for which results are currently lacking. The numerical case studies presented in the following sections suggest that global convergence might be guaranteed under mild regularity conditions nonetheless, and this will be the topic of further research.

Numerical implementation and case studies {#Sec15}
=========================================

This section presents numerical case studies that illustrate the convergence properties of Chebyshev models, and makes comparisons with Taylor models in terms of tightness and computational effort. All the computations that led to these results are performed with the library MC++ (<http://omega-icl.bitbucket.org/mcpp/>), which features classes for the evaluation of factorable functions in Taylor and Chebyshev arithmetics (along with other arithmetics). Verified interval libraries, such as PROFIL (<http://www.ti3.tu-harburg.de/>) or FILIB++ ([http://www2.math.uni-wuppertal.de/\~xsc/](http://www2.math.uni-wuppertal.de/%7exsc/)), are used to perform the remainder interval computations, but any round-off caused by operations between the polynomial coefficients is not accounted for in the current implementation. We note that approaches to account for errors due to floating-point arithmetic are well-documented for Taylor models \[[@CR36], [@CR43]\], and the very same approaches may be used for Chebyshev models in order to arrive at a fully verified implementation.

MC++ implements the Chebyshev model arithmetic and range bounding operations described in Sect. [3](#Sec3){ref-type="sec"}, based on a dense representation of the multivariate polynomials. If not otherwise noted, Chebyshev models for the intrinsic functions are computed based on the Chebyshev interpolation polynomial. Moreover, the exact remainder bounding approach established in Sect. [3.4](#Sec7){ref-type="sec"} is used for the intrinsic functions $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\exp (x)$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\log (x)$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\sqrt{x}$$\end{document}$, and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\frac{1}{x}$$\end{document}$, and the default range bounder uses exact bounding of first- and diagonal second-order terms as described in Sect. [3.5](#Sec8){ref-type="sec"}. Another improvement involves intersecting the polynomial model bounds with those derived from simple interval analysis at each operation as a means of tightening the remainder bound, avoiding division by zero, etc; this extension is similar in essence to the mixed affine-arithmetic/interval-arithmetic model by Stolfi and Figueiredo \[[@CR19]\].

The most critical operation in an efficient implementation is multiplication of Chebyshev models, since it is used extensively in univariate outer-composition operations besides bivariate product operations. The approach used in MC++ involves constructing the index sets $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathbb {P}_q(\kappa )$$\end{document}$ given by ([15](#Equ15){ref-type=""}), for all $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\kappa \in \mathbb {N}^n$$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|\kappa |\le q$$\end{document}$, prior to propagating the Chebyshev models through a DAG of the factorable function. The number of floating-point operations (FLOP) required to compute the product polynomial $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{f_1 f_2,Y}^{q}$$\end{document}$ in ([17](#Equ17){ref-type=""}) is shown on the left plot of Fig. [1](#Fig1){ref-type="fig"} for Chebyshev expansions of order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=1,\ldots ,4$$\end{document}$ and for functions with up to $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$n=19$$\end{document}$ variables. Even for low expansion orders does performing polynomial products in Chebyshev basis result in an increase by 1 or 2 orders of magnitude in the number of FLOP, compared with monomial basis which is shown on the right plot. The better approximation capability of Chebyshev models may thus come at the price of a much higher computational burden compared with Taylor models, at least in a dense implementation. This refinement calls for the development of sparse implementations for both Taylor and Chebyshev model arithmetics in order to overcome this limitation and tackle larger-scale problems, a topic for future research.Fig. 1Number of FLOP for performing polynomial product in monomial basis (*left plot*) and in Chebyshev basis (*right plot*)

Example 1 {#Sec16}
---------

We consider the univariate function *f* given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} f(x) = \exp \left( -x^2+\frac{1}{x}\right) , \end{aligned}$$\end{document}$$for the variable $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x\in X:=[0.3,2]$$\end{document}$. A number of Chebyshev model enclosures for various expansion orders *q* in the range 2-7 are shown on the top-left plot of Fig. [2](#Fig2){ref-type="fig"}, along with the corresponding approximation errors on the opposite plot. The scheme of Chebyshev models appears to be globally convergent as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\rightarrow \infty $$\end{document}$, both with the derivative-based remainder formula and with the exact remainder formula for the univariate outer-composition operations, as shown on the bottom-right plot of Fig. [2](#Fig2){ref-type="fig"}. However, the use of the exact remainder significantly tightens the bounds and improves the rate of convergence. In contrast, classical Taylor models fail to converge to the function as the order *q* is increased, despite the function *f* being analytic on the set *X*. This is because the variable range of interest is partially outside the radius of convergence for this function, which also means that the width of the interval remainder increases with *q*. Finally, the bottom-left plot on Fig. [2](#Fig2){ref-type="fig"} shows that local convergence of *q*th-order Chebyshev model is of order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q+1$$\end{document}$ around $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x=1$$\end{document}$. This observation is consistent with the analysis conducted in Sect. [4](#Sec11){ref-type="sec"} since *f* is analytic.Fig. 2Application of Chebyshev models to Example 1. *Top-left plot* polynomial model enclosures of the function *f* (*black line*) for expansion orders $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=2,3,5,7$$\end{document}$. *Top-right plot* approximation error for expansion orders $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=2,3,5,7$$\end{document}$. *Bottom-left plot* local convergence analysis for a sequence of shrinking interval sets $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Y\subset X$$\end{document}$ around $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x=1$$\end{document}$. *Bottom-right plot* global convergence analysis on *X* as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\rightarrow \infty $$\end{document}$ with the exact remainder estimate compared with the high-order derivative-based remainder formula Fig. 3Application of Chebyshev models to Example 2. *Top-left plot* Polynomial model enclosures of the function *f* (*black line*) for expansion orders $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=2,4,8,16$$\end{document}$. *Top-right plot* approximation error for expansion orders $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=2,4,8,16$$\end{document}$. *Bottom-left plot* local convergence analysis for a sequence of shrinking interval sets $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$Y\subset X$$\end{document}$ around $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x=0$$\end{document}$. *Bottom-right plot* global convergence analysis on *X* as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\rightarrow \infty $$\end{document}$ with the derivative-based remainder formula or the reformulation as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|x|=\sqrt{x^2}$$\end{document}$

Example 2 {#Sec17}
---------

We consider the univariate function *f* given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} f(x) = \sin |x|, \end{aligned}$$\end{document}$$for the variable $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x\in X:=\left[ -\frac{\pi }{2},\frac{\pi }{4}\right] $$\end{document}$. Notice that classical Taylor models may not be used to bound the function *f* on *X* given the nonsmoothness at $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x=0$$\end{document}$. In contrast, Chebyshev models can be constructed since the univariate $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\sin (x)$$\end{document}$ is analytic in $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathbb {R}$$\end{document}$, whereas the derivative of \|*x*\| has total bounded variation on *X*. The top plots in Fig. [3](#Fig3){ref-type="fig"} shows how the bounds improve with a larger expansion order *q*, up to $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=16$$\end{document}$. It is evident that a much higher-order expansion is necessary in comparison with the smooth function in the previous example (Sect. [5.1](#Sec16){ref-type="sec"}), and that the approximation error is the largest close to $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x=0$$\end{document}$. In agreement with the theory in Sect. [4](#Sec11){ref-type="sec"}, the local convergence of any *q*th-order Chebyshev model at this point is of order 1, regardless of the expansion order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q \ge 0$$\end{document}$; a prediction that is confirmed in the bottom-left plot of Fig. [3](#Fig3){ref-type="fig"}. But despite the nonsmoothness at $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x=0$$\end{document}$ and the linear local convergence around this point, Chebyshev models are found to be globally convergent on *X* as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\rightarrow \infty $$\end{document}$, when the remainder formula ([29](#Equ29){ref-type=""}) is used for bounding the remainder of the univariate term \|*x*\| (bottom-right plot on Fig. [2](#Fig2){ref-type="fig"}).

Notice that this function may also be reformulated by substitution of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|x| = \sqrt{x^2}$$\end{document}$, thus providing an alternative route to creating the Chebyshev model. This reformulation makes the use of the exact remainder estimate for the square-root term possible, which explains why tighter bounds for Chebyshev models of order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=5$$\end{document}$ or less are obtained here with the reformulation approach (bottom-right plot of Fig. [3](#Fig3){ref-type="fig"}). However, the reformulation also causes Chebyshev models to diverge on *X* as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\rightarrow \infty $$\end{document}$, a behavior attributed to the fact that the total variation of the square-root term is unbounded on intervals containing zero.Fig. 4Application of Chebyshev models to Example 3. *Top-left plot* polynomial model enclosing the function *f* for $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$n=2$$\end{document}$ with expansion orders $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=8$$\end{document}$. *Bottom-left plot* corresponding approximation error for expansion orders $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=8$$\end{document}$. *Right plot* computational effort and global convergence of the Chebyshev remainder $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$${{\mathrm{rad\,}}}\mathcal R^q_{f,X}$$\end{document}$ with expansion orders $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=2,\ldots ,20$$\end{document}$ for dimensions $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$n=2$$\end{document}$ to 4 (*solid lines*), and comparison with their Taylor models counterparts (*dashed lines*)

Example 3 {#Sec18}
---------

As a final example, we consider the multivariate function *f* given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} f(x) = \exp \left( \sum _{i=1}^n x_i\right) \cos \left( 2\pi \sum _{i=1}^n x_i\right) , \end{aligned}$$\end{document}$$for the variable $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$x\in X := [-0.6,0.6]^n$$\end{document}$, with the dimension *n* taken as 2, 3 or 4. The Chebyshev model enclosure and corresponding approximation error in the 2-variable case are presented on the left part of Fig. [4](#Fig4){ref-type="fig"} for an expansion of order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=8$$\end{document}$, showing a tight approximation of the function. From the right plot of Fig. [4](#Fig4){ref-type="fig"}, it is seen that the computational burden increases with both the expansion order *q* and the dimension *n* of the problem. Moreover, computing Chebyshev models (solid lines) is more demanding than Taylor models (dashed lines), due to the extra complexity of the polynomial product in Chebyshev basis (see Fig. [1](#Fig1){ref-type="fig"} above and related discussion). However, this extra complexity pays off in terms of the tightness of the resulting bounds. For $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$n=2$$\end{document}$, one must consider expansions of order greater than $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=18$$\end{document}$ in order for a Taylor model's remainder bound to be as tight as those given by a second-order Chebyshev models only. For $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$n=3$$\end{document}$ or $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$n=4$$\end{document}$ likewise, although the Taylor models appear to be globally convergent on *X*, an expansion order far greater than $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q=20$$\end{document}$ would be necessary in order to be competitive with just second-order Chebyshev models. Overall, Chebyshev models thus outperform their Taylor counterparts. Finally, it is interesting to note that as the function *f* being even, odd ordered polynomial models prove to provide tighter bounds than even ones, thus explaining why there is sometimes little improvement (or even deterioration) in the bounds even though the computational time increases.

Conclusions {#Sec19}
===========

This paper has formalized an arithmetic for the propagation of Chebyshev models through binary sums, binary products and univariate outer-composition operations in factorable functions. A simple method to calculate the exact remainder bounds for certain univariate terms, including $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\exp (\cdot )$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\log (\cdot )$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\sqrt{\cdot }$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\cdot )^{-1}$$\end{document}$, is established. Adaptations of existing range bounders for Taylor models in order to use them with Chebyshev models are also discussed, which are a vital part of polynomial model arithmetic. The local convergence of Chebyshev models has been analyzed and proven to be equivalent to Taylor models, although the Chebyshev model remainder bounds are often found to be orders of magnitude better than their Taylor model counterparts over larger variable domains. This behavior is supported by the result that Chebyshev models are globally convergent through addition/subtraction and multiplication operations. The global convergence property is also conjectured for composition of Chebyshev models under mild conditions, although sufficient conditions are yet to be formally established. Such convergence properties have been illustrated through several numerical examples, based on an implementation in the authors' in-house library MC++. These examples also illustrate some of the advantages of Chebyshev models over Taylor models, including being able to bound functions with points of non-smoothness or a divergent Taylor expansion. Even though Chebyshev models of equivalent order are computationally more expensive to create than their Taylor counterparts, mainly due to the much greater number of operations required for the binary product operation, they provide benefits by being significantly tighter. This may allow for a much lower order polynomial model to be used and creates a net benefit in terms of computational effort.

On the computational side, the current implementation creates a polynomial which is of the same order for all the variables and is stored using a dense representation. For many applications however, a large number of coefficients may be equal to zero, thereby calling for a sparse implementation of Chebyshev model arithmetic. Such an implementation would reduce both the computational time and the memory requirements. Another area of improvement would concern the binary product operation, whose complexity currently scales exponentially in the number of variables and the expansion order. Several studies have explored approaches to speeding up the multiplication of univariate and bivariate polynomials in Chebyshev basis (e.g., \[[@CR2], [@CR5], [@CR22], [@CR45]\]). Ways to use some of these developments for multivariate Chebyshev polynomials could be explored as part of future work. Binary product operations are also where much of the overestimation occurs, given that the polynomial part of the Chebyshev model product may not match the Chebyshev expansion of the product function. Instead, it might be possible to directly create multivariable Chebyshev approximations, for instance by following a similar approach as in chebfun2 \[[@CR56]\]. This idea is also similar to the recent work on multivariate McCormick relaxations \[[@CR59]\], where instead of decomposing the factorable function down to binary sums, products and univariate compositions, it would be possible to directly bound some multivariate terms.

Appendix 1: Proof of Lemma [1](#FPar2){ref-type="sec"} {#Sec20}
======================================================

Since *F* in Lipschitz continuous on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$, it has a unique representation as an absolutely and uniformly convergent Chebyshev series, and we have$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \xi \in [-1,1], \quad F(\xi )-\mathcal {P}^q_F(\xi ) = \sum _{k=q+1}^\infty a_k T_k(\xi ), \end{aligned}$$\end{document}$$with the Chebyshev coefficients given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} a_k:=\frac{2}{\pi }\int _0^\pi f(\cos (\theta ))\cos (k\theta )d\theta . \end{aligned}$$\end{document}$$Observe that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\Vert F-\mathcal {P}^q_F\Vert = \sum \nolimits _{k=q+1}^\infty |a_k|$$\end{document}$ under the following scenarios:(i)the coefficients $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a_k$$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$k\ge q+1$$\end{document}$ are of the same sign, either non-negative or non-positive, in which case $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|F(1)-\mathcal {P}^q_F(1)| = \sum \nolimits _{k=q+1}^\infty |a_k|$$\end{document}$ by the property that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$T_k(1)=1, \forall k\ge 0$$\end{document}$;(ii)the coefficients $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a_k$$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$k\ge q+1$$\end{document}$ alternate sign, either non-negative for even coefficients and non-positive for odd coefficients or vice-versa, in which case $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|F(-1)-\mathcal {P}^q_F(-1)| = \sum \nolimits _{k=q+1}^\infty |a_k|$$\end{document}$ by the property that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$T_k(-1)=(-1)^k, \forall k\ge 0$$\end{document}$.Without loss of generality, we shall prove that case (i) holds with non-negative Chebyshev coefficients $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a_k$$\end{document}$ when *F* has an analytic extension on the closed unit disk and all its successive derivatives are nonnegative on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$; the other three cases follow readily by symmetry, considering $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$-F(\xi )$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F(-\xi )$$\end{document}$ or $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$-F(-\xi )$$\end{document}$.

Under the assumption that *F* has an analytic extension on the closed unit disk, its Taylor series at any point $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\bar{\xi }\in [-1,1]$$\end{document}$ converges to *F* on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$. At $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\bar{\xi }=0$$\end{document}$ in particular, we have$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \xi \in [-1,1], \quad F(\xi ) = \sum _{k=0}^\infty \frac{F^{(k)}(0)}{k!}\xi ^k = \sum _{k=0}^\infty a_k T_k(\xi ^k). \end{aligned}$$\end{document}$$Recall that any polynomial $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}^q$$\end{document}$ of order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\ge 1$$\end{document}$, given in monomial basis by $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}^q(\xi ):=\sum \nolimits _{k=0}^q \alpha _k\xi ^k$$\end{document}$, has a unique representation in Chebyshev basis as $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}^q(\xi )=\sum \nolimits _{k=0}^q \beta _kT_k(\xi )$$\end{document}$ with$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall k\in \{1,\ldots , q\}, \quad \beta _k := \sum _{i=1}^q \frac{\alpha _i}{2^i} \sum _{j=0}^i \left( {\begin{array}{c}i\\ j\end{array}}\right) \delta _{k,|i-2j|}. \end{aligned}$$\end{document}$$By passing to the limit $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$q\rightarrow \infty $$\end{document}$, we obtain from ([38](#Equ38){ref-type=""}) that$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall k\ge 1, \quad a_k := \sum _{i=1}^\infty \frac{F^{(i)}(0)}{2^i i!} \sum _{j=0}^i \left( {\begin{array}{c}i\\ j\end{array}}\right) \delta _{k,|i-2j|}. \end{aligned}$$\end{document}$$If follows from the nonnegativity of all the successive derivatives $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$F^{(i)}$$\end{document}$ that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a_k\ge 0$$\end{document}$ for every $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$k\ge 1$$\end{document}$.

Regarding Chebyshev interpolation, on the other hand, we have$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \xi \in [-1,1], \quad F(\xi )-\widehat{\mathcal {P}}^q_F(\xi ) = \sum _{k=0}^q \left( a_k-\widehat{a}_k\right) T_k(\xi ) + \sum _{k=q+1}^\infty a_k T_k(\xi ), \end{aligned}$$\end{document}$$with the coefficients $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\widehat{a}_k$$\end{document}$ given by$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \widehat{a}_k := \frac{2}{q+1}\sum _{j=0}^q F(\zeta _j) T_k(\zeta _j) \quad \text {with} \zeta _j := \cos \left( \frac{\pi }{2}\frac{j+1}{q+1}\right) . \end{aligned}$$\end{document}$$Like previously, we shall prove that$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \left\| F -\widehat{\mathcal {P}}^q_F\right\| = \left| F(1)-\widehat{\mathcal {P}}^q_F(1)\right| = \pm \left( a_0-\widehat{a}_0\right) + \sum _{k=1}^q \left( a_k-\widehat{a}_k\right) + \sum _{k=q+1}^\infty a_k \end{aligned}$$\end{document}$$in the case that *F* has an analytic extension on the closed unit disk and all its successive derivatives are nonnegative on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[-1,1]$$\end{document}$---the other cases following by symmetry. Since $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a_k\ge 0$$\end{document}$ for every $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$k\ge 1$$\end{document}$ under these assumptions, it is sufficient to prove that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a_k \ge \widehat{a}_k$$\end{document}$ for each $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$k=1,\ldots , q$$\end{document}$. By the aliasing property of the Chebyshev polynomials \[[@CR37]\], we have$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall k\in \{0,\ldots ,q\},\quad \widehat{a}_k = a_k - a_{2q+2-k} - a_{2q+2+k} + a_{4q+4-k} + a_{4q+4+k} - \cdots \end{aligned}$$\end{document}$$and, therefore, the result follows by showing that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$a_k\ge a_{k+2}$$\end{document}$ for every $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$k\ge 1$$\end{document}$. Considering the case of odd coefficients first, ([39](#Equ39){ref-type=""}) gives$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} a_{2k-1}&= \sum _{i=k}^\infty \frac{F^{(2i-1)}(0)}{2^{2i-2} (2i-1)!} \left( {\begin{array}{c}2i-1\\ i-k\end{array}}\right) \\&= \frac{F^{(2k-1)}(0)}{2^{2k-2} (2k-1)!} + \sum _{i=k+1}^\infty \frac{F^{(2i-1)}(0)}{2^{2i-2} (2i-1)!} \left( {\begin{array}{c}2i-1\\ i-k\end{array}}\right) \\&= \frac{F^{(2k-1)}(0)}{2^{2k-2} (2k-1)!} + \sum _{i=k+1}^\infty \frac{F^{(2i-1)}(0)}{2^{2i-2} (2i-1)!} \left( {\begin{array}{c}2i-1\\ i-k-1\end{array}}\right) \frac{i+k}{i-k}\\&\ge \frac{F^{(2k-1)}(0)}{2^{2k-2} (2k-1)!} + a_{2k+1}, \end{aligned}$$\end{document}$$and a similar result holds for even coefficients. $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\square $$\end{document}$

Appendix 2: Proof of Lemma [3](#FPar8){ref-type="sec"} {#Sec21}
======================================================

Let the Chebyshev coefficients of the function *F* on the subset $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$W\subset X$$\end{document}$ be denoted as$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall k\ge 1, \quad a_{W,k}:=\frac{2}{\pi }\int _{-1}^{1} \frac{F({{\mathrm{mid\,}}}W+\xi {{\mathrm{rad\,}}}W)T_k(\xi )}{\sqrt{1-\xi ^2}}d\xi . \end{aligned}$$\end{document}$$Since *F* is *s*-times continuously differentiable, Taylor's theorem asserts that, for all $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\xi \in [-1,1]$$\end{document}$, there exists $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\zeta _\xi \in [-1,1]$$\end{document}$ such that$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} F({{\mathrm{mid\,}}}W+\xi {{\mathrm{rad\,}}}W) = \sum _{i=0}^{s-1} \frac{F^{(i)}({{\mathrm{mid\,}}}W)}{i!}(\xi {{\mathrm{rad\,}}}W)^i + \frac{F^{(s)}(\zeta _\xi )}{s!}(\xi {{\mathrm{rad\,}}}W)^s. \end{aligned}$$\end{document}$$Therefore, letting $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$M_i := \Vert F^{(i)}\Vert _X$$\end{document}$ for $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$0\le i\le s$$\end{document}$, we have$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} |a_{W,k}|&\le \sum _{i=0}^{s-1} \frac{2({{\mathrm{rad\,}}}W)^i M_i}{\pi \,i!} \left| \int _{-1}^{1} \frac{\xi ^i T_k(\xi )}{\sqrt{1-\xi ^2}}d\xi \right| + \frac{2({{\mathrm{rad\,}}}W)^s M_s}{\pi \,s!} \int _{-1}^{1} \frac{\left| \xi ^s T_k(\xi )\right| }{\sqrt{1-\xi ^2}}d\xi . \end{aligned}$$\end{document}$$Since $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\int _{-1}^{1} \frac{\xi ^i T_k(\xi )}{\sqrt{1-\xi ^2}}d\xi =0$$\end{document}$ for all $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$i<k$$\end{document}$, $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$|\int _{-1}^{1} \frac{\xi ^i T_k(\xi )}{\sqrt{1-\xi ^2}}d\xi |\le \pi $$\end{document}$ for all $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$i\ge 0$$\end{document}$, and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\int _{-1}^{1} \frac{|\xi ^s T_k(\xi )|}{\sqrt{1-\xi ^2}}d\xi \le \pi $$\end{document}$, we obtain$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall k\ge 1, \quad |a_{W,k}|&\le ({{\mathrm{rad\,}}}W)^k \sum _{i=k}^{s} \frac{2({{\mathrm{rad\,}}}W)^{i-k}M_i}{i!} + ({{\mathrm{rad\,}}}W)^s \frac{2M_s}{s!}\\&= \mathcal O\left( ({{\mathrm{rad\,}}}W)^{\min \{k,s\}} \right) . \end{aligned}$$\end{document}$$Using Chebyshev interpolation likewise, a direct consequence of the aliasing property ([41](#Equ41){ref-type=""}) is that$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall k\in \{1,\ldots ,q\}, \quad \widehat{a}_{W,k}&= \mathcal O\left( ({{\mathrm{rad\,}}}W)^{\min \{k,s\}} \right) . \end{aligned}$$\end{document}$$The result of the lemma follows by noting that $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$${{\mathrm{rad\,}}}\mathcal {R}_{F,W}^{q} = \mathcal O(({{\mathrm{rad\,}}}W)^{\min \{s,q+1\}})$$\end{document}$ by ([26](#Equ26){ref-type=""}). $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\square $$\end{document}$

Appendix 3: Proof of Lemma [4](#FPar10){ref-type="sec"} {#Sec22}
=======================================================

We start by showing, using finite recursion on $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$k=q+2,q+1,\ldots , 1$$\end{document}$, that the intermediate schemes $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\beta _k,Y}^{q},\mathcal {R}_{\beta _k,Y}^{q})_{Y\subset Z}$$\end{document}$ all satisfy$$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} \forall \kappa \in \mathbb {N}^n, \quad b^q_{k,Y,\kappa }&= \left\{ \begin{array}{ll} \mathcal O\left( ({{\mathrm{rad\,}}}Y)^{\min \{|\kappa |+k,r,s\}}\right) &{}\quad \text {if }\,0\le |\kappa |\le q-k\\ 0 &{}\quad \text {otherwise}\end{array}\right. \end{aligned}$$\end{document}$$ $$\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\begin{aligned} {{\mathrm{rad\,}}}\mathcal {R}_{\beta _k,Y}^{q}&= \mathcal O\left( ({{\mathrm{rad\,}}}Y)^{\min \{r,q+1\}}\right) , \end{aligned}$$\end{document}$$where $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$b^q_{k,Y,\kappa }$$\end{document}$ denote the coefficients of $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{\beta _k,Y}^{q}$$\end{document}$.This property is trivially satisfied for the schemes $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\beta _k,Y}^{q+1},\mathcal {R}_{\beta _k,Y}^{q+1})_{Y\subset Z}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\beta _k,Y}^{q+2},\mathcal {R}_{\beta _k,Y}^{q+2})_{Y\subset Z}$$\end{document}$.Next, assume that the property holds for both $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\beta _k,Y}^{k+1},\mathcal {R}_{\beta _k,Y}^{k+1})_{Y\subset Z}$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\beta _k,Y}^{k+2},\mathcal {R}_{\beta _k,Y}^{k+2})_{Y\subset Z}$$\end{document}$, and consider the resulting scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\beta _k,Y}^{k},\mathcal {R}_{\beta _k,Y}^{k})_{Y\subset Z}$$\end{document}$ through ([21](#Equ21){ref-type=""}). Since $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[\mathcal {P}_{\varPhi ,Y}^{q}]=\mathcal O(1)$$\end{document}$ and $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$[\mathcal {P}_{\beta _k,Y}^{k+1}]=\mathcal O(({{\mathrm{rad\,}}}Y)^{\min \{k+1,r,s\}})$$\end{document}$, the remainder of the product scheme $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$(\mathcal {P}_{\varPhi ,Y}^{q},\mathcal {R}_{\varPhi ,Y}^{q}) (\mathcal {P}_{\beta _{k+1},Y}^{q},\mathcal {R}_{\beta _{k+1},Y}^{q})$$\end{document}$ has order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal O(({{\mathrm{rad\,}}}Y)^{\min \{r,q+1\}})$$\end{document}$, and so $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {R}_{\beta _k,Y}^{k}$$\end{document}$ satisfies condition ([43](#Equ43){ref-type=""}). Concerning the polynomial part, the effect of multiplying $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{\beta _{k+1},Y}^{q}$$\end{document}$ with $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{\varPhi ,Y}^{q}$$\end{document}$ is a reduction of the order of all the coefficients $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$b^q_{k+1,Y,\kappa }$$\end{document}$ by one order. Moreover, subtracting $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal {P}_{\beta _{k+2},Y}^{q}$$\end{document}$ does not change the order of the coefficients and by Lemma [3](#FPar8){ref-type="sec"} the lead coefficient $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\varphi ^q_{B^q_{Y},k}$$\end{document}$ has order $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\mathcal O(({{\mathrm{rad\,}}}Y)^{\min \{k,s\}})$$\end{document}$. Therefore, the coefficients $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$b^q_{k,Y,\kappa }$$\end{document}$ satisfy condition ([42](#Equ42){ref-type=""}), too.The same reasoning applies in taking the recurrence one step further through ([20](#Equ20){ref-type=""}), thereby giving the result. $\documentclass[12pt]{minimal}
                \usepackage{amsmath}
                \usepackage{wasysym} 
                \usepackage{amsfonts} 
                \usepackage{amssymb} 
                \usepackage{amsbsy}
                \usepackage{mathrsfs}
                \usepackage{upgreek}
                \setlength{\oddsidemargin}{-69pt}
                \begin{document}$$\square $$\end{document}$

Financial support from Marie Curie Career Integration Grant PCIG09-GA-2011-293953 and the Engineering and Physical Sciences Research Council (EPSRC) under Grant EP/J006572/1 is gratefully acknowledged. JR thanks EPSRC and the Department of Chemical Engineering at Imperial College London for doctoral training awards. MEV thanks CONACYT for doctoral scholarship.
