1
$\begingroup$

I have an integral of the form $$ I(n) = \int_0^1 dx f(x) \cos(n \pi x) , $$ where $n$ is an integer. In other words, I calculate the cosine Fourier coefficients of function $f$, which is real and continuous on the interval. I need to calculate this for a large number of large values of $n$. Currently I'm just doing this:

iint = array([ integrate.quad(lambda x: f(x)*cos(n*pi*x), 0, 1, limit=1000)[0] for n in range(0,nlim) ])

where $nlim$ is typically of the order of $\sim 10^3$. I have the feeling that this is not the most optimal way of doing this. How do I make my calculation faster?

Improve this question
$\endgroup$
5
  • 2
    $\begingroup$ what is f(x) ....? $\endgroup$
    – Nasser
    Commented Apr 11, 2014 at 11:31
  • $\begingroup$ Is there a reason you don't want to use one of the Fourier Transform libraries? Most of them can compute Discrete Cosine Transforms, too. $\endgroup$
    – Bill Barth
    Commented Apr 11, 2014 at 11:53
  • $\begingroup$ $f(x)$ can be one of several well-behaved (smooth, no fast oscillations, real valued, positive) functions. The FFT library in scipy requires me to input the function at discrete space points. If I want to get the output for large values of $n$, I'd need to first make a very dense grid, calculate this function on this grid and then give it to the function that calculates FFT. I doubt this can be much faster than what I'm doing now, but correct me if I'm wrong. $\endgroup$
    – Echows
    Commented Apr 11, 2014 at 13:02
  • $\begingroup$ It depends, do you need this for only certain values of $n$ or for all $n$ in a range $[1,N]$ where $N$ is big? $\endgroup$
    – Bill Barth
    Commented Apr 11, 2014 at 14:10
  • $\begingroup$ i'd build my own quadrature $\endgroup$
    – Aksakal almost surely binary
    Commented Apr 11, 2014 at 19:02

3 Answers 3

Reset to default
4
$\begingroup$

This is part of the (complex-valued) Fourier transform for which there is (provably) no more efficient way than the Fast Fourier Transform (FFT) if you want to compute the integrals for at least a significant fraction of all frequencies $n$ between zero and the largest frequency you care about.

In other words, while your intuition may tell you that evaluating the function at a large number of points is probably in efficient, the intuition is wrong in this case.

Improve this answer
$\endgroup$
7
  • 3
    $\begingroup$ To amplify Wolfgang: At the very least, you will need a large number of function evaluations to capture the oscillatory nature of $cos(n\pi x)$ if $N$ is large ($O(N)$ of them) and an equivalent number of arithmetic operations on the results. The only thing that using quadrature buys you is the ability not to store the $f(x)$ values at the integration points. By switching to FFTs, you get the advantage of getting all the cosine transform values between 1 and $N$. $\endgroup$
    – Bill Barth
    Commented Apr 12, 2014 at 14:23
  • $\begingroup$ Thats a really good way to think about this. $\endgroup$
    – meawoppl
    Commented Apr 17, 2014 at 16:15
  • $\begingroup$ @BillBarth: Your comment about large numbers of function evaluations isn't necessarily true, e.g. if the integrals are evaluated using Filon transforms. In that case, the accuracy actually increases as the frequency increases. $\endgroup$
    – Lysistrata
    Commented Mar 13, 2016 at 12:47
  • $\begingroup$ @Lysistrata, maybe you could link to an example? $\endgroup$
    – Bill Barth
    Commented Mar 13, 2016 at 16:30
  • $\begingroup$ @BillBarth: I think the best place for examples is in some of the publications and reports by Arieh Iserles and his group at the University of Cambridge. damtp.cam.ac.uk/user/ai/Arieh_Iserles/Publications.html In particular, see: Iserles, A., "On the numerical quadrature of highly-oscillating'integrals I: Fourier transforms", DAMTP Tech. Report NA2003/05, University of Cambridge, UK, June 6, 2003, pp. 26. damtp.cam.ac.uk/user/na/NA_papers/NA2003_05.ps.gz $\endgroup$
    – Lysistrata
    Commented Mar 13, 2016 at 18:16
2
$\begingroup$

Isn't this the same as the real part of the Fourier transform of f(x) $$ \begin{align} I(n) &= \int^1_0dxf(x)\cos(n\pi x) \\ &=\Re{{\int^1_0dx f(x) e^{-j n \pi x}}} \\ &= \Re\int^1_0dx f(x) (cos(n\pi x) - j\sin(n \pi x)) \\ &= \int^1_0dxf(x)\cos(n\pi x) \end{align} $$

Improve this answer
$\endgroup$
2
$\begingroup$

If $f$ is very nice, then you can approximate $f$ with a sequence of piecewise polynomials and integrate over the resulting intervals exactly. This may be much, much cheaper than using an FFT. This is also true for approximations of $f$ by any functions where you can easily know or determine the exact antiderivative in each interval. If $f$ is a black box, then you'll have more trouble with this method.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.