Asymptotics of Implicit Functions and Computer Algebra
Bruno Salvy
Inria Rocquencourt
Algorithms Seminar
September 9, 1996
[summary by Joris Van der Hoeven]
A properly typeset version of this document is available in
postscript and in
pdf.
If some fonts do not look right on your screen, this might be
fixed by configuring your browser (see the documentation here).
Abstract
We describe several algorithms for the asymptotic inversion of functions,
from the computer algebra viewpoint. We start with some classical results,
such as the inversion theorem of Lagrange. We next consider functions
which can only be expanded in more complicated asymptotic scales.
Finally, we briefly discuss the problem of finding the asymptotic
behaviour of implicit functions.
1 Introduction
Let f be a sufficiently regular function which admits
a functional inverse f^{inv} in the neighbourhood of some point.
One often encounters the problem of determining the asymptotic
behaviour of f^{inv} in terms of the asymptotics of f.
For instance, this problem arises systematically
when computing integral transforms by the saddlepoint method.
If the function admits a power series expansion near
the point of interest, then Lagrange's inversion theorem
can be used, or Brent and Kung's fast algorithm
for power series inversion (see section 2).
However, one often has to deal with logarithmic or exponential singularities,
whence a more general device for asymptotic functional inversion is needed.
A first approach would be to consider more general classes of singularities
of a particular form. For instance, one may consider functions which
admit an asymptotic expansion of the form
f= e^{x} x 



d_{n} x^{n} (x ® ¥); 
see section 3. More generally, one might consider the class
of explog functions, i.e. functions built up from Q and x by
+, , ×, /, exp and log. A first algorithm to obtain
asymptotic information about functional inverses of such
functionsin the form of a socalled nested expansionwas
given
by Salvy and Shackell [6]; see also section 4.
An interesting more general problem than the inversion problem
is to find an algorithm to determine the asymptotic behaviour
of implicit functions, say the solutions of explog functions
in two variables. A partial algorithm to determine
the potential nested expansions of such functions was proposed
in [7] and will be discussed in section 5.
More generally, one might require complete asymptotic expansions of
the inverses of explog functions and solutions of implicit equations.
A nice algorithm to compute complete asymptotic expansions of
explog functions themselves was given in [4]. Notice also that extending the techniques from [4],
solutions for the inversion and implicit function problems for explog
functions were given independently in [10] and [11].
2 Classical results
Let f(y)= f_{0}+ f_{1} y+ f_{2} y^{2}+ ··· be a power series with f_{0} ¹ 0
and assume that y(x) satisfies
y= x f(y).
Then Lagrange's inversion theorem gives a formula for [x^{n}] y:
[x^{n}] y= 

[y^{n1}] f^{n}. 
The formula can be generalized to
[x^{n}] G= 

[y^{n1}] G' f^{n}, 
for an arbitrary power series G.
If f resp. G have an analytic meaning (i.e. the power
series converge, or the power series are the asymptotic expansion
of some functions), then the same holds for the above formulae.
Example 1 The generating function y of Cayley trees satisfies
y= x e^{y}.
Hence [x^{n}] y= 1/n [y^{n1}] e^{ny}= n^{n1}/n!.
Other trees can be treated in a similar way.
Assume now that we do not merely want the nth coefficient of
the inverse of a power series x(y)= y+ a_{2} y^{2}+ a_{3} y^{3}+ ···= F(y),
but the first N coefficients of y(x). This problem can
be reduced to the problem of functional composition by
the Newton method: define the sequence y_{n}(x) by y_{0}=x and
Brent and Kung have shown that the number of correct terms
doubles at each step [1]. Suitably truncating the Newton iterations
then yields an inversion algorithm for power series of
the same complexity as functional composition,
i.e. O(M(N) (N log N)^{1/2})= O((N log N)^{3/2}) [1].
Example 2
How to find the nth positive root x_{n} of tan x= x?
Looking at the graph, we see that x_{n}= (2n+1) p/2 u_{n}
with u_{n} ® 0. Applying tan we get
(2n+1) 

 u_{n}= tan 
æ ç ç è 


 u_{n} 
ö ÷ ÷ ø 
, 
whence
with t_{n} ® 0. Inverting the above power series, we get
x_{n}= p n+ 

 

+ 



æ ç ç è 


+ 


ö ÷ ÷ ø 

+ ··· . 
3 Inverses of more general functions
In this section we are interested in computing the asymptotic
inverse of a function which admits
e^{y} y 

D(y^{1}) (d_{n} ¹ 0, y ® ¥) 
as an asymptotic expansion, where D(y^{1})= å_{n ³ 0} d_{n} y^{n}.
In other words, we want the asymptotic solution to
in y. A first method [2] is to do this is to take logarithms
y= log
x+
a log
y log
D(
y^{1}),
(1)
and to replace y by the right hand side in an iterative manner.
This yields an expansion of the form
y » log 
x+ 


P_{n}(log log x) 

log^{n} x 

, 
where the P_{n} are polynomials.
A faster method to compute the P_{n} was proposed in [5].
Setting z= log log x and t= (log x)^{1}, we set
y log 
x= P(z, t)= 

P_{n}(z) t^{n}. 
Taking logarithms, this yields
log y= log (log x+ P)= z+ log (1+ tP).
From (1), we get
P= a z+ a log (1+ tP)
log D 
æ ç ç è 



ö ÷ ÷ ø 
. 
In this equation z and t are decoupled. We get
P(0,t)= a log (1+ tP(0,t))
log D 
æ ç ç è 



ö ÷ ÷ ø 
and

æ ç ç è 


 t 
ö ÷ ÷ ø 

+
t^{2} 

= 1. 
These two relations enable us to compute the coefficients of P efficiently.
Example 3
The number of prime numbers p(x) smaller than x satisfies
p(x) » 


æ ç ç è 
1+ 

+


+ ··· 
ö ÷ ÷ ø 
. 
Using the above method we find the asymptotic expansion
p(n) » n log n 
æ ç ç è 
1+ 

+


+ ··· 
ö ÷ ÷ ø 
for the nth prime number p(n).
The above result can be extended to get asymptotic expansions
for log y and e^{b y} y^{g}G(y^{1}) [5].
4 Inversion of explog functions
Let f be an explog function (i.e. a function built up from
Q and x by +, , ×, /, exp and log),
which tends to infinity for x ® ¥. We denote by
exp_{i} resp. log_{i} the ith iterate of exp resp. log.
We write f <<< g, if f is of a smaller comparability class
than g, i.e. (log f)/(log g) tends to 0.
In this section, we present an algorithm to compute
a ``nested expansion'' for the asymptotic functional inverse of f.
Such a nested expansion often (although not always) yields an asymptotic
expansion of the function. In any case, the limit behaviour of
a function can be determined from its nested expansion.
One first defines a nested form as being a finite sequence
{ (s_{i}, e_{i}, m_{i}, d_{i}, f_{i}) } for 1 £ i £ n.
Here s_{i}, m_{i} Î N, d_{i} Î R, e_{i}= ± 1,
f_{i} lies in a Hardy field and
f_{i1}(x)= exp 

(log 

(x) f_{i}(x)), 
with f_{i} <<< log_{m}_{i} for 2 £ i £ n.
We also require that
 f_{n} has a finite limit l;
 d_{n} ¹ 1 unless s_{n}=0 or m_{n}=0;
 d_{i}>0 unless s_{i}=0.
Continuing the process of taking nested forms with f_{n}l
(if possible), we obtain a nested expansion.
Example 4
exp [log^{2} x exp [(log log x)^{1/2} (7+ f(x))]]
is a nested form, if f(x) tends to 0.
John Shackell was the first to prove in [8] that any explog
function admits a nested form, and even a nested expansion.
Nested forms and expansions are particularly well suited
to the purpose of functional inversion. The algorithm
proceeds as follows:
Algorithm
An explog function f tending to infinity
is given at input, and we compute a nested expansion of its functional
inverse at infinity.
 Rewrite f as a NF
f(y)= exp_{s} (log_{m}^{d} (y) f_{1}(y))= x.
 Invert
y= exp_{m} (log_{s}^{1/D} x g_{1}(x)). (E)
 Iterate

Compute NF(g_{i})= exp_{s}_{i}^{e}_{i} (log_{m}_{i}^{d}_{i} y G_{i}(y)),
 Substitute (E) in log_{m}_{i}^{d}_{i} y,
 Rewrite to get g_{i+1},
until NF(g_{i})= c+ e (y), with e(y) ® 0, yielding NF(y(x)).
 Repeat the above procedure for e(y), whence NE(y(x)).
Example 5
Step 1: rewrite f as NF
exp 
[ 
log^{2} y exp [ (log_{2} y)^{1/2} (1+ W)] 
] 
, 
with
W= log_{2}^{1/2} y log (1+ log 
^{1} y e 

). 
Step 2: Invert
y= exp 
é ê ê ë 
(log x)^{1/2} exp 
æ ç ç è 
 

(log_{2} y)^{1/2} (1+W) 
ö ÷ ÷ ø 
ù ú ú û 
. 
Step 3: Iterate
(log_{2} y)^{1/2}= 

(log_{2} x)^{1/2}

æ ç ç è 
1+ 


ö ÷ ÷ ø 

, 
whence the nested form for y:
y= exp 
é ê ê ë 
(log x)^{1/2} exp 
æ ç ç è 
 

(log_{2} x)^{1/2} (1+ e (x)) 
ö ÷ ÷ ø 
ù ú ú û 
. 
Step 4: Repeat
y= exp 
é ê ê ê ê ê ê ë 


· e^{1/8} · 
æ ç ç è 
1 

log_{2}^{1/2} x+ 




log_{2}^{1} x+


log_{2}^{3/2} x+ ··· 
ö ÷ ÷ ø 
ù ú ú û 
. 

5 Asymptotic expansion of implicit functions
Assume now that H is an explog function in two variables x and y.
By theorems of Khovanskii and van den Dries [3, 9],
if y(x) is a real solution
of H(x,y)=0 (for x ® ¥), then (the germ of) y(x) belongs
to a Hardy field. In particular, y(x) does not present any oscillatory
phenomena at infinity. As a corollary [7], y(x) has a nested form
and even a nested expansion. A question is how to compute all
possible nested expansions of such solutions.
The idea from the algorithm in [7] is to compute the asymptotic
behaviour for H(x,y) for all possible asymptotic behaviours of x and y.
In order to do so, one often has to distinguish several number of cases,
but always finitely many, depending on the values of x and y.
The strategy is best illustrated on an example.
Example 6
H(x,y)= 
exp (x^{2}+ 2 x log^{2} x+ y) 

exp (x^{2}+x log^{2} x+y) 

1. 
Three cases are distinguished:
H ~ 
ì í î 
1, 
if y ® ¥ and x <<< y; 
?, 
if log y ~ 2 log x; 
exp (x log^{2} x), 
otherwise. 


Since H=0, we find that log y ~ 2 log x,
whence y ~ x^{2}. Continuing the process, we find
y= x^{2} 2xlog^{2} x+ 




+ ··· . 
We remark that it is not claimed that the obtained nested expansions
are indeed nested expansions of actual solutions. For a solution to
this problem, we refer to [11].
References
 [1]

Brent (R. P.) and Kung (H. T.). 
Fast algorithms for manipulating formal power series. Journal of
the ACM, vol. 25, 1978, pp. 581595.
 [2]

De Bruijn (N. G.). 
Asymptotic Methods in Analysis. 
Dover, 1981. A reprint of the third North Holland
edition, 1970 (first edition, 1958).
 [3]

Khovanskii (A. G.). 
Fewnomials and Pfaff manifolds. In Proceedings of the
International Congress of Mathematicians, pp. 549564. 
1983.
 [4]

Richardson (Dan), Salvy (Bruno), Shackell (John), and Van der Hoeven
(Joris). 
Asymptotic expansions of explog functions. In Lakshman (Y. N.)
(editor), Symbolic and Algebraic Computation. pp. 309313. 
New York, 1996. Proceedings ISSAC'96. Zürich.
 [5]

Salvy (Bruno). 
Fast computation of some asymptotic functional inverses. Journal
of Symbolic Computation, vol. 17, 1994, pp. 227236.
 [6]

Salvy (Bruno) and Shackell (John). 
Asymptotic expansions of functional inverses. In Wang (Paul S.)
(editor), Symbolic and Algebraic Computation. pp. 130137. 
New York, 1992. Proceedings of ISSAC'92, Berkeley.
 [7]

Salvy (Bruno) and Shackell (John). 
Symbolic Asymptotics: Functions of Two Variables, Implicit
Functions. 
Research Report n°2883, Institut National de Recherche en
Informatique et en Automatique, 1996. Submitted to the Journal of Symbolic Computation.
 [8]

Shackell (John). 
Growth estimates for explog functions. Journal of Symbolic
Computation, vol. 10, n°6, December 1990, pp. 611632.
 [9]

Van den Dries (Lou). 
Analytic Hardy fields and exponential curves in the real plane.
American Journal of Mathematics, vol. 106, 1984,
pp. 149167.
 [10]

van der Hoeven (J.). 
General algorithms in asymptotics II: Common operations. 
Technical Report n°LIX/RR/94/10, LIX, École polytechnique,
France, 1994.
 [11]

Van der Hoeven (Joris). 
Asymptotique Automatique. 
PhD thesis, École polytechnique, Palaiseau, France,
1997.
This document was translated from L^{A}T_{E}X by
H^{E}V^{E}A.