%%%First we have a character check
%%%
%%% ! exclamation mark    " double quote  
%%% # hash                ` opening quote (grave)
%%% & ampersand           ' closing quote (acute)
%%% $ dollar              % percent       
%%% ( open parenthesis    ) close paren.  
%%% - hyphen              = equals sign
%%% | vertical bar        ~ tilde         
%%% @ at sign             _ underscore
%%% { open curly brace    } close curly   
%%% [ open square         ] close square bracket
%%% + plus sign           ; semi-colon    
%%% * asterisk            : colon
%%% < open angle bracket  > close angle   
%%% , comma               . full stop
%%% ? question mark       / forward slash 
%%% \ backslash           ^ circumflex
%%%
%%% ABCDEFGHIJKLMNOPQRSTUVWXYZ 
%%% abcdefghijklmnopqrstuvwxyz 
%%% 1234567890
%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%
%\documentstyle{ioplppt}                % use this for journal style
\documentstyle[12pt]{ioplppt}          % use this for preprint style
%
\begin{document}
\jl{1}
\title{Summation of asymptotic expansions of multiple-valued
functions using algebraic approximants:\ application
to anharmonic oscillators}[Algebraic approximants]

\author{Alexei V Sergeev\ftnote{1}{Permanent address:
S. I. Vavilov State Optical Institute, Tuchkov
Pereulok 1, St. Petersburg, 199034 Russia}
and David Z Goodson}

\address{Department of Chemistry, Southern Methodist University,
Dallas, TX 75275, USA}


\begin{abstract}
The divergent Rayleigh-Schr\"odinger perturbation expansions for
energy eigenvalues of cubic, quartic, sextic and
octic oscillators are summed using algebraic approximants.  These
approximants are generalized Pad\'e approximants that are obtained
from an algebraic equation of arbitrary degree.  Numerical results
indicate that given enough terms in the asymptotic expansion the
rate of convergence of the diagonal staircase approximant sequence
increases with the degree.  Different branches of the approximants
converge to different branches of the function.  The success of the
high-degree approximants is attributed to their ability to model
the function on multiple sheets of the Riemann surface and to
reproduce the correct singularity structure in the limit of large
perturbation parameter.  An efficient recursive algorithm for
computing the diagonal approximant sequence is presented.
\end{abstract}

%
%  Uncomment out if preprint format required
%
\pacs{03.65.Ge, 02.30, 02.60.Gf}
%% http://www.th.physik.uni-frankfurt.de/~cbest/pacs.numbers.00
%%    02.30.Dk Functions of a complex variable
%%    02.30.Lt Sequences, series, and summability
%%    02.30.Mv Approximations and expansions
%%    02.60.-x Numerical approximation and analysis
%%    02.60.Gf Algorithms for functional approximation
%%    03.65.-w Quantum theory; quantum mechanics
%%    03.65.Db Functional analytical methods
%%    03.65.Ge Solutions of wave equations: bound states
\submitted          
\maketitle

\section{Introduction}

The Schr\"odinger equation $H\psi=E\psi$, where
\begin{equation}
\label{anhq}
H=\case{1}{2}p^2+\case{1}{2}x^2+\lambda x^\beta ,
\end{equation}
gives rise to a well-known example of singular perturbation theory
(Bender and Orszag 1978).  This problem is of physical interest,
as a prototypical quantum field theory and as a model for molecular
vibrations, and of mathematical interest, on account of the rich
singularity structure of the function $E(\lambda)$.  A
characteristic feature of $E(\lambda)$ is an infinite sequence of
branch points approaching a limit point at $\lambda=0$.
(Bender and Wu 1969, Simon 1970, Shanley 1986,
Alvarez 1995, Bender and Orszag 1978 p 350-61).
The asymptotic expansions for the energy,
$E(\lambda)\sim\sum_{n=0}^\infty E_n\lambda^n$, are therefore
divergent for all $\lambda$.  These expansions have become a
standard test case for new summation procedures
(Reid 1967, Graffi \etal 1970, Seznec and Zinn-Justin 1979,
Caswell 1979, Dmitrieva and Plindov 1980, Drummond 1981,
\v{C}\'\ii\v{z}ek and Vrscay 1982, Cohen and Kais 1986,
Weniger \etal 1993),
in part because the $E_n$ can be easily computed even for extremely
large $n$.

The multiple-valued nature of $E(\lambda)$ causes trouble for
summation approximants that are single valued.
Consider for example the quartic oscillator, $\beta =4$.  It can be
proved in that case (Loeffel \etal 1969) that $E(\lambda)$ is
Pad\'e summable as long as $|\arg\lambda|<\pi$.
In practice, the approximants place a sequence
of poles along the negative real axis, simulating a branch cut
(Baker 1975).  The convergence slows as $|\arg\lambda|$
approaches $\pi$ and fails completely at $\pi$.
If $\lambda$ is pure real and negative then the eigenvalue
corresponds to a double-valued complex resonance energy,
$E=E_r\pm\i\Gamma/2$, where $\Gamma$ can be identified, at least
approximately, with the resonance width (Connor and Smith 1981).
The plus and minus signs correspond to the incoming and outgoing
wave boundary conditions respectively.  Since the $E_n$ are real
the Pad\'e approximants for negative real $\lambda$ are also real,
and cannot converge to the correct result.

The Pad\'e approximant $E_{[L,M]}(\lambda)$ is a rational function
$P_L(\lambda)/Q_M(\lambda)$ that is asymptotically equal to
the expansion of $E(\lambda)$.  $P_L$ and $Q_M$ are polynomials
of degree $L$ and $M$ respectively that satisfy the linear equation
\begin{equation}
P(\lambda)-Q(\lambda)E(\lambda)= \Or (\lambda^{L+M+1}).
\label{lineq}
\end{equation}
``$E(\lambda)$'' in \eref{lineq} represents the power series in
$\lambda$.  A natural generalization for double-valued functions
is a quadratic Pad\'e approximant (Shafer 1974). The $[L,M,N]$
quadratic approximant is a function of 3 polynomials $A$, $B$ and
$C$, of degree $L$, $M$ and $N$, respectively, that satisfy
\begin{equation}
A(\lambda)+B(\lambda)E(\lambda)
+C(\lambda)E^2(\lambda)= \Or (\lambda^{L+M+N+2}).
\label{quadeq}
\end{equation}
The approximant $E_{[L,M,N]}(\lambda)$ is obtained by solving
the quadratic equation
\begin{equation}
A(\lambda)+B(\lambda)E_{[L,M,N]}(\lambda)
+C(\lambda)E^2_{[L,M,N]}(\lambda)=0.
\label{quadapprox}
\end{equation}
Convergence theorems for these approximants exist for certain
special cases (Baker and Graves-Morris 1996 pp~544--569),
but in practice applicability to functions of
physical interest has been justified numerically (Short 1979,
Jeziorski \etal 1980,
Liu and Bergersen 1981, Common 1982, Mayer and Tong 1985,
Va\u{\ii}nberg \etal 1986, De'Bell 1992,
Goodson \etal 1992, Hamer \etal 1992, Germann and Kais 1993).
These approximants do not need to simulate the branch cut with
poles; they explicitly contain square-root branch points.  For
the quartic oscillator they converge
at negative real $\lambda$ (Sergeev 1995).

In the same spirit, an algebraic approximant of arbitrary degree $m$
(Short 1979, Brak and Guttmann 1990) can be defined by the equation
\begin{equation}
\sum_{k=0}^m A^{(k)}(\lambda)E^k_{[p_0,p_1,...,p_k]}(\lambda)=0,
\label{Eeq}
\end{equation}
where $A^{(k)}(\lambda)$ are polynomials in $\lambda$ of degree $p_k$
that satisfy
\begin{equation}
\sum_{k=0}^m A^{(k)}(\lambda)E^k(\lambda)= \Or (\lambda^{n+1}),
\qquad n=m-1+\sum_{k=0}^m p_k.
\label{O}
\end{equation}
This is a special case of a class of summation schemes known collectively
as Pad\'e-Hermite approximation (Hermite 1893, Pad\'e 1894, Della Dora
and Di Crescenzo 1979, Baker and Graves-Morris 1996 pp 524-69).
We will use high-degree algebraic approximants to sum the expansion
of the multiple-valued oscillator eigenvalue.

The paper of Short (1979) is the closest prototype of the
present study.  There, similar multiple-valued approximants were
constructed so as to incorporate the known
branch-point structure of Feynman matrix elements.
For the multiple-valued function $\ln(1-z)$,
Short observed that quadratic approximants reduce the
error by roughly two orders of magnitude compared with Pad\'e
approximants.  He found similar improvement in accuracy for
quadratic and especially cubic approximants for certain Feynman
integrals and found that these approximants provide, in effect,
analytic continuations of the asymptotic expansion from the first
Riemann sheet to the second.

Here we present a further demonstration of the power of algebraic
approximants to describe functions with complicated branch-point
structure.  We find that the convergence for the ground-state
energy of anharmonic oscillators improves with
approximant degree, given enough terms in the expansion.  High-degree
approximants yield very high accuracy
for the principal value of $E(\lambda)$ and reasonably good results
on higher sheets.  Certain choices of the approximant indices give
the known large-$\lambda$ singularity structure, and the
corresponding numerical results are especially accurate.
\Sref{calculation} presents an efficient algorithm for
computing algebraic approximants.  In \sref{singularities} we analyze
the singularity structure of the approximants and in \sref{logarithm},
as a simple demonstration, we determine the rate of convergence
for the infinite-valued function $\lambda^{-1}\ln (1+\lambda)$.
\Sref{oscillators} contains results for quartic, cubic, sextic,
and octic oscillators.


\section{Computational algorithm}
\label{calculation}

The most direct method for calculating algebraic approximants is to
solve the system of $n+1$ linear homogeneous equations for the
coefficients of the polynomials $A^{(0)}$, $A^{(1)}$,\dots,
$A^{(m)}$, that follows from \eref{O} after collecting terms by
powers of $\lambda$ (Della Dora and Di Crescenzo 1979).  This
approach is appropriate for low-order analyses, but the
number of arithmetic operations increases very rapidly with
increasing $n$.

For the conventional ``linear'' Pad\'e approximants ($m=1$) the
$[L,M]$ approximants with $L\approx M$ tend in general to be
the most accurate (Baker 1975).  Our experience is that
quadratic and higher-degree approximants with approximately
equal indexes are also the most accurate.
Here we present an algorithm for computing such approximants that is
much faster and needs much less computer memory at large orders than
solving the system of linear equations.  It yields what we will
call the {\it diagonal staircase} sequence of degree-$m$ approximants,
$E_{\{ m,n\}}(\lambda)$, $n=m-1, m, m+1, m+2, \dots$, with
\begin{equation}
{\{ m,n\}\equiv \Big[\atop\,}
{\underbrace{j,j,...,j}\atop i}{,\atop\,}
{\underbrace{j-1,j-1,...,j-1}\atop m+1-i}{\Big],\atop\,}
\label{index}
\end{equation}
where $j$ satisfies the equation $(m+1)j=n-i+2$ with $1\leq i\leq m+1$.
\Tref{appr} lists representative examples of the index sequences,
illustrating the correspondence between $\{ m,n\}$ and
$[p_0,\dots,p_m]$.
\begin{table}
\noindent
\caption{Indices of approximants that comprise the
diagonal staircase sequences.
$m$ is the degree of the approximant and $n$ is the highest order
in the asymptotic expansion that is needed in the calculation.
The entries $[p_0,p_1,\dots ,p_m]$ give the degrees of the
constituent polynomials.}
\begin{indented}
\label{appr}
\item[]\begin{tabular}{@{}ccccc}
\br
$n$&$m=1$&$m=2$&$m=3$&$m=4$\\
\mr
0&$[0,0]$\\
1&$[1,0]$&$[0,0,0]$\\
2&$[1,1]$&$[1,0,0]$&$[0,0,0,0]$\\
3&$[2,1]$&$[1,1,0]$&$[1,0,0,0]$&$[0,0,0,0,0]$\\
4&$[2,2]$&$[1,1,1]$&$[1,1,0,0]$&$[1,0,0,0,0]$\\
5&$[3,2]$&$[2,1,1]$&$[1,1,1,0]$&$[1,1,0,0,0]$\\
$\vdots$&$\vdots$&$\vdots$&$\vdots$&$\vdots$\\
50&$[25,25]$&$[17,16,16]$&$[12,12,12,12]$&$[10,10,9,9,9]$\\
\br
\end{tabular}
\end{indented}
\end{table}

Our algorithm is an extension to arbitrary degree of the
Berlekamp-Massey algorithm (Baker and Graves-Morris pp 153-66).
It was used previously by Mayer and Tong (1985)
for calculating quadratic approximants and is a special case of a
more general algorithm, for Pad\'e-Hermite approximants,
derived by Sergeev (1986).  Let $R_n(\lambda)$ be a sequence of
residual functions,
\begin{equation}
R_n(\lambda)=\sum_{i=0}^\infty r_{n,i}\lambda^i,
\end{equation}
such that
\begin{equation}
\sum_{k=0}^m A_n^{(k)}E^k(\lambda)=\lambda^{n+1}R_n(\lambda).
\label{eqB1}
\end{equation}
Note that we have added the subscript $n$ to $A^{(k)}$ in order
to indicate which $E_{\{ m,n\}}$ it corresponds to in the
diagonal staircase sequence.
We will assume that $r_{n,0}\ne 0$ for all $n$.

The lowest-order approximant of degree $m$
will have each $A_n^{(k)}(\lambda)$ equal to a constant.
Any solution for the set $\{ A_n^{(0)},A_n^{(1)},\dots,A_n^{(m)}\}$
can be multiplied by a common factor.  Thus, one of these constants
is arbitrary.  The remaining $m$ constants are determined from the
accuracy-through-order conditions, \eref{O}.  The lowest-order
approximant corresponds to $n=m-1$, so that there can be $m$
such conditions.  The solution for this approximant is
$A_{m-1}^{(k)}(\lambda)=\left({m\atop k}\right)(-E_0)^{m-k}$, as
can be verified by substitution into \eref{eqB1}.  It follows that
\begin{equation}
\lambda^m R_{m-1}(\lambda)=[E(\lambda)-E_0]^m.
\label{eqB3}
\end{equation}
If $n<m-1$ then the resulting approximant cannot
be of degree $m$.  However, if we define
\begin{equation}
A_n^{(k)}(\lambda)=
\cases{0,&$k>n+1$,\\
       \left({n+1\atop k}\right)(-E_0)^{n+1-k},&$k\le m-j$},
\label{eqB4}
\end{equation}
for $n=-1,0,1,\dots,m-1$, then \eref{eqB1} will be satisfied for
all $n\ge 0$, with
\begin{equation}
\lambda^{n+1}R_n(\lambda)=[E(\lambda)-E_0]^{n+1}.
\label{eqB5}
\end{equation}

Let $\{ c_{n,1}, c_{n,2}, \dots, c_{n,m}\}$ be a set of constants
such that
\begin{equation}
\label{constants}
R_{n-m-1}(\lambda)
+ \sum_{j=1}^m c_{n,j}\lambda^{j-1}R_{n-m-1+j}(\lambda)
= \Or (\lambda^m),
\label{eqB6}
\end{equation}
where the $R_k$ are residuals of degree-$m$ approximants
according to \eref{eqB1}.  Then
the constituent polynomials of
the diagonal approximant sequence satisfy the recursion
\begin{equation}
A_n^{(k)}(\lambda) =
\lambda A_{n-m-1}^{(k)}(\lambda)
+ \sum_{j=1}^m c_{n,j}A_{n-m-1+j}^{(k)}(\lambda)
\label{eqB7}
\end{equation}
for $n\ge m$, with $A_n^{(k)}$ at $n<m$ given by and
\eref{eqB4}, and the corresponding residuals satisfy
\begin{equation}
\lambda^m R_n(\lambda) = R_{n-m-1}(\lambda)
+ \sum_{j=1}^m c_{n,j}\lambda^{j-1}R_{n-m-1+j}(\lambda),
\label{eqB8}
\end{equation}
with $R_n$ at $n<m$ given by \eref{eqB5}.  These equations are easily
verified by substituting \eref{eqB7} into \eref{eqB1}, giving
\begin{equation}
\sum_{k=0}^m A_n^{(k)}(\lambda) E^k(\lambda) =
\lambda^{n-m+1}R_{n-m+1}(\lambda)
+\sum_{j=1}^m c_{n,j}\lambda^{n-m+j}R_{n-m-1+j}(\lambda),
\end{equation}
which according to \eref{eqB6} goes to zero asymptotically
as $\Or (\lambda^{n+1})$.  Thus, \eref{eqB7} satisfies the
accuracy-through-order condition.  Comparison with \eref{eqB1}
establishes \eref{eqB8}.

Collecting terms in \eref{eqB6} according to order in $\lambda$
gives the set of equations
\begin{equation}
r_{n-m-1,k-1}+\sum_{j=1}^k c_{n,j}r_{n-m-1+j,k-j}=0,
\label{ceqs}
\end{equation}
corresponding to $k=1,2,\dots,m$.  These can be solved recursively
for the $c_{n,j}$ starting with $c_{n,1}=-r_{n-m-1,0}/r_{n-m,0}$.
Thus, we can calculate the set of $c_{m,j}$ from the residuals
given by \eref{eqB5}.  Substituting the $c_{m,j}$ into
\eref{eqB7} and \eref{eqB8} gives the $A_m^{(k)}$ and
$R_m$.  The coefficients of $R_m$ can then be used in \eref{ceqs}
to calculate the $c_{m+1,j}$, and so on.

Alternative recursive calculation of $A_n^{(k)}$ and $R_n$
(Sergeev 1986) is as follows.
First, define
subsidiary polynomials $a_{0}^{(k)}$ ($k=0,1,...,m$)
and a subsidiary remainder function $r_0$ as
\begin{eqnarray}
\label{begin}
a_0^{(k)}(\lambda)=\lambda A_{n-m-1}^{(k)}(\lambda),\\
r_{0}(\lambda)=R_{n-m-1}(\lambda).
\end{eqnarray}
Second, define subsequently for $p=1,2,...,m$ another
subsidiary polynomials and remainder functions as
\begin{eqnarray}
\label{reca}
a_{p}^{(k)}(\lambda)=a_{p-1}^{(k)}(\lambda)
+c_{n,p}A_{n-m+p-1}^{(k)}(\lambda),\\
r_{p}(\lambda)=\lambda^{-1}[r_{p-1}(\lambda)
+c_{n,p}R_{n-m+p-1}(\lambda)],
\end{eqnarray}
where
\begin{equation}
c_{n,p}=-r_{p-1}(0)/R_{n-m+p-1}(0).
\label{alternative}
\end{equation}
And third, calculate
\begin{eqnarray}
A_n^{(k)}(\lambda)=a_{m}^{(k)}(\lambda),\\
\label{end}
R_n(\lambda) = r_m(\lambda).
\end{eqnarray}

Equations \eref{begin} -- \eref{end} are very feasible
for computer programming. Note that the constants defined by
\eref{constants} are the same as in \eref{alternative},
\begin{eqnarray}
a_p^{(k)}(\lambda) =
\lambda A_{n-m-1}^{(k)}(\lambda)
+ \sum_{j=1}^p c_{n,j}A_{n-m-1+j}^{(k)}(\lambda),\\
\lambda^p r_p(\lambda) = R_{n-m-1}(\lambda)
+ \sum_{j=1}^p c_{n,j}\lambda^{j-1}R_{n-m-1+j}(\lambda).
\end{eqnarray}
Note also, that
polynomials $A_n^{(k)}$ generated by this algorithm
are normalized so that the leading-order coefficient
in a polynomial $A_n^{(i-1)}$ (i. e. a coefficient before
$\lambda^j$) equals to one where $i$ and $j$ are
defined by \eref{index}.



\section{Branches of the approximants}
\label{singularities}

Consider the diagonal approximant sequence $E_{\{ m,n\}}(\lambda)$,
as defined in \sref{calculation}.  These approximants can have as
many as $m$ branches, corresponding to the $m$ roots of \eref{Eeq}.
Let us determine the asymptotic behavior of these branches
at $\lambda \rightarrow 0$ and at $\lambda \rightarrow \infty$.

It follows from \eref{Eeq} that $E_{\{ m,n\} }(0)$ is a root of a
polynomial
${\cal P}(E_{\{ m,n\} })=\sum_{k=0}^m A_n^{(k)}(0)E_{\{ m,n\} }^k$.
(We do not consider here the case of a multiple root, when
$E_{\{ m,n\} }(0)$ is also a root of $d{\cal P}/dE_{\{ m,n\} }$,
which may occur accidentally but is rare in practice.)
According to a theorem of Baker (Baker and Graves-Morris
1995, pp~534--5) the root of equation \eref{Eeq}
for which $E_{\{ m,n\}}(0)$ is equal to $E(0)$ differs from
$\sum_{i=0}^\infty E_i\lambda^i$ by an error that is at worst
$\Or (\lambda^{n+1})$.  We call this root the {\it principal branch}
and call the sheet of the Riemann surface on which it approximates
the function the {\it principal sheet}.
In general, for all other branches the $\lambda\to 0$ limit
of the approximant will not be equal to the $\lambda\to 0$
limit of the asymptotic expansion $\sum_{n=0}^\infty E_n\lambda^n$.
However, these approximants can, at least in principle, describe
the function $E(\lambda)$ on other sheets.

The large-$\lambda$ behavior is described by the following theorem.

\noindent
{\it Theorem:\ } Let $i$ be the index defined in \eref{index}
that describes the pattern of the degrees of the $A_n^{(k)}$.
In the limit $\lambda\to\infty$, $i-1$ branches of
$E_{\{ m,n\}}(\lambda)$ tend to constants, while the remaining
$m+1-i$ branches behave as $\lambda^{1/(m+1-i)}$.

\noindent
{\it Proof:\ }
In the limit $\lambda \rightarrow \infty$ the asymptotic
behavior of $E_{\{ m,n\}}(\lambda)$ is given by the equation
\begin{equation}
\label{infla}
\sum_{k=0}^{i-1} A_{n,0}^{(k)}E^k_{\{ m,n\} }(\lambda) +
{1\over \lambda}\sum_{k=i}^m A_{n,0}^{(k)}E^k_{\{ m,n\} }(\lambda)=0,
\end{equation}
where $A_{n,0}^{(k)}$ are the leading-order
(i.e. $\lambda^j$ or $\lambda^{j-1}$ according to \eref{index})
coefficients of $A_{n}^{(k)}(\lambda)$.
Assume that $E_{\{ m,n\} }(\lambda)\sim c\lambda^\alpha$, where
$c$ and $\alpha$ are constants.
We consider 3 possibilities: $\alpha=0$, $\alpha>0$, and $\alpha<0$. 
If $\alpha=0$ then at leading order we have
$\sum_{k=0}^{i-1} A_{n,0}^{(k)}c^k=0$, which has $i-1$
solutions for $c$ if $i>1$ and no solutions if $i=1$.
If $\alpha>0$ then
$A_{n,0}^{(m)}c^{m+1-i}\lambda^{(m+1-i)\alpha-1}-A_{n,0}^{(i-1)}=0$.
This has a solution only if $\alpha=m+1-i$.
If $\alpha<0$, then at leading order in $\lambda$
we have $A_{n,0}^{(0)}=0$, which has no solution.


\section{A simple example of a multiple-valued function}
\label{logarithm}

Consider the infinite-valued function
$F(\lambda)=\lambda^{-1}[\ln(1+\lambda)+2\pi K\i]$, where $K$ is
an integer indicating the branch.  For the principal branch,
$K=0$, $F$ has the asymptotic expansion
\begin{equation}
F(\lambda)\sim 1-\case12\lambda+\case13\lambda^2-\case14\lambda^3+
\case15\lambda^4-\case16\lambda^5+\cdots .
\label{log}
\end{equation}
It has been proved (Bender and Orszag 1978 pp 402-3) that the
convergence of linear approximants for \eref{log} to the
principal branch is geometric.
We know of no such theoretical estimates for
higher-degree approximants but find numerically that the diagonal
sequences of approximants of any degree also converge
geometrically, with 
\begin{equation}
\Big| F_{\{ m,n\}}(\lambda)-F(\lambda)\Big|\propto
\left| { 1-(1+\lambda)^{(m+1)^{-1}}
   \exp\left({2\pi\i\over m+1}\bigg[{m+1\over 2}\bigg]\right)\over
        1-(1+\lambda)^{(m+1)^{-1}}    
   \exp\left({2\pi\i\over m+1}K\right) }\right|^{\, n} ,
\label{Fconv}
\end{equation}
where $\left[{m+1\over 2}\right]$ is the greatest integer less than
or equal to ${m+1\over 2}$.  It follows that the approximant
sequence will converge on those branches for which 
$|K|\le\left[{m-1\over 2}\right]$.  On such branches, \eref{Fconv}
in the limit of large $m$ reduces to
\begin{equation}
\Big| F_{\{ m,n\}}-F\Big|\propto
\left[{\scriptstyle{1\over 4}}\ln^{\! 2}\! (1+\lambda)
+\pi^2K^2\right]^{n/2}m^{-n}.
\end{equation}
Thus, increasing the degree $m$ always
increases the rate of convergence in the limit of large order $n$.
However, the convergence rate is slower for larger $K$,
corresponding to branches that are more distant from the principal
branch.

Of particular interest is the determination of the optimal $m$
for given $n$.  \Fref{logacc} shows the accuracy vs $m$ at given
values of $n$ for the principal branch and for the $|K|=2$
branches.  In the Appendix we develop the following expression for
the optimal $m$:
\begin{equation}
m\approx (n+2)^{1/2}-1.
\label{moptimal}
\end{equation}
As shown in \fref{logacc}, this expression does in practice give
very nearly the highest accuracy.


\section{Anharmonic oscillators}
\label{oscillators}

\subsection{Quartic oscillator}

We have computed the exact asymptotic expansion coefficients for the
ground-state energy of the quartic oscillator through 600th order
using a linear algebraic method (Va\u{\ii}nberg \etal 1988,
Dunn \etal 1994).  The coefficients are rational numbers.
Calculations of diagonal staircase approximant sequences
were carried out with Mathematica (Wolfram 1991) in
multiple-precision arithmetic (5000 digits), because the recursive
algorithm is numerically unstable.
The accuracy of $E_{\{m,n\}}$ for various $\lambda$ is shown in
figures~\ref{anhq0}--\ref{panel2}.

\Fref{anhq0} shows the convergence at $\lambda=1/2$ for the
principal branch of $E(\lambda)$, which corresponds to the
ground-state energy.  The 20th-degree approximant sequence appears to
converge to 80 decimal digits, which far surpasses the accuracy
reported previously (Mei\ss ner and Steinborn 1997).
At large $n$ we find that
$\ln |E_{\{m,n\}}-E|\sim -n^\alpha$, with the parameter $\alpha$
increasing with $m$.  For $m=1$ we find numerically that
$\alpha=0.50$, which is the same convergence as linear approximants
for the simple Stieltjes series
$\sum(-\lambda)^n n!$ (Bender and Orszag 1978 pp. 404-5).
For $m=2$, 3, 4, and 20, respectively,
we find $\alpha=0.59$, 0.67, 0.68, and 0.69.

\Fref{panel1} compares the convergence for different $\lambda$
on the circle $|\lambda|=1/2$.  The curves here are
polynomial fits, which supress the relatively small fluctuations
around regular trends.  Convergence at $\lambda=\i/2$ is similar
to that at $\lambda=1/2$, but the accuracy is slightly poorer.
(The physical meaning of imaginary coupling constants will be
discussed below.)
At $\lambda=-1/2$ the potential does not support bound
states.  Linear approximants no longer converge,
but quadratic and, especially, cubic and higher-degree
approximants converge fairly well to a complex energy corresponding
to a quasi-stationary state.  Using the scaling transformation of
the variable $x=x'\exp(-\pi\i/4)$ in the Hamiltonian \eref{anhq},
with $\beta =4$, one can
prove (Crutchfield 1978, Seznec and Zinn-Justin 1978)
that $\lambda=\exp(\frac32\i\pi)\lambda'$ corresponds to a
double-well problem
\begin{equation}
H'=\case12 p^2 -\case12 x'^2 +\lambda' x'^4,
\end{equation}
with eigenvalues $E'(\lambda')=-\i E(\lambda)$.
These lie on the second sheet of Riemann surface.
The bottom panel of \fref{panel1} corresponds to this branch,
with $\lambda'=1/2$.

\Fref{panel2} shows the accuracy of
results for $\lambda'=1/10$ and $\lambda'=3/100$.
Convergence improves significantly with increasing degree
$m$ of the approximant sequences, especially for smaller $\lambda'$.
We attribute this to the presence on this sheet of the infinite sequence
of square-root branch-point pairs, identified by Bender and Wu (1969).
The positions of these branch points are shown in \fref{brpts}.
The closer $\lambda$ is to these points, the greater is the
advantage of increasing the degree of the approximants.
$\lambda'=3/100$ corresponds to a $\lambda$ deep in the heart of
the branch-point region.  We find that approximants with degree
less than 4 do not converge at all at this point while the 20th
degree approximants show slow but steady convergence to a pure
imaginary energy,
\begin{equation}
E\big(\case{3}{100}\e^{3\i\pi/2}\big) =\i E'(\lambda')=
-1.411\,819\,732\,54\i,
\end{equation}
which corresponds to the ground-state energy in the double well.
Moreover, another branch of the degree-20 approximants converges
(at very high order) to $-0.312\,162\,1\i$, which corresponds
to the energy of the second excited state in the double well.
These two branches meet at the branch cut between the first
branch points of the sequence, $\pm 0.0319934-0.0367596\i$.

The principal branch of the function $E(\lambda)$ at
$\lambda=-\i\lambda'$ corresponds to the complex energy of the
barrier resonance in the double well,
$E_{\rm DW}^{\rm r}(\lambda')=-\i E(-\i\lambda')$.
The small-coupling expansion
\begin{equation}
E_{\rm DW}^{\rm r}(\lambda')={-\i\over 2}+{3\over 4}\lambda'
-{21\over 8}\i\lambda'^2+{333\over 16}\i\lambda'^3+\cdots
\end{equation}
represents a formal Rayleigh--Schr\"odinger perturbation series for
the anharmonic oscillator $\frac12\omega^2 x^2+\lambda'x^4$
with an imaginary frequency $\omega=-\i$.
A similar perturbation theory for resonances was recently used by
Fern\'andez (1996).
Such broad resonances with the real part of the energy near the
potential maximum are associated with chemical reaction thresholds
(Friedman \etal 1995).

The case $\lambda=-1/1000$ shown in \fref{1000th}
corresponds to a quasistationary state with extremely small
width, $\Im E\approx\pm 4.319\cdot 10^{-144}$.  The linear approximants
are all pure real.  Their error decreases with $n$ until it becomes
approximately equal to $|\Im E|$ and then it holds steady at that level.
This level of accuracy is eventually reached also by partial
summation, just before the divergence sets in.
Approximants with $m\ge 2$ are real at low $n$, and their accuracy
stalls at the same level as for $m=1$, but at large $n$, once an
imaginary part appears, the accuracy increases rapidly.  This
behavior is qualitatively similar to that observed in a study of
molecular resonances with two degrees of freedom using quadratic
approximants (Suvernev and Goodson, 1997).

The cases $\lambda=100$ and $\lambda=10^6$ displayed at the bottom of
\fref{panel2} correspond to a strong-coupling region, in
which $E(\lambda)\sim b_0 \lambda^{1/3}$ (Turbiner and Ushveridze
1988, Guardiola \etal 1992).  Since linear and quadratic
approximants cannot accurately model
cube-root singularities, their convergence is very slow.
Convergence of approximants with $m>2$ is much better.
The 20th degree approximants of the form $[j,j,...,j,j-1,j-1,j-1]$
($n=21j+16$) are marked by crosses.
Their accuracy is signicantly higher than the overall accuracy of the
20th degree approximants (solid curve) because they always have
correct $\lambda^{1/3}$ behavior at large $\lambda$,
according to the theorem in \sref{singularities}.

Apart from dependence of the accuracy on $n$, we have also studied
the dependence on $m$, that is, the covergence along rows in
\tref{appr}.  The behavior is qualitatively the same as that in
\fref{logacc}.  The condition \eref{moptimal} gives nearly optimal
convergence.


\subsection{Cubic oscillator}

The harmonic oscillator with cubic perturbation,
$H=\case12 p^2+\case12 x^2+gx^3$, is a prototypical system
exhibiting resonances.  Its complex eigenvalues have been studied
numerically (Drummond 1981) and analytically (Alvarez 1988, 1995).
One can expect in general that harmonic oscillators with polynomial
perturbations of any degree greater than 2 will have an infinite
sequence of square-root branch points approaching the origin, and
therefore should benefit from the use of high-degree approximants.

Since odd-order terms of the energy series in $g$ are zero
(because the energy is an even function of $g$), we define the
perturbation parameter as $\lambda=g^2$ and analyze the series
$E(\lambda)\sim \case12 -\case{11}{8}\lambda
-\case{465}{32}\lambda^2-\cdots$,
which has nonzero terms at every order.
Convergence of the algebraic approximants for $\lambda=1/4$
(i.e. $g=1/2$), $\lambda=\case14\exp(5\i\pi/2)$, and $\lambda=100$
is shown in \fref{panel3}.  The convergence behavior is quite
similar to that for the quartic oscillator.  However, for the
quartic oscillator, which has a cube-root branch structure,
there was significant improvement from increasing $m$ up to 3 and
more modest improvement for $m>3$.  For the cubic oscillator,
with a fifth-root structure (Alvarez 1995), there is a greater
advantage from increasing $m$ above 3.  This is especially true for
large-$\lambda$ case, shown in the bottom panel, where the asymptotic
$\lambda^{1/5}$ behavior becomes dominant in $E(\lambda)$.
For the $m=20$ case the accuracy of approximants of the form
$[j,\dots,j,j-1,j-1,j-1,j-1,j-1]$ is marked by crosses.  Their
accuracy is consistently higher than the average accuracy of
the $E_{\{20,n\}}$, as expected from the theorem in
\sref{singularities}.

Using the scaling transformation $x=\omega^{1/2} x'$, one can prove
that $\omega E(\omega^{-5}\lambda)$ is an eigenvalue in a
potential $\omega^2 x^2/2+\lambda^{1/2}x^3$.
In particular, for $\omega=\exp(-\pi\i/2)$, the value
$-\i E(\e^{\frac52\i\pi}\lambda)$ is an eigenvalue in a potential
$-x^2/2+\lambda^{1/2}x^3$.
A shift transformation transforms this modified potential back to the
original potential,
\begin{equation}
-\case12 x^2+\lambda^{1/2}x^3=\case12 x'^2+\lambda^{1/2}x'^3
-1/(54\lambda),
\end{equation}
where $x'=x+1/(3\lambda^{1/2})$, which implies that the eigenvalues
$E(\lambda)$ in the original potential can be expressed in terms of
the new eigenvalues according to
\begin{equation}
E(\lambda)=-\i E(\lambda')+1/(54\lambda),
\label{upturn}
\end{equation}
where $\lambda'=\exp(\frac52\i\pi)\lambda$.
The point $\lambda'$ lies
on the second sheet of Riemann surface under the cut $(0,\infty)$.
$E(\lambda')$ can be expressed as $E'(\i\lambda)$
where $E'$ represents the second branch of the function $E$.
Thus, the eigenvalues can be calculated either by
direct summation of the series $E(\lambda)$ on the principal sheet
or by summing the expansion for $E'$ on its second sheet.
The latter approach is equivalent to expanding the potential
$-x^2/2+\lambda^{1/2}x^3$ at its local maximum
and then developing a complex perturbation theory for an upturned
oscillator with pure imaginary frequency with summation
of the energy expansion on the second sheet.  As shown in the
second panel of \fref{panel3}, this indirect approach does
indeed converge to the same result as the direct approach, and
it benefits even more strongly from the use of high-degree
approximants, although the rate of convergence appears to always
be less than that of the direct analysis.


\subsection{Sextic and octic oscillators}

Perturbation theory for sextic ($\lambda x^6$) and octic
($\lambda x^8$) anharmonic oscillators is very strongly
divergent---the $E_n$ grow as $(2n)!$ and $(3n)!$, respectively.
Linear approximants converge very slowly for the sextic oscillator
even at small $\lambda$ and fail to converge at all for the octic
oscillator (Graffi and Grecchi 1978).
\Fref{anhs} shows that increasing the approximant degree for the
sextic oscillator considerably improves the convergence rate.

The problem of the octic oscillator is particularly interesting because
the $[j,j]$ and $[j+1,j]$ sequences of linear approximants
converge to different limits, giving lower and upper bounds to the
true energy (Austin 1984).  We find that higher-degree approximants
of a given index pattern also converge to an
incorrect result but the number of digits of agreement with the
correct result is considerably greater than that
for linear approximants.  This behavior is shown in \fref{anho}
for diagonal quadratic and cubic approximants.
The accuracy obtained both for the sextic and for the
octic oscillator exceeds the accuracy obtained by Weniger \etal
(1993) using a Levin transform of a renormalized series.


\section{Conclusions}
\label{conclusions}

We have demonstrated that algebraic approximants of degree
$m\ge 3$ can be very effective for summing perturbation series
for quantum oscillators, both on the principal sheet and on
nearby sheets of the Riemann surface.  These approximants
can reproduce several sheets of a multiple-valued function starting
from the Taylor expansion of the function on the principal sheet.
The eigenvalues of a given symmetry are branches of a
single multiple-valued function and the branch points form a sequence
with limit point at the origin.  Similar singularity structure
has also been identified for other kinds of potentials,
including the angular spheroidal wave equation (Hunter and
Guerrieri 1982), the two-center Coulomb problem
(Grozdanov and Solov'ev 1990), and analytically solvable
models (Bender \etal 1974, Ushveridze 1988).  To the extent
that such structure is typical of quantum mechanical eigenvalues,
we expect that algebraic approximants will be useful as a
general summation method for perturbation theories of the
Schr\"odinger equation.  A limitation of the method is that
in practice the number of expansion coefficients needed for
a given degree is approximately equal to the square of the degree.  
Therefore, approximants of very high degree will be most useful
for problems with few degrees of freedom, for which the
perturbation theory can be computed to very high order.


\ack
We thank Professor F M Fern\'andez for bringing to our
attention the phenomenon of broad barrier resonances for
inverted wells.  AVS is grateful to Southern Methodist
University for its hospitality and support.
This work was funded by the U. S. National Science Foundation and
the Robert A. Welch Foundation.

{\appendix
\section*{Appendix: Predicting the optimal degree}

The defining relation for
algebraic approximants of the function $E(\lambda)$, which we have
written in \eref{O} as a polynomial in $E$, can also be thought
of as a polynomial equation in terms of $\lambda$.  Thus, if
we substitute the explicit expression
$\sum_{i=0}^{p_k}a_{k,i}\lambda^i$ for $A^{(k)}$, then we
can write \eref{O} in the form
\begin{equation}
\sum_{i=0}^q B^{(i)}(E)\,\lambda^i=\Or (\lambda^{n+1}),
\label{lambdaeq}
\end{equation}
in terms of a set of polynomials $B^{(i)}$, with $q=[(n+1)/(m+1)]$.
It follows that \eref{O} simultaneously defines an algebraic
approximant of degree $m$ for $E(\lambda)$ and algebraic
approximant of degree $q$ for the inverse function $\lambda(E)$.

In principle, if $E(\lambda)$ has an infinite number of branches
then we can expect that the accuracy of the approximant $E_{\{m,n\}}$
will increase with $m$.  However, the error in the approximant
$E_{\{m,n\}}$ is related to the error in the approximant
$\lambda_{\{q,n\}}$.  Let
\begin{equation}
\delta E=E(\lambda_s)-E_{\{m,n\}}(\lambda_s),\quad
\delta\lambda=\lambda(E_s)-\lambda_{\{q,n\}}(E_s),
\end{equation}
where $\lambda_s$ is the point at which $E$ is being summed and
$E_s=E(\lambda_s)$.  Then
$E_{\{m,n\}}(\lambda_s)=E_{\{m,n\}}(\lambda(E_s)+\delta\lambda)$.
Assuming that $\delta\lambda$ is small, it follows that
\begin{equation}
E_{\{m,n\}}(\lambda_s)\approx
E_s+\delta\lambda\,{dE_{\{m,n\}}\over d\lambda}\approx
E_s+\delta\lambda\,{dE\over d\lambda},
\end{equation}
which implies that $\delta E$ is proportional to $\delta\lambda$.
Let us assume that $\lambda(E)$ is a multiple-valued function.
(This is true for our model function $x^{-1}\ln(1+x)$.)  Then
we can expect that the accuracy of the approximants
$\lambda_{\{q,n\}}$ will increase with $q$.  However, $q$ decreases
with $m$.  This implies that the error in $E_{\{m,n\}}$ will
{\it decrease} with $m$.

Thus, we want $m$ somewhat large, to model the singularity
structure of $E(\lambda)$, but we also want $q$ somewhat large,
to model the singularity structure of $\lambda(E)$.
Based on these arguments, we conjecture that the optimal
approximant degree will correspond to $m\approx q$, from which
\eref{moptimal} follows.}
\vfil\eject

\References
\begin{harvard}
\item[] Alvarez G 1988 \PR A {\bf 37} 4079--83
\item[] Alvarez G 1995 \JPA {\bf 27} 4589-4598
\item[] Austin E J 1984 \JPA {\bf 17} 367--74
\item[] Baker G A Jr 1975 {\it The Essentials of Pad\'e
Approximants} (New York: Academic Press).
\item[] Baker G A Jr and Graves-Morris P 1996 {\it Pad\'e Approximants}
(Cambridge: Cambridge University Press)
\item[] Bender C M, Happ H J and Svetitsky B 1974 \PR D {\bf 9},
2324-9
\item[] Bender C M and Orszag S A 1978 {\it Advanced Mathematical
Methods for Scientists and Engineers} (New York: McGraw-Hill)
\item[] Bender C M and Wu T T 1969 \PR {\bf 184} 1231--60
\item[] Brak R and Guttmann A J 1990 \JPA {\bf 23} L1331--7
\item[] Caswell W E 1979 \APNY {\bf 123} 153--184
\item[] \v{C}\'\ii\v{z}ek J and Vrscay E R 1982 {\it Int. J. Quantum
Chem.} {\bf 21} 27--68
\item[] Cohen M and Kais S 1986 \JPA {\bf 19} 683--690
\item[] Common A K 1982 \JPA {\bf 15} 3665--77
\item[] Connor J N L and Smith A D 1981 {\it Molecular Physics}
{\bf 43} 397-414
\item[] Crutchfield W Y II 1978 \PL {\bf 77B} 109--13
\item[] De'Bell K 1992 \JPA {\bf 25} 1815--20
\item[] Della Dora J and Di Crescenzo C 1979 {\it Pad\'e
Approximation and its Applications (Springer Lecture Notes in
Mathematics, no. 765)} (Berlin: Springer) pp~85--115
\item[] Dmitrieva I K and Plindov G I 1980 \PL {\bf 79A} 47--50
\item[] Drummond J E 1981 \JPA {\bf 14} 1651--61
\item[] Dunn M, Germann T C, Goodson D Z, Traynor C A,
Morgan J D III, Watson D K, Herschbach D R 1994 \JCP {\bf 101}
5987--6004
\item[] Fern\'andez F M 1996 \JCP {\bf 105} 10444--8
\item[] Friedman R S, Hullinger V D and Truhlar D G 1995
{\it J. Phys. Chem.} {\bf 99} 3184--94
\item[] Germann T C and Kais S 1993 \JCP {\bf 99}, 7739--7747.
\item[] Goodson D Z, L\'opez-Cabrera M, Herschbach D R and
Morgan J D III 1992 \JCP {\bf 97} 8481--96
\item[] Graffi S and Grecchi V 1978 \JMP {\bf 19} 1002--6
\item[] Graffi S, Grecchi V and Simon B 1970 \PL {\bf 32B} 631--4
\item[] Grozdanov T P and Solov'ev E A 1990 \PR A {\bf 42} 2703--18
\item[] Guardiola R, Sol\'{\ii}s M A and Ros J 1992 {\it Il Nuovo
Cimento} {\bf 107B} 713--24
\item[] Hamer C J, Oitmaa J and Zheng Weihong 1992 \PR D {\bf 45}
4652--8
\item[] Hermite C 1893 {\it Annali di Mathematica Pura e Applicata,
S\'er. 2} {\bf 21} 289--308
\item[] Hunter C and Guerrieri B 1982 {\it Stud. Appl. Math.}
{\bf 66}, 217-40
\item[] Jeziorski B, Schwalm W A and Szalewicz K 1980
\JCP {\bf 73} 6215--24
\item[] Liu K L and Bergersen B 1981 \CJP {\bf 59} 141--9
\item[] Loeffel J J, Martin A, Simon B and Wightman A S 1969 \PL
{\bf 30B} 656--8
\item[] Mayer I L and Tong B Y 1985 \JPC {\bf 18} 3297--318
\item[] Mei\ss ner H and Steinborn E O 1997 \PR A {\bf 56},
1189-200
\item[] Pad\'e H 1894 {\it J. Math. Pures et Appl. 4} {\bf 10}
291--329
\item[] Reid C E 1967 {\it Int. J. Quantum Chem.} {\bf 1} 521--34
\item[] Sergeev A V 1986 {\it Zh. vychisl. Mat. mat. Fiz.} {\bf 26}
348--56 [{\it U.S.S.R. Comput. Math. Math. Phys.} {\bf 26} 17--22]
\item[] Sergeev A V 1995 \JPA {\bf 28} 4157--62
\item[] Seznec R and Zinn-Justin J 1979 \JMP {\bf 20} 1398--408
\item[] Shafer R E 1974 {\it SIAM J. Numer. Analys.} {\bf 11} 447--60
\item[] Shanley P E 1986 \PL {\bf 117A} 161--5
\item[] Short L 1979 \JPG {\bf 5} 167--98
\item[] Simon B 1970 {\it Ann. Phys., NY} {\bf 58} 76--136
\item[] Suvernev A A and Goodson D Z 1997 \JCP {\bf 106}, 2681-4
\item[] Turbiner A V and Ushveridze A G 1988 \JMP {\bf 29} 2053--63
\item[] Ushveridze A G 1988 \JPA {\bf 21} 955-
\item[] Va\u{\ii}nberg V M, Mur V D, Popov V S, and Sergeev A V 1986
{\it Pis'ma Zh. Eksp. Teor. Fiz.} {\bf 44} 9--12 [{\it JETP
Lett.} {\bf 44} 9--13]
\item[] Va\u{\ii}nberg V M, Mur V D, Popov V S, Sergeev A V, and
Shcheblykin 1988 {\it Teor. Mat. Fiz.} {\bf 74} 399--411 [{\it Theo.
Math. Phys.} {\bf 74} 269--278].
\item[] Weniger E J, \v C\'\ii\v zek J and Vinette F 1993 \JMP
{\bf 34} 571--609
\item[] Wolfram S 1991 {\it Mathematica: a system for doing
mathematics by computer} (Redwood City, California: Addison-Wesley)
\end{harvard}
\endrefs

\Figures
\Figure{
Dependence of accuracy on the degree $m$ of the algebraic
approximant of order $n$ in the diagonal staircase sequence
for the function $F(\lambda)=\lambda^{-1}[\ln (1+\lambda)+2\pi\i K]$
on the branches $K=0$ and $K=\pm 2$ at $\lambda=1$.  The optimal
$m$ for given $n$ is indicated by a circle.  The predicted optimal
$m$, according to \protect\eref{moptimal} is indicated by a star.  
The measure of accuracy is $-\lg|F_{\{m,n\}}-F|$, which is roughly
equal to the number of correct digits after the decimal point.
\label{logacc} }

\Figure{Accuracy of diagonal staircase approximant sequences
$E_{\{m,n\}}$ for the ground-state energy of the quartic
oscillator at $\lambda=1/2$.  $n$ is the order of
the perturbation expansion while $m$ is the approximant degree.
The measure of accuracy is $-\lg|E_{\{m,n\}}-E|$.
\label{anhq0} }

\Figure{Accuracy of diagonal staircase approximant sequences
$E_{\{m,n\}}$ for the ground-state energy of the quartic
anharmonic oscillator at different $\lambda$ on the complex
circle $|\lambda|=1/2$.  The curves have been smoothed by
fitting with a polynomial.  The approximant degree is indicated
as follows: $m=1$, \dotted$\;$; $m=2$, \broken{\;}; $m=3$,
\chain$\;$; $m=4$, \dashddot$\;$; $m=20$, \full$\;$.
The measure of accuracy is $-\lg|E_{\{m,n\}}-E|$.
(The vertical scale is different for different $\lambda$.)
\label{panel1}}

\Figure{Accuracy of diagonal staircase approximant sequences
$E_{\{m,n\}}$ for the ground-state energy of the quartic
anharmonic oscillator at various values of $\lambda$.
The curves are smoothed by polynomial fits.  The approximant
degree is indicated by curve type as in
\protect{\fref{panel1}}.  In the bottom two panels approximants
that have correct large-$\lambda$ cube-root behavior are marked
by crosses.
\label{panel2}}

\Figure{Diagram of analytic structure of $E(\lambda)$ for the
quartic oscillator on the second sheet of the Riemann surface
corresponding to the branch cut $(0,\infty)$, showing pairs
of square-root branch points with limit point at the origin
(Shanley 1986).
\label{brpts}}

\Figure{Accuracy of diagonal staircase approximant sequences
$E_{\{m,n\}}$ and of partial sums for a very long-lived
quasistationary state of the quartic oscillator.
Accuracy of approximants with $m>3$ is indistinguishable from
that of cubic approximants within the scale of the figure.
\label{1000th}}

\Figure{Accuracy of diagonal staircase approximant sequences
$E_{\{m,n\}}$ for the ground-state energy of the cubic
anharmonic oscillator at different values of $\lambda$.
The curves are smoothed by polynomial fits.  The approximant
degree is indicated by curve type as in
\protect{\fref{panel1}}.
The energy at $\lambda=\frac14\exp(\frac52\i\pi)$
is closely related to the energy at $\lambda=\frac14$
according to \protect{\eref{upturn}}.
In the bottom panel approximants that have correct
large-$\lambda$ fifth-root behavior are marked by crosses.
\label{panel3} }

\Figure{Accuracy of diagonal staircase approximant sequences
$E_{\{m,n\}}$ for the ground-state energy of the sextic
oscillator at $\lambda=1/10$.
\label{anhs}}

\Figure{Accuracy of diagonal approximants $E_{[j,j,\dots,j]}$
for the ground-state energy of the octic oscillator at
$\lambda=1/100$, with approximant degree $m$ as indicated.
The linear approximants converge to 0.5272
\ (99.1\% of the exact energy),
the quadratic approximants converge to 0.532\,105\ (100.0002\%),
and the cubic approximants converge
to 0.532\,103\,926\ (99.999\,999\,7\%).  The ``exact'' energy
for this system was calculated by numerical integration of
the differential equation.
\label{anho}}


\end{document}