Compute a generalized series expansion
This functionality does not run in MATLAB.
series(f
,x
, <order
>, <Left  Right  Real  Undirected
>, <NoWarning>, <UseGseries>) series(f
,x = x_{0}
, <order
>, <Left  Right  Real  Undirected
>, <options
>)
series(f, x = x_{0})
computes
the first terms of a series expansion of f
with
respect to the variable x
around the point x_{0}
.
series
tries to compute either the Taylor
series, the Laurent series, the Puiseux series, or a generalized series
expansion of f
around x = x_{0}
.
See Series::gseries
for
details on generalized series expansions.
The mathematical type of the series returned by series
can
be queried using the type expression Type::Series
.
If series
cannot compute a series expansion
of f
, a symbolic function call is returned. This
is an expression of type "series"
. Cf. Example 11.
Mathematically, the expansion computed by series
is
valid in some neighborhood of the expansion point in the complex plane.
Usually, this is an open disc centered at x_{0}
.
However, if the expansion point is a branch point, then the returned
expansion may not approximate the function f
for
values of x
close to the branch cut. Cf. Example 12.
Using the options Left
or Right
,
one can compute directed expansions that are valid along the real
axis. With the option Real
, a twosided expansion
along the real axis is computed. See Example 5 and Example 6.
If x_{0}
is infinity
or infinity
,
then a directed series expansion along the real axis from the left
to the positive real infinity or from the right to the negative real
infinity, respectively, is computed. If x_{0}
is complexInfinity
and dir
is
not specified or Undirected
, then an undirected
series expansion around the complex infinity, i.e., the north pole
of the Riemann sphere, is computed. Specifying x_{0}=
infinity
is equivalent to x_{0}=
complexInfinity
and dir = Left
. Similarly, x_{0}=
infinity
is equivalent to x_{0}=
complexInfinity
and dir = Right
. Cf. Example 7.
Such a series expansion is computed as follows: The series variable x
in f
is
replaced by
(or
for x_{0}=
infinity
). Then, a series expansion of f
around u =
0 is computed. Finally,
(or
,
respectively) is substituted in the result.
Mathematically, the result of such a series expansion is a series in . However, it may happen that the coefficients of the returned series depend on the series variable. See the corresponding paragraph below.
The number of requested terms for the expansion is the argument order
if
specified. Otherwise, the value of the environment variable ORDER
is
used. One can change the default value 6 by
assigning a new value to ORDER
.
The number of terms is counted from the lowest degree term on
for finite expansion points, and from the highest degree term on for
expansions around infinity, i.e., "order
"
has to be regarded as a "relative truncation order".
series
implements a limited amount of precision
management to circumvent cancellation. If the number of terms of the
computed expansion is less than order
, a second
series computation with a higher value of order
is
tried automatically, and the result of the latter is returned.
Note: Nevertheless, the actual number of terms in the resulting series expansion may differ from the requested number of terms. See Example 13 and Example 15. 
Taylor/Laurent/Puiseux expansions (all of domain type Series::Puiseux
) can
be restricted easily to an absolute order term by adding an appropriate O
term. Cf. Example 14.
Expansions of symbolic integrals can be computed. Cf. Example 16.
If f
is an expression of type RootOf
,
then series
returns the set of all nonzero series
solutions of the corresponding algebraic equation. Cf. Example 9.
If order
has the value infinity
, then the system tries to convert
the first argument into a formal infinite series, i.e., it computes
a general formula for the nth
coefficient in the Taylor expansion of f
. The result
is an inactive symbolic sum or a polynomial expression.
Cf. Example 10.
If series
returns a series expansion of domain
type Series::Puiseux
,
it may happen that the "coefficients" of the returned
series depend on the series variable. In this case, the expansion
is not a proper Puiseux series in the mathematical sense. See Example 7 and Example 8. However, if the series variable
is x and
the expansion point is x_{0},
then the following is valid for each coefficient function c(x) and
every positive ε: c(x) (x  x_{0})^{ε} converges
to zero and
is
unbounded when x approaches x_{0}.
Similarly, if the expansion point is infinity,
then, for every positive ε,
converges
to zero and c(x) x^{ε} is
unbounded when x approaches infinity.
The function returns a domain object that can be manipulated
by the standard arithmetical operations. Moreover, the following methods
are available: ldegree
returns
the exponent of the leading term; Series::Puiseux::order
returns
the exponent of the error term; expr
converts to an arithmetical expression,
removing the error term; coeff(s, n)
returns the
coefficient of the term of s
with exponent n
; lcoeff
returns the leading
coefficient; revert
computes
the inverse with respect to composition; diff
and int
differentiate and integrate a series
expansion, respectively; map
applies
a function to all coefficients. See the help pages for Series::Puiseux
and Series::gseries
for
further details.
Note:

The function is sensitive to the environment variable ORDER
,
which determines the default number of terms in series computations.
We compute a series expansion of sin(x) around x = 0. The result is a Taylor series:
s := series(sin(x), x)
Syntactically, the result is an object of domain type Series::Puiseux
:
domtype(s)
The mathematical type of the series expansion can be queried
using the type expression Type::Series
:
testtype(s, Type::Series(Taylor))
Various system functions are overloaded to operate on series objects. E.g.,
the function coeff
can
be used to extract the coefficients of a series expansion:
coeff(s, 5)
The standard arithmetical operators can be used to add or multiply series expansions:
s + 2*s, s*s
delete s:
This example computes the composition of s
by
itself, i.e. the series expansion of sin(sin(x)).
s := series(sin(x), x): s @ s = series(sin(sin(x)), x)
delete s:
We compute the series expansion of the tangent function around the origin in two ways:
series(sin(x), x) / series(cos(x), x) = series(tan(x), x)
bool(%)
We compute a Laurent expansion around the point 1:
s := series(1/(x^2  1), x = 1)
testtype(s, Type::Series(Taylor)), testtype(s, Type::Series(Laurent))
Without an optional argument or with the option Undirected
,
the sign
function
is not expanded:
series(x*sign(x^2 + x), x) = series(x*sign(x^2 + x), x, Undirected)
Some simplification occurs if one requests an expansion that is valid along the real axis only:
series(x*sign(x^2 + x), x, Real)
The sign
vanishes
from the result if one requests a onesided expansion along the real
axis:
series(x*sign(x^2 + x), x, Right), series(x*sign(x^2 + x), x, Left)
In MuPAD^{®}, the heaviside
function
is defined only on the real axis. Thus an undirected expansion in
the complex plane does not make sense:
series(x*heaviside(x + 1), x)
Warning: Cannot find an undirected series expansion. Try the 'Left', 'Right', or 'Real' option. [Series::main]
After specifying corresponding options, the system computes an expansion along the real axis:
series(x*heaviside(x + 1), x, Real), series(x*heaviside(x + 1), x, Right)
At the point I
in the complex plane, the
function heaviside
is
not defined, and neither is a series expansion:
series(heaviside(x), x = I, Real)
We compute series expansions around infinity:
s1 := series((x + 1)/(x  1), x = complexInfinity)
s2 := series(psi(x), x = infinity)
domtype(s1), domtype(s2)
Although both expansions are of domain type Series::Puiseux
, s2
is
not a Puiseux series in the mathematical sense, since the first term
contains a logarithm, which has an essential singularity at infinity:
testtype(s1, Type::Series(Puiseux)), testtype(s2, Type::Series(Puiseux))
coeff(s2)
The following expansion is of domain type Series::gseries
:
s3 := series(exp(x)/(1  x), x = infinity, 4)
domtype(s3)
delete s1, s2, s3:
Oscillating but bounded functions may appear in the "coefficients" of a series expansion as well:
s := series(sin(x + 1/x), x = infinity)
domtype(s), testtype(s, Type::Series(Puiseux))
coeff(s, 1)
The algebraic equation y^{5}  y  x = 0 cannot be resolved in terms of radicals:
solve(y^5  y  x, y)
However, series
can compute all series solutions
of this equation around x =
0:
series(%, x = 0)
It may happen that the series solutions themselves are expressed
in terms of RootOf
s:
series(RootOf(y^5 (x + 2*x^2)*y^3  x^3*y^2 + (x^3 + x^4)*y + x^4 + x^5, y), x)
The coefficients of the algebraic equation are allowed to be
transcendental. They are internally converted into Puiseux series
by series
:
series(RootOf(y^3  y  exp(x  1) + 1, y), x = 1, 4)
An error occurs if some coefficient cannot be expanded into a Puiseux series:
series(RootOf(y^3  y  exp(x), y), x = infinity)
Error: Cannot expand the coefficients of 'RootOf(y^3  y  exp(1/x), y)' into a series. [Series::algebraic]
In this example, we compute a formula for the nth
coefficient a_{n} in
the Taylor expansion of the function
around
zero, by specifying infinity
as order
.
The result is a symbolic sum:
series(exp(x), x, infinity)
If the input is a polynomial expression, then so is the output:
series(x^5  1, x = 1, infinity)
No asymptotic expansion exists for the Bessel
J function of unspecified index, and series
returns
a symbolic function call:
series(besselJ(k, x), x=infinity)
domtype(%), type(%)
The branch cut of the logarithm and
the square root is the negative real
axis. For a series expansion on the branch cut, series
uses
the function signIm
to
return an expansion that is valid in an open disc around the expansion
point:
series(ln(x), x = 1, 3)
series(sqrt(x), x = 1, 3)
The situation is more intricate when the expansion point is
a branch point. The following expansion of the function arcsin(x
+ 1)
is valid in an open disc around the branch point 0:
series(arcsin(x + 1), x, 4)
However, the expansion of f = ln(x + I*x^3)
around
the branch point 0 that is
returned by series
does not approximate f
for
values of x
that are close to the negative real
axis:
f := ln(x + I*x^3); g := series(f, x, 4);
DIGITS := 20: float(subs([f, expr(g)], x = 0.01 + 0.0000001*I)); delete DIGITS:
The situation is similar for algebraic branch points:
f := sqrt(x + I*x^3); g := series(f, x, 4);
DIGITS := 20: float(subs([f, expr(g)], x = 0.01 + 0.0000001*I)); delete DIGITS:
delete f, g:
The first six terms, including zeroes, of the following two series expansions agree:
series(sin(tan(x)), x, 12); series(tan(sin(x)), x, 12);
If we want to compute the series expansion of the difference sin(tan(x))
 tan(sin(x))
, cancellation happens and produces too few
terms in the result. series
detects this automatically
and performs a second series computation with increased precision:
series(sin(tan(x))  tan(sin(x)), x, 6)
It may nevertheless happen that the result has too few terms; cf. Example 15.
If rational exponents occur in the series expansion, then it may even happen that the result has more than the number of terms requested by the third argument:
series(x^2*exp(x) + x*sqrt(sin(x)), x, 3)
series
's control of the order term is based
on the concept of `relative order', counting the number of terms beginning
with the lowest order that is present in the expansion. An `absolute
order' control can be achieved by simply adding an appropriate order term to restrict a result returned
by series
:
series(exp(x) + x*sqrt(sin(x)), x, 7)
series(exp(x) + x*sqrt(sin(x)), x, 7) + O(x^4)
Note, however, that the series must have enough terms for the added order term to have any effect:
series(exp(x^2), x, 4)
series(exp(x^2), x, 4) + O(x^8)
If the specified order for the expansion is too small to compute
the reciprocal (due to cancellation), series
returns
a symbolic call:
series(exp(x), x, 4)
series(1/(exp(x)  1  x  x^2/2  x^3/6), x, 2)
After increasing the order, an expansion is computed, but possibly with fewer terms:
series(1/(exp(x)  1  x  x^2/2  x^3/6), x, 3); series(1/(exp(x)  1  x  x^2/2  x^3/6), x, 4)
series
and int
support each other. On the one hand,
series expansions can be integrated:
int(series(1/(2  x), x), x = 0..1)
On the other hand, series
knows how to handle
symbolic integrals:
int(x^x, x)
series(%, x = 0, 3)
int(exp(x*sin(t)), t = 0.. x)
series(%, x = 0)
int(cos((x*t^2 + x^2*t))^(1/3), t = 0..2)
series(%, x)
Users can extend the power of series
by implementing series
attributes
(slots) for their own special mathematical
functions.
We illustrate how to write such a series attribute, using the
case of the exponential function. (Of course, this function already
has a series
attribute in MuPAD, which you
can inspect via expose(exp::series)
.) In order
not to overwrite the already existing attribute, we work on a copy
of the exponential function called Exp
.
The series
attribute must be a procedure
with four arguments. This procedure is called whenever a series expansion
of Exp
with an arbitrary argument is to be computed.
The first argument is the argument of Exp
in the series
call.
The second argument is the series variable; the expansion point is
always the origin 0; other
expansion points are internally moved to the origin by a change of
variables. The third and the fourth argument are identical with the order
and
the dir
argument of series
,
respectively.
For example, the command series(Exp(x^2 + 2), x, 5)
is
internally converted into the call Exp::series(x^2 + x, x,
5, Undirected)
. Here is an example of a series
attribute
for Exp
.
// The series attribute for Exp. It handles the call // series(Exp(f), x = 0, order, dir) ExpSeries := proc(f, x, order, dir) local t, x0, s, r, i; begin // Expand the argument into a series. t := series(f, x, order, dir); // Determine the order k of the lowest term in t, so that // t = c*x^k + higher order terms, for some nonzero c. k := ldegree(t); if k = FAIL then // t consists only of an error term O(..) return(FAIL); elif k < 0 then // This corresponds to an expansion of exp around infinity, // which does not exist for the exponential // function, since it has an essential singularity. Thus we // return FAIL, which makes series return unevaluatedly. For // other special functions, you may add an asymptotic // expansion here. return(FAIL); else // k >= 0 // This corresponds to an expansion of exp around a // finite point x0. We write t = x0 + y, where all // terms in y have positive order, use the // formula exp(x0 + y) = exp(x0)*exp(y) and compute // the series expansion of exp(y) as the functional // composition of the Taylor series of exp(x) around // x = 0 with t  x0. If your special function has // any finite singularities, then they should be // treated here. x0 := coeff(t, x, 0); s := Series::Puiseux::create(1, 0, order, [1/i! $ i = 0..(order  1)], x, 0, dir); return(Series::Puiseux::scalmult(s @ (t  x0), Exp(x0), 0)) end_if end_proc:
This special function must be embedded in a function
environment. The following command defines Exp
as
a function environment and lets the system function exp
do the evaluation.
The subs
command
applied on the result achieves that Exp
with symbolic
arguments is returned as Exp
and not as exp
.
Exp := funcenv(x > subs(exp(x), hold(exp)=hold(Exp))): Exp(1), Exp(1.0), Exp(x^2 + x)
series
can already handle this "new"
function, but it can only compute a Taylor expansion with symbolic
derivatives:
ORDER := 3: series(Exp(x), x = 0)
One can define the series
attribute of Exp
by
assigning the procedure above to its series
slot:
Exp::series := ExpSeries:
Now we can test the new attribute:
series(Exp(x^2 + x), x = 0) = series(exp(x^2 + x), x = 0)
series(Exp(x^2 + x), x = 2) = series(exp(x^2 + x), x = 2)
series(Exp(x^2 + x), x = 0, 1)
series(Exp(x^2 + x), x = infinity)
Another possibility to obtain series expansions of userdefined
functions is to define the diff
attribute
of the corresponding function environment. This is used by series
to
compute a Taylor expansion when no series
attribute
exists. However, this only works when a Taylor expansion exists, whilst
a series
attribute can handle more general types
of series expansions as well.
delete ExpSeries, Exp:

An arithmetical
expression representing a function in 

An identifier 

The expansion point: an arithmetical expression. If not specified, the default expansion point 0 is used. 

The number of terms to be computed: a nonnegative integer or 

If no expansion exists that is valid in the complex plane, this
argument can be used to request expansions that only need to be valid
along the real line. The default is 

Supresses warning messages printed during the series computation.
This can be useful if 

Option, specified as Use 
If order
is a nonnegative integer, then series
returns
either an object of the domain type Series::Puiseux
or Series::gseries
, an expression of type "series"
,
or, if f
is a RootOf
expression, a set of type Type::Set
. If order
= infinity
, then series
returns an arithmetical expression.
f
asympt
 limit
 mtaylor
 O
 ORDER
 RootOf
 Series::gseries
 Series::Puiseux
 solve
 taylor
 Type::Series