Numerical Integration

Calling Sequence


evalf(Int(f, x=a..b, ...))


evalf(Int(f, a..b, ...))


evalf(Int(f, listofequations, ...))


evalf(Int(f, listofranges, ...))


evalf(int(f, x=a..b))






Parameters


f



algebraic expression or procedure; the integrand

x



name; the variable of integration

a, b



endpoints of the interval of integration

listofequations



list of equations [x1=a1..b1, ..., xn=an..bn]

listofranges



list of ranges [a1..b1, ..., an..bn]

...



(optional) zero or more options, as described below





Description


•

The most common command for numerical integration is evalf(Int(f, x=a..b)) where the integration command is expressed in inert form to avoid first invoking the symbolic integration routines. It is also possible to invoke evalf on an unevaluated integral returned by the symbolic int command, as in evalf(int(f, x=a..b)), if it happens that symbolic int fails (returns an unevaluated integral).

•

All numerical integration calling sequences can also be accessed directly from the int command by using the numeric option.

•

You can enter the command evalf/Int using either the 1D or 2D calling sequence. For example, evalf(Int(1/(x^2+1), x=0..infinity)) is equivalent to .

•

The integrand f may be another unevaluated integral, that is, multiple integrals are supported. A special list syntax (see below) can be used to specify multiple integrals, rather than using nested integrals. Integrals expressed in the standard nonlist notation are referred to as 1D (onedimensional) integrals including the case of nested 1D integrals.

•

If the integrand f is specified as a procedure or a Maple operator, then the second argument must be a range a..b and not an equation, that is, a variable of integration must not be specified.

•

Various levels of user information are displayed during the computation if infolevel[`evalf/int`] is assigned values between 1 and 4.



Optional Arguments


•

Additional options may be specified as equations. (For backward compatibility some options are accepted as values rather than equations, as specified below.) An option is one of the following forms:

method = <name> or <name>

digits = <posint> or <posint>

epsilon = <numeric>

methodoptions = <list>

maxintervals = <posint>



•

The specification method = <name> (or simply <name>) indicates a particular numerical integration method to be applied. The methods that can be specified are described below. By default, a hybrid symbolicnumeric strategy is applied.

•

The specification digits = <posint> (or simply <posint>) indicates the number of digits of precision for the computation. Some additional guard digits are carried during the computation to attempt to achieve a result with <posint> correct digits (although a larger tolerance can be specified by using the 'epsilon' option). By default, the Maple environment variable Digits specifies the precision for the computation.

•

The specification epsilon = <numeric> specifies the relative error tolerance for the computed result. The routines attempt to achieve a final result with a relative error less than this value. By default, the relative error tolerance which the routines attempt to achieve for the final result is


where digits is the precision specified for the computation. In attempting to achieve this accuracy, the working value of Digits is increased as deemed necessary. It is an error to specify 'epsilon' smaller than the default value above, and for any value larger than 1e3 the value 1e3 is used instead if the method in use is deterministic (i.e. not the MonteCarlo or Cuba methods).


Note: For some integrands, the numerical accuracy attained when computing values of the integrand may be insufficient to allow the value of the integral to be computed to the default tolerance (even though the computation is using some number of guard digits). In such cases, specifying a larger tolerance (relative to the setting of digits) via the 'epsilon' option may be helpful. Alternatively, increasing Digits and fixing 'epsilon' may provide the desired answer (see the end of the examples section).

•

The specification methodoptions = <list> specifies a list of zero or more options that are specific to a method selected with the method option. In particular, if the _d01ajc or _d01akc method is selected, one can supply an option of the form methodoptions=[maxintervals = <posint>] to specify a maximal number of subintervals that can be used internally by those methods.

•

For backward compatibility, the option maxintervals = <posint> for the _d01ajc and _d01akc methods can also be specified as a separate option, as an argument to Int directly, rather than in the methodoptions option.



Outline of the Numerical Integration Polyalgorithm (1D Integrals)


•

In the default case (no particular method specified), the problem is first passed to NAG integration routines if Digits is not too large (that is, if Digits <= evalhf(Digits)). The NAG routines are in a compiled C library and hence operate at hardware floatingpoint speed. If the NAG routines cannot perform the integration, then some singularity handling may be performed and control may pass back to the NAG routines with a modified problem. Native Maple routines are invoked if the NAG routines cannot solve the problem (for example, if Digits is too large or if the integrand involves functions for which hardware floatingpoint evaluation is not supported).

•

The native Maple hybrid symbolicnumeric solution strategy is as follows. The default numerical method applied is ClenshawCurtis quadrature (_CCquad). If slow convergence is detected then there must be singularities in or near the interval of integration (perhaps in the complex plane). Some techniques of symbolic analysis are used to deal with the singularities. For problems with nonremovable endpoint singularities, an adaptive doubleexponential quadrature method (_Dexp) is applied.

•

If singularities interior to the interval are suspected, then an attempt is made to locate the singularities in order to break up the interval of integration. Finally, if still unsuccessful, then the interval is subdivided and the _Dexp method is applied, or if the method was already _Dexp or _Sinc then an adaptive Gaussian quadrature method (_Gquad) is applied.

•

For the limits of integration, the values infinity and/or infinity are valid, and a symbolicnumeric strategy attempts to deal with singularities. Techniques employed include variable transformations, subtracting out the singularity, and integration of a truncated generalized series near the singularity.

•

No singularity handling is attempted in the case where the integrand f is specified as a procedure or a Maple operator.



Special (List) Syntax for Multiple Integrals


•

A numerical multiple integration problem may be specified in a natural way using nested onedimensional integrals, for example:

evalf( Int(...(Int(Int(f, x1=a1..b1), x2=a2..b2), ...), xn=an..bn) )




where the integrand f depends on x1, x2, ..., xn. Such a problem may also be specified using the following special multiple integration notation with a list as the second argument:

evalf( Int(f, [x1=a1..b1, x2=a2..b2, ..., xn=an..bn]) ) .



•

Additional optional arguments may be stated just as in the case of 1D integration. Also as in 1D integration, the integrand f may be specified as a procedure in which case the second argument must be a list of ranges: [a1..b1, a2..b2, ..., an..bn].

•

Whether a multiple integration problem is stated using nested integrals or using the list notation, the arguments will be extracted so as to invoke the same numerical multiple integration routines.



The Method Names


•

The optional argument method = <name> (or simply <name>) accepts the following method names.

method = _DEFAULT



equivalent to not specifying a method; the solution strategy outlined above is applied for 1D integrals; for multiple integrals, the problem is passed to the _cuhre method and if it fails, then the problem is treated via nested 1D integration.




method = _NoNAG



indicates to avoid calling NAG routines; otherwise follow the _DEFAULT strategy.




method = _NoMultiple



indicates to avoid calling numerical multiple integration routines; compute multiple integrals via nested 1D integration.




Maple Methods


•

Specifying a method indicates to try only that method (in particular, no NAG methods and no singularity handling).

method = _CCquad



ClenshawCurtis quadrature method.




method = _Dexp



adaptive doubleexponential method.




method = _Gquad



adaptive Gaussian quadrature method.




method = _Sinc



adaptive sinc quadrature method.




method = _NCrule



adaptive NewtonCotes method "quanc8". Note that in contrast to the other Maple methods listed here, "quanc8" (method = _NCrule) is a fixedorder method and hence it is not recommended for very high precisions (e.g. Digits > 15).





NAG Methods


•

Specifying a method indicates to try only that method (in particular, no singularity handling and no Maple methods).

method = _d01ajc



for finite interval of integration; allows for badly behaved integrands; uses adaptive Gauss 10point and Kronrod 21point rules.




method = _d01akc



for finite interval of integration, oscillating integrands; uses adaptive Gauss 30point and Kronrod 61point rules.




method = _d01amc



for semiinfinite/infinite interval of integration.





Multiple Integration Methods


•

Specifying a method indicates to try only that method (in particular, do not revert to nested 1D integration).

method = _cuhre



multiple integrals over a hyperrectangle; i.e., the limits of integration are finite constants; dimensions 2 to 15; ACM TOMS Algorithm 698.




method = _MonteCarlo



Monte Carlo method over a hyperrectangle; i.e., the limits of integration are finite constants; for low accuracy only (less than 5 digits of accuracy); NAG routine 'd01gbc'.




method = _CubaVegas



Monte Carlo method over a hyperrectangle; i.e., the limits of integration are finite constants; for low accuracy. For methodspecific options, see evalf/Int/cuba.




method = _CubaSuave



Monte Carlo method over a hyperrectangle; i.e., the limits of integration are finite constants; for low accuracy. For methodspecific options, see evalf/Int/cuba.




method = _CubaDivonne



Monte Carlo method over a hyperrectangle; i.e., the limits of integration are finite constants; for low accuracy. For methodspecific options, see evalf/Int/cuba.




method = _CubaCuhre



Cuhre method over a hyperrectangle; i.e., the limits of integration are finite constants. For methodspecific options, see evalf/Int/cuba.






Examples


>


 (1) 
>


 (2) 
>


 (3) 
The following integrals are computed to higher precision.
>


>


 (4) 
>


>


 (5) 
>


>


 (6) 
>


 (7) 
>


 (8) 
>


 (9) 
The following command returns an error because procedure is invoked with argument , a symbolic name.
>

f := proc(x) if x < 2 then 2*x else x^2 end if; end proc;

 (10) 
>


When the integrand is a procedure, the following syntax should be used.
>


 (11) 
Note that the following command also works by delaying the evaluation of via unevaluation quotes.
>


 (12) 
Multiple integrals may be expressed as nested onedimensional integrals.
>


 (13) 
>


 (14) 
Numerical multiple integration may also be invoked using a list syntax.
>


>


 (15) 
>


 (16) 
When low accuracy is sufficient, the Monte Carlo method may be used.
>


 (17) 
>


 (18) 
Only trust about 3 digits when epsilon = 0.5e2.
>


 (19) 
The following integrand has a region near x=0.5 where evaluation incurs catastrophic cancellation to the extent that the function cannot even be evaluated to 1 significant Digit at standard precision.
>


 (20) 
>


 (21) 
Note that evalf fails to compute this integral with default settings.
>


 (22) 
So to compute the value of this integral to 10 digits, we need to add significant guard digits:
>


 (23) 
>


 (24) 
>


 (25) 
So we add 15 digits to assure we get the answer to 10 digits:
>


 (26) 
In the following example, the default setting for the maximum number of subintervals, in this case, is not enough for successful integration using the NAG method for oscillatory integrands, , which would be suitable for this integrand. By supplying a higher upper bound, we can get successful completion with this method.
>


 (27) 
>


>


 (28) 

