### 1.1 What is a chebfun?

A chebfun is a function of one variable defined on an interval $[a,b]$. The syntax for chebfuns is almost exactly the same as the usual MATLAB syntax for vectors, with the familiar MATLAB commands for vectors overloaded in natural ways. Thus, for example, whereas sum(f) returns the sum of the entries when f is a vector, it returns a definite integral when f is a chebfun.

Chebfun with a capital C is the name of the software system.

The aim of Chebfun is to "feel symbolic but run at the speed of numerics". More precisely, our vision is to achieve for functions what floating-point arithmetic achieves for numbers: rapid computation in which each successive operation is carried out exactly apart from a rounding error that is very small in relative terms [Trefethen 2007].

The implementation of Chebfun is based on the mathematical fact that smooth functions can be represented very efficiently by polynomial interpolation in Chebyshev points, or equivalently, thanks to the Fast Fourier Transform, by expansions in Chebyshev polynomials. For a simple function, 20 or 30 points often suffice, but the process is stable and effective even for functions complicated enough to require 1000 or 1,000,000 points. Chebfun makes use of adaptive procedures that aim to find the right number of points automatically so as to represent each function to roughly machine precision, that is, about 15 or 16 digits of relative accuracy. (Originally Chebfun stored function values at Chebyshev points; in Version 5 it switched to storing Chebyshev expansion coefficients.)

The mathematical foundations of Chebfun are for the most part well established by results scattered throughout the 20th century. A key early figure, for example, was Bernstein in the 1910s. Much of the relevant material can be found collected in the Chebfun-based book Approximation Theory and Approximation Practice [Trefethen 2013].

Chebfun was originally created by Zachary Battles and Nick Trefethen at Oxford during 2002-2005 [Battles & Trefethen 2004]. Battles left the project in 2005, and soon four new members were added to the team: Ricardo Pachon (from 2006), Rodrigo Platte (from 2007), and Toby Driscoll and Nick Hale (from 2008). In 2009, Asgeir Birkisson and Mark Richardson also became involved, and other contributors included Pedro Gonnet, Joris Van Deun, and Georges Klein. Nick Hale served as Director of the project during 2010-2014. The Chebfun Version 5 rewrite was directed by Nick Hale during 2013-2014, and the team included Anthony Austin, Asgeir Birkisson, Toby Driscoll, Hrothgar, Mohsin Javed, Hadrien Montanelli, Alex Townsend, Nick Trefethen, Grady Wright, and Kuan Xu. October 2014 brough new arrivals Jared Aurentz, and Behnam Hashemi. In 2019 the team includes also Nicolas Boulle, Abi Gopal, Yuji Nakatsukasa, and Ryan Sherbo. Further information about Chebfun history is available at the Chebfun web site, http://www.chebfun.org, where one can also find a discussion of other software projects related to Chebfun. This Guide is based on Chebfun Version 5.7.0, released in June 2017.

### 1.2 Constructing simple chebfuns

The chebfun command constructs a chebfun from a specification such as a string or an anonymous function. If you don't specify an interval, then the default interval $[-1,1]$ is used. For example, the following command makes a chebfun corresponding to $\cos(20x)$ on $[-1,1]$ and plots it.

  f = chebfun('cos(20*x)');
plot(f), ylim([-1.2,1.2]) From this little experiment, you cannot see that f is represented by a polynomial. One way to see this is to find the length of f:

  length(f)
ans =
51


Another is to remove the semicolon that suppresses output:

  f
f =
chebfun column (1 smooth piece)
interval       length     endpoint values
[      -1,       1]       51      0.41     0.41
vertical scale =   1


These results tell us that f is represented by a polynomial interpolant through 51 Chebyshev points, i.e., a polynomial of degree 50. These numbers have been determined by an adaptive process. We can see the data points by plotting f with the '.-' option:

  plot(f,'.-'), ylim([-1.2 1.2]) The formula for $N+1$ Chebyshev points in $[-1,1]$ is $$x(j) = -\cos(j \pi/N), \quad j = 0:N,$$ and in the figure we can see that the points are clustered accordingly near $1$ and $-1$. Note that in the middle of the grid, there are about 5 points per wavelength, which is evidently what it takes to represent this cosine to 15 digits of accuracy. For intervals other than $[-1,1]$, appropriate Chebyshev points are obtained by a linear scaling.

The curve between the data points is the polynomial interpolant, which can be evaluated by the barycentric formula introduced by Salzer [Berrut & Trefethen 2004, Salzer 1972]. This method of evaluating polynomial interpolants is stable and efficient even if the degree is in the millions [Higham 2004]. Chebfun actually evaluates polynomials from their Chebyshev series rather than by barycentric interpolation; the difference in the two methods is little.

What is the integral of $f$ from $-1$ to $1$? Here it is:

  sum(f)
ans =
0.091294525072763


This number was computed by integrating the polynomial (Clenshaw-Curtis quadrature -- see Section 2.1), and it is interesting to compare it to the exact answer from calculus:

  exact = sin(20)/10
exact =
0.091294525072763


Here is another example, now with the chebfun defined by an anonymous function instead of a string. In this case the interval is specified as $[0,100]$.

  g = chebfun(@(t) besselj(0,t),[0,100]);
plot(g), ylim([-.5 1]) The function looks complicated, but it is actually a polynomial of surprisingly small degree:

  length(g)
ans =
89


Is it accurate? Well, here are three random points in $[0,100]$:

  format long
x = 100*rand(3,1)
x =
81.472368639317892
90.579193707561927
12.698681629350606


Let's compare the chebfun to the true Bessel function at these points:

  exact = besselj(0,x);
error = g(x) - exact;
[g(x) exact error]
ans =
0.048059538377880   0.048059538377880   0.000000000000000
-0.021311086924089  -0.021311086924090   0.000000000000000
0.176415464952875   0.176415464952875  -0.000000000000000


If you want to know the first 5 zeros of the Bessel function, here they are:

  r = roots(g); r = r(1:5)
r =
2.404825557695756
5.520078110286287
8.653727912911013
11.791534439014278
14.930917708487790


Notice that we have just done something nontrivial and potentially useful. How else would you find zeros of the Bessel function so readily? As always with numerical computation, we cannot expect the answers to be exactly correct, but they will usually be very close. In fact, these computed zeros are accurate to close to machine precision:

  besselj(0,r)
ans =
1.0e-14 *
0.891222737425086
-0.809539285445395
0.010653221425832
-0.058908346770306
-0.092296985262750


Most often we get a chebfun by operating on other chebfuns. For example, here is a sequence that uses plus, times, divide, and power operations on an initial chebfun x to produce a famous function of Runge:

  x = chebfun('x');
f = 1/(1+25*x^2);
length(f)
clf, plot(f)
ans =
187 ### 1.3 Operations on chebfuns

There are more than 300 commands that can be applied to a chebfun. For a list of many of them you can type methods:

  methods chebfun
Methods for class chebfun:

abs              csc              jaccoeffs        range
acos             cscd             join             rank
acosd            csch             jump             rdivide
acosh            ctranspose       kron             real
acot             cummax           ldivide          reallog
acotd            cummin           le               realpow
acoth            cumprod          legcoeffs        realsqrt
acsc             cumsum           length           rem
acscd            cylinder         log              remez
acsch            deriv            log10            removeDeltas
airy             diff             log1p            repmat
all              dirac            log2             residue
and              disp             logical          restrict
angle            display          loglog           roots
any              domain           lt               round
arcLength        eigs             lu               sec
area             ellipj           mat2cell         secd
arrowplot        ellipke          max              sech
asec             end              mean             semilogx
asecd            eq               measure          semilogy
asech            erf              merge            sign
asin             erfc             mesh             simplify
asind            erfcinv          min              sin
asinh            erfcx            minandmax        sinc
atan             erfinv           minus            sind
atan2            exp              mldivide         sinh
atan2d           expm             mod              size
atand            expm1            movie            sound
atanh            feval            mpower           spy
besselh          fill             mrdivide         sqrt
besseli          find             mtimes           std
besselj          fix              nchoosek         subsasgn
besselk          fliplr           ne               subspace
bessely          flipud           newDomain        subsref
bvp4c            floor            nextpow2         sum
bvp5c            fracDiff         norm             surf
cat              fracInt          normal           surface
ceil             fred             normest          surfc
cf               gamma            not              svd
cheb2cell        ge               null             tan
cheb2quasi       get              num2cell         tand
chebcoeffs       gmres            or               tanh
chebellipseplot  gt               orth             times
chebfun          heaviside        overlap          transpose
chebtune         hypot            permute          trigremez
circconv         imag             pinv             truncate
comet            innerProduct     plot             uminus
comet3           integral         plot3            unwrap
complex          inv              plotcoeffs       uplus
compose          isPeriodicTech   plotregion       vander
cond             isdelta          plus             var
conj             isempty          poly             vertcat
conv             isequal          polyfit          volt
cos              isfinite         polyval          vscale
cosd             ishappy          pow2             waterfall
cosh             isinf            power            why
cot              isnan            prod             xor
cotd             isreal           qr
coth             issing           quantumstates
cov              iszero           quasi2cheb

Static methods:

dct              idst             nufft            pchip
dlt              interp1          nufft2           spline
dst              inufft           ode113           update
idct             lagrange         ode15s
idlt             ndct             ode45



To find out what a command does, you can use help.

  help chebfun/sum
 SUM   Definite integral of a CHEBFUN.
SUM(F) is the integral of a column CHEBFUN F over its domain of definition.

SUM(F, A, B), where A and B are scalars, integrates a column CHEBFUN F over
[A, B], which must be a subdomain of F.domain:

B
/
SUM(F) =  | F(t) dt.
/
A

SUM(F, A, B), where A and B are CHEBFUN objects, returns a CHEBFUN S which
satisfies

B(s)
/
S(s) =  | F(t) dt.
/
A(s)

SUM(F, DIM), where DIM is one of 1, 2, sums F over the dimension DIM. If F
is a column CHEBFUN and DIM = 1 or if F is a row CHEBFUN and DIM = 2 then
this integrates in the continuous dimension of F, as described above.
Otherwise, SUM(F, DIM) sums across the columns (rows) of the column (row)
CHEBFUN F.



Most of the commands in the list above exist in ordinary MATLAB; some exceptions are domain, restrict, and chebcoeffs. We have already seen length and sum in action. In fact we have already seen subsref too, since that is the MATLAB command for (among other things) evaluating arguments in parentheses. Here is another example of its use:

  f(0.5)
ans =
0.137931034482759


Here for comparison is the true result:

  1/(1+25/4)
ans =
0.137931034482759


In this Runge function example, we have also implicitly seen times, plus, power, and rdivide, all of which have been overloaded from their usual MATLAB uses to apply to chebfuns.

In the next part of this tour we shall explore many of these commands systematically. First, however, we should see that chebfuns are not restricted to smooth functions.

### 1.4 Piecewise smooth chebfuns

Many functions of interest are not smooth but piecewise smooth. In this case a chebfun may consist of a concatenation of smooth pieces, each with its own polynomial representation. Each of the smooth pieces is called a "fun". This enhancement of Chebfun was developed initially by Ricardo Pachon during 2006-2007, then also by Rodrigo Platte starting in 2007 [Pachon, Platte and Trefethen 2010]. Essentially funs are the "classic chebfuns" for smooth functions on $[-1,1]$ originally implemented by Zachary Battles in Chebfun Version 1.

Later we shall describe the options in greater detail, but for the moment let us see some examples. One way to get a piecewise smooth function is directly from the Chebfun constructor, taking advantage of its capability of automatic edge detection. For example, in the default "splitting off" mode a function with a jump in its derivative produces a warning message,

  f = chebfun('abs(x-.3)');
Warning: Function not resolved using 65537 pts. Have you tried 'splitting on'?


The same function can be successfully captured with splitting on:

  f = chebfun('abs(x-.3)','splitting','on');

The length command reveals that f is defined by four data points, two for each linear interval:

  length(f)
ans =
4


We can see the structure of f in more detail by typing f without a semicolon:

  f
f =
chebfun column (2 smooth pieces)
interval       length     endpoint values
[      -1,     0.3]        2       1.3        0
[     0.3,       1]        2   5.6e-17      0.7
vertical scale = 1.3    Total length = 4


This output confirms that f consists of two funs, each defined by two points and two corresponding function values. The functions live on intervals defined by breakpoints at $-1$, $1$, and a number very close to $0.3$.

Another way to make a piecewise smooth chebfun is to construct it explicitly from various pieces. For example, the following command specifies three functions $x^2$, $1$, and $4-x$, together with a vector of endpoints indicating that the first function applies on $[-1,1]$, the second on $[1,2]$, and the third on $[2,4]$:

  f = chebfun({@(x) x^2, 1, @(x) 4-x},[-1 1 2 4]);
plot(f) We expect f to consist of three pieces of lengths 3, 1, and 2, and this is indeed the case:

  f
f =
chebfun column (3 smooth pieces)
interval       length     endpoint values
[      -1,       1]        3         1        1
[       1,       2]        1         1        1
[       2,       4]        2         2        0
vertical scale =   2    Total length = 6


Our eyes see pieces, but to Chebfun, f is just another function. For example, here is its integral.

  sum(f)
ans =
3.666666666666667


Here is an algebraic transformation of f, which we plot in another color for variety.

  plot(1/(1+f),'r') Some Chebfun commands naturally introduce breakpoints in a chebfun. For example, the abs command first finds zeros of a function and introduces breakpoints there. Here is a chebfun consisting of 6 funs:

  f = abs(exp(x).*sin(8*x));
plot(f) And here is an example where breakpoints are introduced by the max command, leading to a chebfun with 13 pieces:

  f = sin(20*x);
g = exp(x-1);
h = max(f,g);
plot(h), ylim([0 1.2]) As always, h may look complicated to a human, but to Chebfun it is just a function. Here are its mean, standard deviation, minimum, and maximum:

  mean(h)
ans =
0.578242020778010

  std(h)
ans =
0.280937455806246

  min(h)
ans =
0.135335283236613

  max(h)
ans =
1.000000000000000


A final note about piecewise smooth chebfuns is that the automatic edge detection or "splitting" feature, when it is turned on, may subdivide functions even though they do not have clean point singularities, and this may be desirable or undesirable depending on the application. For example, considering $\sin(x)$ over $[0,1000\pi]$ with splitting on, we end up with a chebfun with many pieces:

  tic, f = chebfun('sin(x)',[0 1000*pi],'splitting','on'), toc
f =
chebfun column (32 smooth pieces)
interval       length     endpoint values
[       0,      98]       91   2.3e-14    -0.71
[      98,   2e+02]       88     -0.71        1
[   2e+02, 2.9e+02]       88         1    -0.71
[ 2.9e+02, 3.9e+02]       87     -0.71  1.5e-13
[ 3.9e+02, 4.9e+02]       87   2.5e-13     0.71
[ 4.9e+02, 5.9e+02]       87      0.71       -1
[ 5.9e+02, 6.9e+02]       86        -1     0.71
[ 6.9e+02, 7.9e+02]       85      0.71   -6e-13
[ 7.9e+02, 8.8e+02]       86  -9.8e-13    -0.71
[ 8.8e+02, 9.8e+02]       86     -0.71        1
[ 9.8e+02, 1.1e+03]       86         1    -0.71
[ 1.1e+03, 1.2e+03]       86     -0.71 -8.4e-14
[ 1.2e+03, 1.3e+03]       85  -4.3e-13     0.71
[ 1.3e+03, 1.4e+03]       86      0.71       -1
[ 1.4e+03, 1.5e+03]       86        -1     0.71
[ 1.5e+03, 1.6e+03]       85      0.71 -1.3e-12
[ 1.6e+03, 1.7e+03]       86  -9.9e-13    -0.71
[ 1.7e+03, 1.8e+03]       86     -0.71        1
[ 1.8e+03, 1.9e+03]       86         1    -0.71
[ 1.9e+03,   2e+03]       85     -0.71  6.7e-13
[   2e+03, 2.1e+03]       85   1.7e-12     0.71
[ 2.1e+03, 2.2e+03]       86      0.71       -1
[ 2.2e+03, 2.3e+03]       86        -1     0.71
[ 2.3e+03, 2.4e+03]       87      0.71 -1.5e-12
[ 2.4e+03, 2.5e+03]       85     1e-12    -0.71
[ 2.5e+03, 2.6e+03]       86     -0.71        1
[ 2.6e+03, 2.7e+03]       84         1    -0.71
[ 2.7e+03, 2.7e+03]       85     -0.71   -2e-12
[ 2.7e+03, 2.8e+03]       85  -1.2e-14     0.71
[ 2.8e+03, 2.9e+03]       84      0.71       -1
[ 2.9e+03,   3e+03]       86        -1     0.71
[   3e+03, 3.1e+03]       85      0.71 -2.8e-12
vertical scale =   1    Total length = 2752
Elapsed time is 0.265972 seconds.


In this case it is more efficient -- and more interesting mathematically -- to omit the splitting and construct one global chebfun:

  tic, f2 = chebfun('sin(x)',[0 1000*pi]), toc
f2 =
chebfun column (1 smooth piece)
interval       length     endpoint values
[       0, 3.1e+03]     1684   1.1e-14 -2.8e-13
vertical scale =   1
Elapsed time is 0.035212 seconds.


Splitting on and off are discussed further in Section 8.3.

### 1.5 Infinite intervals and infinite function values

A major change from Chebfun Version 2 to Version 3 was the generalization of chebfuns to allow certain functions on infinite intervals or which diverge to infinity; the initial credit for these innovations belongs to Nick Hale, Rodrigo Platte, and Mark Richardson, and many later developments are due to Kuan Xu. For example, here is a function on the whole real axis,

  f = chebfun('exp(-x^2/16).*(1+.2*cos(10*x))',[-inf,inf]);
plot(f) and here is its integral:

  sum(f)
Warning: Result may not be accurate as the function decays slowly at infinity.
ans =
7.089815403622063


Here's the integral of a function on $[1,\infty)$:

  sum(chebfun('1/x^4',[1 inf]))
ans =
0.333333333328451


Notice that several digits of accuracy have been lost here. Be careful! -- operations involving infinities in Chebfun are not always as accurate and robust as their finite counterparts.

Here is an example of a function that diverges to infinity, which we can capture with the 'exps' flag; see Chapter 7 for details:

  h = chebfun('(1/pi)/sqrt(1-x^2)','exps',[-.5 -.5]);
plot(h) In this case the integral comes out just right:

  sum(h)
ans =
1.000000000000000


For more on the treatment of infinities in Chebfun, see Chapter 9.

### 1.6 Periodic functions

Until 2014, Chebfun used only nonperiodic representations, based on Chebyshev polynomials. Beginning with Version 5, a new capability was introduced for representing sufficiently smooth periodic functions by trigonometric polynomials instead, that is, Fourier series. Such an object is still called a chebfun, but it is a periodic one, and the signal to invoke such capabilities is the string trig. For abbreviation, we call a periodic chebfun a "trigfun". This section gives a quick introduction, and more details can be found in Chapter 11.

Trigfuns were initiated by Grady Wright in the first half of 2014 [Wright et al. 2015]. Another project along the same lines was carried out independently by Kristyn McLeod and Rodrigo Platte at Arizona State University [McLeod 2014].

For example, here is a periodic function on $[-\pi,\pi]$ represented in the usual way by a Chebyshev series.

ff = @(t) sin(t) + cos(2*t) - cos(t)/3 + cos(100*t)/6;
f = chebfun(ff,[-pi,pi]);
max(f)
plot(f)
ans =
1.614526099978750 Its length, very roughly, is $100 \pi$,

length(f)
ans =
383


Here is the same function represented by a Fourier series:

f2 = chebfun(ff,[-pi,pi],'trig');
max(f2)
plot(f2,'m')
ans =
1.614526099978745 Its length is now only about $200$ (exactly 201). This improvement by a factor of about $\pi/2$ is typical.

length(f2)
ans =
201


Sampling at a few arbitrary points confirms that the two functions agree closely:

xx = [1/3 sqrt(2) exp(1)];
f(xx) - f2(xx)
ans =
1.0e-14 *
0.355271367880050   0.222044604925031   0.111022302462516


Readers may be interested to compare plotcoeffs applied to the first and second versions of $f$. Rather than display that here we shall turn to a simpler example involving a shorter Fourier series. Consider the function

f = chebfun('7 + sin(t) + exp(1)*cos(2*t)',[-pi,pi],'trig')
f =
chebfun column (1 smooth piece)
interval       length     endpoint values trig
[    -3.1,     3.1]        5       9.7      9.7
vertical scale = 9.7


Here are the coefficients of $f$ as an expansion in sines and cosines:

[a,b] = trigcoeffs(f)
a =
7.000000000000000
0.000000000000001
2.718281828459045
b =
1.000000000000000
0


Here they are as an expansion in complex exponentials:

c = trigcoeffs(f)
c =
1.359140914229522 + 0.000000000000000i
0.000000000000000 + 0.500000000000000i
7.000000000000000 + 0.000000000000000i
0.000000000000000 - 0.500000000000000i
1.359140914229522 + 0.000000000000000i


Bookkeeping of Fourier coefficients can often be a headache. If these examples don't make the patterns clear, details can be found with help trigcoeffs.

For a mathematically less trivial example, here is the cosine expansion of a function whose Fourier series coefficients are known to be values of a Bessel function:

f = chebfun('exp(cos(t))',[-pi pi],'trig');
[a,b] = trigcoeffs(f);
n = floor(length(f)/2);
exact = 2*besseli(0:n,1); exact(1) = exact(1)/2;
disp('        computed             exact')
disp([a exact'])
        computed             exact
1.266065877752008   1.266065877752008
1.130318207984970   1.130318207984970
0.271495339534077   0.271495339534077
0.044336849848664   0.044336849848664
0.005474240442094   0.005474240442094
0.000542926311914   0.000542926311914
0.000044977322954   0.000044977322954
0.000003198436462   0.000003198436462
0.000000199212481   0.000000199212481
0.000000011036772   0.000000011036772
0.000000000550590   0.000000000550590
0.000000000024980   0.000000000024980
0.000000000001039   0.000000000001039
0.000000000000040   0.000000000000040
0.000000000000001   0.000000000000001


### 1.7 Rows, columns, and quasimatrices

MATLAB doesn't only deal with column vectors: there are also row vectors and matrices. The same is true of Chebfun. The chebfuns shown so far have all been in column orientation, which is the default, but one can also take the transpose, compute inner products, and so on:

  x = chebfun(@(x) x)
x =
chebfun column (1 smooth piece)
interval       length     endpoint values
[      -1,       1]        2        -1        1
vertical scale =   1

  x'
ans =
chebfun row (1 smooth piece)
interval       length     endpoint values
[      -1,       1]        2        -1        1
vertical scale =   1

  x'*x
ans =
0.666666666666667


One can also make matrices whose columns are chebfuns or whose rows are chebfuns, like this:

  A = [1 x x.^2]
A =
chebfun column1 (1 smooth piece)
interval       length     endpoint values
[      -1,       1]        3         1        1
vertical scale =   1
chebfun column2 (1 smooth piece)
interval       length     endpoint values
[      -1,       1]        3        -1        1
vertical scale =   1
chebfun column3 (1 smooth piece)
interval       length     endpoint values
[      -1,       1]        3         1        1
vertical scale =   1

  A'*A
ans =
2.000000000000000  -0.000000000000000   0.666666666666667
0   0.666666666666667   0.000000000000000
0.666666666666667                   0   0.400000000000000


These are called quasimatrices, and they are discussed in Chapter 6.

### 1.8 Chebfun features not in this Guide

Some of Chebfun's most remarkable features haven't made it into this edition of the Guide. Here are some of our favorites:

o leg2cheb and cheb2leg for fast Legendre-Chebyshev conversions (also legvals2chebcoeffs, chebcoeffs2legpts, and ten more)

o conv for convolution,

o The 'equi' flag to the Chebfun constructor for equispaced data,

o polyfit for least-squares fitting in the continuous context,

o inv for computing the inverse of a chebfun,

o pde15s for PDEs in one space and one time variable.

To learn about any of these options, try the appropriate help command. Just as a taster, here's a hint of how fast Chebfun can convert a ten-thousand coefficient Chebyshev expansion to Legendre coefficients and back again using an algorithm from [Hale & Townsend 2013]:

tic
ccheb = randn(10000,1);
cleg = cheb2leg(ccheb);
ccheb2 = leg2cheb(cleg);
norm(ccheb-ccheb2,inf)
toc
ans =
1.905142710256769e-13
Elapsed time is 0.131316 seconds.


### 1.9 Chebfun example galleries

MATLAB has long had a gallery command to generate interesting matrices, and Chebfun has an analogous gallery command to generate interesting functions.

Here is what is currently available:

help cheb.gallery
 CHEB.GALLERY   Chebfun example functions.
F = CHEB.GALLERY(NAME) returns a chebfun or a quasimatrix corresponding to
NAME.  See the listing below for available names.

For example,  plot(cheb.gallery('zigzag'))  plots a degree 10000 polynomial
that doesn't look like a polynomial, and  plot(cheb.gallery('gamma'))  shows
a chebfun with poles. For details of how each function is constructed, try
type +cheb/gallery  or  edit cheb.gallery.

[F,FA] = CHEB.GALLERY(NAME) also returns the anonymous function FA used to
define the function. Some gallery functions are generated by operations
beyond the usual Chebfun constructor (e.g. by solving ODEs), so FA in those
cases simply evaluates the chebfun.

CHEB.GALLERY with no input argument returns a random function from the
gallery.

CHEB.GALLERY with no output argument creates a plot of the selected
function.

airy         Airy Ai function on [-40,40]
bessel       Bessel function J_0 on [-100,100]
bump         C-infinity function with compact support
blasius      Blasius function on [0,10]
chirp        Sine with exponentially increasing frequency
daubechies   Approximation to Daubechies phi_2 wavelet scaling function
erf          Error function on [-10,10]
fishfillet   Wild oscillations from Extreme Extrema example
gamma        Gamma function on [-4,4]
gaussian     Gaussian function on [-Inf,Inf]
jitter       A piecewise constant function generated by ROUND
kahaner      Challenging integrand with four spikes
motto        Chebfun motto (Gilbert Strang)
random       Polynomial interpolant through random data in Chebyshev points
rose         A complex-valued sinusoid
runge        Runge function
seismograph  Tanh plus growing oscillation
Si           Sine integral on [-50,50]
sinefun1     As smooth as it looks
sinefun2     Not as smooth as it looks
spikycomb    25 peaks, each sharper than the last
stegosaurus  max(wiggly, x/10)
vandercheb   Chebyshev-Vandermonde quasimatrix
vandermonde  Vandermonde quasimatrix
wiggly       One of the Chebfun team's favorites
wild         An iteratively defined function on [-1 1]
zigzag       Degree 10000 polynomial that looks piecewise linear

Gallery functions are subject to change in future releases of Chebfun.



For example, here is a chebfun representing the Airy function,

plot(cheb.gallery('airy')), ylim([-.8 .8])
title('Airy function') In this instance the underlying code fits in a line,

fa = @airy; p = chebfun(fa, [-40,40]);

Some examples make use of more complicated code, like this approximation to a Daubechies wavelet scaling function (accurate to about 3 digits of accuracy; the underlying function is a fractal):

plot(cheb.gallery('daubechies')), ylim([-0.5 1.5])
title('Daubechies scaling function') To find out how a gallery example was generated, take a look at the code with type +cheb/gallery or edit +cheb/gallery.

Like the MATLAB gallery command, cheb.gallery produces a plot if you call it without specifying output variables. To illustrate, let us finish with an example the Chebfun team enjoys from the appendix to [Trefethen 2013], "Six myths of polynomial interpolation and quadrature":

cheb.gallery('zigzag') This function looks piecewise linear, but in fact, it is a polynomial of degree 10000. This serves no purpose from an approximation point of view -- one would never represent this function in this manner -- but it illustrates the robustness of high-degree polynomial approximation.

If you call cheb.gallery without any input arguments, it selects a gallery function at random.

Other collections worth exploring include cheb.gallerytrig for periodic functions and cheb.gallery2 for 2D functions.

### 1.10 How this Guide is produced

This guide is produced in MATLAB using the publish command with a style sheet somewhat different from the usual; the output of publish is then processed by Markdown. To publish a chapter for yourself, make sure the chebfun guide directory is in your path and then type, for example, open(publish('guide1')). The formatting may not be exactly right but it should certainly be intelligible.

### 1.11 References

[Battles & Trefethen 2004] Z. Battles and L. N. Trefethen, "An extension of MATLAB to continuous functions and operators", SIAM Journal on Scientific Computing, 25 (2004), 1743-1770.

[Berrut & Trefethen 2005] J.-P. Berrut and L. N. Trefethen, "Barycentric Lagrange interpolation", SIAM Review 46, (2004), 501-517.

[Hale & Townsend 2013] N. Hale and A. Townsend, A fast, simple, and stable Chebyshev--Legendre transform using an asymptotic formula, SIAM Journal on Scientific Computing, 36 (2014), A148-A167.

[Higham 2004] N. J. Higham, "The numerical stability of barycentric Lagrange interpolation", IMA Journal of Numerical Analysis, 24 (2004), 547-556.

[McLeod 2014] K. N. McLeod, "Fourfun: A new system for automatic computations using Fourier expansions," SIAM Undergraduate Research Online, 7 (2014), http://dx.doi.org/10.1137/14S013238.

[Pachon, Platte & Trefethen 2010] R. Pachon, R. B. Platte and L. N. Trefethen, "Piecewise-smooth chebfuns", IMA J. Numer. Anal., 30 (2010), 898-916.

[Salzer 1972] H. E. Salzer, "Lagrangian interpolation at the Chebyshev points cos(nu pi/n), nu = 0(1)n; some unnoted advantages", Computer Journal 15 (1972), 156-159.

[Trefethen 2007] L. N. Trefethen, "Computing numerically with functions instead of numbers", Mathematics in Computer Science 1 (2007), 9-19. Revised and reprinted in Communications of the ACM 58 (2014), 91-97.

[Trefethen 2013] L. N. Trefethen, Approximation Theory and Approximation Practice, SIAM, 2013.

[Wright et al. 2015] G. B. Wright, M. Javed, H. Montanelli, and L. N. Trefethen, Extension of Chebfun to periodic functions, SIAM J. Sci. Comp. 37 (2015), C554-C573.