multidimensional integration

multidimensional integration

[Harald Wiedemann]

Hello everybody,

I am using gsl since recently and now I was looking for
a routine for integration on a multidimensional square
(in my case the dimension was just 2), but in the docs
I found only routines for 1d-integration. Do I have to
integrate my function in two steps, then? Or did I miss
the right routines in the documentation?

With best regards,

Harald Wiedemann

Re: multidimensional integration

[Wartan Hachaturow]

On Mon, Jul 15, 2002 at 02:06:58PM +0200, Harald Wiedemann wrote:

> a routine for integration on a multidimensional square
> (in my case the dimension was just 2), but in the docs
> I found only routines for 1d-integration. Do I have to

Try using Monte Carlo integration, gsl_monte_*_integration().

-- 
Regards, Wartan.
"Computers are not intelligent. They only think they are."

Re: multidimensional integration

[Lukas Dobrek]

On Mon, Jul 15, 2002 at 02:06:58PM +0200, Harald Wiedemann wrote:
> Hello everybody,
> 
> I am using gsl since recently and now I was looking for
> a routine for integration on a multidimensional square
> (in my case the dimension was just 2), but in the docs
> I found only routines for 1d-integration. Do I have to
> integrate my function in two steps, then? Or did I miss
> the right routines in the documentation?

There is no. You can do it by compose integration or using monte 
carlo (This is supported).  The secound methiod of course has its
limitations.   
> With best regards,
Lukasz
-- 
£ukasz Dobrek
   An optimist believes that we live in the best of all possible worlds.
   A pessimist is sure that this must be so.
http://www.pld-linux.org

multidimensional integration

[Robert G. Brown]

[-- Attachment #1: Type: text/plain, Size: 2333 bytes --]

Hi fellow GSL'ers.

We have a postdoc in our department who is preparing to integrate
something. In his previous position at another place, he used NAG to do
this, and has the requisite code already in place.  He requested that we
buy and install a single copy of NAG just for him and a student to be
able to use this one routine to do this one integral on just one
computer, at a cost of many hundreds of dollars.

I suggested that he look into using the GSL instead, since it is a very
high-quality library to my own direct experience and of course is both
free and universally installed in our department.  GSL and NAG both use
QUADPACK as the basis for their 1D integrals (and have nearly identical
call structure) so I figured that the transition would actually be
painless.

However, the integrand he has to integrate is actually defined and
integrated over somewhere between 5 to 7 dimensions (with rectangular
limits).  The routine he used from NAG was actually d01fcc, which is NOT
from QUADPACK but rather implements the multidimensional adaptive
routine HALF with a custom interval rule.  When I looked at GSL's online
manual (version 1.6 as of this last December) I didn't see a
multidimensional integration routine equivalent to d01fcc.

SO, questions:

  a) Is a multidimensional integration routine equivalent to d01fcc
implemented or under development, and if so, where is it and/or how do I
get a version that has it?  I looked at the CVS tree and didn't
immediately see one.  In principle I could probably use e.g. a
multidimentional ODE solver but I'd think that having a d01fcc
equivalent would be much more efficient.

  b) If not, does anybody have any suggestions on the "best" way to
attack this sort of integral using existing tools?  At five dimensions I
suspect that just calling 1 dim integrations five levels deep would
result in an awful lot of wasted energy and time.  Framing it as an ODE
set also seems like it would work but likely not be terribly efficient
or terribly easy to control error-wise.

  c) On a related note, has anybody done a head-to-head performance
comparison of GSL with NAG -- either time/efficiency performance or
numerical accuracy type performance?  This isn't a significant issue on
this particular project but is an issue that I expect to see come up in
the future.

   rgb

[-- Attachment #2: Type: application/pgp-signature, Size: 189 bytes --]

Re: multidimensional integration

[Manoj Warrier]

Is not Monte-Carlo integration the recommended method for
multi-dimensional integration?
The GSL specific routines for this are described at:
http://www.gnu.org/software/gsl/manual/gsl-ref_23.html#SEC371

I think your point (b) below hits the nail on the head.
I have not used these techniques, and would like to hear other,
more expert opinion.

Best Regards.
Manoj

On Wed, 29 Jun 2005, Robert G. Brown wrote:

> Hi fellow GSL'ers.
> 
> We have a postdoc in our department who is preparing to integrate
> something. In his previous position at another place, he used NAG to do
> this, and has the requisite code already in place.  He requested that we
> buy and install a single copy of NAG just for him and a student to be
> able to use this one routine to do this one integral on just one
> computer, at a cost of many hundreds of dollars.
> 
> I suggested that he look into using the GSL instead, since it is a very
> high-quality library to my own direct experience and of course is both
> free and universally installed in our department.  GSL and NAG both use
> QUADPACK as the basis for their 1D integrals (and have nearly identical
> call structure) so I figured that the transition would actually be
> painless.
> 
> However, the integrand he has to integrate is actually defined and
> integrated over somewhere between 5 to 7 dimensions (with rectangular
> limits).  The routine he used from NAG was actually d01fcc, which is NOT
> from QUADPACK but rather implements the multidimensional adaptive
> routine HALF with a custom interval rule.  When I looked at GSL's online
> manual (version 1.6 as of this last December) I didn't see a
> multidimensional integration routine equivalent to d01fcc.
> 
> SO, questions:
> 
>   a) Is a multidimensional integration routine equivalent to d01fcc
> implemented or under development, and if so, where is it and/or how do I
> get a version that has it?  I looked at the CVS tree and didn't
> immediately see one.  In principle I could probably use e.g. a
> multidimentional ODE solver but I'd think that having a d01fcc
> equivalent would be much more efficient.
> 
>   b) If not, does anybody have any suggestions on the "best" way to
> attack this sort of integral using existing tools?  At five dimensions I
> suspect that just calling 1 dim integrations five levels deep would
> result in an awful lot of wasted energy and time.  Framing it as an ODE
> set also seems like it would work but likely not be terribly efficient
> or terribly easy to control error-wise.
> 
>   c) On a related note, has anybody done a head-to-head performance
> comparison of GSL with NAG -- either time/efficiency performance or
> numerical accuracy type performance?  This isn't a significant issue on
> this particular project but is an issue that I expect to see come up in
> the future.
> 
>    rgb
> 

-- 
Manoj Warrier (manoj.warrier@ipp.mpg.de)

Stellaratortheorie, Max-Planck Institut Fur Plasmaphysik
TeilInstitut Greifswald Wendelsteinstrasse 1
D-17491 Greifswald Germany Tel: +49-3834-882434

--------- History of Computing 10-11-3003 ---------------
Then there used to be this great user friendly OS which
overwrote your MBR whenever you installed it.
---------------------------------------------------------

Re: multidimensional integration

[Robert G. Brown]

[-- Attachment #1: Type: text/plain, Size: 5753 bytes --]

Manoj Warrier writes:

> Is not Monte-Carlo integration the recommended method for
> multi-dimensional integration?
> The GSL specific routines for this are described at:
> http://www.gnu.org/software/gsl/manual/gsl-ref_23.html#SEC371

I >>think<< that this is in an edge case of sorts.  MC is definitely
right for very high dimensionalities (say, 100 dimensions).  In those
cases it is the only game in down since even bisecting each dimension a
single time leaves one with 2^100 cells to sum over, which is "bad".  In
1 dimension, though, it is very slow to converge compared to more
analytic and adaptive methods, especially for smooth integrands with
support in only certain parts of the integration region.

In between there is clearly a crossover.  In (say) five dimensions
bisection introduces 2^5 or 32 cells -- covering any one of those cells
more finely in an adaptive mesh isn't THAT expensive in its scaling if
the function has reasonable support and smoothness.  Adaptive methods in
this case give you the benefit of good control over errors provided only
that significant "features" if the integrand are large enough to be
picked up by the initial meshes.

> I think your point (b) below hits the nail on the head.
> I have not used these techniques, and would like to hear other,
> more expert opinion.

The person who wants to do this problem has already indicated that MC
methods are very slow and not horribly accurate for this problem,
suggesting that it is on the multidimensional mesh side of the
crossover.  In any event, even in 2 dimensions it is likely desirable to
have a real multidimensional integration routine instead of just calling
a 1 d routine twice, once inside the other.  In 2d the refinement mesh
would have to cover the regions where the integrand is supported in
stripes that reach all the way across the long dimensions, resulting in
lots of calls at points where the integrand has little support -- this
gets to be LOTS of calls in signficant subvolumes in d = 3 or 4 or 5 as
one has to refine the mesh across all dimensions where the support has a
signficant projection, even if the support itself has relatively little
volume.

I'm guessing that HALF does something like mesh refinement by halving
the scale (doubling the mesh) successively in subvolumes only to the
extent required to get each successive subvolume to converge, so that a
coarse mesh is used over all regions where the integrand is smooth and
small but a fine one covers regions where it varies rapidly and/or is
large.  How smart it is with respect to interpolating functions or the
like, I don't know.

It sounds like there is nothing like this in the GSL, though, for the
intermediate/low dimensions where it would solidly beat MC.

    rgb

> 
> Best Regards.
> Manoj
> 
> On Wed, 29 Jun 2005, Robert G. Brown wrote:
> 
>> Hi fellow GSL'ers.
>> 
>> We have a postdoc in our department who is preparing to integrate
>> something. In his previous position at another place, he used NAG to do
>> this, and has the requisite code already in place.  He requested that we
>> buy and install a single copy of NAG just for him and a student to be
>> able to use this one routine to do this one integral on just one
>> computer, at a cost of many hundreds of dollars.
>> 
>> I suggested that he look into using the GSL instead, since it is a very
>> high-quality library to my own direct experience and of course is both
>> free and universally installed in our department.  GSL and NAG both use
>> QUADPACK as the basis for their 1D integrals (and have nearly identical
>> call structure) so I figured that the transition would actually be
>> painless.
>> 
>> However, the integrand he has to integrate is actually defined and
>> integrated over somewhere between 5 to 7 dimensions (with rectangular
>> limits).  The routine he used from NAG was actually d01fcc, which is NOT
>> from QUADPACK but rather implements the multidimensional adaptive
>> routine HALF with a custom interval rule.  When I looked at GSL's online
>> manual (version 1.6 as of this last December) I didn't see a
>> multidimensional integration routine equivalent to d01fcc.
>> 
>> SO, questions:
>> 
>>   a) Is a multidimensional integration routine equivalent to d01fcc
>> implemented or under development, and if so, where is it and/or how do I
>> get a version that has it?  I looked at the CVS tree and didn't
>> immediately see one.  In principle I could probably use e.g. a
>> multidimentional ODE solver but I'd think that having a d01fcc
>> equivalent would be much more efficient.
>> 
>>   b) If not, does anybody have any suggestions on the "best" way to
>> attack this sort of integral using existing tools?  At five dimensions I
>> suspect that just calling 1 dim integrations five levels deep would
>> result in an awful lot of wasted energy and time.  Framing it as an ODE
>> set also seems like it would work but likely not be terribly efficient
>> or terribly easy to control error-wise.
>> 
>>   c) On a related note, has anybody done a head-to-head performance
>> comparison of GSL with NAG -- either time/efficiency performance or
>> numerical accuracy type performance?  This isn't a significant issue on
>> this particular project but is an issue that I expect to see come up in
>> the future.
>> 
>>    rgb
>> 
> 
> -- 
> Manoj Warrier (manoj.warrier@ipp.mpg.de)
> 
> Stellaratortheorie, Max-Planck Institut Fur Plasmaphysik
> TeilInstitut Greifswald Wendelsteinstrasse 1
> D-17491 Greifswald Germany Tel: +49-3834-882434
> 
> --------- History of Computing 10-11-3003 ---------------
> Then there used to be this great user friendly OS which
> overwrote your MBR whenever you installed it.
> ---------------------------------------------------------

[-- Attachment #2: Type: application/pgp-signature, Size: 189 bytes --]

Re: multidimensional integration

[Brian Gough]

Robert G. Brown writes:
 > However, the integrand he has to integrate is actually defined and
 > integrated over somewhere between 5 to 7 dimensions (with rectangular
 > limits).  The routine he used from NAG was actually d01fcc, which is NOT
 > from QUADPACK but rather implements the multidimensional adaptive
 > routine HALF with a custom interval rule.  When I looked at GSL's online
 > manual (version 1.6 as of this last December) I didn't see a
 > multidimensional integration routine equivalent to d01fcc.

See Steven Johnson's recent post for adaptive multidimensional integration.
http://sources.redhat.com/ml/gsl-discuss/2005-q2/msg00020.html

Otherwise MISER or VEGAS are the only existing options.

Or take out a support contract and I'll write a replacement for d01fcc ;-)

-- 
Brian Gough

Network Theory Ltd,
Commercial support for GSL --- http://www.network-theory.co.uk/gsl/

Re: multidimensional integration

[Giulio Bottazzi]

[-- Attachment #1: Type: text/plain, Size: 224 bytes --]

Hello,
with last CVS

FAIL:   SV_decomp_jacobi (3x3) A=[ 0, 1, 1; 0, 1, 1; 0, 1, 1] [67731]
...
FAIL:   SV_decomp_jacobi (5x5) case=33489456 [167015]
FAIL: test

basically all tests on SV_decomp_jacobi FAIL.

	G.

[-- Attachment #2: Type: application/pgp-signature, Size: 189 bytes --]