how can I evaluate a four fold integral with four variables numerically ?

i used the functions quad, dblquad, and triplequad to evaluate up to three integrals, but i need to evaluate 4 integrals, how can i do that?

[Muthu Annamalai]

Can you separate the kernel of the integral? This is usally the better way to solve the problem. It also provides a computational speeds up to the solution as well from O(n^4) to say 2*O(n^2) if you can split kernel into a function (product?) of 2 double integrals.

HTH.

[ebtehal]

No, I can't separate the integrals to make them multiplied functions. I was asking because I thought there is a way in Matlab to solve higher order integrals that I didn't know.

Answer 1

[Mike Hosea]

There's one thing that people often find difficult: the requirement that the integrand accepts arrays and returns arrays, operating element-wise.

integral(@(x)integral3(@(y,z,w)f(x,y,z,w),ymin,ymax,zmin,zmax,wmin,wmax),xmin,xmax,'ArrayValued',true)

For smooth integrands over finite regions, a double integral of a double integral is usually a lot faster

integral2(@(x,y)arrayfun(@(x,y)integral2(@(z,w)f(x,y,z,w),zmin,zmax,wmin,wmax),x,y),xmin,xmax,ymin,ymax)

[Mike Hosea]

I got a little tired of giving the same answers over the years, and I'm sure the syntax of these calls is confusing for people, so I finally just wrote something and put it on the file exchange. It's called integralN, and it will set up and call INTEGRAL2 and INTEGRAL3 iteratively to compute order 4, 5, and 6 integrals. It's sort of mechanical to build it to do even or odd orders, but I couldn't quite work out how to generalize setting up the function handles, and I figured it didn't matter anyway, because order 6 is pretty excruciating as it is. There are better ways... But this should work pretty well if you've got the time.

[ebtehal]

Thank you Mike, I will try to understand how to apply your method and use it. The problem is i have a very complex integral where the variables can't be separated, they are multiplied together inside exponential functions and fractions.

[Mike Hosea]

You don't need to separate the variables at all. When you create a function handle like so @(z,w)f(x,y,z,w), it takes a snapshot of the values of x and y at the moment of creation. When that happens nested into the function calls, the creation happens over and over again during the integration as the values of x and y are changed. The function itself can be arbitrarily complicated and inseparable. However, just download integralN. Calling it is just a natural extension of the way interal2 and integral3 are called.

http://www.mathworks.com/matlabcentral/fileexchange/47919-integraln-m

[ebtehal]

Hello Mike, Are the functions integral and integral2 and integral3 already built in the newest version of Matlab, I am working with a bit old version of matlab, I am afraid it doesn't have integral2 or integral3. can i use your method with quad and dblquad or triplequad? Thank you

[Mike Hosea]

You can't use integralN, but you should be able to do this:

quad2d(@(x,y)arrayfun(@(x,y)quad2d(@(z,w)f(x,y,z,w),zmin,zmax,wmin,wmax),x,y),xmin,xmax,ymin,ymax)

If your version is so old that you don't have QUAD2D, you can try this with DBLQUAD, but I don't recommend it. BTW, depending on how x and y are used in f(x,y,z,w), you might need to expand those arguments manually like so:

quad2d(@(x,y)arrayfun(@(x,y)quad2d(@(z,w)f(x*ones(size(z)),y*ones(size(z)),z,w),zmin,zmax,wmin,wmax),x,y),xmin,xmax,ymin,ymax)

Answer 2

[Roger Stafford]

Set up a function which, when given the value of the outermost variable of integration, calculates the triple integral of the inner three iterated integrals for the particular value of that fourth variable. Then take the integral of the value of this newly-defined function over the appropriate limits for that fourth variable. There's nothing very difficult about that. Let the computer do all the hard work. You can expect it to take a fairly long time at it. That's inherent in doing integration in four dimensions.

[ebtehal]

Thank you Roger, I'll consider how to apply your method. The idea is my integration is very complex, and I experienced a another one with triple integration which really took me long time too, and all the integrals are from 0 to infinity, which makes me perform many trials to reach a steady value.