htswanalysis/MACS/lib/gsl/gsl-1.11/doc/sum.texi

   1 @cindex acceleration of series
   2 @cindex summation, acceleration
   3 @cindex series, acceleration
   4 @cindex u-transform for series
   5 @cindex Levin u-transform
   6 @cindex convergence, accelerating a series
   7
   8 The functions described in this chapter accelerate the convergence of a
   9 series using the Levin @math{u}-transform.  This method takes a small number of
  10 terms from the start of a series and uses a systematic approximation to
  11 compute an extrapolated value and an estimate of its error.  The
  12 @math{u}-transform works for both convergent and divergent series, including
  13 asymptotic series.
  14
  15 These functions are declared in the header file @file{gsl_sum.h}.
  16
  17 @menu
  18 * Acceleration functions::
  19 * Acceleration functions without error estimation::
  20 * Example of accelerating a series::
  21 * Series Acceleration References::
  22 @end menu
  23
  24 @node Acceleration functions
  25 @section Acceleration functions
  26
  27 The following functions compute the full Levin @math{u}-transform of a series
  28 with its error estimate.  The error estimate is computed by propagating
  29 rounding errors from each term through to the final extrapolation.
  30
  31 These functions are intended for summing analytic series where each term
  32 is known to high accuracy, and the rounding errors are assumed to
  33 originate from finite precision. They are taken to be relative errors of
  34 order @code{GSL_DBL_EPSILON} for each term.
  35
  36 The calculation of the error in the extrapolated value is an
  37 @math{O(N^2)} process, which is expensive in time and memory.  A faster
  38 but less reliable method which estimates the error from the convergence
  39 of the extrapolated value is described in the next section.  For the
  40 method described here a full table of intermediate values and
  41 derivatives through to @math{O(N)} must be computed and stored, but this
  42 does give a reliable error estimate.
  43
  44 @deftypefun {gsl_sum_levin_u_workspace *} gsl_sum_levin_u_alloc (size_t @var{n})
  45 This function allocates a workspace for a Levin @math{u}-transform of @var{n}
  46 terms.  The size of the workspace is @math{O(2n^2 + 3n)}.
  47 @end deftypefun
  48
  49 @deftypefun void gsl_sum_levin_u_free (gsl_sum_levin_u_workspace * @var{w})
  50 This function frees the memory associated with the workspace @var{w}.
  51 @end deftypefun
  52
  53 @deftypefun int gsl_sum_levin_u_accel (const double * @var{array}, size_t @var{array_size}, gsl_sum_levin_u_workspace * @var{w}, double * @var{sum_accel}, double * @var{abserr})
  54 This function takes the terms of a series in @var{array} of size
  55 @var{array_size} and computes the extrapolated limit of the series using
  56 a Levin @math{u}-transform.  Additional working space must be provided in
  57 @var{w}.  The extrapolated sum is stored in @var{sum_accel}, with an
  58 estimate of the absolute error stored in @var{abserr}.  The actual
  59 term-by-term sum is returned in @code{w->sum_plain}. The algorithm
  60 calculates the truncation error (the difference between two successive
  61 extrapolations) and round-off error (propagated from the individual
  62 terms) to choose an optimal number of terms for the extrapolation.
  63 All the terms of the series passed in through @var{array} should be non-zero.
  64 @end deftypefun
  65
  66
  67 @node Acceleration functions without error estimation
  68 @section Acceleration functions without error estimation
  69
  70 The functions described in this section compute the Levin @math{u}-transform of
  71 series and attempt to estimate the error from the ``truncation error'' in
  72 the extrapolation, the difference between the final two approximations.
  73 Using this method avoids the need to compute an intermediate table of
  74 derivatives because the error is estimated from the behavior of the
  75 extrapolated value itself. Consequently this algorithm is an @math{O(N)}
  76 process and only requires @math{O(N)} terms of storage.  If the series
  77 converges sufficiently fast then this procedure can be acceptable.  It
  78 is appropriate to use this method when there is a need to compute many
  79 extrapolations of series with similar convergence properties at high-speed.
  80 For example, when numerically integrating a function defined by a
  81 parameterized series where the parameter varies only slightly. A
  82 reliable error estimate should be computed first using the full
  83 algorithm described above in order to verify the consistency of the
  84 results.
  85
  86 @deftypefun {gsl_sum_levin_utrunc_workspace *} gsl_sum_levin_utrunc_alloc (size_t @var{n})
  87 This function allocates a workspace for a Levin @math{u}-transform of @var{n}
  88 terms, without error estimation.  The size of the workspace is
  89 @math{O(3n)}.
  90 @end deftypefun
  91
  92 @deftypefun void gsl_sum_levin_utrunc_free (gsl_sum_levin_utrunc_workspace * @var{w})
  93 This function frees the memory associated with the workspace @var{w}.
  94 @end deftypefun
  95
  96 @deftypefun int gsl_sum_levin_utrunc_accel (const double * @var{array}, size_t @var{array_size}, gsl_sum_levin_utrunc_workspace * @var{w}, double * @var{sum_accel}, double * @var{abserr_trunc})
  97 This function takes the terms of a series in @var{array} of size
  98 @var{array_size} and computes the extrapolated limit of the series using
  99 a Levin @math{u}-transform.  Additional working space must be provided in
 100 @var{w}.  The extrapolated sum is stored in @var{sum_accel}.  The actual
 101 term-by-term sum is returned in @code{w->sum_plain}. The algorithm
 102 terminates when the difference between two successive extrapolations
 103 reaches a minimum or is sufficiently small. The difference between these
 104 two values is used as estimate of the error and is stored in
 105 @var{abserr_trunc}.  To improve the reliability of the algorithm the
 106 extrapolated values are replaced by moving averages when calculating the
 107 truncation error, smoothing out any fluctuations.
 108 @end deftypefun
 109
 110
 111 @node Example of accelerating a series
 112 @section Examples
 113
 114 The following code calculates an estimate of @math{\zeta(2) = \pi^2 / 6}
 115 using the series,
 116 @tex
 117 \beforedisplay
 118 $$
 119 \zeta(2) = 1 + 1/2^2 + 1/3^2 + 1/4^2 + \dots
 120 $$
 121 \afterdisplay
 122 @end tex
 123 @ifinfo
 124
 125 @example
 126 \zeta(2) = 1 + 1/2^2 + 1/3^2 + 1/4^2 + ...
 127 @end example
 128
 129 @end ifinfo
 130 @noindent
 131 After @var{N} terms the error in the sum is @math{O(1/N)}, making direct
 132 summation of the series converge slowly.
 133
 134 @example
 135 @verbatiminclude examples/sum.c
 136 @end example
 137
 138 @noindent
 139 The output below shows that the Levin @math{u}-transform is able to obtain an
 140 estimate of the sum to 1 part in
 141 @c{$10^{10}$}
 142 @math{10^10} using the first eleven terms of the series.  The
 143 error estimate returned by the function is also accurate, giving
 144 the correct number of significant digits.
 145
 146 @example
 147 $ ./a.out
 148 @verbatiminclude examples/sum.out
 149 @end example
 150
 151 @noindent
 152 Note that a direct summation of this series would require
 153 @c{$10^{10}$}
 154 @math{10^10} terms to achieve the same precision as the accelerated
 155 sum does in 13 terms.
 156
 157 @node Series Acceleration References
 158 @section References and Further Reading
 159
 160 The algorithms used by these functions are described in the following papers,
 161
 162 @itemize @asis
 163 @item
 164 T. Fessler, W.F. Ford, D.A. Smith,
 165 @sc{hurry}: An acceleration algorithm for scalar sequences and series
 166 @cite{ACM Transactions on Mathematical Software}, 9(3):346--354, 1983.
 167 and Algorithm 602 9(3):355--357, 1983.
 168 @end itemize
 169
 170 @noindent
 171 The theory of the @math{u}-transform was presented by Levin,
 172
 173 @itemize @asis
 174 @item
 175 D. Levin,
 176 Development of Non-Linear Transformations for Improving Convergence of
 177 Sequences, @cite{Intern.@: J.@: Computer Math.} B3:371--388, 1973.
 178 @end itemize
 179
 180 @noindent
 181 A review paper on the Levin Transform is available online,
 182 @itemize @asis
 183 @item
 184 Herbert H. H. Homeier, Scalar Levin-Type Sequence Transformations,
 185 @uref{http://arxiv.org/abs/math/0005209}.
 186 @end itemize