%
% LaTeX source file
%
\documentstyle[A4,maple,12pt]{article}
\begin{document}
\font\mengft=cmss9 scaled \magstep1
\def \nzo{\hbox{\mengft{I \hskip -1.7mm N}$_0$}}
\newenvironment{listvar}[3]%
{\begin{list}{}
{\settowidth{\labelwidth}{#2 #1}
\setlength{\leftmargin}{\labelwidth}
\addtolength{\leftmargin}{\labelsep}
\parsep 0.2ex plus 0.1ex minus 0.1ex
\itemsep 0.0ex #3
\renewcommand{\makelabel}[1]{#2 ##1\hfill}}}%
{\end{list}}
\newcounter{saveeqn}
\newcommand{\alpheqn}{\setcounter{saveeqn}
{\value{equation}}
\addtocounter{saveeqn}{1}
\setcounter{equation}{0}
\renewcommand{\theequation}{\mbox{\arabic{saveeqn}\alph{equation}}}}
\newcommand {\reseteqn}{\setcounter{equation}{\value{saveeqn}}
\renewcommand{\theequation}{\arabic{equation}}}
%
% title
%
\title{The {\tt trans} package for transforming series, sequences and
functions}
\author{J. Grotendorst\\
        Zentralinstitut f\"ur Angewandte Mathematik\\
        Forschungszentrum J\"ulich GmbH\\
        Postfach 1913\\
        D-5170 J\"ulich\\
        Fed. Rep. Germany\\
        E-Mail: j.grotendorst@kfa.juelich.de}
\date{}
\maketitle
\begin{abstract}
The {\tt trans} package is a collection of Maple procedures for
transforming symbolic and numerical series (sequences) as well as
functions having formal series expansions.  It
includes most of the currently known transformations, mainly
nonlinear methods, for increasing the rate and expanding the
domain of the convergence of series (sequences). For numerical 
purposes the generation of FORTRAN functions for the resulting rational 
expressions is also supported.
\end {abstract}
\thispagestyle{empty}
\pagebreak
\section{Program description}
The {\tt trans} package is loaded by executing the following Maple statements:

\begin{maple}

> read `trans.m`:
> with(trans);

  [aitken, aitkenhf, eps, epshf, gb, gbhf, lev, levhf, ratgen, ratser, rho,

   rhohf,  rhoit, rhoithf, rich, richhf, sidi, sidihf, theta, thetahf,
 
   thetait, thetaithf]


\end{maple}
\noindent
The program {\bf ratser} is designed for the symbolic transformation
of series expansions with data type ``series''. The program
{\bf ratgen} transforms generalized series expansions, e.g.
asymptotic expansions, with data type \`{}+\`{}.
{\bf ratser} has the following input parameters:
\begin{listvar}{Obligatory:}{\it}{\rm}
\item[Obligatory:]
\item[s]
series expansion with data type ``series''
(has to be the first parameter)
\item[Optional:]
\item[typ]
type of transformation
\item[norder]
order of transformation (default: number of terms of the series)
\item[expts]
list with expansion points of data type ``real constant''
for segmented approximation
(default: [ ], empty list)
\item[switch]
switch for the generation of a FORTRAN function: fort or
nofort (default: nofort)
\item[precision]
type of the FORTRAN function: single or double (default:
single)
\item[evalorder]
computational form for the rational function:
horner or confrac (default: horner)
\item[fn]
filename for the FORTRAN function program. fn has to be of data type
``string'' (default: \`{}\`{}).
\end{listvar}
The formal parameters of procedure {\bf ratgen} are:
\begin{listvar}{Obligatory:}{\it}{\rm}
\item[Obligatory:]
\item[s]
generalized series with data type \`{}+\`{}
(has to be the first parameter)
\item[Optional:]
\item[typ]
type of transformation
\item[norder]
order of transformation (default: number of terms of the series)
\item[switch]
switch for the generation of a FORTRAN function: fort or
nofort (default: nofort)
\item[precision]
type of the FORTRAN function: single or double (default:
single)
\item[evalorder]
computational form for the rational function:
horner (default: \`{}\`{})
\item[fn]
filename for the FORTRAN function program (default:
\`{}\`{}).
\end{listvar}
\pagebreak
The order of the optional parameters in {\bf ratser} and {\bf ratgen} is
arbitrary.
In both procedures the formal input parameter {\it typ}
for the type of transformation can have one of the following values:
\begin{listvar}{thetit}{\tt}{\rm}
\item[u]      {\em u} transformation
\item[t]      {\em t} transformation
\item[d]      {\em d} transformation
\item[v]      {\em v} transformation
\item[su]     {\em su} transformation
\item[st]     {\em st} transformation
\item[sd]     {\em sd} transformation
\item[sv]     {\em sv} transformation
\item[sh]     Shanks transformation
\item[gb]     Germain-Bonne transformation
\item[rh]     $\varrho$ transformation
\item[rhit]   iterated ${\varrho}_{2}$ transformation
\item[r]      Richardson extrapolation process
\item[th]     $\vartheta$ algorithm
\item[thit]   iterated $\vartheta_{2}$ algorithm
\item[at]     iterated Aitken transformation
\end{listvar}
In case of the Richardson extrapolation (r), $\varrho$ transformation
(rh) or iterated $\varrho_{2}$ transformation (rhit)
 the user can specify a list of auxiliary parameters
$x_{n},\ n\in\nzo$
(default: $x_{n} = 1/(n+1),\ n\in\nzo$,
for the Richardson extrapolation process,
and $x_{n} = n+1,\ n\in\nzo$,
for the $\varrho$ transformation and the iterated $\varrho_{2}$
transformation).
Then the first element in the list characterizes the type of the
transformation, i.e. it must be equal to \`{}r\`{}, \`{}rh\`{}
or \`{}rhit\`{}.
The main structure of {\bf ratser}
and {\bf ratgen} is as follows: First it is checked whether forbidden
actual parameters were used as input parameters.
If an error occurs at this level the procedure terminates after an
error message. Then the sequence transformation specified via the
input parameter {\it typ} is performed
by calling the appropriate procedure.
Finally, the resulting expression returned by this subprocedure
is normalized and, if specified by the user,
translated into a FORTRAN function program. The FORTRAN functions
can be generated as single (7 significant decimal digits)
or double precision  (16 significant decimal digits) functions.
In order to optimize the numerical evaluation of
rational functions two computational methods are supported:
(i) numerator and denominator polynomials are evaluated using nested
Horner rules, (ii) the associated continued fraction of the
rational function is evaluated.
Furthermore, procedure {\bf ratser} supports segmented rational
approximation
to improve the quality of the approximation over the whole interval.
With both procedures
{\bf ratser} and {\bf ratgen} there is associated a
remember table, i.e. whenever the value of
one of these procedures is computed, the Maple system makes an entry
in a table ``remembering'' the value of the procedure for the
given arguments. When the procedure is called again with the same
arguments, the Maple system retrieves the result from the procedure's
remember table. In this way, it is possible to avoid redundant
computations.
\pagebreak

\noindent
Examples (see also the help files ?ratser and ?ratgen):
\begin{maple}


> a1 := series(GAMMA(x), x = x0, 3);

  a1 := GAMMA(x0) + Psi(x0) GAMMA(x0) (x - x0)
   
                                               2                    2
      + (1/2 Psi(1, x0) GAMMA(x0) + 1/2 Psi(x0)  GAMMA(x0)) (x - x0)
   
                  3
      + O((x - x0) )
   
> a2 := ratser (a1, sh, horner);
   
                                                              2
  a2 := (- 2 Psi(x0) GAMMA(x0)+(- Psi(1, x0) GAMMA(x0)+Psi(x0)  GAMMA(x0)) x0
   
                                         2
        + (Psi(1, x0) GAMMA(x0) - Psi(x0)  GAMMA(x0)) x)
   
         /                                     2                          2
        /  (- 2 Psi(x0)+(- Psi(1, x0) - Psi(x0) ) x0 + (Psi(1, x0)+Psi(x0) ) x)
       /

   
> ratser (a1, sh, fort, horner, double, [1,2,3]);
   
      function shtr (x)
      double precision shtr, x
      if (x .le. 150.D-2) then
      shtr = (-0.157324813237442D0+0.1311756143040508D1*x)/(-0.823680660
     #8528794D0+0.1978111990655945D1*x)
      elseif (x .gt. 150.D-2 .and. x .le. 250.D-2) then
      shtr = (-0.1777943615884081D1+0.4661874728435735D0*x)/(-0.24929299
     #91902693D1+0.8236806608528794D0*x)
      elseif (x .gt. 250.D-2) then
      shtr = (-0.9515561668645066D0-0.9131937245097873D0*x)/(-0.55849636
     #58050974D1+0.1246464995951347D1*x)
      end if
      end 

> a3 := series (sqrt(sin(x)), x, 8);
   
               1/2         5/2           9/2            13/2      15/2
        a3 := x    - 1/12 x    + 1/1440 x    - 1/24192 x     + O(x    )
   
> ratgen (a3, [rh, 1, 6, 9, 16]);
   
                       1/2          5/2         9/2       13/2
               846720 x    - 34560 x    - 2412 x    - 35 x     
               -----------------------------------------------
                                              2 
                              846720 + 36000 x 

\end{maple} 

\pagebreak
  
\begin{maple}


> a4 := asympt (Psi (2*exp(x))-x, x, 4);
   
                                  1           1            1
                a4 := ln(2) - -------- - ---------- + O(-------)
                              4 exp(x)            2           4
                                         48 exp(x)      exp(x)
   
> ratgen (a4, theta);
   
                                         2
                          48 ln(2) exp(x)  - 12 exp(x) - 1
                     1/48 --------------------------------
                                             2



\end{maple}



\noindent
The programs
\begin{maple}

 aitken, eps, gb, lev, rho, rhoit, rich, sidi, theta, thetait

\end{maple}
represent the Maple implementations of the recursive algorithms
described in section 2.
These functions can be used for the rational approximation of numerical 
series (sequences), i.e. for the purpose of accelerating the convergence of 
slowly convergent series or the summation of divergent series.
The Maple functions 
\begin{maple}

 aitkenhf, epshf, gbhf, levhf, rhohf, rhoithf, richhf, sidihf, thetahf, 

 thetaithf
\end{maple}
are identical to the functions {\tt aitken, ..., thetait} but 
use the hardware floating point system.
The following list describes the obligatory input parameters of the
procedures:
\begin{listvar}{aitken}{\bf}{\it}
\item[eps]    (s\_n, n, arr1)
\item[aitken] (s\_n, n, arr1)
\item[thetit] (s\_n, n, arr1)
\item[theta]  (s\_n, n, arr1, arr2)
\item[rho]    (s\_n, x\_n, n, arr1, arr2)
\item[rhoit]  (s\_n, x\_n, n, arr1, arr2)
\item[rich]   (s\_n, x\_n, n, arr1, arr2)
\item[gb]     (s\_n, a\_n, n, arr1, arr2)
\item[lev]    (s\_n, a\_n, n, arr1, arr2, typ)
\item[sidi]   (s\_n, a\_n, n, arr1, arr2, typ)
\end{listvar}
The meaning of the formal parameters is as follows
\begin{listvar}{array}{\it}{\rm}
\item[s\_n] element $s_{n}$ of the sequence $\{s_n\}_{n \geq 0}$
\item [x\_n] element $x_{n}$ of an auxiliary sequence $\{x_n\}_{n \geq 0}$
\item [a\_n] term $a_{n}={s}_{n}-{s}_{n-1}$ (of a series)
\item [n]    number of the last sequence element
\item [vec1]  a vector to store intermediate
             elements of the transformation table
\item [vec2] another vector to store intermediate
             elements of the transformation table
\item [typ]  type of transformation: u, t, d, v 
\end{listvar}
 
\noindent
These programs are input-directed, i.e.
the elements $s_{n},\ n\in\nzo$,  of the
sequence to be transformed have to be computed and then read in
successively starting with $s_{0}$.
Whenever a new sequence element
$s_{n}$ is computed the Maple procedure has to be called again and
a new approximation to the limit of the sequence is returned
according to the rules defined in section 2. 
Example (see also the help file ?sidi):
\begin{maple}


# Summation of the divergent Euler series,
# alternating n! series
#
Digits := 16:
vec1 := array (-1 .. 30):
vec2 := array (-1 .. 30):
vec3 := array (-1 .. 30):
vec4 := array (-1 .. 30):
sum0 := 0:
printlevel := 0:
for i from 0 to 29 do
   term0   := (-1) ^ i * i!:
   sum0    := sum0 + term0:
   app [1] := i:
   app [2] := evalf (sum0):
   app [3] := evalf (sidi (sum0, term0, i, vec1, vec2, d)):
   app [4] := sidihf (sum0, term0, i, vec3, vec4, d):
   if i >= 23 then
      printlevel := 1
   fi:
   app [m] $ m = 1 .. 4;
   printlevel := 0:
od;

                                  23
         23, -.2477678962998852*10  , .5963473623231885, .5963473624247936
  
                                  24
         24,  .5956716121032509*10  , .5963473623231899, .5963473625137528
   
                                  26
         25, -.1491553843122774*10  , .5963473623231924, .5963473622882467
   
                                  27
         26,  .3883759226953779*10  , .5963473623231936, .5963473597424967
   
                                  29
         27, -.1050049352772297*10  , .5963473623231940, .5963473472547150
   
                                  30
         28,  .2943878510839909*10  , .5963473623231941, .5963473053673139
   
                                  31
         29, -.8547374142655711*10  , .5963473623231941, .5963472184698268
   
#Comparison:
   
> int(exp(-t)/(1+t), t=0 .. infinity);
   
                                  - exp(1) Ei(-1)
   
> evalf (");
                                 .5963473623231941

\end{maple}
\pagebreak
\section{Recursive algorithms}
In this section we give a short description of the algorithms for
computing the various sequence transformations included in the 
package trans.
%
\subsection{Wynn's $\epsilon$ algorithm}
Wynn [1] showed that the Shanks transforms $e_{k}(s_{n})$ [2] can
be computed by an economic computational scheme usually being
called the $\epsilon$ algorithm. The $\epsilon$ algorithm
proceeds as follows. Starting from the boundary values
\alpheqn
\begin{equation}
{\epsilon}_{-1}^{(n)} =  0 \ {\rm and}\  {\epsilon}_{0}^{(n)} = s_{n},
 \ n\in\nzo,
\end{equation}
the nonlinear 4-term recursion
\begin{equation}
\epsilon_{k+1}^{(n)} = \epsilon_{k-1}^{(n+1)} + \frac{1}
                       {\epsilon_{k}^{(n+1)}-\epsilon_{k}^{(n)}},
                       \ k,n\in\nzo,
\end{equation}
\reseteqn
is used to construct the 2-dimensional array
$\{\epsilon_{k}^{(n)}\}_{k,n \geq 0}$,
called $\epsilon$ table, with
\begin{equation}
\epsilon_{2k}^{(n)} = e_{k}(s_{n}), \ k,n\in\nzo,
\end{equation}
and
\begin{equation}
\epsilon_{2k+1}^{(n)} = 1/e_{k}(\Delta s_{n}),\ k,n\in\nzo.
\end{equation}
\par\noindent
The computational tableau of algorithm (10)
can be displayed as follows:
\stepcounter{equation}
\leftline{}
$$\begin{array}{ccccccccccc}
& &        \epsilon_{0}^{(0)}   &  &      &          & & & & &\\
&                 &     & \searrow &  & &        &  & & &\\
\epsilon_{-1}^{(1)} &  & \longrightarrow &  & \epsilon_{1}^{(0)}
& & & & & &\\
&       &  & \nearrow  &        & \searrow  &       & & & &\\
& & \epsilon_{0}^{(1)} & & \longrightarrow &  & \epsilon_{2}^{(0)}
& & & &\\
&    &   & \searrow &  & \nearrow   &  &\searrow & & & \\
\epsilon_{-1}^{(2)} & & \longrightarrow & & \epsilon_{1}^{(1)}&
&\longrightarrow & &\epsilon_{3}^{(0)}        & & \\
& & & \nearrow & &   \searrow &    &   \nearrow &       & &\\
& &\epsilon_{0}^{(2)} & &  \longrightarrow & & \epsilon_{2}^{(1)}&
 & & &\epsilon_{4}^{(0)}  \\
&  & &          \searrow  & & \nearrow  &  &  & & & \vdots \\
\epsilon_{-1}^{(3)} & & \longrightarrow & & \epsilon_{1}^{(2)}&
& & &\epsilon_{3}^{(1)}     & &\\
&  & &     \nearrow  &  &  & &   & \vdots     &  & \\
& & \epsilon_{0}^{(3)} &  & & & \epsilon_{2}^{(2)}& & & & \\
&     & &            &     & &  \vdots     &  & & &\\
\epsilon_{-1}^{(4)} &  & & & \epsilon_{1}^{(3)}&  & & & & &\\
\vdots   &  & &     & \vdots &          &  & & & & \\
& & \epsilon_{0}^{(4)}&  & & & & &\\
& & \vdots &  &  & & & & \\
\end{array}$$
Because of the relations in (3) and (4)
only the elements $\epsilon_{2k}^{(n)}$ with even subscripts serve as
approximations to the limit whereas the elements $\epsilon_{2k+1}^
{(n)}$ with odd subscripts are only auxiliary quantities.
Therefore, if $n$, the index of the last sequence element $s_{n}$,
is even, $n=2\mu$, our approximation to the limit of the sequence
is the transformation
\begin{equation}
\{s_{0},s_{1}, ..., s_{2\mu}\} \longrightarrow
\epsilon_{2\mu}^{(0)}=e_{\mu}(s_{0}),
\end{equation}
and if $n$ is odd, $n=2\mu + 1$, we use the transformation
\begin{equation}
\{s_{1},s_{2}, ... ,s_{2\mu +1}\} \longrightarrow
\epsilon_{2\mu}^{(1)}=e_{\mu}(s_{1}).
\end{equation}
%
\subsection{Aitken's iterated $\Delta^{2}$ process}
The special Shanks transformation
\begin{equation}
e_{1}(s_{n})=s_{n}-\frac{{(\Delta s_{n})}^{2}}{{\Delta}^{2}s_{n}}
            ={\epsilon}_{2}^{(n)}, \ n\in\nzo,
\end{equation}
represents Aitken's famous $\Delta^{2}$ process [3].
Identifying the sequence elements $s_{n}$ with the
starting values ${\cal A}_{0}^{(n)}$ and then iterating
eq. (8) yields the following nonlinear recursive scheme:
\alpheqn
\begin{eqnarray}
{\cal A}_{0}^{(n)}   & = & s_{n}, \ n\in\nzo,\\
{\cal A}_{k+1}^{(n)} & = & {\cal A}_{k}^{(n)}-
                         \frac{{[\Delta {\cal A}_{k}^{(n)}]}^{2}}
                         {\Delta^{2}{\cal A}_{k}^{(n)}}, \ k,n\in\nzo.
\end{eqnarray}
\reseteqn
If we arrange the elements ${\cal A}_{k}^{(n)}$ of the 2-dimensional
array $\{{\cal A}_{k}^{(n)}\}_{k,n \geq 0}$
in a rectangular scheme,
\stepcounter{equation}
\leftline{}
$$\begin{array}{ccccc}
{\cal A}_{0}^{(0)} & {\cal A}_{1}^{(0)} & {\cal A}_{2}^{(0)} &
{\cal A}_{3}^{(0)} & \ldots \\
 &      &   &       &   \\
{\cal A}_{0}^{(1)} & {\cal A}_{1}^{(1)} & {\cal A}_{2}^{(1)} &
{\cal A}_{3}^{(1)} & \ldots \\
&     &           &     &    \\
{\cal A}_{0}^{(2)} & {\cal A}_{1}^{(2)} & {\cal A}_{2}^{(2)} &
{\cal A}_{3}^{(2)} & \ldots \\
&    &            &     &    \\
{\cal A}_{0}^{(3)} & {\cal A}_{1}^{(3)} & {\cal A}_{2}^{(3)} &
{\cal A}_{3}^{(3)} & \ldots \\
\vdots   &  \vdots & \vdots &   \vdots & \ddots
\end{array}\eqno(\theequation)$$
then the 4 elements, which are connected by the recurrence relation
(9b), form the following pattern:
\stepcounter{equation}
\leftline{}
$$\begin{array}{cc}
{\cal A}_{k}^{(n)}   & {\cal A}_{k+1}^{(n)} \\
 &           \\
{\cal A}_{k}^{(n+1)}   &  \\
&         \\
{\cal A}_{k}^{(n+2)}   & \\
\end{array}\eqno(\theequation)$$
It follows from the recursive algorithm (9) that the computation
of ${\cal A}_{k}^{(0)}$ requires the sequence elements $s_{0},s_{1},
...s_{2k}$. Therefore, as in Wynn's $\epsilon$ algorithm,
we distinguish two distinct approximations for the limit
of the sequence depending upon the index $n$ of the last
sequence element $s_{n}$ which was used by the algorithm. For
$n=2\mu$ we use
\begin{equation}
\{s_{0},s_{1}, ..., s_{2\mu}\} \longrightarrow {\cal A}_{2\mu}^{(0)}
\end{equation}
and for $n=2\mu+1$ we use
\begin{equation}
\{s_{1},s_{2}, ... ,s_{2\mu +1}\} \longrightarrow
{\cal A}_{2\mu}^{(1)}.
\end{equation}
%
\subsection{The $\varrho$ algorithm}
The $\varrho$ algorithm looks very much like
Wynn's $\epsilon$ algorithm:
With the initial values
\alpheqn
\begin{equation}
{\varrho}_{-1}^{(n)} =  0 \ {\rm and}\  {\varrho}_{0}^{(n)} = s_{n},
 \ n\in\nzo,
\end{equation}
the nonlinear 4-term recursion
\begin{equation}
\varrho_{k+1}^{(n)} = \varrho_{k-1}^{(n+1)} + \frac{x_{n+k+1}-x_{n}}
                       {\varrho_{k}^{(n+1)}-\varrho_{k}^{(n)}},
                       \ k,n\in\nzo,
\end{equation}
\reseteqn
is used to construct the $\varrho$ table, a 2-dimensional array
$\{\varrho_{k}^{(n)}\}_{k,n \geq 0}$ of transforms.
Algorithm (14) is a rhombus rule
and yields exactly the same computational
scheme as for the $\epsilon$ algorithm.
The only difference between the $\epsilon$ algorithm and the
$\varrho$ algorithm is that the $\varrho$ algorithm also involves
a sequence of interpolation points ${\{x_{n}\}}_{n \geq 0}$.
Here also, only the even order transforms $\varrho_{2k}^{(n)}$
are meaningful as approximations to the limit. Thus, we distinguish
the two transformations:
\begin{equation}
\{s_{0}, x_{0}; s_{1}, x_{1};...; s_{2\mu}, x_{2\mu}\} \longrightarrow
\varrho_{2\mu}^{(0)}
\end{equation}
for $n$ even, $n=2\mu$, and
\begin{equation}
\{s_{1},x_{1}; s_{2},x_{2};... ;s_{2\mu +1}, x_{2\mu +1}\}\longrightarrow
\varrho_{2\mu}^{(1)}
\end{equation}
for $n$ odd, $n=2\mu+1$.
With the choice $x_{n}=n+1,\ n\in\nzo$, the $\varrho$ transformation
assumes its original form  due to Wynn [4]:
\alpheqn
\begin{equation}
{\varrho}_{-1}^{(n)}  =  0, {\varrho}_{0}^{(n)} = s_{n},\ n\in\nzo,
\end{equation}
\begin{equation}
\varrho_{k+1}^{(n)}  =  \varrho_{k-1}^{(n+1)} + \frac{k+1}
                       {\varrho_{k}^{(n+1)}-\varrho_{k}^{(n)}},
                       \ k,n\in\nzo ,
\end{equation}
\reseteqn
which is very powerful
for logarithmically convergent sequences, i.e. sequences
${\{s_{n}\}}_{n=0}^{\infty}$ satisfying
\begin{equation}
\lim_{n \rightarrow \infty} \frac{s_{n+1}-s} {s_{n}-s} = 1.
\end{equation}
%
\subsection{The iterated $\varrho_{2}$ algorithm}
First we have to express the transform $\varrho_{2}^{(n)}$
in terms of some sequence elements $s_{n}$ and interpolation
points $x_{n}$. The resulting expression for $\varrho_{2}^{(n)}$
will then be iterated. For $k=2$ we get from the $\varrho$ algorithm
(14):
\begin{equation}
\varrho_{2}^{(n)}=s_{n+1} +
\frac{(x_{n+2}-x_{n})[\Delta s_{n+1}] [\Delta s_{n}]}
{[\Delta x_{n+1}][\Delta s_{n}]-[\Delta x_{n}][\Delta s_{n+1}]},
 \ n\in\nzo.
\end{equation}
Iterating this expression yields the following recursive algorithm
(see ref. [5], p. 234):
\alpheqn
\begin{eqnarray}
{\cal W}_{0}^{(n)} & = & s_{n},  \ n\in\nzo,\\
{\cal W}_{k+1}^{(n)} & = &
\frac{(x_{n+2k+2}-x_{n})[\Delta {\cal W}_{k}^{(n+1)}]
                        [\Delta {\cal W}_{k}^{(n)}]}
{(x_{n+2k+2}-x{n})[\Delta {\cal W}_{k}^{(n)}]-
(x_{n+2k+1}-x_{n})[\Delta {\cal W}_{k}^{(n+1)}]},\nonumber\\
 &   & \ k,n\in\nzo.
\end{eqnarray}
\reseteqn
With the choice $x_{n}=n+1,\ n\in\nzo$, the ${\cal W}$ transformation
assumes its standard form:
\alpheqn
\begin{eqnarray}
{\cal W}_{0}^{(n)} & = & s_{n},  \ n\in\nzo,\\
{\cal W}_{k+1}^{(n)} & = &
\frac{(2k+2)[\Delta {\cal W}_{k}^{(n+1)}] [\Delta {\cal W}_{k}^{(n)}]}
{(2k+1){\Delta}^{2} {\cal W}_{k}^{(n)}}, \ k,n\in\nzo.
\end{eqnarray}
\reseteqn
The recurrence formulas for Aitken's iterated ${\Delta}^{2}$ process
and the ${\cal W}_{k}^{(n)}$ transforms are structually identical.
Hence, the computational scheme for the transforms ${\cal W}_{k}^{(n)}$
is the same as for the transforms ${\cal A}_{k}^{(n)}$.
The approximation
to the limit depends upon the index $n$ of the last sequence element
$s_{n}$ which was used in the recursion. We distinguish
the two transformations:
\begin{equation}
\{s_{0}, x_{0}; s_{1}, x_{1};...; s_{2\mu}, x_{2\mu}\} \longrightarrow
{\cal W}_{2\mu}^{(0)}
\end{equation}
for $n$ even, $n=2\mu$, and
\begin{equation}
\{s_{1},x_{1}; s_{2},x_{2};... ;s_{2\mu +1}, x_{2\mu +1}\}\longrightarrow
{\cal W}_{2\mu}^{(1)}
\end{equation}
for $n$ odd, $n=2\mu+1$.
%
\subsection{Brezinski's $\vartheta$ algorithm}
It is well known that Wynn's $\epsilon$ algorithm accelerates linear
convergence quite efficiently and in addition is able to sum divergent
series. However, the $\epsilon$ algorithm is not able to accelerate
logarithmic convergence. In order to have a sequence transformation
capable to accelerate  linearly as well as logarithmically convergent
sequences (series) Brezinski [6] constructed a new algorithm, called
$\vartheta$ algorithm, by modifying the $\epsilon$ algorithm in a special
way. The defining equations of the recursive $\vartheta$ scheme are:
\alpheqn
\begin{eqnarray}
{\vartheta}_{-1}^{(n)} & = &  0, {\vartheta}_{0}^{(n)} = s_{n},\\
\vartheta_{2k+1}^{(n)} & = & \vartheta_{2k-1}^{(n+1)} + 1/
                         [\Delta \vartheta_{2k}^{(n)}],\\
\vartheta_{2k+2}^{(n)} & = & \vartheta_{2k}^{(n+1)} +
                             \frac{[\Delta \vartheta_{2k}^{(n+1)}]
                            [\Delta \vartheta_{2k+1}^{(n+1)}]}
                            {{\Delta}^{2} \vartheta_{2k+1}^{(n)}},
                            \ k,n\in\nzo.
\end{eqnarray}
\reseteqn
If we arrange the elements $\vartheta_{k}^{(n)}$ in a rectangular
scheme
\stepcounter{equation}
\leftline{}
$$\begin{array}{ccccc}
{\vartheta}_{0}^{(0)} & {\vartheta}_{1}^{(0)} &
{\vartheta}_{2}^{(0)} & {\vartheta}_{3}^{(0)} & \ldots \\
 &      &   &       &   \\
{\vartheta}_{0}^{(1)} & {\vartheta}_{1}^{(1)} &
{\vartheta}_{2}^{(1)} & {\vartheta}_{3}^{(1)} & \ldots \\
&     &           &     &    \\
{\vartheta}_{0}^{(2)} & {\vartheta}_{1}^{(2)} &
{\vartheta}_{2}^{(2)} & {\vartheta}_{3}^{(2)} & \ldots \\
&    &            &     &    \\
{\vartheta}_{0}^{(3)} & {\vartheta}_{1}^{(3)} &
{\vartheta}_{2}^{(3)} & {\vartheta}_{3}^{(3)} & \ldots \\
\vdots   &  \vdots & \vdots &   \vdots & \ddots
\end{array}\eqno(\theequation)$$
then the geometric pattern, which represents the recurrence
relation (24b), has the form
\leftline{}
\stepcounter{equation}
$$\begin{array}{ccc}
& {\vartheta}_{2k}^{(n)}   & {\vartheta}_{2k+1}^{(n)} \\
 & &           \\
{\vartheta}_{2k-1}^{(n+1)}   & {\vartheta}_{2k}^{(n+1)}  \\
\end{array}\eqno(\theequation)$$
and eq. (24c) is represented by
\stepcounter{equation}
\leftline{}
$$\begin{array}{ccc}
& {\vartheta}_{2k+1}^{(n)}   & {\vartheta}_{2k+2}^{(n)} \\
 & &           \\
{\vartheta}_{2k}^{(n+1)}   & {\vartheta}_{2k+1}^{(n+1)}  \\
 & &           \\
{\vartheta}_{2k}^{(n+2)}   & {\vartheta}_{2k+1}^{(n+2)}  \\
\end{array}\eqno(\theequation)$$
As with the $\epsilon$ algorithm, only transforms
$\vartheta_{2k}^{(n)}$ with even numbered subscripts serve as
approximations. The two geometric pattern (26) and (27) and the
corresponding eqs. (24b) and (24c) show that for the
computation of the transform $\vartheta_{2k}^{(0)}$ the sequence
elements $s_{0},s_{1}, ...s_{3k}$ have to be known. Hence,
our choice of approximants has to depend upon the index $n$ of the
last sequence element $s_{n}$ which was used in the recursion.
If $n$ is a multiple of 3, $n=3\mu$, we use the transformation
\begin{equation}
\{s_{0},s_{1}, ..., s_{3\mu}\} \longrightarrow
\vartheta_{2\mu}^{(0)},
\end{equation}
if $n=3\mu + 1$, we use the transformation
\begin{equation}
\{s_{1},s_{2}, ... ,s_{3\mu +1}\} \longrightarrow
\vartheta_{2\mu}^{(1)},
\end{equation}
and if $n=3\mu + 2$, we use
\begin{equation}
\{s_{2},s_{3}, ... ,s_{3\mu +2}\} \longrightarrow
\vartheta_{2\mu}^{(2)}.
\end{equation}
%
\subsection{The iterated $\vartheta_{2}$ algorithm}
>From the recursive scheme (24) we obtain
the following expression for $\vartheta_{2}^{(n)}$:
\begin{equation}
\vartheta_{2}^{(n)}=s_{n+1} -
\frac{[\Delta s_{n}] [\Delta s_{n+1}][{\Delta}^{2} s_{n+1}]}
{[\Delta s_{n+2}][{\Delta}^{2} s_{n}]-[\Delta s_{n}]
[{\Delta}^{2} s_{n+1}]},  \ n\in\nzo.
\end{equation}
Iterating eq. (31) we get the following recursive scheme (see ref. [5],
p. 282):
\alpheqn
\begin{eqnarray}
{\cal T}_{0}^{(n)}  & = &  s_{n},  \ n\in\nzo,\\
{\cal T}_{k+1}^{(n)} & = & {\cal T}_{k}^{(n+1)} -
\frac{[\Delta {\cal T}_{k}^{(n)}] [\Delta {\cal T}_{k}^{(n+1)}]
[{\Delta}^{2} {\cal T}_{k}^{(n+1)}]}
{[\Delta {\cal T}_{k}^{(n+2)}][{\Delta}^{2} {\cal T}_{k}^{(n)}]-
[\Delta {\cal T}_{k}^{(n)}][{\Delta}^{2} {\cal T}_{k}^{(n+1)}]},
\nonumber\\
 &   & \ k,n\in\nzo.
\end{eqnarray}
\reseteqn
The 5 elements, which are connected by the recurrence relation
(32b), form the following pattern in a rectangular scheme:
\stepcounter{equation}
\leftline{}
$$\begin{array}{cc}
{\cal T}_{k}^{(n)}   & {\cal T}_{k+1}^{(n)} \\
 &           \\
{\cal T}_{k}^{(n+1)}   &  \\
&         \\
{\cal T}_{k}^{(n+2)}   & \\
&         \\
{\cal T}_{k}^{(n+3)}   & \\
\end{array}\eqno(\theequation)$$
As in the case of the $\vartheta$ algorithm the computation
of the ${\cal T}_{k}^{(n)}$ transforms require the knowlegde of
the sequence elements
$s_{n},s_{n+1}, ...s_{n+3k}$. Therefore,
we distinguish the following approximations:
\begin{equation}
\{s_{0},s_{1}, ..., s_{3\mu}\} \longrightarrow
{\cal T}_{2\mu}^{(0)},
\end{equation}
for $n=3\mu$,
\begin{equation}
\{s_{1},s_{2}, ... ,s_{3\mu +1}\} \longrightarrow
{\cal T}_{2\mu}^{(1)},
\end{equation}
for $n=3\mu + 1$ and
\begin{equation}
\{s_{2},s_{3}, ... ,s_{3\mu +2}\} \longrightarrow
{\cal T}_{2\mu}^{(2)}.
\end{equation}
for $n=3\mu + 2$.
%
\subsection{The Richardson extrapolation process}
The linear Richardson extrapolation for sequences (series) [7]
can be computed by the following variant of Neville's recursive
scheme (see p. 6 of ref. [8]):
\alpheqn
\begin{eqnarray}
{\cal N}_{0}^{(n)} & = &  s_{n},  \ n\in\nzo, \\
{\cal N}_{k+1}^{(n)} & = & \frac{x_{n} {\cal N}_{k}^{(n+1)}-
                        x_{n+k+1} {\cal N}_{k}^{(n)}}
                       {x_{n} - x_{n+k+1}}, \ k,n\in\nzo
\end{eqnarray}
\reseteqn
For the special auxiliary sequence $x_{n} = \Delta s_{n}, \ n\in\nzo$,
we obtain the sequence transformation of Germain-Bonne [9].
In this case the recursive scheme (37) assumes the form
\alpheqn
\begin{eqnarray}
{\cal B}_{0}^{(n)} & = & s_{n},  \ n\in\nzo, \\
{\cal B}_{k+1}^{(n)}& = &\frac{\Delta s_{n} {\cal B}_{k}^{(n+1)}-
                        \Delta s_{n+k+1} {\cal B}_{k}^{(n)}}
                       {\Delta s_{n} - \Delta s_{n+k+1}},
                       \ k,n\in\nzo .
\end{eqnarray}
\reseteqn
The computational tableau of the deltoid scheme (38) can be
displayed as follows:
\stepcounter{equation}
\leftline{}
$$\begin{array}{ccccccccc}
{\cal N}_{0}^{(0)}   &  &   &  & & & & &\\
& \searrow &  & &  &  & & &\\
&  & {\cal N}_{1}^{(0)} & & & & & &\\
& \nearrow  &        & \searrow  &       & & & &\\
{\cal N}_{0}^{(1)} & &   &  & {\cal N}_{2}^{(0)} & & & &\\
& \searrow &  & \nearrow   &  &\searrow & & & \\
& & {\cal N}_{1}^{(1)}& &  & &{\cal N}_{3}^{(0)}& & \\
& \nearrow & &   \searrow &    &   \nearrow & \vdots   & &\\
{\cal N}_{0}^{(2)} & &    & & {\cal N}_{2}^{(1)} &  & & & \\
&          \searrow  & & \nearrow  & \vdots &  & & & \\
& & {\cal N}_{1}^{(2)}& & & &    & &\\
&     \nearrow  & \vdots &  & &   &  &  & \\
{\cal N}_{0}^{(3)} &  & & & & & & & \\
\vdots  &  &  & &   &  &   &  & \\
\end{array}\eqno(\theequation)$$
As approximation to the limit of the sequence we use
\begin{equation}
\{s_{1},x_{1}; s_{2},x_{2};... ;s_{n}, x_{n}\} \longrightarrow
{\cal N}_{n}^{(0)},
\end{equation}
and in case of Germain-Bonne's transformation we use
\begin{equation}
\{s_{1},\Delta s_{1}; s_{2},\Delta s_{2};...;
s_{n}, \Delta s_{n}\} \longrightarrow  {\cal B}_{n}^{(0)}.
\end{equation}
%
\subsection{Recursive evaluation of the Levin transformation}
The general sequence transformation of Levin [10] is given by
the following ratio of finite sums:
\begin{equation}
{\cal L}_{k}(s_{n}) =
\frac
{\sum\limits_{j=0}^{k} (-1)^{j} {k \choose j}
(n+j+1)^{k-1}
\frac{s_{n+j}}{R_{n+j}}}
{\sum\limits_{j=0}^{k} (-1)^{j} {k \choose j}
(n+j+1)^{k-1}
\frac{1}{R_{n+j}}}.
\end{equation}
Here, ${\{R_{n}\}}_{n \geq 0}$ is a sequence of remainder estimates
of a given sequence ${\{s_{n}\}}_{n \geq 0}$.
The numerator and denominator sums in eq. (42) can be computed
recursively with the help of a linear 3-term recurrence relation
of Fessler $et \ al$. [11]
\begin{equation}
L_{k+1}^{(n)} = L_{k}^{(n+1)}-
                \frac {(n+1) {(n+k+1)}^{k-1}}{(n+k+2)^{k}}
                L_{k}^{(n)}, \ k,n\in\nzo .
\end{equation}
The starting values for this recursive formula (43) are
\begin{equation}
L_{0}^{(n)} =  s_{n}/R_{n},  \ n\in\nzo,
\end{equation}
for the numerator in eq. (42) and
\begin{equation}
L_{0}^{(n)} =  1/R_{n},  \ n\in\nzo,
\end{equation}
for the denominator in eq. (42).
The computational tableau produced by the deltoid algorithm (43),
which either represent numerators or denominators of the general
Levin transformation (42), is exactly the same as in the case
of the Richardson extrapolation process. Our approximation
to the limit of the sequence ${\{s_{n}\}}_{n \geq 0}$
will be:
\begin{equation}
\{s_{0}, R_{0}; s_{1}, R_{1};...; s_{n}, R_{n}\} \longrightarrow
{\cal L}_{n}(s_{0}).
\end{equation}
%
\subsection{Recursive evaluation of the Sidi transformation}
Formally the transformation of Sidi can be obtained
by replacing the powers $(n+j+1)^{k-1}$ in the explicit expression
(42) for Levin's general sequence transformation by the Pochhammer
symbols $(n+j+1)_{k-1}$
\begin{equation}
{\cal S}_{k}(s_{n}) =
\frac
{\sum\limits_{j=0}^{k} (-1)^{j} {k \choose j}
(n+j+1)_{k-1}
\frac{s_{n+j}}{R_{n+j}}}
{\sum\limits_{j=0}^{k} (-1)^{j} {k \choose j}
(n+j+1)_{k-1}
\frac{1}{R_{n+j}}}.
\end{equation}
Here also, ${\{R_{n}\}}_{n \geq 0}$ is a sequence of remainder
estimates of a given sequence ${\{s_{n}\}}_{n \geq 0}$.
Transformation (47) was first used by Sidi [12] for deriving Pad\'{e}
approximants of some hypergeometric functions. Later Weniger [5]
introduced the transformation (47) as a nonlinear sequence transformation
for the acceleration of convergence and the summation of divergent
series.
The recursive computation of the numerator
and denominator sums in eq. (47) is performed
using the following linear 3-term recurrence relation
of Weniger (see ref. [5], p. 259):
\begin{equation}
S_{k+1}^{(n)} = S_{k}^{(n+1)}-
                \frac {(n+k)(n+k+1)}{(n+2k+1)(n+2k)}
                S_{k}^{(n)}, \ k,n\in\nzo .
\end{equation}
The starting values for the recursive formula in eq. (48) are
\begin{equation}
S_{0}^{(n)} =  s_{n}/R_{n},  \ n\in\nzo,
\end{equation}
in case of the numerator and
\begin{equation}
S_{0}^{(n)} =  1/R_{n},  \ n\in\nzo,
\end{equation}
in case of the denominator in eq. (47).
The computaional scheme for the deltoid algorithm (51)
is the same as for the general Levin transformation.
We use the following approximation to the limit of
the sequence ${\{s_{n}\}}_{n \geq 0}$:
\begin{equation}
\{s_{0}, R_{0}; s_{1}, R_{1};...; s_{n}, R_{n}\} \longrightarrow
{\cal S}_{n}(s_{0}).
\end{equation}
\pagebreak
\parskip=0pt plus 1pt
\begin{thebibliography}{wimp}
%\bibitem{jg1} J. Grotendorst, Comput. Phys. Commun. 59 (1990) 289.
%\bibitem{cb5} C. Brezinski, Numer. Math. 35 (1980) 176.
%\bibitem{th} T. H\aa vie, BIT 19 (1979) 204.
\bibitem{pw2} P. Wynn, Math. Tables Aids Comput. 10 (1956) 91.
\bibitem{ds} D. Shanks, J. Math. and Phys.
 (Cambridge, Mass.) 34 (1955) 1.
\bibitem{aca} A. C. Aitken, Proc. Roy. Soc. Edinburgh 46 (1926) 289.
\bibitem{pw1} P. Wynn, Proc. Camb. Phil. Soc. 52 (1956) 633.
\bibitem{ejw1} E. J. Weniger, Nonlinear Sequence Transformations for the
 Acceleration of Convergence and the Summation of Divergent Series,
 Comput. Phys. Rep. 10 (1989) 189.
\bibitem{cb4} C. Brezinski, C. R. Acad. Sc. Paris 273 (1971) 727.
\bibitem{lfr} L. F. Richardson, Phil. Trans. Roy. Soc. London A 226
 (1927) 229.
\bibitem{cb2} C. Brezinski, Algorithmes d'Acc\'{e}l\'{e}ration
 de la Convergence - \'{E}tude Num\'{e}rique (\'{E}ditions
 Technip, Paris, 1978).
\bibitem{bgb} B. Germain-Bonne, Estimation de la Limite des Suites et
 Formalisation de Pr\'{e}c\'{e}des d'Acc\'{e}l\'{e}ration
 de Convergence (Thes\`{e}, Universit\'{e} de Lille, 1978).
\bibitem{dl} D. Levin, Int. J. Comput. Math. B3 (1973) 371.
\bibitem{tf} T. Fessler, W. F. Ford and D. A. Smith, ACM Trans.
 Math. Software 9 (1983) 346.
\bibitem{as} A. Sidi, J. Comput. Appl. Math. 7 (1981) 37.
%\bibitem{bwc} B. W. Char, K. O. Geddes, G. H. Gonnet, M. B. Monagan
% and S. M. Watt, Maple Reference Manual 5th ed. (Watcom Pub.,
% Waterloo, Ontario, 1988).
\end{thebibliography}
\end{document}

