This section contains some functions for numerical analysis, taken from Pari.

Subsections

• Summation of Infinite Series
• Integration

There are three functions for evaluating infinite sums of real numbers. The sum should be specified as a map m from the integers to the free real field, such that m(n) is the n^( th) term of the sum. The summation begins at term i. The precision of the result will be the default precision of the free real field.

InfiniteSum(m, i) : Map, RngIntElt -> FldPrElt

An approximation to the infinite sum m(i) + m(i + 1) + m(i + 2) + ... . This function also works for maps to the free complex field.

An approximation to the infinite sum m(i) + m(i + 1) + m(i + 2) + ... . Designed for series in which every term is positive, it uses van Wijngaarden's trick for converting the series into an alternating one. Due to the stopping criterion, terms equal to 0 will create problems and should be removed.

    Al: MonStgElt                       Default: "Villegas"

An approximation to the infinite sum m(i) + m(i + 1) + m(i + 2) + ... . Designed for series in which the terms alternate in sign. The optional argument Al can be used to specify the algorithm used. The possible values are "Villegas" (the default), and "EulerVanWijngaarden". Due to the stopping criterion, terms equal to 0 will create problems and should be removed.

A number of `Romberg-like' integration methods have been taken from Pari. The precision should not be made too large for this, and singularities are not allowed in the interval of integration (including its boundaries).

Integral(m, a, b) : Map, FldPrElt, FldPRElt -> FldPrElt

    Al: MonStgElt                       Default: Open/Infinite

Given a function f, (at least) continuous everywhere on the open interval (a, b) and defined on [a, b], for free real variables a and b, this function numerically integrates f from a to b. The function f must be defined as a map with as its domain the free real field and as its codomain either the free real or the free complex field. The precision should be not much higher than 28. For improper integrals, where one or both of the limits of integration are plus or minus infinity, the integrand must decrease sufficiently rapidly at infinity (this can often be achieved by integration by parts).

Using the optional argument Al it is possible to choose a special algorithm. With Al := "Infinite" a method tailored for use whe either a or b is `infinite' is used; in that case -a or b could be chosen `large' (how large depends on how rapidly the function tends to 0 at +- Infinity), but it must be the case that a and b have the same sign (and it is suggested that a and b have magnitude at least 1). When neither -a nor b are very large, f is smooth and it can be evaluated everywhere on the closed interval, the method chosen with Al := "Simple" is the fastest. The case Al := "Open" is best suited for cases in which f is undefined (but continuous) at one of the end-points.

The default algorithm depends on the limits of integration. It always uses the "Open" or "Infinite" algorithm (or some combination of these) and so the function need not be defined at a and b.

> R := RealField();
> f := map< R -> R | x :-> Exp(-x^2) >;
> Integral(f, 0, 10);
0.88622692545275801364908374188
> f := map< R -> R | x :-> x^2 >;
> Integral(f, -1, 2);
2.99999999999999999999999999888

> f := map< R -> R | x :-> 1/(x+4)^2 >;
> Integral(f, 1, 10^4000: Al := "Infinite");
0.199999999999999999999999989688
> Integral(f, 1, 10^40: Al := "Infinite");
0.199999999999999999999999989688

> R := RealField();
> C<i> := ComplexField();
> f := map< R -> C | x :-> i*x^2 + 5*x>;
> Integral(f, 0, 10: Al := "Simple");
250.000000000000000000000000 + 333.333333333333333333333333335487*i

> f := map< R -> R | x :-> Sin(x)/x >;
> Integral(f, 0, 1: Al := "Open");
0.94608307036718301494135330891