## `series` -- compute a (generalized) series expansion

#### Introduction

`series``(f, x = x0)` computes the first terms of a series expansion of `f` with respect to the variable `x` around the point `x0`.

#### Call(s)

`series(f, x < = x0> <, order> <, dir> <, NoWarning>)`

#### Parameters

`f` | - | an arithmetical expression representing a function in `x`
`x` | - | an identifier
`x0` | - | the expansion point: an arithmetical expression. If not specified, the default expansion point 0 is used.
`order` | - | the number of terms to be computed: a nonnegative integer or `infinity`. The default order is given by the environment variable `ORDER` (default value 6).

#### Options

`dir` | - | either Left, Right, or Real. 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.
`NoWarning` | - | supresses warning messages printed during the series computation. This can be useful if `series` is called within user-defined procedures.

#### Returns

If `order` is a nonnegative integer, then `series` returns either an object of the domain type `Series::Puiseux` or `Series::gseries`, or an expression of type `"series"`. If `order = infinity`, then `series` returns an arithmetical expression.

#### Side Effects

The function is sensitive to the environment variable `ORDER`, which determines the default number of terms in series computations.

### Overloadable:

`f`

#### Related Functions

`asympt`, `limit`, `O`, `ORDER`, `Series::gseries`, `Series::Puiseux`, `taylor`, `Type::Series`

#### Details

- `series` tries to compute either the Taylor series, the Laurent series, the Puiseux series, or a generalized series expansion of `f` around `x = x0`. 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 9.

- Mathematically, the expansion computed by `series` is valid in some neighborhood of the expansion point in the complex plane. Using the options Left or Right, one can compute directional expansions that are valid along the real axis. With the option Real, a two-sided expansion along the real axis is computed. Cf. examples 4 and 5.

- If `x0` is `infinity` or `-infinity`, then a directional 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. Cf. example 7.

  Such a series expansion is computed as follows: The series variable `x` in `f` is replaced by x = 1/u. Then, a directional series expansion of `f` around u = 0+ is computed. Finally, u = 1/x is substituted in the result.

  Mathematically, the result of such a series expansion is a power series in 1/x. 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 `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, i.e., ```order`'' has to be regarded as a ``relative truncation order''.

  The actual number of terms in the resulting series expansion may differ from the requested number of terms. Cf. examples 10 and 11.

- In some cases, when cancellation occurs, it may happen that the requested order is too small to compute a series expansion. In such a case, the computation is aborted with an error message. Cf. example 10.

- 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 n-th coefficient in the Taylor expansion of `f`. The result is a symbolic sum. Cf. example 8.

- 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. Cf. example 7. However, if the series variable is x and the expansion point is x0, then the following is valid for each coefficient function c(x) and every positive e: c(x)*(x-x0)^e converges to zero and c(x)*(x-x0)^(-e) is unbounded when x approaches x0. Similarly, if the expansion point is infinity, then, for every positive e, c(x)*x^(-e) converges to zero and c(x)*x^e 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` differentiates a series expansion; `map` applies a function to all coefficients. See the help pages for `Series::Puiseux` and `Series::gseries` for further details.

#### Example 1

We compute a series expansion of sin(x) around x = 0. The result is a Taylor series:

    >> s := series(sin(x), x)

                                    3    5
                                   x    x        6
                               x - -- + --- + O(x )
                                   6    120

Syntactically, the result is an object of domain type `Series::Puiseux`:

    >> domtype(s)

                                  Series::Puiseux

The mathematical type of the series expansion can be queried using the type expression `Type::Series`:

    >> testtype(s, Type::Series(Taylor))

                                       TRUE

Various system functions were 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)

                                       1/120

The standard arithmetical operators can be used to add or multiply series expansions:

    >> s + 2*s, s*s

                          3    5                4      6
                         x    x       6    2   x    2 x       7
                   3 x - -- + -- + O(x ), x  - -- + ---- + O(x )
                         2    40               3     45

    >> delete s:

#### Example 2

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)

                          3    5                3    5
                         x    x       6        x    x       6
                     x - -- + -- + O(x ) = x - -- + -- + O(x )
                         3    10               3    10

    >> delete s:

#### Example 3

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)

                        3      5                3      5
                       x    2 x       6        x    2 x       6
                   x + -- + ---- + O(x ) = x + -- + ---- + O(x )
                       3     15                3     15

    >> bool(%)

                                       TRUE

#### Example 4

Without an optional argument, the `sign` function is not expanded:

    >> series(x*sign(x^2 + x), x)

                                          2       7
                              x sign(x + x ) + O(x )

Some simplification occurs if one requests an expansion that is valid along the real axis only:

    >> series(x*sign(x^2 + x), x, Real)

                                                6
                                 x sign(x) + O(x )

The `sign` vanishes from the result if one requests an one-sided expansion along the real axis:

    >> series(x*sign(x^2 + x), x, Right),
       series(x*sign(x^2 + x), x, Left)

                                     6            6
                              x + O(x ), - x + O(x )

#### Example 5

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: Could not find undirected series expansion; try optio\
          n `Left', `Right', or `Real' [Series::main]

                           series(x heaviside(x + 1), x)

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)

                                      7          7
                               x + O(x ), x + O(x )

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)

          Error: heaviside is not defined for non-real expansion points \
          [heaviside::series]

#### Example 6

We compute a Laurent expansion around the point 1:

    >> series(1/(x^2 - 1), x = 1)

                                                 2          3
              1             / x       \   (x - 1)    (x - 1)
          --------- - 1/4 + | - - 1/8 | - -------- + -------- +
          2 (x - 1)         \ 8       /      16         32

                      4
             O((x - 1) )

#### Example 7

We compute series expansions around infinity:

    >> s1 := series((x + 1)/(x - 1), x = infinity)

                              2   2    2    2     / 1  \
                          1 + - + -- + -- + -- + O| -- |
                              x    2    3    4    |  5 |
                                  x    x    x     \ x  /

    >> s2 := series(psi(x), x = infinity)

                               1      1       1       / 1  \
                      ln(x) - --- - ----- + ------ + O| -- |
                              2 x       2        4    |  5 |
                                    12 x    120 x     \ x  /

    >> domtype(s1), domtype(s2)

                         Series::Puiseux, Series::Puiseux

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:

    >> coeff(s2)

                           ln(x), -1/2, -1/12, 0, 1/120

The following expansion is of domain type `Series::gseries`:

    >> s3 := series(exp(x)/(1 - x), x = infinity, 4)

                       exp(x)   exp(x)   exp(x)    / exp(x) \
                     - ------ - ------ - ------ + O| ------ |
                         x         2        3      |    4   |
                                  x        x       \   x    /

    >> domtype(s3)

                                  Series::gseries

    >> delete s1, s2, s3:

#### Example 8

In this example, we compute a formula for the n-th coefficient a[n] in the Taylor expansion of the function exp(-x) = sum a[n]*x^n around zero, by specifying `infinity` as `order`. The result is a symbolic sum:

    >> series(exp(-x), x, infinity)

                          /   n1     n1                    \
                          |  x   (-1)                      |
                       sum| ------------, n1 = 0..infinity |
                          \ n1 gamma(n1)                   /

#### Example 9

The sine function has an essential singularity at infinity. `series` cannot compute a series expansion and returns a symbolic function call:

    >> series(sin(x), x = infinity)

                           series(sin(x), x = infinity)

    >> domtype(%), type(%)

                                DOM_EXPR, "series"

#### Example 10

In the following example, the specified order for the expansion is too small to compute the reciprocal, due to cancellation:

    >> series(exp(x), x, 3)

                                         2
                                        x       3
                                1 + x + -- + O(x )
                                        2

    >> series(1/(exp(x) - 1 - x - x^2/2), x, 3)

          Error: order too small [Series::Puiseux::_invert]

After increasing the order, an expansion is computed, but with fewer terms:

    >> series(1/(exp(x) - 1 - x - x^2/2), x, 5)

                                6     3      / 1 \
                                -- - ---- + O| - |
                                 3      2    \ x /
                                x    2 x

#### Example 11

Here are some examples where the actual number of computed terms differs from the requested number:

    >> series(sin(x^2), x, 5)

                                     2      5
                                    x  + O(x )

    >> series((sin(x^4) - tan(x^4)) / x^10, x, 15)

                                      2
                                     x       5
                                   - -- + O(x )
                                     2

#### Example 12

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, FAIL)`. In this case, the fourth argument `FAIL` indicates that an undirected series expansion is to be computed. 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(..)
           error(örder too small");

         elif k < 0 then
           // This corresponds to an expansion of exp around infinity
           // or -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);
           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 copies the code for evaluating the system function `exp`. The `subsop` command achieves that `Exp` with symbolic arguments is returned as `Exp` and not as `exp`, see the help page for `DOM_PROC`.

    >> Exp := funcenv(subsop(op(exp, 1), 6 = hold(Exp)), NIL, NIL):
       Exp(1), Exp(-1.0), Exp(x^2 + x)

                                                        2
                         Exp(1), 0.3678794412, Exp(x + 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)

                                        2
                                       x  D(D(Exp))(0)      3
                     1 + x D(Exp)(0) + --------------- + O(x )
                                              2

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)

                               2                      2
                            3 x       3            3 x       3
                    1 + x + ---- + O(x ) = 1 + x + ---- + O(x )
                             2                      2

    >> series(Exp(x^2 + x), x = 2) = series(exp(x^2 + x), x = 2)

                                                       2
                                      27 Exp(6) (x - 2)             3
          Exp(6) + 5 Exp(6) (x - 2) + ------------------ + O((x - 2) ) =
                                              2

                                                          2
                                         27 exp(6) (x - 2)             3
             exp(6) + 5 exp(6) (x - 2) + ------------------ + O((x - 2) )
                                                 2

    >> series(Exp(x^2 + x), x = 0, 0)

          Error: order too small [ExpSeries]

    >> series(Exp(x^2 + x), x = infinity)

                                         2
                         series(Exp(x + x ), x = infinity)

Another possibility to obtain series expansions of user-defined 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:

#### Changes

- The new options Real and NoWarning were introduced.

Do you have questions or comments?
Copyright © SciFace Software GmbH & Co. KG 2000