%% Exact Chebyshev expansion coefficients of a function
% Mark Richardson, 13th June 2012
%%
% (Chebfun example approx/ExactChebCoeffs.m)
% [Tags: #Chebyshevcoefficients, #residue]
%% 1. Introduction
% In this example, we shall compare the results of the Chebfun construction
% process to a known closed-form formula for the Chebyshev expansion
% coefficients of a function with a pole.
%
% This gives us an excellent excuse to visit some interesting approximation
% theory from the 1960s!
%% 2. The residue method of Elliott
% For certain functions, explicit formulas for the coefficients in the
% Chebsyhev series expansion may be obtained using a contour-integral
% technique described by Elliott [1]. Here is how it works:
%
% The Chebsyhev coefficients of a Lipschitz-continuous function f can be
% determined by the integral
% $$ a_n = \frac{2}{\pi} \int_{-1}^{1} \frac{f(x)T_n(x)}{\sqrt{1-x^2}}
% {\rm d} x.$$
%%
% If f(x) is analytic within a particular contour C in the complex-plane,
% then by Cauchy's integral formula, we can also write
% $$ f(x) = \frac{1}{2 \pi i} \int_{C} \frac{f(z)}{z-x} {\rm d} z.$$
%%
% If C is large enough to enclose the unit interval, then the second of
% these two formulas can be subsituted into the first and the orders of
% integration interchanged to give
%
% $$ a_n = \frac{1}{\pi^2 i} \int_C f(z) \int_{-1}^{1}
% \frac{T_n(x) {\rm d} x}{(z-x)\sqrt{1-x^2}} {\rm d} z.$$
%%
% The integral with respect to x can be computed exactly so that we end up
% with
% $$a_n = \frac{1}{\pi i} \int_{C} \frac{f(z)}{\sqrt{z^2-1}(z \pm
% \sqrt{z^2-1})^n} {\rm d} z.$$
%%
% Here, we note that $\rho = |z \pm \sqrt{z^2-1}|$ is the parameter of the
% usual Bernstein Ellipse $E_\rho$ with the point $z = x + iy$ on its
% boundary.
%
% So, if $f$ is a function with a pole at $z_0$ whose integral around
% $E_\rho$ tends to zero as $\rho \to \infty$, then by the Residue theorem
% we have
% $$ a_n = \frac{-2r_0}{\sqrt{z_0^2-1}(z_0 \pm \sqrt{z_0^2-1})^n},$$
% where $r_0$ is the residue of the pole at $z_0$.
%% 3. A function with pole
% As an example, consider the function $$ f(x) = \frac{1}{5 + x} .$$
%%
% This function can be represented in Chebfun by an interpolant in
% 17 points.
f = @(x) 1./(5+x);
fc = chebfun(f);
k = 1:length(fc);
%%
% The function has a pole at -5 and corresponding residue 1. Substituing
% these values into the above formula then gives the us following exact
% expression for the Chebyshev expansion coefficients
% $$ a_n = \frac{1}{\sqrt{6}} \frac{(-1)^n}{(5+\sqrt{24})^n} $$
%%
% The theoretical coefficients match the ones computed by Chebfun, apart
% from floating point representation and aliasing effects. The $a_0$
% coefficient is out by the usual factor of 2.
exact_coeffs = flipud( (1/sqrt(6)*(-1).^(k-1)./(5+sqrt(24)).^(k-1))' );
cheb_coeffs = chebpoly(fc)';
display([exact_coeffs cheb_coeffs exact_coeffs-cheb_coeffs])
chebpolyplot(fc)
title('Chebyshev coefficients of 1/(5+x)')
xlabel('n'), ylabel('log(|a_n|)'), grid on
%%
% References:
%
% [1] D. Elliott, The evaluation and estimation of the coefficients in the
% Chebyshev series expansion of a function, Math. Comp. 18 (1964).