Numerical Summation of a Series 

Produce a table of the values of the series

\begin{displaymath}\psi (x) = \sum_{k=1}^{\infty} {1 \over {k(k+x)}}
\end{displaymath} (1)

for the 3001 values of x, $x= 0.0, 0.1, 0.2, \dots, 300.0$. All entries of the table must have a relative error less than 1.0e-10 (10 digits of precision). This problem is based on a problem from Hamming (1962), when mainframes were very slow by today's microcomputer standards.

Input 

This problem has no input file.

Output 

The output is to be written into a file. The output is to be formatted as two columns with the values of x and y(x) printed as in the C printf or the Pascal writeln.

printf("%6.2f %16.12f\n", x, psix )

writeln(x:6:2, psix:16:2)

As an example, the sample output below shows 4 acceptable lines out of 3001, which might appear in the output file.

The values of x should start at 0.00 and increase by 0.1 until the line with x=300.00 is output.

Process 

The problem with summing the sequence in equation 1 is that too many terms may be required to complete the summation in the given time. Additionally, if enough terms were to be summed, roundoff would render any typical double precision computation useless for the desired precision.


To improve the convergence of the summation process note that


\begin{displaymath}{1 \over {k(k+1)}} = {1 \over k} - {1 \over {k+1}}
\end{displaymath} (2)

which implies y(1)=1.0. One can then produce a series for y(x) - y(1) which converges faster than the original series. This series not only converges much faster, it also reduces roundoff loss.


This process of finding a faster converging series may be repeated again on the second series to produce a third sequence, which converges even more rapidly than the second.


The following equation is helpful in determining how may items are required in summing the series above.


\begin{displaymath}\sum_{k=n}^{\infty} {1 \over {k^r}} < \int_{n-1}^{\infty} {1 \over
{x^r}} dx \quad \mbox{ for } k>1 \mbox{ and } r \ge 1
\end{displaymath} (3)

Sample Output 

  0.00   1.644934066848
  0.10   1.534607244904
...
  1.00   1.000000000000
...
  2.00   0.750000000000
...
300.00   0.020942212934



Miguel Revilla
2000-08-31