[R] bivariate vector numerical integration with infinite range

[Hans W Borchers]

baptiste auguie <baptiste.auguie <at> googlemail.com> writes:

> 
> Dear list,
> 
> I'm seeking some advice regarding a particular numerical integration I
> wish to perform.
> 
> The integrand f takes two real arguments x and y and returns a vector
> of constant length N. The range of integration is [0, infty) for x and
> [a,b] (finite) for y. Since the integrand has values in R^N I did not
> find a built-in function to perform numerical quadrature, so I wrote
> my own after some inspiration from a post in R-help,

The function adaptIntegral() in the 'cubature' package integrates 
multi-valued functions over n-dimensional finite hypercubes, as do the
functions in 'R2Cuba'.
If the hypercube is (partly) infinite, a transformation such as x --> 1/x
per infinite axis (as in NR) has to be applied.

For two dimensions, another approach could be to apply the integrate()
function twice, because this 1-dimensional integration function can handle
infinite intervals.
Hint: First inegrate over the finite dimension(s).

Example: Integrate sin(x)/exp(y) for 0 <= x <= pi, 0 <= y <= Inf (value: 2)
----
    f1 <- function(y) 1/exp(y)
    f2 <- function(x) sin(x) * integrate(f1, 0, Inf)$value
    integrate(f2, 0, pi)
    # 2 with absolute error < 2.2e-14
----
Note that the absolute error is not correct, as the first integration has
its own error term. You have to do your own error estimation.

Hans Werner

> 
> [...]
> 
> So it seems to work. I wonder though if I may have missed an easier
> (and more reliable) way to perform such integration using base
> functions or an add-on package that I may have overlooked.
> 
> Best regards,
> 
> baptiste
> 
>