I hope that numerical questions are also permitted here.

I want to expand a smooth functions $f \in C^{\infty}$in terms of Legendre polynomials. Thus I need to calculate integrals of the form $\int_{-1}^{1} f(x) P_n(x)$, where $n$ becomes sufficiently large (between 40 and 80). In that regime, the Legendre polynomials oscillate heavily, so my question is whether you are aware of a particular method that is good for integrating such things, cause the standard MATLAB method ( probably some low-level Newton Cotes method) cannot do it sufficiently accurate.

Although this is not a pure math question, I think that this question is of some particular value in applications, so please be not to hard with me ;-)

  • 1
    $\begingroup$ Are you aware of anything clever for Fourier series (which suffer from the same oscillation problem)? $\endgroup$ – Igor Rivin Sep 22 '14 at 11:39
  • $\begingroup$ @IgorRivin unfortunately not, sorry. $\endgroup$ – user54300 Sep 22 '14 at 11:46
  • $\begingroup$ See if you can borrow some ideas in chebfun which is centered around doing numerics with Chebyshev polynomials... $\endgroup$ – Suvrit Sep 22 '14 at 22:36
  • 1
    $\begingroup$ I am quite unfamiliar with approximations and the like. But, let me add my "one rupee worth": the integral may be converted to one of the form $\int _0 ^1 f(x) \frac{d^n}{dx^n}(\frac{x^n(1-x)^n}{n!})$. Up to $\pm 1$ this is (by integration by parts), the same as $\int _0 ^1 f^{(n)}(x)\frac{x^n(1-x)^n}{n!}$, which may be more tractable. $\endgroup$ – Venkataramana Sep 23 '14 at 4:51
  • $\begingroup$ @Venkataramana I upvoted, since depending on what $f$ is, this might be helpful to somebody, although this is by no mean a 'general method'. $\endgroup$ – user54300 Sep 23 '14 at 9:27

It seems fundamentally ill-conditioned.

Since $\int_{-1}^{+1} x^rP_n(x),dx=0$ for $r=0,1,\dots,n-1$, your integral is unchanged if you subtract a polynomial of degree $n-1$ from $f(x)$. I'm guessing that if you can do that very accurately for some polynomial that approximates $f(x)$, the resulting integral won't be so ill-conditioned as the original.

| cite | improve this answer | |


@article {MR1078802,
    AUTHOR = {Alpert, Bradley K. and Rokhlin, Vladimir},
     TITLE = {A fast algorithm for the evaluation of {L}egendre expansions},
   JOURNAL = {SIAM J. Sci. Statist. Comput.},
  FJOURNAL = {Society for Industrial and Applied Mathematics. Journal on
              Scientific and Statistical Computing},
    VOLUME = {12},
      YEAR = {1991},
    NUMBER = {1},
     PAGES = {158--179},
      ISSN = {0196-5204},
     CODEN = {SIJCD4},
   MRCLASS = {65D20 (41A50)},
  MRNUMBER = {1078802 (91i:65042)},
       DOI = {10.1137/0912009},
       URL = {http://dx.doi.org/10.1137/0912009},
| cite | improve this answer | |
  • $\begingroup$ so as far as I see this, this paper enables me to evaluate the Legendre polynomials efficiently at the nodes of Chebyshev polynomials. Thus, one could some Gauss quadrature method for example. I will try this. $\endgroup$ – user54300 Sep 22 '14 at 13:09
  • $\begingroup$ @TobiasHurth I find it hard to believe that this algorithm is not already implemented somewhere (in particular, in MATLAB). However, the only thing I could find was this: mathworks.com/matlabcentral/fileexchange/… $\endgroup$ – Igor Rivin Sep 22 '14 at 13:18
  • $\begingroup$ yes, I already looked at this code and as somebody said in the comments, it is highly(!) unaccurate for higher Legendre polynomials near $\pm 1$. $\endgroup$ – user54300 Sep 22 '14 at 13:39

I will make the following assumptions

  • you can evaluate your function wherever your like.
  • you are unconcerned with machine precision (this should not be an issue for the number of coefficients you are after).

Then the solution is to use an n-point Gauss-Legendre quadrature in order to optimally extract n Legendre coefficients.

To obtain this result, the key point is that the function space spanned by the first $n$ Legendre polynomials is simply the space of the zero-th to n-th polynomial moments. More specifically:

$$ \mathbb{P}_n\equiv\mathrm{span}(P_0(x),P_1(x),P_2(x),\ldots,P_n(x)) \equiv \mathrm{span}(1,x,x^2,\ldots,x^n).$$

Hence to perfectly integrate the $\{0,\ldots,n-1\}$ Legendre moments, we are really looking to integrate an order-$2n-2$ polynomial exactly. We are looking for an answer to the following.

Question. Given that $f$ is a polynomial of at most order $n$, is there a numerical integration scheme that will integrate the product $f$ exactly?

In this form, the question has a known, provable solution.

Theorem 1 (Gauss-Legendre Quadrature). The $n$-node Gaussian quadrature scheme, whose nodes are defined at the roots of the order $n$ Legendre polynomial, and whose weights are defined via Lagrange interpolants, will peform the integration exactly$$\int_0^1f(x)\,dx=\sum_{i=1}^m w_i f(x_i)$$for all $f\in\mathbb{P}_{2n-1}$, where $\mathbb{P}_{2n-1}$ denotes the space of order $2n-1$ polynomials.

In terms of where to find code for Gauss-Legendre Quadrature, I would look for it on MATLAB central file exchange.

| cite | improve this answer | |

The LegendreSeries command of the OrthogonalExpansions package of Maple does the job.

| cite | improve this answer | |

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy