13
$\begingroup$

Suppose I have the following interesting function: $$ f(x) = \sum_{k\geq1} \frac{\cos k x}{k^2(2-\cos kx)}. $$ It has some unpleasant properties, like its derivative not being continous at rational multiples of $\pi$. I suspect a closed form does not exist.

I can compute it by computing partial sums and using Richardson extrapolation, but the problem is that it is too slow to compute the function to a good number of decimal digits (100 would be nice, for example).

Is there a method that can handle this function better?

Here's a plot of $f'(\pi x)$ with some artifacts:

Function's derivative, $f'(\pi x)$

Improve this question
$\endgroup$
5
  • 1
    $\begingroup$ Maybe you can use the fact that $cos(kx) = T_k(x)$, where $T_k(x)$ is a Chebyshev polynomial. Then the summation starts to look like a series of rational polynomials. Then if you can turn the series in to a rational polynomial in a Chebyshev basis, it would allow a very efficient way to sum it up. If you aren't familiar with Chebyshev polynomials and basis, the Numerical Recipes in C has a good primer, as well as this: www2.maths.ox.ac.uk/chebfun/ATAP/ATAPfirst6chapters.pdf $\endgroup$ – Jay Lemmon Jan 2 '14 at 23:56
  • 1
    $\begingroup$ er, that should say $cos(kx) = T_k(cos(x))$ $\endgroup$ – Jay Lemmon Jan 3 '14 at 0:03
  • $\begingroup$ @JayLemmon Thanks you for that link. I'll have a look and see if it helps. $\endgroup$ – Kirill Jan 5 '14 at 3:51
  • $\begingroup$ I'm joining this party a bit late, but have you tried using Padé approximants, i.e. the $\varepsilon$-Algorithm instead of Richardson extrapolation? $\endgroup$ – Pedro Feb 4 '14 at 19:57
  • $\begingroup$ By analogy with the case of highly oscillatory integrals, I don't think you'll be able to do a good job without some knowledge of the separation between oscillatory and nonoscillatory parts. If you have such a separation, the Fourier series answer gives you easy exponential convergence. $\endgroup$ – Geoffrey Irving Feb 5 '14 at 2:31
7
$\begingroup$

If analytic techniques are disallowed but the periodic structure is known, here is one approach. Let $$g(x) = \frac{\cos x}{2-\cos x}$$ be periodic with period $2 \pi$, so that $$g(x) = \sum_j w_j e^{ijx}$$ where $$w_j = \frac{1}{2\pi}\int_0^{2\pi} g(x) e^{-ijx} dx$$ Thus, $$\begin{aligned} f(x) &= \sum_{k \ge 1} \frac{g(kx)}{k^p} \\ &= \sum_{k \ge 1} \frac{1}{k^p} \sum_j w_j e^{ijkx} \\ &= \sum_j w_j \sum_{k \ge 1} \frac{\left(e^{ijx}\right)^k}{k^p} \\ &= \sum_j w_j \operatorname{Li}_p(e^{ijx}) \end{aligned}$$ You can either approximate the integrals $w_j$ directly or compute a bunch of $f(x)$ values and use a DFT. In either case, you can potentially apply Richardson extrapolation to the result. Since in your case $g(x)$ is analytic within a neighborhood of $\mathbb{R}$, the final series converges exponentially even without Richardson.

Improve this answer
$\endgroup$
1
  • $\begingroup$ I presume you mean $g(x) = \cos(x) / (2 - \cos(x))$? $\endgroup$ – Geoff Oxberry Feb 4 '14 at 0:51
3
$\begingroup$

For $x = 2\pi a/b$ with $a,b$ integer, we have $$\begin{aligned} f(x) &= \sum_{k \ge 1} \frac{\cos {kx}}{k^2 (2 - \cos {kx})} \\ &= \sum_{k=1}^b \frac{\cos {kx}}{2-\cos {kx}} \sum_{n \ge 0} \frac{1}{(k+bn)^2} \\ &= \sum_{k=1}^b \frac{\cos {kx}}{2-\cos {kx}} \frac{\psi^1(k/b)}{b^2} \end{aligned}$$ where $\psi^1(z)$ is the trigamma function (http://en.wikipedia.org/wiki/Polygamma). Here are plots of the function and its derivative with the artifacts removed: Values and derivatives for the series

Improve this answer
$\endgroup$
2
  • $\begingroup$ Thank you. The problem is that I selected this specific function as a model for another more complicated function that I actually wanted to evaluate, having similar features, but not actually the same. I'm aware of the closed form from this question on MSE. I meant this as a question about summing an infinite series numerically without closed form. $\endgroup$ – Kirill Feb 3 '14 at 20:04
  • $\begingroup$ Maybe my other answer is better then? $\endgroup$ – Geoffrey Irving Feb 4 '14 at 0:29
0
$\begingroup$

How about the Levin u-transform? In addition to Fortan codes, there are several versions in the GSL: `gsl_sum_levin_u*'. Matlab's MuPAD and Maple use this scheme.

Improve this answer
$\endgroup$
1
  • $\begingroup$ I tried it, but I found Richardson extrapolation to be more accurate. $\endgroup$ – Kirill Nov 3 '13 at 1:26

Your Answer

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

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