Chebyshev Coefficients for 1/(1+x^2); interpretation of equations in paper

Specific Domains Numerics
1

I need to calculate C2k from the paper, Modification of of Owen’s Method for Computing the Bivariate Normal Integral (DM Borth, 1973).

The relevant bit of math is here:

borthC2k
borthC2k676×528 70.1 KB

Note the sentence,

Letting\ C_{2k} \ be\ the\ sum\ of\ the\ coefficients\ of\ x^{2k} \ from\ the\ various\ Chebyshev\ polynomials ...

Later, a table of the coefficients for m=6 is given:

BorthC2kTable

So I’ve tried to calculate C2k for m=6.

Using ApproxFun 
cf = Fun(x->1/(1+x^2),Chebyshev())
cx = coefficients(cf)

julia> cx
43-element Vector{Float64}:
  0.7071067811865475
  0.0
 -0.24264068711928516
  0.0
  0.04163056034261582
  0.0
 -0.007142674936409816
  ⋮
  0.0
 -4.036230797029855e-15
  0.0
  6.467742735441352e-16
  0.0
 -9.003474632934613e-17

If I try to reproduce the table, I get (note translation from 0 indexing to 1 indexing):

julia> for k in 0:6;println(k,"\t",cx[2k+1]);end
0       0.7071067811865475
1       -0.24264068711928516
2       0.04163056034261582
3       -0.007142674936409816
4       0.0012254892758431872
5       -0.00021026071864915557
6       3.6075036051548e-5

As you see, it does not match!

Any ideas on how I should be computing C2k?

2

The coefficients of the Chebyshev polynomials (what ApproxFun gives you, and the coefficients of S in your equation 9) are not the same as the coefficients of the monic polynomials (your C2k).

Furthermore, it looks like that paper may be using the exact minimax polynomial fit, whereas ApproxFun uses a (very good) approximation thereof. In any case, I would suggest comparing ApproxFun’s coefficients to those of S_2k in equation (9).

In general, there is no need to compute the coefficients of a polynomial in a monic basis. You can evaluate a polynomial directly from its coefficients in a Chebyshev basis using a Clenshaw recurrence.

4 Likes
3

It is possible to compute the monomial expansion from the Chebyshev one if you use high precision (BigFloat). Unfortunately I don’t think this is currently implemented in an easy-to-use way, though it wouldn’t be too hard.

4

I am not sure, but is there a chance you’re mixing the polynomials and the coefficients?
Namely, and I vaguely remember, Chebyshev Polynomials are the basis and the coefficients define the linear combination.

So in @blackeneth’s problem {S}_{2K} \left( x \right) are the Chebyshev Polynomials but the code he calculated are the coefficients of {S}_{2k} \left( x \right) but not itself (As implied from your answer: " and the coefficients of S in your equation 9").

So ApproxFun gives the coefficients for Chebyshev Polynomials while his problem his the coefficients of the polynomial of the form {x}^{2k} . In his problem they derive it by expanding the coefficients of {x}^{2k} from the equation 9.

So, in order to get {C}_{2k} @blackeneth needs to put the definition of Chebyshev Polynomials (Not the coefficients) in (9) and create the form in (10) with some basic algebra (Hopefully).

5

I don’t “get it” :frowning:

6

To illustrate simple example, think you have the polynomial basis:

{s}_{1} \left( x \right) = 2 x, \quad {s}_{2} \left( x \right) = 4 {x}^{2}

Now you have a polynomial which is:

P \left( x \right) = \sum_{i = 1}^{2} {\left( ii + 1 \right)}^{2} {s}_{i} \left( x \right)

You could write:

P \left( x \right) = \sum_{i = 1}^{2} {c}_{i} {x}^{i}

Where {c}_{1} = 8, \; {c}_{2} = 36 .

This is what happens in your problem essentially.

@mbauman , Do you know why the display mode math doesn’t work?

P. S.
I figured the math thing. Discourse require Display Mode to be written:

$$
...LaTeX Code...
$$
7

I’ll leave it to others to discuss the merits of doing this, but these coefficients can found using Polynomials:

using Polynomials
P = ChebyshevT # or Chebyshev in SpecialPolynomials
q = 3 - 2*sqrt(2)
P2m(m) = 8q/(1-q^2) * (1/2 +
                       sum((-1)^k * q^k * Polynomials.basis(P, 2k) for k ∈ 1:m-1) +
                       (-1)^m * q^m/(1-q^2)*Polynomials.basis(P, 2m))

Q2m(m) = convert(Polynomial, P2m(m))

julia> coeffs(Q2m(6))[1:2:end]
7-element Vector{Float64}:
  0.9999936227743424
 -0.9992988806338389
  0.9872975750776766
 -0.9109972523246309
  0.6829098229737781
 -0.33602102198833345
  0.07612251134666158
1 Like