# On Exponential Convergence of Gegenbauer Interpolation and Spectral
Differentiation

Ziqing Xie1, Li-Lian Wang2 and Xiaodan Zhao2

###### Abstract.

This paper is devoted to a rigorous analysis of exponential convergence of
polynomial interpolation and spectral differentiation based on the Gegenbauer-
Gauss and Gegenbauer-Gauss-Lobatto points, when the underlying function is
analytic on and within an ellipse. Sharp error estimates in the maximum norm
are derived.

###### Key words and phrases:

Bernstein ellipse, exponential accuracy, Gegenbauer polynomials, Gegenbauer
Gauss-type interpolation and quadrature, spectral differentiation, maximum
error estimates

###### 1991 Mathematics Subject Classification:

65N35, 65E05, 65M70, 41A05, 41A10, 41A25

1 School of Mathematics and Computer Science, Guizhou Normal University,
Guiyang, Guizhou 550001, China; College of Mathematics and Computer Science,
Hunan Normal University, Changsha, Hunan 410081, China. The research of the
author is partially supported by the NSFC (10871066) and the Science and
Technology Grant of Guizhou Province (LKS[2010]05).
2 Division of Mathematical Sciences, School of Physical and Mathematical
Sciences, Nanyang Technological University, 637371, Singapore. The research of
the authors is partially supported by Singapore AcRF Tier 1 Grant RG58/08,
Singapore MOE Grant T207B2202 and Singapore NRF2007IDM-IDM002-010.

##  1\. Introduction

Perhaps the most significant advantage of the spectral method is its high-
order of accuracy. The typical convergence rate of the spectral method is
$O(n^{-m})$ for every $m,$ provided that the underlying function is
sufficiently smooth [22, 4, 6, 28, 8]. If the function is suitably analytic,
the expected rate is $O(q^{n})$ with $0<q<1.$ This is the so-called
exponential convergence, which is well accepted among the community. There has
been much investigation on exponential decay of spectral expansions of
analytic functions. For instance, the justification for Fourier and/or
Chebyshev series can be found in [38, 10, 36, 5, 33, 41]. In the seminal work
of Gottlieb and Shu [24, 23], on the resolution of the Gibbs phenomenon (the
interested readers are referred to Gustafsson [29] for a review of this
significant contribution), the exponential convergence, in the maximum norm
(termed as the so-called regularization error), of Gegenbauer polynomial
expansions was derived, when the index (denoted by $\lambda$ below) grows
linearly with the degree $n$. Boyd [7] provided an insightful study of the
“diagonal limit” (i.e., $\lambda=\beta n$ for some constant $\beta>0$)
convergence of the Gegenbauer reconstruction algorithm in [24]. We remark that
under the assumption of analyticity in [24], the exponential accuracy of
Gegenbauer expansions is actually valid for fixed $\lambda$ (see Appendix A
for the justification).

It is known that the heart of a collocation/pseudospectral method is the
spectral differentiation process. That is, given a set of Gauss-type points
$\\{x_{j}\\}_{j=0}^{n},$ e.g., on $[-1,1],$ the derivative values
$\\{u^{\prime}(x_{j})\\}$ can be approximated by an exact differentiation of
the polynomial interpolant $\\{(I_{n}u)^{\prime}(x_{j})\\}.$ Such a direct
differentiation technique is also called a differencing method in the
literature (see, e.g., [39, 19, 35]). For the first time, Tadmor [39] showed
the exponential accuracy of differencing analytic functions on Chebyshev-Gauss
points, where the main argument was based on analyzing the (continuous)
Chebyshev coefficients and the aliasing error, and where the intimate relation
between Fourier and Chebyshev basis functions played an essential role in the
analysis. Reddy and Weideman [35] took a different approach and improved the
estimate in [39] for Chebyshev differencing of functions analytic on and
within an ellipse with foci $\pm 1.$ As pointed out in [35], although the
exponential convergence of spectral differentiation of analytic functions is
appreciated and mentioned in the literature (see, e.g., [19, 40, 6]), the
rigorous proofs (merely for the Fourier and Chebyshev methods) can only be
found in [39, 35]. Indeed, to the best of our knowledge, the theoretical
justification even for the Legendre method is lacking. It is worthwhile to
point out that under the regularity condition (M):
$\|u^{(k)}\|_{L^{\infty}}\leq cM^{k},$ the super-geometric convergence of
Legendre spectral differentiation was proved by Zhang [46]:

$\max_{0\leq j\leq n}|(u-I_{n}u)^{\prime}(x_{j})|\leq
C\Big{(}\frac{eM}{2(n+1)}\Big{)}^{n+2},$ (1.1)

where $\\{x_{j}\\}_{j=0}^{n}$ are the Legendre-Gauss points. Similar estimate
was nontrivially extended to the Chebyshev collocation method in [47]. The
condition (M) covers a large class of functions, but it is even more
restrictive than analyticity. On the other hand, the regularity index $k$
could be infinite, while the dependence of (1.1) on $k$ is not clear.

The main concern of this paper is to show exponential convergence, in the
maximum norm, of Gegenbauer interpolation and spectral differentiation on
Gegenbauer-Gauss and Gegenbauer-Gauss-Lobatto points, provided that the
underlying function is analytic on and within a sizable ellipse. The essential
argument is based on the classical Hermite’s contour integral (see (2.19)
below), and a delicate estimate of the asymptotics of the Gegenbauer
polynomial on the ellipse of interest. It is important to remark that the
Chebyshev polynomial on the ellipse takes a very simple explicit form (see,
e.g., [10] or (3.5) below), but the Gegenbauer polynomial has a complicated
expression. Accordingly, compared with the Chebyshev case in [35], the
analysis in this paper is much more involved. The Chebyshev and Legendre
methods are commonly used in spectral approximations, but we have also
witnessed renewed applications of the Gegenbauer (or more general Jacobi)
polynomial based methods in, e.g., defeating Gibbs phenomenon (see, e.g., [24,
21]), $hp$-elements (see, e.g., [17, 2, 31, 25]), and numerical solutions of
differential equations (see, e.g., [3, 26, 27, 12, 13, 14, 15, 16, 44, 37])
and integral equations (see, e.g., [9]). The results in this paper might be
useful for a better understanding of the methods and have implications in
other applications.

The rest of the paper is organized as follows. As some preliminaries, we
briefly review basic properties of Gegenbauer polynomials, Gamma functions and
analytic functions in Section 2. We study the asymptotics of the Gegenbauer
polynomials in Section 3, and present the main results on exponential
convergence of interpolation and spectral differentiation, together with some
numerical results and extensions in the last section.

##  2\. Preliminaries

In this section, we collect some relevant properties of Gegenbauer polynomials
and assorted facts to be used throughout the paper.

###  2.1. Gegenbauer polynomials

The analysis heavily relies on the normalization of [38], so we define the
Gegenbauer polynomials111Historically, they are sometimes called
“ultraspherical polynomials” (see, e.g., the footnote on Page 80 of [38] and
Page 302 of [1]). by the three-term recurrence:

$\begin{split}&nC_{n}^{\lambda}(x)=2\left(n+\lambda-1\right)xC_{n-1}^{\lambda}(x)-(n+2\lambda-2)C_{n-2}^{\lambda}(x),\quad
n\geq 2,\\\ &C_{0}^{\lambda}(x)=1,\quad C_{1}^{\lambda}(x)=2\lambda
x,\quad\lambda>-1/2,\;\;x\in[-1,1].\end{split}$ (2.1)

Notice that if $\lambda=0,$ $C_{n}^{\lambda}(x)$ vanishes identically for
$n\geq 1.$ This corresponds to the Chebyshev polynomial, and there holds

$\lim_{\lambda\to
0}\lambda^{-1}C_{n}^{\lambda}(x)=\frac{2}{n}T_{n}(x)=\frac{2}{n}\cos(n\ {\rm
arccos}(x)),\quad n\geq 1.$ (2.2)

Hereafter, if not specified explicitly, we assume $\lambda\not=0,$ and refer
to [35] for the analysis of the Chebyshev case. Notice that for $\lambda=1/2,$
$C^{\lambda}_{n}(x)=L_{n}(x),$ i.e., the usual Legendre polynomial of degree
$n.$

The Gegenbauer polynomials are orthogonal with respect to the weight function
$(1-x^{2})^{\lambda-1/2},$ namely,

$\int_{-1}^{1}C_{n}^{\lambda}(x)C_{m}^{\lambda}(x)(1-x^{2})^{\lambda-1/2}dx=h_{n}^{\lambda}\delta_{mn},$
(2.3)

where $\delta_{mn}$ is the Kronecker symbol, and

$h_{n}^{\lambda}=\frac{2^{1-2\lambda}\pi}{\Gamma^{2}(\lambda)}\frac{\Gamma(n+2\lambda)}{n!(n+\lambda)}.$
(2.4)

Moreover, we have

$C_{n}^{\lambda}(-x)=(-1)^{n}C_{n}^{\lambda}(x),\quad
C_{n}^{\lambda}(1)=\frac{\Gamma{(n+2\lambda)}}{n!\Gamma{(2\lambda)}},$ (2.5)

and

$\frac{d}{dx}C_{n}^{\lambda}(x)=2\lambda C_{n-1}^{\lambda+1}(x).$ (2.6)

By Formula (4.7.1) and Theorems 7.32.1 and 7.33.1 of Szeg$\ddot{o}$ [38], we
have

$|C_{n}^{\lambda}(x)|\leq C_{n}^{\lambda}(1)\;\;\;{\rm
if}\;\lambda>0;\quad|C_{n}^{\lambda}(x)|\leq
D_{\lambda}\,n^{\lambda-1}\;\;\;{\rm if}\;-\frac{1}{2}<\lambda<0\;{\rm
and}\;n\gg 1,$ (2.7)

where $D_{\lambda}$ is a positive constant independent of $n.$ A tight bound
can be found in [34] (also see [32]):

$\max_{|x|\leq
1}\Big{\\{}(1-x^{2})^{\lambda}\big{(}C_{n}^{\lambda}(x)\big{)}^{2}\Big{\\}}\leq\frac{2e(2+\sqrt{2}\lambda)}{\pi}h_{n}^{\lambda},\quad\lambda>0,\;n\geq
0.$ (2.8)

###  2.2. Gamma and incomplete Gamma functions

The following properties of the Gamma and incomplete Gamma functions (cf.
[45]) are found useful. The Gamma function satisfies

$\Gamma(x)\Gamma(x+1/2)=2^{1-2x}\sqrt{\pi}\,\Gamma(2x),\quad\forall x\geq 0,$
(2.9)

and

$\Gamma{(x)}\Gamma{(-x)}=-\frac{\pi}{x\sin(\pi x)},\quad\forall x>0.$ (2.10)

We would like to quote the Stirling’s formula (see, e.g., [23]):

$\sqrt{2\pi}x^{x+1/2}e^{-x}\leq\Gamma(x+1)\leq\sqrt{2\pi}x^{x+1/2}e^{-x}e^{\frac{1}{12x}},\quad\forall
x\geq 1.$ (2.11)

We also need to use the incomplete Gamma function defined by

$\Gamma(\alpha,x)=\int_{x}^{\infty}t^{\alpha-1}e^{-t}dt,\quad\alpha>0,\;\;x\geq
0,$ (2.12)

which satisfies (see P. 899 of [30])

$\Gamma(n+1,x)=n!e^{-x}\sum_{k=0}^{n}\frac{x^{k}}{k!},\quad n=0,1,\cdots.$
(2.13)

###  2.3. Basics of analytic functions

Suppose that $u(x)$ is analytic on $[-1,1].$ Based on the notion of analytic
continuation, there always exists a simple connected region $R$ in the complex
plane containing $[-1,1]$ into which $f(x)$ can be continued analytically. The
analyticity may be characterized by the growth of the derivatives of $u$. More
precisely, let ${\mathcal{C}}$ be a simple positively oriented closed contour
surrounding $[-1,1]$ and lying in $R.$ Then we have (see, e.g., [10]):

$\frac{|u^{(m)}(x)|}{m!}\leq\frac{\max_{z\in{\mathcal{C}}}|u(z)|L({\mathcal{C}})}{2\pi\delta^{m+1}},\quad\forall
x\in[-1,1],$ (2.14)

where $L({\mathcal{C}})$ is the length of ${\mathcal{C}},$ and $\delta$ is the
distance from ${\mathcal{C}}$ to $[-1,1]$ (can be viewed as the distance from
$[-1,1]$ to the nearest singularity of $u$ in the complex plane).
Mathematically, an appropriate contour to characterize the analyticity is the
so-called Bernstein ellipse:

${\mathcal{E}}_{\rho}:=\Big{\\{}z\in{\mathbb{C}}~{}:~{}z=\frac{1}{2}(w+w^{-1})\;\;\text{with}\;\;w=\rho
e^{{\rm i}\theta},\;\theta\in[0,2\pi]\Big{\\}},\quad\rho>1,$ (2.15)

where ${\mathbb{C}}$ is the set of all complex numbers, and ${\rm
i}=\sqrt{-1}$ is the complex unit. The ellipse ${\mathcal{E}}_{\rho}$ has the
foci at $\pm 1$ and the major and minor semi-axes are, respectively

$a=\frac{1}{2}\big{(}\rho+\rho^{-1}\big{)},\quad
b=\frac{1}{2}(\rho-\rho^{-1}),$ (2.16)

so the sum of two axes is $\rho.$ As illustrated in Figure 2.1, $\rho$ is the
radius of the circle $w=\rho e^{{\rm i}\theta}$ that is mapped to the ellipse
${\mathcal{E}}_{\rho}$ under the conformal mapping: $z=\frac{1}{2}(w+w^{-1}).$

Figure 2.1. Circle (left): $|w|=\rho=1.5$, and Bernstein ellipse (right):
${\mathcal{E}}_{\rho}$ with foci $\pm 1$ linked by the conformal mapping:
$z=\frac{1}{2}(w+w^{-1})$.

According to [35], the perimeter of ${\mathcal{E}}_{\rho}$ satisfies

$L({\mathcal{E}}_{\rho})\leq\pi\sqrt{\rho^{2}+\rho^{-2}},$ (2.17)

which overestimates the perimeter by less than 12 percent. The distance from
${\mathcal{E}}_{\rho}$ to the interval $[-1,1]$ is

$\delta_{\rho}=\frac{1}{2}(\rho+\rho^{-1})-1.$ (2.18)

We are concerned with the interpolation and spectral differentiation of
analytic functions on the Gegenbauer-Gauss-type points. Let
$\\{x_{j}:=x_{j}(\lambda,n)\\}_{j=0}^{n}$ be the Gegenbauer-Gauss points
(i.e., the zeros of $C_{n+1}^{\lambda}(x)$) or the Gegenbauer-Gauss-Lobatto
points (i.e., the zeros of $(1-x^{2})C_{n-1}^{\lambda+1}(x)$). The associated
Lagrange interpolation polynomial of $u$ is given by
$I_{n}u\in{\mathbb{P}}_{n}$ (the set of all polynomials of degree $\leq n$)
such that $(I_{n}u)(x_{j})=u(x_{j})$ for $0\leq j\leq n.$ Our starting point
is the Hermite’s contour integral (see, e.g., [11]):

$(u-I_{n}u)(x)=\frac{1}{2\pi{\rm
i}}\oint_{{\mathcal{E}}_{\rho}}\frac{Q_{n+1}(x)}{Q_{n+1}(z)}\frac{u(z)}{z-x}dz,\quad\forall
x\in[-1,1],$ (2.19)

where $Q_{n+1}(x)=C_{n+1}^{\lambda}(x)$ or $(1-x^{2})C_{n-1}^{\lambda+1}(x).$
Consequently, we have

$(u-I_{n}u)^{\prime}(x_{j})=\frac{1}{2\pi{\rm
i}}\oint_{\mathcal{E}_{\rho}}\frac{Q_{n+1}^{\prime}(x_{j})}{Q_{n+1}(z)}\frac{u(z)}{z-x_{j}}dz,\quad
0\leq j\leq n.$ (2.20)

A crucial component of the error analysis is to obtain a sharp asymptotic
estimate of $Q_{n+1}(z)$ on ${\mathcal{E}}_{\rho}$ with large $n.$ This will
be the main concern of the forthcoming section.

##  3\. Asymptotic estimate of Gegenbauer polynomials on
${\mathcal{E}}_{\rho}$

Much of our analysis relies on the following representation of the Gegenbauer
polynomial.

######  Lemma 3.1.

Let $z=\frac{1}{2}(w+w^{-1}).$ We have

$C_{n}^{\lambda}(z)=\sum_{k=0}^{n}g_{k}^{\lambda}g_{n-k}^{\lambda}w^{n-2k},\quad
n\geq 0,\;\;\lambda>-1/2,$ (3.1)

where

$g_{0}^{\lambda}=1,\quad g_{k}^{\lambda}={{k+\lambda-1}\choose
k}=\frac{\Gamma(k+\lambda)}{k!\Gamma(\lambda)},\quad 1\leq k\leq n.$ (3.2)

This formula is derived from the three-term recurrence formula (2.1) and the
mathematical induction. Its proof is provided in Appendix B.

######  Remark 3.1.

Some consequences of Lemma 3.1 are in order.

  * (a)

Comparing the coefficients $w^{n}$ on both sides of (3.1), we find that the
leading coefficient of $C_{n}^{\lambda}$ is $2^{n}g_{n}^{\lambda}.$ This can
be also verified from (2.1) by the mathematical induction.

  * (b)

If $\lambda>0,$ then $g_{k}^{\lambda}>0$ for all $0\leq k\leq n.$ On the other
hand, if $\lambda<0,$ we find from (2.10) and (3.2) that

$g_{k}^{\lambda}=\frac{\sin(\pi\lambda)}{\pi}\frac{\Gamma{(k+\lambda)}\Gamma{(1-\lambda)}}{k!}<0,\quad
1\leq k\leq n.$ (3.3)

  * (c)

If $\lambda=1/2,$ it follows from (2.9) and (3.2) that

$g_{k}^{\lambda}=\frac{(2k)!}{(k!)^{2}2^{2k}},\quad 0\leq k\leq n.$ (3.4)

Such a representation for the Legendre polynomial can be found in, e.g., [11]
and [38], but the derivation is quite different.

  * (d)

If $\lambda=0,$ then by (2.2),

$T_{n}(z)=\frac{1}{2}\big{(}w^{n}+w^{-n}\big{)},\quad n\geq 1.$ (3.5)

  * (e)

If $\lambda=1,$ then $g_{k}^{\lambda}\equiv 1$ for $0\leq k\leq n.$ Therefore,
the Chebyshev polynomial of the second kind has the representation

$U_{n}(z)=\frac{w^{n+1}-w^{-(n+1)}}{w-w^{-1}}=w^{n}\sum_{k=0}^{n}w^{-2k}={C_{n}^{1}(z)},\quad
n\geq 0,$ (3.6)

which can be found in [33]. ∎

It is interesting to observe from (3.6) that for $\lambda=1,$
$C_{n}^{\lambda}(z)/w^{n}$ converges to $(1-w^{-2})^{-\lambda}$ uniformly for
all $|w|>1,$ that is,

$\sum_{k=0}^{\infty}w^{-2k}=\frac{1}{1-w^{-2}},\quad|w|>1.$

In what follows, we show a similar property holds for general $\lambda>-1/2$
and $\lambda\not=0.$ More precisely, we estimate the upper bound of remainder:

$\begin{split}\Big{|}\big{(}1-w^{-2}\big{)}^{-\lambda}-\frac{C_{n}^{\lambda}(z)}{g_{n}^{\lambda}w^{n}}\Big{|}&\leq\sum_{k=1}^{n}|d_{n,k}^{\lambda}||g^{\lambda}_{k}|\rho^{-2k}+\sum_{k=n+1}^{\infty}|g_{k}^{\lambda}|\rho^{-2k}:=R_{n}(\rho,\lambda),\end{split}$
(3.7)

where $z\in{\mathcal{E}}_{\rho}$ with $|w|=\rho>1,$ and

$d_{n,k}^{\lambda}=1-\frac{g_{n-k}^{\lambda}}{g_{n}^{\lambda}},\;\;\;1\leq
k\leq n.$ (3.8)

To obtain a sharp estimate of $R_{n}(\rho,\lambda),$ it is necessary to
understand the behavior of the coefficients $\\{g_{k}^{\lambda}\\}_{k=1}^{n}$
and $\\{d_{n,k}^{\lambda}\\}_{k=1}^{n},$ which are summarized in the following
two lemmas.

######  Lemma 3.2.

For $\lambda>-1/2,\,k\geq 1$ and $k+\lambda\geq 1,$

$c_{1}\,\Big{(}1+\frac{\lambda}{k}\Big{)}^{k+1/2}e^{-\lambda}\leq\frac{\Gamma(\lambda)\,g_{k}^{\lambda}}{(k+\lambda)^{\lambda-1}}\leq
c_{2}\,\Big{(}1+\frac{\lambda}{k}\Big{)}^{k+1/2}e^{-\lambda},$ (3.9)

where $c_{1}=e^{-\frac{1}{12k}}$ and $c_{2}=e^{\frac{1}{12(k+\lambda)}}.$

###### Proof.

Applying the Stirling’s formula (2.11) to

$(k+\lambda)\Gamma(\lambda)\,g_{k}^{\lambda}=\frac{\Gamma(k+\lambda+1)}{\Gamma(k+1)}$

leads to (3.9). ∎

######  Lemma 3.3.

Let $\\{g_{k}^{\lambda}\\}_{k=1}^{n}$ and $\\{d_{n,k}^{\lambda}\\}_{k=1}^{n}$
be the sequences as defined in (3.2) and (3.8), respectively.

  * (i)

If $\lambda>1$, then there holds

$0<d_{n,1}^{\lambda}<d_{n,2}^{\lambda}<\cdots<d_{n,n}^{\lambda}<1.$ (3.10)

  * (ii)

If $-1/2<\lambda<1$ and $\lambda\not=0,$ then

$\cdots<|g_{k+1}^{\lambda}|<|g_{k}^{\lambda}|<\cdots<|g_{1}^{\lambda}|<g_{0}^{\lambda}=1,$
(3.11)

and we have

$0<-d_{n,1}^{\lambda}<-d_{n,2}^{\lambda}<\cdots<-d_{n,n-1}^{\lambda},$ (3.12)

and

$|d_{n,k}^{\lambda}g_{k}^{\lambda}|<1\;\;\;{\rm for}\;\;\;1\leq k\leq
n-1,\;n\geq 3.$ (3.13)

###### Proof.

By (3.2),

$\frac{g_{k+1}^{\lambda}}{g_{k}^{\lambda}}=\frac{k+\lambda}{k+1}\;\;\;{\rm
for}\;\;\lambda\not=0.$ (3.14)

Thus for $\lambda>1,$ $\\{g_{k}^{\lambda}\\}$ is strictly increasing with
respect to $k,$ which, together with the fact $g_{k}^{\lambda}>0,$ implies

$0<d_{n,k}^{\lambda}=1-\frac{g_{n-k}^{\lambda}}{g_{n}^{\lambda}}<1,$

and

$d_{n,k+1}^{\lambda}-d_{n,k}^{\lambda}=\frac{g_{n-k}^{\lambda}-g_{n-k-1}^{\lambda}}{g_{n}^{\lambda}}>0.$

This completes the proof of (i).

The property (3.11) is a direct consequence of (3.14), and (3.12) can be
proved in a fashion similar to (3.10). It remains to verify (3.13). If $k=1,$
a direct calculation by using (3.2) yields

$|d_{n,1}^{\lambda}g_{1}^{\lambda}|=|\lambda|\frac{1-\lambda}{n-(1-\lambda)}<1,\quad\forall
n\geq 3.$

For $2\leq k\leq n-1,$ it follows from (3.2) and (3.12) that

$\begin{split}|d_{n,k}^{\lambda}g_{k}^{\lambda}|&=\Big{(}\frac{g_{n-k}^{\lambda}}{g_{n}^{\lambda}}-1\Big{)}|g_{k}^{\lambda}|<\frac{g_{n-k}^{\lambda}}{g_{n}^{\lambda}}|g_{k}^{\lambda}|=\left(\prod_{j=0}^{k-2}\frac{1-\frac{1-\lambda}{k-j}}{1-\frac{1-\lambda}{n-j}}\right)\frac{|\lambda|}{1-\frac{1-\lambda}{n-k+1}}<1.\end{split}$

This ends the proof. ∎

With the above preparation, we are ready to present the main result on the
upper bound of $R_{n}(\rho,\lambda)$ in (3.7).

######  Theorem 3.1.

Let $R_{n}(\rho,\lambda)$ with $\rho>1$ be the remainder as defined in (3.7).

  * (i)

If $\lambda>1,$ then for all $n\geq m\geq 1$ and

$m+2\geq(\lambda-1)\Big{(}\frac{1}{2\ln\rho}-1\Big{)},$ (3.15)

we have

$R_{n}(\rho,\lambda)\leq
d_{n,m}^{\lambda}\Big{(}(1-\rho^{-2})^{-\lambda}-1\Big{)}+A\frac{[\lambda]!}{(2\ln\rho)^{\lambda}}\frac{(m+\lambda)^{[\lambda]}}{\rho^{2(m-1)}},$
(3.16)

where $[\lambda]$ is the largest integer $\leq\lambda,$ and

$A=\frac{1}{\Gamma(\lambda)}{\rm
exp}\Big{(}{\frac{1}{12(m+1+\lambda)}+\frac{\lambda}{2(m+1)}}\Big{)}.$ (3.17)

  * (ii)

If $-1/2<\lambda<1$ and $\lambda\not=0,$ then for all $n\geq m\geq 1$ and
$n\geq 3,$

$R_{n}(\rho,\lambda)\leq|d_{n,m}^{\lambda}|\big{|}(1-\rho^{-2})^{-\lambda}-1\big{|}+\frac{\rho^{-2m}}{\rho^{2}-1}+2\rho^{-2n}.$
(3.18)

Here, the factor $d_{n,m}^{\lambda}$ is given by (3.8).

###### Proof.

(i) For $\lambda>1,$ we obtain from (3.10) that

$\begin{split}R_{n}(\rho,\lambda)&=\sum_{k=1}^{m}d_{n,k}^{\lambda}g_{k}^{\lambda}\rho^{-2k}+\sum_{k=m+1}^{n}d_{n,k}^{\lambda}g_{k}^{\lambda}\rho^{-2k}+\sum_{k=n+1}^{\infty}g_{k}^{\lambda}\rho^{-2k}\\\
&\overset{(\ref{dnkcase1})}{\leq}d_{n,m}^{\lambda}\sum_{k=1}^{m}g_{k}^{\lambda}\rho^{-2k}+\sum_{k=m+1}^{n}g_{k}^{\lambda}\rho^{-2k}+\sum_{k=n+1}^{\infty}g_{k}^{\lambda}\rho^{-2k}\\\
&\leq
d_{n,m}^{\lambda}\Big{(}(1-\rho^{-2})^{-\lambda}-1\Big{)}+\sum_{k=m+1}^{\infty}g_{k}^{\lambda}\rho^{-2k}.\end{split}$
(3.19)

By Lemma 3.2,

$\begin{split}\sum_{k=m+1}^{\infty}g_{k}^{\lambda}\rho^{-2k}&\leq{(\Gamma(\lambda))^{-1}}{e^{-\lambda}e^{\frac{1}{12(m+1+\lambda)}}}\sum_{k=m+1}^{\infty}\frac{(k+\lambda)^{\lambda-1}}{\rho^{2k}}\Big{(}1+\frac{\lambda}{k}\Big{)}^{k+1/2}\\\
&\leq{(\Gamma(\lambda))^{-1}}e^{-\lambda}e^{\frac{1}{12(m+1+\lambda)}}\sum_{k=m+1}^{\infty}\frac{(k+\lambda)^{\lambda-1}}{\rho^{2k}}e^{\lambda+\frac{\lambda}{2k}},\end{split}$

where we used the inequality $1+x<e^{x}$ for $x>0.$ Hence,

$\sum_{k=m+1}^{\infty}g_{k}^{\lambda}\rho^{-2k}\leq
A\sum_{k=m+1}^{\infty}\frac{(k+\lambda)^{\lambda-1}}{\rho^{2k}},$ (3.20)

where $A$ is given by (3.17). One verifies that under the condition (3.15),
${(k+\lambda)^{\lambda-1}}/{\rho^{2k}}$ is decreasing with respect to $k.$
Therefore, by (2.12) and (2.13),

$\begin{split}\sum_{k=m+1}^{\infty}\frac{(k+\lambda)^{\lambda-1}}{\rho^{2k}}&\leq\int_{m}^{\infty}{(x+\lambda)^{\lambda-1}}\rho^{-2x}dx=\frac{\rho^{2\lambda}}{(2\ln\rho)^{\lambda}}\int_{2(m+\lambda)\ln\rho}^{\infty}x^{\lambda-1}e^{-x}dx\\\
&=\frac{\rho^{2\lambda}}{(2\ln\rho)^{\lambda}}\Gamma\big{(}\lambda,2(m+\lambda)\ln\rho\big{)}\leq\frac{\rho^{2\lambda}}{(2\ln\rho)^{\lambda}}\Gamma\big{(}[\lambda]+1,2(m+\lambda)\ln\rho\big{)}\\\
&=\frac{[\lambda]!\rho^{-2m}}{(2\ln\rho)^{\lambda}}\sum_{k=0}^{[\lambda]}\frac{(m+\lambda)^{k}(2\ln\rho)^{k}}{k!}\leq\frac{[\lambda]!}{(2\ln\rho)^{\lambda}}\frac{(m+\lambda)^{[\lambda]}}{\rho^{2m}}\sum_{k=0}^{\infty}\frac{(2\ln\rho)^{k}}{k!}\\\
&=\frac{[\lambda]!}{(2\ln\rho)^{\lambda}}\frac{(m+\lambda)^{[\lambda]}}{\rho^{2(m-1)}}.\end{split}$

A combination of the above estimates leads to (3.16).

(ii) Now, we turn to the proof of the second case: $-1/2<\lambda<1$ and
$\lambda\not=0.$ By Lemma 3.3,

$\begin{split}R_{n}(\rho,\lambda)&=\sum_{k=1}^{m}|d_{n,k}^{\lambda}||g_{k}^{\lambda}|\rho^{-2k}+\sum_{k=m+1}^{n}|d_{n,k}^{\lambda}||g_{k}^{\lambda}|\rho^{-2k}+\sum_{k=n+1}^{\infty}|g_{k}^{\lambda}|\rho^{-2k}\\\
&\overset{(\ref{dnkcase2})}{\leq}|d_{n,m}^{\lambda}|\sum_{k=1}^{m}|g_{k}^{\lambda}|\rho^{-2k}\overset{(\ref{dnkcase3})}{+}\sum_{k=m+1}^{n-1}\rho^{-2k}+|d_{n,n}^{\lambda}g_{n}^{\lambda}|\rho^{-2n}\overset{(\ref{dnkcase2aa})}{+}\sum_{k=n+1}^{\infty}\rho^{-2k}\\\
&\quad\leq|d_{n,m}^{\lambda}|\big{|}(1-\rho^{-2})^{-\lambda}-1\big{|}+\frac{\rho^{-2m}}{\rho^{2}-1}+2\rho^{-2n},\end{split}$

where in the last step, we used the following facts:

$\begin{split}\sum_{k=1}^{m}|g_{k}^{\lambda}|\rho^{-2k}&\leq{\rm
sign}(\lambda)\sum_{k=1}^{\infty}g_{k}^{\lambda}\rho^{-2k}={\rm
sign}(\lambda)\big{(}(1-\rho^{-2})^{-\lambda}-1\big{)},\quad\rho>1,\end{split}$

(note: ${\rm sign}(\lambda)$ is the sign of $\lambda$), and
$|d_{n,n}^{\lambda}g_{n}^{\lambda}|=|g_{n}^{\lambda}-1|<2,$ thanks to (3.11).
∎

The estimate in Theorem 3.1 is quite tight and is valid even for small $n.$ By
choosing a suitable $m$ to balance the two error terms in the upper bound, we
are able to derive the anticipated asymptotic estimate.

######  Theorem 3.2.

For any $z\in{\mathcal{E}}_{\rho}$ with $|w|=\rho>1,$ and any $\lambda>-1/2$
and $\lambda\not=0,$ there exists $0<\varepsilon\leq 1/2$ such that

$\Big{|}\big{(}1-w^{-2}\big{)}^{-\lambda}-\frac{C_{n}^{\lambda}(z)}{g_{n}^{\lambda}w^{n}}\Big{|}\leq
A(\rho,\lambda)n^{\varepsilon-1}+O(n^{-1}),$ (3.21)

where

$A(\rho,\lambda)=|1-\lambda|\big{|}(1-\rho^{-2})^{-\lambda}-1\big{|}.$ (3.22)

###### Proof.

We first estimate $|d_{n,m}^{\lambda}|$ in Theorem 3.1, when $n-m$ is large.
Using the Stirling’s formula (2.11) and (3.2) that

$\begin{split}\frac{g_{n-m}^{\lambda}}{g_{n}^{\lambda}}&=\Big{(}1+\frac{1-\lambda}{n+\lambda-1}\Big{)}^{n+\frac{1}{2}}\Big{(}1-\frac{1-\lambda}{n-m}\Big{)}^{n-m+\frac{1}{2}}\Big{(}1-\frac{m}{n+\lambda-1}\Big{)}^{\lambda-1}\Big{\\{}1+O\Big{(}\frac{1}{n-m}\Big{)}\Big{\\}}\\\
&=\Big{(}1-\frac{m}{n+\lambda-1}\Big{)}^{\lambda-1}\Big{\\{}1+O\Big{(}\frac{1}{n-m}\Big{)}\Big{\\}}\\\
&=\Big{\\{}1+\frac{(1-\lambda)m}{n+\lambda-1}+O\Big{(}\frac{m^{2}}{n^{2}}\Big{)}\Big{\\}}\Big{\\{}1+O\Big{(}\frac{1}{n-m}\Big{)}\Big{\\}}.\end{split}$

Hereafter, taking $m=[n^{\varepsilon}]$ with $0<\varepsilon\leq 1/2$ yields

$\begin{split}\frac{g_{n-m}^{\lambda}}{g_{n}^{\lambda}}&=1+(1-\lambda)n^{\varepsilon-1}+O\Big{(}\frac{1}{n-n^{\varepsilon}}\Big{)}\\\
&\Longrightarrow\;\;d_{n,m}^{\lambda}=(\lambda-1)n^{\varepsilon-1}+O(n^{-1}).\end{split}$
(3.23)

One verifies readily that for $\lambda>1$ and any $0<\varepsilon\leq 1/2,$

$\frac{m^{[\lambda]}}{\rho^{2m}}\leq\frac{1}{n}\;\;\Longleftrightarrow\;\;\frac{\ln
n}{n^{\varepsilon}}\leq\frac{2\ln\rho}{1+\varepsilon[\lambda]},$ (3.24)

which, together with (3.16) and (3.23), implies (3.21) with $\lambda>1.$

If $-1/2<\lambda<1$ and $\lambda\not=0,$ it follows from (3.24) that
$\rho^{-2m}\leq n^{-1}$ for any $0<\varepsilon\leq 1/2.$ This validates the
desired estimate. ∎

A direct consequence of Theorem 3.2 is that

$\lim_{n\to\infty}\frac{C_{n}^{\lambda}(z)}{g_{n}^{\lambda}}=\lim_{n\rightarrow\infty}\sum_{k=0}^{n}\frac{g_{n-k}^{\lambda}}{g_{n}^{\lambda}}g_{k}^{\lambda}w^{-2k}=(1-w^{-2})^{-\lambda},$
(3.25)

for all $z\in{\mathcal{E}}_{\rho}$ with $|w|=\rho>1,$ and any $\lambda>-1/2$
and $\lambda\not=0.$

######  Remark 3.2.

Based on a completely different argument, Elliott [18] derived an asymptotic
expansion for large $n$ near $z=1$ (but not near $z=-1$):
$C_{n}^{\lambda}(z)\sim\frac{B(n,\lambda)}{(z^{2}-1)^{\lambda/2}},$ where $B$
is a series involving modified Bessel functions, and some other asymptotic
expansions for $|z|$ large and $n$ fixed. Although they are valid for general
$z$ off the interval $[-1,1],$ our results in Theorems 3.1 and 3.2 provide
tighter and sharper bounds when $z$ is sitting on ${\mathcal{E}}_{\rho}.$ ∎

As the end of this section, we provide some numerical results to illustrate
the tightness of the upper bound in (3.21). Denote by

$E_{n}(\rho;\lambda):=\frac{1}{A(\rho,\lambda)}\max_{z\in{\mathcal{E}}_{\rho}}\Big{|}\big{(}1-w^{-2}\big{)}^{-\lambda}-\frac{C_{n}^{\lambda}(z)}{g_{n}^{\lambda}w^{n}}\Big{|}.$
(3.26)

To approximate the maximum value, we sample a set of points dense on the
ellipse ${\mathcal{E}}_{\rho}$ based on the conformal mapping
$z=\frac{1}{2}(w+w^{-1})$ of the Fourier points on the circle $w=\rho e^{{\rm
i}\theta}.$ We plot in Figure 3.1 in (Matlab) log-log scale of
$E_{n}(\rho;\lambda),n^{-1}$ and $n^{\varepsilon-1}$ (with $\varepsilon=0.1$)
for several sets of parameters $\lambda$ and $\rho,$ and for large $n.$
According to Theorem 3.2, $E_{n}$ should be bounded by $n^{\varepsilon-1}$
from above, and it is anticipated to be bounded below by $n^{-1},$ if the
estimate is tight. Indeed, we observe from Figure 3.1 such a behavior when $n$
is large.

Figure 3.1. $E_{n}$ against $n^{\varepsilon-1}$ (with $\varepsilon=0.1$) and
$n^{-1}$ for large $n$.

##  4\. Error estimates of interpolation and spectral differentiation

After collecting all the necessary results, we are ready to estimate
exponential convergence of interpolation and spectral differentiation of
analytic functions.

Hereafter, the notation $a_{n}\cong b_{n}$ means that $a_{n}/b_{n}\to 1$ as
$n\to\infty,$ for any two sequences $\\{a_{n}\\}$ and $\\{b_{n}\\}$ (with
$b_{n}\not=0$) of complex numbers.

###  4.1. Gegenbauer-Gauss interpolation and differentiation

We start with the analysis of interpolation and spectral differentiation on
zeros of the Gegenbauer polynomial $C_{n+1}^{\lambda}(x).$

######  Theorem 4.1.

Let $u$ be analytic on and within the ellipse ${\mathcal{E}_{\rho}}$ with foci
$\pm 1$ and $\rho>1$ as defined in (2.15), and let $(I_{n}u)(x)$ be the
interpolant of $u(x)$ at the set of $(n+1)$ Gegenbauer-Gauss points.

  * (i)

If $\lambda>0,$ we have

$\max_{|x|\leq
1}\big{|}(u-I_{n}u)(x)\big{|}\leq\frac{c\Gamma(\lambda)M_{\rho}\sqrt{\rho^{2}+\rho^{-2}}}{\Gamma(2\lambda)(\rho-1)^{2}(1+\rho^{-2})^{-\lambda}}\frac{n^{\lambda}}{\rho^{n}}.$
(4.1)

  * (ii)

If $-1/2<\lambda<0,$ we have

$\max_{|x|\leq
1}\big{|}(u-I_{n}u)(x)\big{|}\leq\frac{cD_{\lambda}|\Gamma(\lambda)|M_{\rho}\sqrt{\rho^{2}+\rho^{-2}}}{(\rho-1)^{2}(1-\rho^{-2})^{-\lambda}}\frac{1}{\rho^{n}}.$
(4.2)

Here, $M_{\rho}=\max_{z\in{\mathcal{E}_{\rho}}}|u(z)|$, $D_{\lambda}$ is
defined in (2.7), and $c\cong 1$ is a generic positive constant.

###### Proof.

By the formula (2.19) with $Q_{n+1}=C_{n+1}^{\lambda}$ and (2.17)-(2.18), we
have the bound of the point-wise error:

$\begin{split}\big{|}(u-I_{n}u)(x)\big{|}&\leq\frac{|C_{n+1}^{\lambda}(x)|}{2\pi}\frac{\max_{z\in\mathcal{E}_{\rho}}|u(z)|}{\min_{z\in{\mathcal{E}_{\rho}}}|C_{n+1}^{\lambda}(z)|}\oint_{\mathcal{E}_{\rho}}\frac{|dz|}{|z-x|}\\\
&\leq\frac{M_{\rho}L({\mathcal{E}_{\rho}})}{2\pi\delta_{\rho}}\frac{|C_{n+1}^{\lambda}(x)|}{\min_{z\in{\mathcal{E}_{\rho}}}|C_{n+1}^{\lambda}(z)|}\\\
&\leq\frac{M_{\rho}\sqrt{\rho^{2}+\rho^{-2}}}{\rho+\rho^{-1}-2}\frac{|C_{n+1}^{\lambda}(x)|}{\min_{z\in{\mathcal{E}_{\rho}}}|C_{n+1}^{\lambda}(z)|},\quad
x\in[-1,1],\;n\geq 0.\end{split}$ (4.3)

Therefore, it is essential to obtain the lower bound of
$|C_{n+1}^{\lambda}(z)|.$ Recall that for any two complex numbers $z_{1}$ and
$z_{2},$ we have $\big{|}|z_{1}|-|z_{2}|\big{|}\leq|z_{1}-z_{2}|.$ It follows
from Theorem 3.2 that

$\bigg{|}|1-w^{-2}|^{-\lambda}-\frac{|C_{n+1}^{\lambda}(z)|}{|g_{n+1}^{\lambda}|\rho^{n+1}}\bigg{|}\leq
A(\rho,\lambda)n^{\varepsilon-1}+O(n^{-1}),$

which implies

$\begin{split}\big{|}1-w^{-2}\big{|}^{-\lambda}&-A(\rho,\lambda)n^{\varepsilon-1}-O(n^{-1})\leq\frac{|C_{n+1}^{\lambda}(z)|}{|g_{n+1}^{\lambda}|\rho^{n+1}}\\\
&\leq\big{|}1-w^{-2}\big{|}^{-\lambda}+A(\rho,\lambda)n^{\varepsilon-1}+O(n^{-1}).\end{split}$
(4.4)

Notice that

$1-\rho^{-2}\leq|1-w^{-2}|\leq 1+\rho^{-2}.$ (4.5)

Consequently,

$|C_{n+1}^{\lambda}(z)|\geq
c\frac{n^{\lambda-1}\rho^{n+1}}{|\Gamma(\lambda)|}\begin{cases}(1+\rho^{-2})^{-\lambda},\quad&{\rm
if}\;\;\lambda>0,\\\ (1-\rho^{-2})^{-\lambda},\quad&{\rm
if}\;\;\lambda<0,\end{cases}$ (4.6)

where we used (3.9), and the constant $c\cong 1.$

On the other hand, we derive from (2.5), (2.7) and (2.11) that if $\lambda>0,$

$\max_{|x|\leq
1}|C_{n+1}^{\lambda}(x)|=C_{n+1}^{\lambda}(1)\cong\frac{n^{2\lambda-1}}{\Gamma(2\lambda)}.$
(4.7)

Hence, a combination of (4.3), (4.6) and (4.7) leads to (4.1). Similarly, for
$-1/2<\lambda<0,$ we use (2.7) to derive (4.2). ∎

######  Remark 4.1.

For $\lambda>0,$ we obtain from (2.4), (2.8) and (2.11) that

$|C_{n+1}^{\lambda}(x)|\leq\frac{c2^{1-\lambda}\sqrt{e(2+\sqrt{2}\lambda)}}{\Gamma(\lambda)}n^{\lambda-1}(1-x^{2})^{-\lambda/2},\;\;|x|<1.$
(4.8)

Replacing (4.7) by this bound in the above proof, we can derive the point-wise
estimate for $\lambda>0:$

$\big{|}(u-I_{n}u)(x)\big{|}\leq
D(\rho,\lambda)\frac{(1-x^{2})^{-\lambda/2}}{\rho^{n}},\quad|x|<1,$ (4.9)

where the positive constant $D(\rho,\lambda)$ can be worked out as well. It
appears to be sharper than (4.1) at the points which are not too close to the
endpoints $x=\pm 1.$ A similar remark also applies to the Gegenbauer-Gauss-
Lobatto interpolation to be addressed in a minute. ∎

Now, we turn to the estimate of spectral differentiation.

######  Theorem 4.2.

Let $u$ be analytic on and within the ellipse ${\mathcal{E}_{\rho}}$ with foci
$\pm 1$ and $\rho>1$ as defined in (2.15), and let $(I_{n}u)(x)$ be the
interpolant of $u(x)$ at $(n+1)$ Gegenbauer-Gauss points
$\\{x_{j}\\}_{j=0}^{n}$. Then we have

$\max_{0\leq j\leq
n}\big{|}(u-I_{n}u)^{\prime}(x_{j})\big{|}\leq\Lambda(\rho,\lambda)\frac{n^{\lambda+2}}{\rho^{n}},$
(4.10)

where the constant

$\Lambda(\rho,\lambda)=\frac{2c\Gamma(\lambda+1)M_{\rho}\sqrt{\rho^{2}+\rho^{-2}}}{\Gamma(2\lambda+2)(\rho-1)^{2}}\begin{cases}(1+\rho^{-2})^{\lambda},\quad&{\rm
if}\;\;\lambda>0,\\\ (1-\rho^{-2})^{\lambda},\quad&{\rm
if}\;\;\lambda<0,\end{cases}$ (4.11)

and $c,M_{\rho}$ are the same as in Theorem 4.1.

###### Proof.

In view of (2.19) and (2.20), it is enough to replace $x$ and
$C_{n+1}^{\lambda}(x)$ by $x_{j}$ and $\frac{d}{dx}C_{n+1}^{\lambda}(x),$
respectively, in (4.3). Thus, we have

$\begin{split}\big{|}(u-I_{n}u)^{\prime}(x_{j})\big{|}&\leq\frac{M_{\rho}}{2\pi}\frac{|(C_{n+1}^{\lambda})^{\prime}(x_{j})|}{\min_{z\in\mathcal{E}_{\rho}}|C_{n+1}^{\lambda}(z)|}\oint_{\mathcal{E}_{\rho}}\frac{|dz|}{|z-x_{j}|}\\\
&\leq\frac{M_{\rho}\sqrt{\rho^{2}+\rho^{-2}}}{\rho+\rho^{-1}-2}\frac{|(C_{n+1}^{\lambda})^{\prime}(x_{j})|}{\min_{z\in{\mathcal{E}_{\rho}}}|C_{n+1}^{\lambda}(z)|},\quad
0\leq j\leq n.\end{split}$ (4.12)

By (2.6) and (4.7),

$\max_{|x|\leq
1}\big{|}(C_{n+1}^{\lambda})^{\prime}(x)\big{|}=2|\lambda||C_{n}^{\lambda+1}(1)|\cong\frac{2|\lambda|}{\Gamma(2\lambda+2)}n^{2\lambda+1},$
(4.13)

which, together with (4.6) and (4.12), leads to the desired estimate. ∎

######  Remark 4.2.

Obviously, by (2.19),

$(u-I_{n}u)^{\prime}(x)=\frac{1}{2\pi{\rm
i}}\oint_{{\mathcal{E}}_{\rho}}\Big{(}\frac{(C_{n+1}^{\lambda})^{\prime}(x)}{z-x}+\frac{C_{n+1}^{\lambda}(x)}{(z-x)^{2}}\Big{)}\frac{u(z)}{C_{n+1}^{\lambda}(z)}dz.$
(4.14)

If $x\not=x_{j},$ we need to estimate the second term in the summation, which
can be done in the same fashion in the proof of Theorem 4.1. The first term is
actually estimated above. Consequently, we have

$\max_{|x|\leq
1}\big{|}(u-I_{n}u)^{\prime}(x)\big{|}\leq\Lambda(\rho,\lambda)\frac{n^{\lambda+2}}{\rho^{n}}+\frac{1}{\delta_{\rho}}\max_{|x|\leq
1}\big{|}(u-I_{n}u)(x)\big{|},$

where $\delta_{\rho}$ is given by (2.18). Hence, by Theorem 4.1,

$\max_{|x|\leq
1}\big{|}(u-I_{n}u)^{\prime}(x)\big{|}=O\Big{(}\frac{n^{\lambda+2}}{\rho^{n}}\Big{)}.$
(4.15)

In fact, the results for higher-order derivatives can be derived recursively,
and it is anticipated that

$\max_{|x|\leq
1}\big{|}(u-I_{n}u)^{(k)}(x)\big{|}=O\Big{(}\frac{n^{\lambda+2k}}{\rho^{n}}\Big{)},\quad
k\geq 1.$ (4.16)

A similar remark applies to the Gegenbauer-Gauss-Lobatto case below. ∎

###  4.2. Gegenbauer-Gauss-Lobatto interpolation and differentiation

We are now in a position to estimate the Gegenbauer-Gauss-Lobatto
interpolation and spectral differentiation. In this case,
$Q_{n+1}(x)=(1-x^{2})C_{n-1}^{\lambda+1}(x)$ in (2.19)-(2.20). The main result
is stated as follows.

######  Theorem 4.3.

Let $u$ be analytic on and within the ellipse ${\mathcal{E}_{\rho}}$ with foci
$\pm 1$ and $\rho>1$ as defined in (2.15), and let $(I_{n}u)(x)$ be the
interpolant of $u(x)$ at the set of $(n+1)$ Gegenbauer-Gauss-Lobatto points.

  * (a)

We have the interpolation error:

$\max_{|x|\leq
1}\big{|}(u-I_{n}u)(x)\big{|}\leq\frac{4cM_{\rho}\sqrt{\rho^{2}+\rho^{-2}}(1+\rho^{-2})^{\lambda+1}}{(1-\rho^{-1})^{2}(\rho-\rho^{-1})^{2}}\frac{\Gamma(\lambda+1)}{\Gamma(2\lambda+2)}\frac{n^{\lambda+1}}{\rho^{n}}.$
(4.17)

  * (b)

We have the estimate:

$\begin{split}&\max_{0\leq j\leq
n}\big{|}(u-I_{n}u)^{\prime}(x_{j})\big{|}\leq\frac{8cM_{\rho}\sqrt{\rho^{2}+\rho^{-2}}(1+\rho^{-2})^{\lambda+1}}{(1-\rho^{-1})^{2}(\rho-\rho^{-1})^{2}}\frac{\Gamma(\lambda+2)}{\Gamma(2\lambda+4)}\frac{n^{\lambda+3}}{\rho^{n}}.\end{split}$
(4.18)

Here, $c\cong 1$ and $M_{\rho}=\max_{z\in{\mathcal{E}}_{\rho}}|u(z)|.$

###### Proof.

For any $z\in{\mathcal{E}}_{\rho},$ one verifies that

$\frac{1}{4}(\rho-\rho^{-1})^{2}\leq|z^{2}-1|\leq\frac{1}{4}(\rho+\rho^{-1})^{2},$

and

$\min_{z\in{\mathcal{E}_{\rho}}}\big{|}(1-z^{2})C_{n-1}^{\lambda+1}(z)\big{|}\geq\frac{1}{4}(\rho-\rho^{-1})^{2}\min_{z\in{\mathcal{E}_{\rho}}}\big{|}C_{n-1}^{\lambda+1}(z)\big{|}.$
(4.19)

(a) By (2.19) with $Q_{n+1}(x)=(1-x^{2})C_{n-1}^{\lambda+1}(x)$ and
(2.17)-(2.18), we have the bound of the point-wise error:

$\begin{split}\big{|}(u-I_{n}u)(x)\big{|}&\leq\frac{\big{|}(1-x^{2})C_{n-1}^{\lambda+1}(x)\big{|}}{2\pi}\frac{\max_{z\in\mathcal{E}_{\rho}}|u(z)|}{\min_{z\in{\mathcal{E}_{\rho}}}\big{|}(1-z^{2})C_{n-1}^{\lambda+1}(z)\big{|}}\oint_{\mathcal{E}_{\rho}}\frac{|dz|}{|z-x|}\\\
&\overset{(\ref{minCestGL})}{\leq}\frac{M_{\rho}\sqrt{\rho^{2}+\rho^{-2}}}{\rho+\rho^{-1}-2}\frac{4(\rho-\rho^{-1})^{-2}|C_{n-1}^{\lambda+1}(x)|}{\min_{z\in{\mathcal{E}_{\rho}}}|C_{n-1}^{\lambda+1}(z)|},\quad
x\in[-1,1],\;n\geq 0.\end{split}$ (4.20)

Thus, the estimate (4.17) follows from (4.6) and (4.7).

(b) Similarly, we have

$\begin{split}\big{|}(u-I_{n}u)^{\prime}(x_{j})\big{|}&\leq\frac{M_{\rho}}{2\pi}\frac{\big{|}[(1-x^{2})C_{n-1}^{\lambda+1}(x)]^{\prime}(x_{j})\big{|}}{\min_{z\in\mathcal{E}_{\rho}}\big{|}(1-z^{2})C_{n-1}^{\lambda+1}(z)\big{|}}\oint_{\mathcal{E}_{\rho}}\frac{|dz|}{|z-x_{j}|}\\\
&\leq\frac{4M_{\rho}\sqrt{\rho^{2}+\rho^{-2}}}{\rho+\rho^{-1}-2}\frac{\big{|}[(1-x^{2})C_{n-1}^{\lambda+1}(x)]^{\prime}(x_{j})\big{|}}{(\rho-\rho^{-1})^{2}\min_{z\in{\mathcal{E}_{\rho}}}\big{|}C_{n-1}^{\lambda+1}(z)\big{|}},\quad
0\leq j\leq n.\end{split}$ (4.21)

A direct calculation leads to

$\big{[}(1-x^{2})C_{n-1}^{\lambda+1}(x)\big{]}^{\prime}=-2xC_{n-1}^{\lambda+1}(x)+(1-x^{2})\big{[}C_{n-1}^{\lambda+1}(x)\big{]}^{\prime},$
(4.22)

which, together with (4.7) and (4.13), gives

$\max\limits_{|x|\leq
1}\big{|}[(1-x^{2})C_{n-1}^{\lambda+1}(x)]^{\prime}\big{|}\cong\frac{2(\lambda+1)}{\Gamma(2\lambda+4)}n^{2\lambda+3}.$
(4.23)

A combination of (4.6), (4.21) and (4.23) yields the desired estimate. ∎

######  Remark 4.3.

Similar to Remarks 4.1 and 4.2, we can derive a sharper point-wise estimate
and analyze higher-order derivatives of interpolation errors. ∎

###  4.3. Analysis of quadrature errors

Recall the interpolatory Gegenbauer-Gauss-type quadrature formula:

$\int_{-1}^{1}u(x)(1-x^{2})^{\lambda-1/2}dx\approx\sum_{j=0}^{n}u(x_{j})\omega_{j}=\int_{-1}^{1}(I_{n}u)(x)(1-x^{2})^{\lambda-1/2}dx,$
(4.24)

where the quadrature weights $\\{\omega_{j}\\}_{j=0}^{n}$ are expressed by the
Lagrange basis polynomials (see, e.g., [38, 20]). Observe that

$\Big{|}\int_{-1}^{1}(u-I_{n}u)(x)(1-x^{2})^{\lambda-1/2}dx\Big{|}\leq
h_{0}^{\lambda}\max_{|x|\leq 1}|(u-I_{n}u)(x)|,$ (4.25)

where $h_{0}^{\lambda}$ is given by (2.4). With the aid of interpolation error
estimates in Theorem 4.1 and Theorem 4.3, we are able to derive the quadrature
errors immediately.

###  4.4. Numerical results

In what follows, we provide two numerical examples to demonstrate the
sharpness of the estimates established in Theorems 4.2 and 4.3.

####  4.4.1. Example 1

We take

$u(x)=\frac{1}{x^{2}+1},$ (4.26)

which has two simple poles at $\pm{\rm i}.$ By (2.16), we are free to choose
${\mathcal{E}}_{\rho}$ with $\rho$ in the range

$1<\rho<1+\sqrt{2}\approx 2.414,$ (4.27)

such that $u$ is analytic on and within ${\mathcal{E}}_{\rho}.$ To compare the
(discrete) maximum error of spectral differentiation with the upper bound, we
sample about $2000$ values of $\rho$ equally from $(1,1+\sqrt{2}),$ and find a
tighter upper bound (which is usually attained when $\rho$ is close to
$1+\sqrt{2}$ ). In Figure 4.1 (a)-(b), we plot the (discrete) maximum errors
of spectral differentiation against the upper bounds. We visualize the
exponential decay of the errors, and the upper bounds and the errors decay at
almost the same rate. Moreover, it seems that the bounds are slightly sharper
in the Gegenbauer-Gauss case.

####  4.4.2. Example 2

The estimates indicates that the errors essentially depend on the location of
singularity (although it affects the constant $M_{\rho}$) rather than the
behavior of the singularity. To show this, we test the function with poles at
$\pm{\rm i}$ of order $2:$

$u(x)=\frac{1}{(x^{2}+1)^{2}}.$ (4.28)

We plot in Figure 4.1 (c)-(d) the errors and upper bounds as in (a)-(b).
Indeed, a similar convergence behavior is observed. Indeed, Boyd [6] pointed
out that the type of singularity might change the rate of convergence by a
power of $n,$ but not an exponential function of $n.$

(a) Upper bound vs Max.-error (Example 1)

(b) Upper bound vs Max.-error (Example 1)

(c) Upper bound vs Max.-error (Example 2)

(d) Upper bound vs Max.-error (Example 2)

Figure 4.1. The (discrete) maximum errors of Gegenbauer-Gauss and Gegenbauer-
Gauss-Lobatto spectral differentiation against the upper bounds in Theorems
4.2 and 4.3 with $\lambda=1/2,3/2.$ Example 1 (a)-(b), and Example 2 (c)-(d).
In the legend, GG SP err (resp. GGL SP err) represents the Gegenbauer-Gauss
(resp. Gegenbauer-Gauss-Lobatto) spectral differentiation error.

Concluding remarks In this paper, the exponential convergence of Gegenbauer
interpolation, spectral differentiation and quadrature of functions analytic
on and within a sizable ellipse is analyzed. Sharp estimates in the maximum
norm with explicit dependence on the important parameters are obtained.
Illustrative numerical results are provided to support the analysis. For
clarity of presentation, it is assumed that $\lambda$ is fixed in our
analysis, but the dependence of the error on this parameter can also be
tracked if it is necessary.

Acknowledgement

The first author would like to thank the Division of Mathematical Sciences at
Nanyang Technological University in Singapore for the hospitality during the
visit. The second author is grateful to Prof. Jie Shen at Purdue University
for bringing this topic to attention. The authors would also like to thank the
anonymous referees for their valuable comments and suggestions.

##  Appendix A Exponential convergence of Gegenbauer expansions

To this end, we discuss the exponential convergence of Gegenbauer polynomial
(with fixed $\lambda$) expansions of functions satisfying the assumption of
analyticity in [24], that is, (2.14) in the form:

$\max_{|x|\leq 1}\Big{|}\frac{d^{k}u}{dx^{k}}(x)\Big{|}\leq
C(\delta)\frac{k!}{\delta^{k}},\quad\delta\geq 1,\;\;\forall k\geq 0,$ (A.1)

where $C(\delta)$ is a positive constant depending only on $\delta.$ We write

$u(x)=\sum_{l=0}^{\infty}\hat{u}_{l}^{\lambda}C_{l}^{\lambda}(x)\;\;\;{\rm
with}\;\;\;\hat{u}_{l}^{\lambda}=\frac{1}{h_{l}^{\lambda}}\int_{-1}^{1}u(x)C_{l}^{\lambda}(x)(1-x^{2})^{\lambda-1/2}dx,$

where $h_{l}^{\lambda}$ is defined in (2.4).

######  Proposition A.1.

If $u(x)$ is an analytic function on $[-1,1]$ satisfying the assumption (A.1),
then for fixed $\lambda>-1/2,$ we have

$\max_{|x|\leq 1}\big{|}(\pi_{n}^{\lambda}u-u)(x)\big{|}\leq
C(\delta,\lambda)\frac{n^{\lambda}}{(2\delta)^{n}},$ (A.2)

where
$(\pi_{n}^{\lambda}u)(x)=\sum_{l=0}^{n}\hat{u}_{l}^{\lambda}C_{l}^{\lambda}(x),$
and $C(\delta,\lambda)$ is a positive constant depending on $\delta$ and
$\lambda.$

###### Proof.

The estimate follows directly from Theorem 4.3 in [24]. Indeed, by Theorem 4.3
in [24], we have

$\begin{split}\max_{|x|\leq 1}\big{|}(\pi_{n}^{\lambda}u-u)(x)\big{|}\leq
A\frac{C(\delta)\Gamma(\lambda+1/2)\Gamma(n+2\lambda+1)}{n\sqrt{\lambda}(2\delta)^{n}\Gamma(2\lambda)\Gamma(n+\lambda)},\quad
n\geq 1,\;\lambda>-1/2,\end{split}$

where $A$ is a generic positive constant independent of $\lambda,n$ and $u.$
For fixed $\lambda>-1/2,$ using (2.11) leads to

$\begin{split}\frac{\Gamma(n+2\lambda+1)}{\Gamma(n+\lambda+1)}&=\frac{n+\lambda+1}{n+2\lambda+1}\frac{\Gamma(n+2\lambda+2)}{\Gamma(n+\lambda+2)}\\\
&\leq
A(n+2\lambda+1)^{\lambda}\Big{(}1+\frac{\lambda}{n+\lambda+1}\Big{)}^{n+\lambda+3/2}e^{-\lambda}\\\
&\leq A(n+2\lambda+1)^{\lambda}e^{\frac{\lambda}{2(n+\lambda+1)}},\end{split}$

where in the last step, we used the inequality $1+x\leq e^{x}$ for all real
$x.$ Therefore, we obtain the estimate

$\begin{split}\max_{|x|\leq 1}\big{|}(\pi_{n}^{\lambda}u-u)(x)\big{|}\leq
C(\delta,\lambda)\frac{n+\lambda}{n}\frac{(n+2\lambda+1)^{\lambda}}{(2\delta)^{n}}.\end{split}$

This implies (A.2). ∎

######  Remark A.1.

The error $\max_{|x|\leq 1}\big{|}(\pi_{n}^{\lambda}u-u)(x)\big{|}$ is termed
as the regularization error in [24, 23]. This, together with the so-called
truncation error, contributes to the error of the Gegenbauer reconstruction.
As shown in [24], an exponential convergence could be recovered from Fourier
partial sum of a nonperiodic analytic function, when $\lambda$ grows linearly
with $n.$ However, the estimate (A.2) shows that the regularization error
alone decays exponentially for fixed $\lambda$. ∎

######  Remark A.2.

It is also interesting to study the exponential convergence of Gegenbauer
expansions of analytic functions, which are analytic on and within the
Bernstein ellipse ${\mathcal{E}}_{\rho}$ in (2.15). We refer to [10, 36, 33]
for the results on Chebyshev expansions. However, the analysis for the
Legendre expansions is much more involved. Davis [10] stated an estimate due
to K. Neumann (see Page 312 of [10]), that is, if $u$ is analytic on and
inside ${\mathcal{E}}_{\rho}$ with $\rho>1,$ then the Legendre expansion
coefficient satisfies

$\underset{{n\to\infty}}{\lim\sup}\,|\hat{u}_{n}|^{1/n}=\rho^{-1}.$ (A.3)

A more informative estimate:

$|\hat{u}_{n}|\leq C(\rho)\frac{n}{\rho^{n+1}},\quad\forall n\geq 0,$ (A.4)

was presented in the very recent paper [42]. Some discussions on the
Gegenbauer expansions can be found in [23]. We shall report the analysis for
the general Jacobi expansions in [43]. ∎

##  Appendix B Proof of Lemma 3.1

We carry out the proof by induction.

Apparently, by (2.1), $C_{0}^{\lambda}(z)=1,$ so (3.1)-(3.2) holds for $n=0.$
Similarly, we can verify the case with $n=1.$

Assume that the formula holds for $C_{n-2}^{\lambda}(z)$ and
$C_{n-1}^{\lambda}(z)$ with $n\geq 2$. It follows from the three-term
recurrence relation (2.1) that

$\begin{split}nC_{n}^{\lambda}(z)&=2(n+\lambda-1)z\sum_{k=0}^{n-1}g_{k}^{\lambda}g_{n-1-k}^{\lambda}w^{n-2k-1}\\\
&\quad-(n+2\lambda-2)\sum_{k=0}^{n-2}g_{k}^{\lambda}g_{n-2-k}^{\lambda}w^{n-2k-2}\\\
&=(n+\lambda-1)\sum_{k=0}^{n-1}g_{k}^{\lambda}g_{n-1-k}^{\lambda}w^{n-2k}+(n+\lambda-1)\sum_{k=0}^{n-1}g_{k}^{\lambda}g_{n-1-k}^{\lambda}w^{n-2k-2}\\\
&\quad-(n+2\lambda-2)\sum_{k=0}^{n-2}g_{k}^{\lambda}g_{n-2-k}^{\lambda}w^{n-2k-2}\\\
&=(n+\lambda-1)g_{n-1}^{\lambda}(w^{n}+w^{-n})+\sum_{k=1}^{n-1}D_{n,k}^{\lambda}\,g_{k}^{\lambda}g_{n-k}^{\lambda}w^{n-2k},\end{split}$

where

$D_{n,k}^{\lambda}=(n+\lambda-1)\left(\frac{g_{n-k-1}^{\lambda}}{g_{n-k}^{\lambda}}+\frac{g_{k-1}^{\lambda}}{g_{k}^{\lambda}}\right)-(n+2\lambda-2)\frac{g_{n-k-1}^{\lambda}}{g_{n-k}^{\lambda}}\frac{g_{k-1}^{\lambda}}{g_{k}^{\lambda}}.$

One verifies from (3.2) that

${g_{n}^{\lambda}}=\frac{n+\lambda-1}{n}{g_{n-1}^{\lambda}},\quad
D_{n,k}^{\lambda}=n.$ (B.1)

This completes the induction. ∎

## References

  * [1] G.E. Andrews, R. Askey, and R. Roy.  Special Functions, volume 71 of Encyclopedia of Mathematics and its Applications.  Cambridge University Press, Cambridge, 1999.
  * [2] I. Babuška and B.Q. Guo.  Direct and inverse approximation theorems for the $p$-version of the finite element method in the framework of weighted Besov spaces. I. Approximability of functions in the weighted Besov spaces.  SIAM J. Numer. Anal., 39(5):1512–1538 (electronic), 2001/02.
  * [3] C. Bernardi, M. Dauge, and Y. Maday.  Spectral Methods for Axisymmetric Domains, volume 3 of Series in Applied Mathematics (Paris).  Gauthier-Villars, Éditions Scientifiques et Médicales Elsevier, Paris, 1999.  Numerical algorithms and tests due to Mejdi Azaïez.
  * [4] C. Bernardi and Y. Maday.  Spectral Methods.  In Handbook of Numerical Analysis, Vol. V, Handb. Numer. Anal., V, pages 209–485. North-Holland, Amsterdam, 1997.
  * [5] J.P. Boyd.  The rate of convergence of Fourier coefficients for entire functions of infinite order with application to the Weideman-Cloot sinh-mapping for pseudospectral computations on an infinite interval.  J. Comput. Phys., 110(2):360–372, 1994.
  * [6] J.P. Boyd.  Chebyshev and Fourier Spectral Methods.  Dover Publications Inc., 2001.
  * [7] J.P. Boyd.  Trouble with Gegenbauer reconstruction for defeating Gibbs’ phenomenon: Runge phenomenon in the diagonal limit of Gegenbauer polynomial approximations.  J. Comput. Phys., 204(1):253–264, 2005.
  * [8] C. Canuto, M.Y. Hussaini, A. Quarteroni, and T.A. Zang.  Spectral Methods: Fundamentals in Single Domains.  Springer, Berlin, 2006.
  * [9] Y.P. Chen and T. Tang.  Convergence analysis of the Jacobi spectral-collocation methods for Volterra integral equations with a weakly singular kernel.  Math. Comp., 79(269):147–167, 2010.
  * [10] P.J. Davis.  Interpolation and Approximation.  Dover Publications, Inc, New Year, 1975.
  * [11] P.J. Davis and P. Rabinowitz.  Methods of Numerical Integration.  Computer Science and Applied Mathematics. Academic Press Inc., Orlando, FL, second edition, 1984.
  * [12] E.H. Doha.  The coefficients of differentiated expansions and derivatives of ultraspherical polynomials.  Comput. Math. Appl., 21(2-3):115–122, 1991.
  * [13] E.H. Doha.  The ultraspherical coefficients of the moments of a general-order derivative of an infinitely differentiable function.  J. Comput. Appl. Math., 89(1):53–72, 1998.
  * [14] E.H. Doha and A.H. Bhrawy.  Efficient spectral-Galerkin algorithms for direct solution for second-order differential equations using Jacobi polynomials.  Numer. Algorithms, 42(2):137–164, 2006.
  * [15] E.H. Doha and A.H. Bhrawy.  A Jacobi spectral Galerkin method for the integrated forms of fourth-order elliptic differential equations.  Numer. Methods Partial Differential Equations, 25(3):712–739, 2009\.
  * [16] E.H. Doha, A.H. Bhrawy, and W.M. Abd-Elhameed.  Jacobi spectral Galerkin method for elliptic Neumann problems.  Numer. Algorithms, 50(1):67–91, 2009.
  * [17] M. Dubiner.  Spectral methods on triangles and other domains.  J. Sci. Comput., 6(4):345–390, 1991.
  * [18] David Elliott.  Uniform asymptotic expansions of the Jacobi polynomials and an associated function.  Math. Comp., 25:309–315, 1971.
  * [19] B. Fornberg.  A Practical Guide to Pseudospectral Methods.  Cambridge Univ Pr, 1998.
  * [20] D. Funaro.  Polynomial Approxiamtions of Differential Equations.  Springer-Verlag, 1992.
  * [21] A. Gelb and Z. Jackiewicz.  Determining analyticity for parameter optimization of the Gegenbauer reconstruction method.  SIAM J. Sci. Comput., 27(3):1014–1031, 2006.
  * [22] D. Gottlieb and S.A. Orszag.  Numerical Analysis of Spectral Methods: Theory and Applications.  Society for Industrial Mathematics, 1977.
  * [23] D. Gottlieb and C.W. Shu.  On the Gibbs phenomenon and its resolution.  SIAM Rev., 39(4):644–668, 1997.
  * [24] D. Gottlieb, C.W. Shu, A. Solomonoff, and H. Vandeven.  On the Gibbs phenomenon. I. Recovering exponential accuracy from the Fourier partial sum of a nonperiodic analytic function.  J. Comput. Appl. Math., 43(1-2):81–98, 1992.  Orthogonal polynomials and numerical methods.
  * [25] B.Q. Guo.  Mathematical Foundation and Applications of the $p$ and $h$-$p$ Finite Element Methods, volume 4 of Series in Analysis.  World Scientific, 2011.
  * [26] B.Y. Guo.  Gegenbauer approximation and its applications to differential equations with rough asymptotic behaviors at infinity.  Appl. Numer. Math., 38(4):403–425, 2001.
  * [27] B.Y. Guo and L.L. Wang.  Jacobi interpolation approximations and their applications to singular differential equations.  Adv. Comput. Math., 14(3):227–276, 2001.
  * [28] B.Y. Guo and L.L. Wang.  Jacobi approximations in non-uniformly Jacobi-weighted Sobolev spaces.  J. Approx. Theory, 128(1):1–41, 2004.
  * [29] B. Gustafsson.  The work of David Gottlieb: A success story.  Commun. Comput. Phys., 9(3):481–496, 2011.
  * [30] A. Jeffrey.  Table of Integrals, Series, and Products (Seventh Edition).  Elsevier, San Diego, CA, 2007.
  * [31] G.E. Karniadakis and S.J. Sherwin.  Spectral/$hp$ Element Methods for Computational Fluid Dynamics.  Numerical Mathematics and Scientific Computation. Oxford University Press, New York, second edition, 2005.
  * [32] I. Krasikov.  An upper bound on Jacobi polynomials.  J. Approx. Theory, 149(2):116–130, 2007.
  * [33] J.C. Mason and D.C. Handscomb.  Chebyshev Polynomials.  Chapman & Hall/CRC, Boca Raton, FL, 2003.
  * [34] P. Nevai, T. Erdélyi, and A.P. Magnus.  Generalized Jacobi weights, Christoffel functions, and Jacobi polynomials.  SIAM J. Math. Anal., 25(2):602–614, 1994.
  * [35] S.C. Reddy and J.A.C. Weideman.  The accuracy of the Chebyshev differencing method for analytic functions.  SIAM J. Numer. Anal., 42(5):2176–2187 (electronic), 2005.
  * [36] T.J. Rivlin and V. Kalashnikov.  Chebyshev Polynomials: From Approximation Theory to Algebra and Number Theory.  Wiley New York, 1990.
  * [37] J. Shen, T. Tang, and L.L. Wang.  Spectral Methods : Algorithms, Analysis and Applications, volume 41 of Series in Computational Mathematics.  Springer, 2011.
  * [38] G. Szegö.  Orthogonal Polynomials (Fourth Edition).  AMS Coll. Publ., 1975.
  * [39] E. Tadmor.  The exponential accuracy of Fourier and Chebyshev differencing methods.  SIAM J. Numer. Anal., 23(1):1–10, 1986.
  * [40] L.N. Trefethen.  Spectral Methods in Matlab.  Software, Environments, and Tools. Society for Industrial and Applied Mathematics (SIAM), Philadelphia, PA, 2000.
  * [41] L.N. Trefethen.  Is Gauss quadrature better than Clenshaw-Curtis?  SIAM Rev., 50(1):67–87, 2008.
  * [42] H.Y. Wang and S.H. Xiang.  On the convergence rates of Legendre approximation.  Math. Comput., In press, 2011.
  * [43] L.L. Wang, Z.Q. Xie, and X.D. Zhao.  A note on exponential convergence of Jacobi polynomial expansions.  In preparation, 2011.
  * [44] Z.Q. Wang and B.Y. Guo.  Jacobi rational approximation and spectral method for differential equations of degenerate type.  Math. Comp., 77(262):883–907, 2008.
  * [45] G.N. Watson.  A Treatise on the Theory of Bessel Functions.  Cambridge University Press, 1995.
  * [46] Z. Zhang.  Superconvergence of spectral collocation and $p$-version methods in one dimensional problems.  Math. Comp., 74(252):1621–1636 (electronic), 2005.
  * [47] Z. Zhang.  Superconvergence of a Chebyshev spectral collocation method.  J. Sci. Comput., 34(3):237–246, 2008.


[a]Oscar Garcia-Montero

# Charge and energy deposition in the McDipper framework

Sören Schlichting  Hannah Elfner

###### Abstract

In this short note we present aspects of the energy and charge deposition
within the McDipper, a novel 3D resolved model for the initial state of
ultrarelativistic Heavy-Ion collisions based on the $k_{\perp}$-factorized
Color Glass Condensate hybrid approach. This framework is a initial-state
Monte Carlo event generator which deposits the relevant conserved charges
(energy, charge and baryon densities) both in the midrapidity and
forward/backward regions of the collision. The event-by-event generator
computes the gluon and (anti-) quark phase-space densities using the IP-Sat
model, from where the conserved charges can be extracted directly. In this
work we present the centrality and collision energy dependence for the
deposited conserved quantities at midrapidity and the full event, the so-
called $4\pi$ solid angle range.

##  1 Introduction

It is widely accepted that during a Heavy Ion Collision (HICs) the formation
of a complex, hot and dense system comprised of quarks and gluons is formed,
the so-called Quark-Gluon Plasma. In recent years, an interest has risen over
observations of a rich longitudinal structure of the long-range correlations
of the system [1, 2, 3, 4]. Such structure is thought to be created through
diverse mechanisms during the initial instants, and translated to momentum
space through interactions. Therefore, new experimental and theoretical
insights towards the forward/backward rapidity windows of heavy ion collisions
call for a fast and phenomenologically available initial conditions model with
non-trivial rapidity dependence and strongly inspired by first principles
computations. Motivated by these ideas, as well as by the need for a framework
for comprehensive comparison of saturation physics and models in HICs and the
upcoming Electron-Ion Collider, we recently presented the McDipper framework
[5]. This saturation based frameworks provides a flexible and systematically
improvable 3D initial conditions for HICs, firmly rooted in high energy
Quantum Chromodynamics.

Expanding on our recent work, we present a study on the initial charge and
energy deposition, valuable for many aspects in high energy nuclear physics,
ranging from the initial dynamics of the pre-equilibrium stage, the underlying
understanding of net-particle yields and correlations and the production of
electromagnetic, among others. We will focus on the deposition at midrapidity
of energy, $\mathrm{d}E/\mathrm{d}\eta_{s}$, electric charge,
$\mathrm{d}Q/\mathrm{d}\eta_{s}$, and baryon charge ,
$\mathrm{d}B/\mathrm{d}\eta_{s}$. We will compare this to the event total,
i.e. the rapidity integrated-charges. In what follows, we will call this
rapidity integrated charges the full solid angle or $4\pi$ charges, e.g.
$Q_{4\pi}\equiv\int\mathrm{d}\eta_{s}\mathrm{d}Q/\mathrm{d}\eta_{s}$.

##  2 The McDipper framework

The McDipper framework computes the initial energy and charge deposition in
high-energy HICs within the dilute-dense approximation of the $k_{\perp}$
factorized Color Glass Condensate Effective Field Theory [6]. Using the single
particle production formulas, one can compute inclusive gluon
$(\mathrm{d}N_{g}/{\mathrm{d}^{2}\mathbf{x}\mathrm{d}^{2}\mathbf{p}\mathrm{d}y})$
and (net-) quark distributions
$\mathrm{d}N_{\bar{q}-q}/{\mathrm{d}^{2}\mathbf{x}\mathrm{d}^{2}\mathbf{p}\mathrm{d}y}$
as a function of (momentum) rapidity $y$ transverse momentum $\mathbf{p}$ and
transverse position $\mathbf{x}$. This is done by evaluating the leading order
cross-sections for the transverse momentum dependent dipole gluon
distributions from the IPSat model [7, 8] and collinear parton distributions
(PDFs) from the LHAPDF library [10]. For details on the computation of the
gluon and quark single particle distributions please see the main reference,
ref. [5]. The relevant macroscopic conserved quantities can be calculated by
taking moments of the single particle distributions and noting that at the LO
in the high energy limit, the indentification of the momentum and spacetime
rapidities can be taken, $y=\eta_{s}$. Energy deposition in this framework
obeys the relation

$(e\tau)_{0}=\int\mathrm{d}^{2}\mathbf{p}~{}|\mathbf{p}|~{}\left[K_{g}\frac{\mathrm{d}N_{g}}{\mathrm{d}^{2}\mathbf{x}\mathrm{d}^{2}\mathbf{p}\mathrm{d}y}+\sum_{f,\bar{f}}\frac{\mathrm{d}N_{q_{f}}}{\mathrm{d}^{2}\mathbf{x}\mathrm{d}^{2}\mathbf{p}\mathrm{d}y}\right]_{y=\eta_{s}}\;,$
(1)

where the prefactor $K_{g}$ is a free paramter which accounts for the
uncertainties coming from higher order perturbative corrections to the total
cross-sections. Notice that since the IP-Sat model have been fitted to HERA
data [9], the $K_{g}$ prefactor in the gluon energy is the only free parameter
in the framework and is fixed by fitting to the deposited energy of $pp$ and
$pPb$ collisions at minimum bias for some collisional center-of-mass energy
per nucleon,$\sqrt{s_{\rm NN}}=5.02$. We choose to fix at $\sqrt{s_{\rm
NN}}=5.02$ TeV, leading to $K_{g}=1.85$ for the IP-Sat model. On the other
hand, the quark number densities $n_{f}$ with $f=(u,d,s)$, can be computed
using

$(n_{f}\tau)_{0}=\int\mathrm{d}^{2}\mathbf{p}~{}\left[\frac{\mathrm{d}N_{q_{f}}}{\mathrm{d}^{2}\mathbf{x}\mathrm{d}^{2}\mathbf{p}\mathrm{d}y}-\frac{\mathrm{d}N_{\bar{q}_{f}}}{\mathrm{d}^{2}\mathbf{x}\mathrm{d}^{2}\mathbf{p}\mathrm{d}y}\right]_{y=\eta_{s}}$
(2)

The electric charge density is then computed as $(\tau
Q)_{0}=\sum_{f}Q_{f}\,n_{f}$ where $Q_{f}=-1/3$ for $f=d,s$ quarks and
$Q_{f}=2/3$ for the $u$ quark. while the baryon charge is given by $(\tau
B)_{0}=Q_{B}\sum_{f}\,n_{f}$, with $Q_{B}=1/3$ the baryon charge of a valence
quark. The framework does not yet include spatial charge fluctuation coming
from fluctations in the PDFs, and therefore the strangeness density vanishes
trivially. With these conditions, a typical event in the McDipper framework
can be pictured in Fig. 1, where energy (left) and baryon charge (right) are
pictured for a single event at impact parameter $b=2\,$fm for a Pb-Pb
collision at $\sqrt{s_{\rm NN}}=5.02$ TeV. As one can expect, the energy
density is symmetric and higher at midrapidity. As we will see in the
following section this is thanks to the gluon domination of the energy density
at high collisional energies close to $\eta_{s}=0$. The interesting
visualization is the right hand panel where we find the baryon density, which,
as expected for a high energy event, is very small at midrapidity and peaks in
the forward/backward regions. The charge distributions are in general
asymmetric, thanks to their dependence on the nuclear density. As has been
more extensively explained in ref. [5], this dependence is not trivial, but
one can schematically think of it roughly as $T_{1}T_{2}^{a}+T_{2}T_{1}^{a}$,
where$T_{1/2}$ denotes the nuclear thickness, and $a\sim 1/2$. Additionally,
if one wants to find the charge and energy yields, e.g. the baryon charge
yield, one needs to use $(\tau
B)_{0}=\mathrm{d}B/\mathrm{d}^{2}\mathbf{x}\mathrm{d}\eta_{s}$, which means
that $B_{4\pi}=\int\mathrm{d}^{2}\mathbf{x}\,\mathrm{d}\eta_{s}(\tau B)_{0}$.
The same is true for the energy and electric charge densities.

Figure 1: Full spatial resolution for the deposited energy (left) and baryon
charge (right) of a single event of a Pb-Pb collision at$\sqrt{s_{\rm
NN}}=5.02$ TeV and impact parameter $b=2$ fm. The computation is done for the
IP-Sat model with the CT18NNLO parton distribution set [7, 8, 10].

##  3 Results

In Fig. 2 the reader can find several properties of the deposited energy in
the McDipper framework. In the left panel, the excitation function (the
$\sqrt{s_{\rm NN}}$ dependence) of both the midrapidity and total energy are
plotted for different centrality classes. Additionally, the reader can observe
that the midrapidity slice and the total energy deposited scale roughly with
$\sqrt{s_{\rm NN}}$ as a power law sufficiently slower than linear increase.
The center panel shows the deposited energy at midrapidity as normalized by
the full $4\pi$ integrated energy. The ratio exhibits a very minor residual
dependence on $\sqrt{s_{\rm NN}}$, from which we can conclude that even though
midrapidity energy increases with $\sqrt{s_{\rm NN}}$, more energy is flowing
towards the fragmentation regions at higher collisional energies.

Figure 2:  (Left) Energy deposition for several centralities at midrapidity
($\eta_{s}=0$) and in $4\pi$ as a function collision energy $\sqrt{s_{\rm
NN}}$. (Center) Ratio of the total energy deposited at midrapidity normalized
with the total deposited energy. (Right) Dependence of the fraction of the
gluon and quark energy contributions as a function of $\sqrt{s_{\rm NN}}$.
Here, the fraction is obtained by normalization to the respective energy at
midrapidity (full symbols), as well as the total energy (open symbols).

In the right panel one can find the relative contribution of quarks (diamonds)
and gluons (circles) to the energy deposition. In this panel, contributions
taken with respect only to the midrapidity are marked as full symbols, while
contributions taken for the whole event are shown as empty markers. For
simplicity, only two centralities are shown, $0-5\%$ (line) and $40-50\%$
(dashed line). In this figure, one can observe that, as intuitively expected
from the model, the gluon contribution becomes more important at higher
energies. This dependence is more intense for the case of the $4\pi$
integrated energy, which can be understood as the addition of the energy
deposited in the fragmentation regions, where quarks with larger $x$ can be
easily deflected and deposited into a medium at lower collisional energies.
Additionally, the relative contributions exhibit a quite soft dependence on
centrality, where the quark contribution becomes more significant for more
peripheral events. In the advent of full 3+1D hybrid models, these last
findings become quite relevant, especially in the early stages after the
initial deposition. Since quark deposition comprises a minimum of 20% on the
current energies available, one can question how much will this inflict the
evolution of the thermalizing plasma. These effects may be particularly
important in observables such as photon and dilepton production, which have
been shown to be sensitive to the process of chemical equilibration of the QGP
[11], and such effects may be addressed in the future in a quantitative manner
using kinetic theory [12, 13].

Figure 3: Deposited electric (left) and baryon (center) charges for different
centrality classes as a function of $\sqrt{s_{\rm NN}}$. Shown in full markers
are the midrapidity slices, while the empty markers correspond to the total
charge deposited with $4\pi$ acceptance. The slim black line is the total
incoming electric and baryon charge, $Q_{0}=2\,Z$ and $B_{0}=2\,A$. (Right)
Ratio of electric and baryon charges at midrapidity to their $4\pi$ integrated
corresponding charge.

In Fig. 3 the reader can find several properties of the deposited electric and
baryon charges, here denoted as $Q$ and $B$, respectively. In the left and
center panels of Fig. 3 one can find $Q$ and $B$, respectively, for the
midrapidity slice as well as for the total integrated case. The bahavior of
these quantitites is qualitatively similar. It is interesting to see that the
total integrated charges are virtually constant with respect to $\sqrt{s_{\rm
NN}}$ , and their deposition depends uniquely on the geometrical constraints,
coming from the sampling of the nucleons from a Woods-Saxon’s profile, and
results from a strict constraining of the conserved charges throught the
parton PDFs. On the other hand the charge deposited at midrapidity is
monotonically decreasing with $\sqrt{s_{\rm NN}}$. As expected, nuclei become
more transparent at midrapidity. The deposition of charge is shifted with
collisional energy more and more towards the fragmentation regions. In the
right panel of Fig. 3 the ratio of the midrapidity and total deposited charge
is plotted for both Q and B. The behavior of both charges is identical when
comparing the ratios, indicating that the behavior is purely geometrical, as
it exhibits exactly the same behavior for electric and baryon charges,
regardless of system (Au-Au/ Pb-Pb) or collisional energy.

##  4 Conclusions and Outlook

In this short note we have shown several aspects of the initial energy and
charge deposition of a 3D resolved heavy-ion collision in the McDipper
framework. For this we have used the language of the so-called excitation
functions to show that the energy deposition presents a soft power-law
dependence on $\sqrt{s_{\rm NN}}$, very similar for both midrapidity and
$4\pi$ deposition. Additionally, we have shown the behavior of the quark and
gluon relative contribution with respect to $\sqrt{s_{\rm NN}}$, which
exhibits a sharper dependence for the $4\pi$ energy as for the midrapidity
slice.

Additionally, we have shown that the full $4\pi$ deposition of the electric
and baryon charges is only geometrical, as it has virtually no dependence on
$\sqrt{s_{\rm NN}}$. On the midrapidity slice both charges exhibit a steep
decrease with rising energy, which is expected from nuclei becoming more
transparent, and the baryon stopping is shifted towards higher values of
$\eta_{s}$. We find that once normalized by their respective full $4\pi$
values of each charge, these ratios exhibit the same values, allowing for a
purely geometrical interpretation. For more detailed on the general aspects of
the framework and the relevant observables, the reader can find more results
on the McDipper framework in the main paper, see ref. [5].

Future endeavors for the framework will contain localized fluctuation of
charge deposition, which will give not only sub-nucleonic hotspots of electric
and baryon charge, but also a locally non-vanishing strangeness distribution.
That said, such fluctations should should not affect the results of this work,
as the variations should vanish for the expectation values of the charge
observables, remaining only at the level of higher order correlations.

## References

  * [1] V. Khachatryan et al. [CMS], Phys. Rev. C 92 (2015) no.3, 034911 doi:10.1103/PhysRevC.92.034911 [arXiv:1503.01692 [nucl-ex]].
  * [2] M. Nie [STAR], Nucl. Phys. A 982 (2019), 403-406 doi:10.1016/j.nuclphysa.2018.09.068
  * [3] H. Petersen, V. Bhattacharya, S. A. Bass and C. Greiner, Phys. Rev. C 84 (2011), 054908 doi:10.1103/PhysRevC.84.054908 [arXiv:1105.0340 [nucl-th]].
  * [4] L. G. Pang, H. Petersen and X. N. Wang, Phys. Rev. C 97 (2018) no.6, 064918 doi:10.1103/PhysRevC.97.064918 [arXiv:1802.04449 [nucl-th]].
  * [5] O. Garcia-Montero, H. Elfner and S. Schlichting, [arXiv:2308.11713 [hep-ph]].
  * [6] L. V. Gribov, E. M. Levin and M. G. Ryskin, Phys. Rept. 100 (1983), 1-150 doi:10.1016/0370-1573(83)90022-4
  * [7] J. Bartels, K. J. Golec-Biernat and H. Kowalski, Phys. Rev. D 66 (2002), 014001 doi:10.1103/PhysRevD.66.014001 [arXiv:hep-ph/0203258 [hep-ph]].
  * [8] H. Kowalski and D. Teaney, Phys. Rev. D 68 (2003), 114005 doi:10.1103/PhysRevD.68.114005 [arXiv:hep-ph/0304189 [hep-ph]].
  * [9] A. H. Rezaeian, M. Siddikov, M. Van de Klundert and R. Venugopalan, Phys. Rev. D 87 (2013) no.3, 034002 doi:10.1103/PhysRevD.87.034002 [arXiv:1212.2974 [hep-ph]].
  * [10] A. Buckley, J. Ferrando, S. Lloyd, K. Nordström, B. Page, M. Rüfenacht, M. Schönherr and G. Watt, Eur. Phys. J. C 75 (2015), 132 doi:10.1140/epjc/s10052-015-3318-8 [arXiv:1412.7420 [hep-ph]].
  * [11] C. Gale, J. F. Paquet, B. Schenke and C. Shen, Phys. Rev. C 105 (2022) no.1, 014909 doi:10.1103/PhysRevC.105.014909 [arXiv:2106.11216 [nucl-th]].
  * [12] O. Garcia-Montero, A. Mazeliauskas, P. Plaschke and S. Schlichting, [arXiv:2308.09747 [hep-ph]].
  * [13] X. Du and S. Schlichting, Phys. Rev. D 104 (2021) no.5, 054011 doi:10.1103/PhysRevD.104.054011 [arXiv:2012.09079 [hep-ph]].