How is the Chebyshev method used by JPL?

Question

In the DE packages, NASA gives us a series of coefficients for the Chebyshev approximation. As I understand it, those approximations are given by a series of polynomials $T_n$ for the interval $[-1,1]$. Then we consider $f(t)$ to be $\sum_i a_i T_i(\tau)$ (where $\tau$ is the normalised time). What is $a_i$ and what are its values?

Then we consider an approximation using a Taylor series defined as $f(t)=\sum_i b_i\tau^i$. What is $b_i$? Is it a Bessel function?

The last two relations represent the Clenshaw algorithm. Here I don't exactly understand how it works. Also does the number of 13 coefficients have something to do with Runge's Phenomenon?

Below are Eqs 3.52 through 3.59, from section 3.3.3 (Chebychev Approximation) in the book Satellite Orbits: Models, Methods and Applications* by De Oliver Montenbruck and Eberhard Gill Springer 2000, ISBN 978-3-642-58351-3 which is viewable in books.google.com

Answer 1

The $a_i$ are the coefficients for the Chebyshev approximation. As you say, NASA gives us those. That is what you find in the DE files, e.g. de430.bsp. (Don't click on that unless you want to download a >100MB file.) NASA/JPL needed a way to provide high resolution and high accuracy functions of time for the positions of the planets, and the most compact way was coefficients of Chebyshev polynomials. The functions are broken up into small intervals, within which the coefficients and resulting functions are valid.

No, the $b_i$ are not Bessel functions. They are the coefficients of a Taylor series that estimates $f$ over the $\tau$ for that interval. That equation was taken out of context, where that context was discussing why Chebyshev polynomials, $T_n\!\left(\tau\right)$ are used for the approximation instead of a Taylor series, which uses the polynomials $\tau^n$. The reason is that a Chebyshev approximation will require fewer terms for the same accuracy.

The Clenshaw algorithm is simply a way to both generate the Chebyshev polynomials and to multiply by the $a_i$ and add them up, all at the same time to minimize the number of operations required. It is well described in the Wikipedia page linked, and can be easily derived from the recurrence relation for the Chebyshev polynomials near the top of your image.

Answer 2

The article Format of the JPL Ephemeris Files has a detailed breakdown of how to use the Development Ephemeris including an example walkthrough and example source code.

I think you have a lot of superfluous information, all you really need is:

$$ \sum_{i=1}^{n}a_i T_i(x) $$

Where $a_i$ are the coefficients, $x$ is the time variable normalized to the interval $[-1,1]$, and $T_i(x)$ is defined as follows:

$T_0(x)=0$

$T_1(x)=x$

$T_{n+1}(x)=2xT_n(x)-T_{n-1}(x)$

Code in Javascript to perform this computation:

function computePolynomial(x,coefficients){
   let T=new Array();

   T[0]=1;
   T[1]=x;
   for(let n=2;n<coefficients.length;n++)  {
     T[n]=2*x*T[n-1] - T[n-2];
   }

  let v=0;
  for(let i=coefficients.length-1;i>=0;i--){
    v+=T[i]*coefficients[i];
  }
  return v;
}   

Note the summation of all of the variables is done in reverse order, from smallest to largest to avoid floating point rounding.

The Github repository gmiller123456/jpl-development-ephemeris has non-optimized source code in several languages which implement the whole process.