Message Boards Message Boards

1
|
2465 Views
|
2 Replies
|
3 Total Likes
View groups...
Share
Share this post:

Integrate recusive function.

Posted 9 years ago

I'm trying to understand how to calculate the volume of a recursive function containing an integral. To be precise; it is the volume V[n,r] of an n-dimensional sphere with radius r expressed as an integral of the expression for an (n-1)-dimensional sphere with radius y, i.e. V[n-1,y]. My first attemt was simply the following:

V[n_, r_] := If[n == 1, 2 r,
Assuming[r > 0, 2*Integrate[V[n - 1, Sqrt[r^2 - x^2]]
, {x, 0, r}]]
]  

But it only works for n =1 and n=2, and returns undefined for n=3. After some searching on the web I found out that the problem could be that the radius from differens steps in the recursion gets mixed up. A fix for that is to introduce a module, like this:

V[n_, r_] := V[n, r] = If[n == 1, 2 r,
   Assuming[r > 0,
    Module[{x}, 2*Integrate[V[n - 1, Sqrt[r^2 - x^2]], {x, 0, r}]
     ]
    ]
   ]

Now the problem is the speed, V[2,r] is calulated in a blink but V[3,r] takes several seconds despite my attemt to use memoization. Anyone who can explain how to get this to work?

regards, Robert

POSTED BY: Robert Johansson
2 Replies

In addition to isolating the integration variables from one another, it might help to force an evaluation of the integrand. Maybe.

vV2[n_, r_] := 
 If[n == 1, 2 r, 
  2*Integrate[Evaluate[vV2[n - 1, Sqrt[r^2 - x[n]^2]]], {x[n], 0, r}, 
    Assumptions -> r > 0]]

AbsoluteTiming[vV2[4, r]]                                               

                    2  4
                  Pi  r
Out[2]= {19.8016, ------}
                    2

This might or might not be a small improvement. And of course it is still slower as we increase n.

Better still is to remove the complexity from the radial parameter. For this one can use a new radius and replace after integrating.

vV3[n_, r_] := 
 If[n == 1, 2 r, 
  2*Module[{x, r2}, 
    Integrate[
     Evaluate[vV3[n - 1, r2] /. r2 -> Sqrt[r^2 - x[n]^2]], {x[n], 0, 
      r}, Assumptions -> r > 0]]]

We get a nice improvement now.

AbsoluteTiming[vV3[6, r]]

(* Out[368]= {7.30655, (\[Pi]^3 r^6)/6} *)
POSTED BY: Daniel Lichtblau

Thanks Daniel, it´s a nice solution. I like the way you separate the integration variables without the using a module. Your solution is faster than my change of variables by hand. I guess one should leave that to Mathematica :-)

V[n_, r_] := If[n == 1, 2 r,
  Assuming[r > 0,
   2*Integrate[
     Evaluate[x[n]*V[n - 1, x[n]]/Sqrt[r^2 - x[n]^2]], {x[n], 0, r}]
   ]
  ]

AbsoluteTiming[V[6, r]]
{15.856, (\[Pi]^3 r^6)/6}

AbsoluteTiming[vV3[6, r]]
{7.24475, (\[Pi]^3 r^6)/6}
POSTED BY: Robert Johansson
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract