4
$\begingroup$

To find the volume of the following region:

fn[x_, y_, z_]:= Abs[0.7*x*Exp[I*y] + 0.3*Sqrt[x^2 + 8*10^-5] 
 + Sqrt[x^2 + 3*10^-3]*0.02*Exp[I*z]]

R = ImplicitRegion[fn[x, y, z]<=3*10^-3, {{x, 0, 0.015}, {y, 0, 2*Pi}, {z, 0, 2*Pi}}]

RegionPlot3D[
 fn[x, y, z] <= 3*10^-3, {x, 0, 0.015}, {y, 0, 2*Pi}, {z, 0, 2*Pi}, 
 PlotPoints -> 50, AxesLabel -> {x, y, z}, 
 PlotStyle -> Directive[Yellow, Opacity[0.5]], Mesh -> None]

Plot of region of interest

Using Volume:

Volume[R]

Volume::nmet: Unable to compute the volume of region

Using NIntegrate:

NIntegrate[
  Boole[fn[x, y, z] <= 3*10^-3], {x, 0, 0.015}, {y, 0, 2*Pi}, {z, 0, 
   2*Pi}, Method -> {"GlobalAdaptive", "MaxErrorIncreases" -> 50000}, 
  PrecisionGoal -> 3] // Timing

NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small.

NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 50000 times. The global error is expected to decrease monotonically after a number of integrand evaluations... NIntegrate obtained 0.12219642052793653 and 0.007515502934285486 for the integral and error estimates.

{19.25, 0.122196}

The result does not converge.

Using MonteCarlo:

NIntegrate[
  Boole[fn[x, y, z] <= 3*10^-3], {x, 0, 0.015}, {y, 0, 2*Pi}, {z, 0, 
   2*Pi}, Method -> {"MonteCarlo", "MaxPoints" -> 10^12, 
    "RandomSeed" -> 9}, PrecisionGoal -> 3] // Timing

{8.39063, 0.121479}

This is slow, and if PrecisionGoal -> 5 is set,

NIntegrate[
 Boole[fn[x, y, z] <= 3*10^-3], {x, 0, 0.015}, {y, 0, 2*Pi}, {z, 0, 
  2*Pi}, Method -> {"MonteCarlo", "MaxPoints" -> 10^12, 
   "RandomSeed" -> 9}, PrecisionGoal -> 5]

It does not give a result even after running for more than 2 minutes.

I hope to find a better way to calculate the volume of this region with:

  1. High accuracy: at least PrecisionGoal -> 5
  2. Good error estimate
  3. Fast
|improve this question
$\endgroup$
  • $\begingroup$ Usually, I would try Volume[DiscretizeRegion[R,MaxCellMeasure->0.01]] with different values for MaxCellMeasure. However, it fails in this case for some reason... $\endgroup$ – Henrik Schumacher Feb 10 '18 at 8:35
  • $\begingroup$ Try to Rationalize the expression you use in fn and use symbolic bounds over numerical constants whenever possible (e.g. 15/1000 instead of 0.015). This often helps, because it enables arbitrary precision results without loss of precision/accuracy, which can cause problems (slow convergence warnings) in numerical methods like NIntegrate. It also helps avoid ugly results when a symbolic integration is possible. $\endgroup$ – Thies Heidecke Feb 10 '18 at 9:20
  • $\begingroup$ @Henrik Schumacher. Yeah, I also found that DiscretizeRegion[R] failed in this case. $\endgroup$ – K.D Feb 10 '18 at 12:26
  • $\begingroup$ @Thies Heidecke Thanks for your advice. I rationalized the expression and found the result is same as the above. $\endgroup$ – K.D Feb 10 '18 at 12:46
5
$\begingroup$

I think the main issue is the axes bounds are quite disproportionate and that's effecting the sampling.

Here's your region scaled to the unit cube:

R2 = ImplicitRegion[fn[3 x/200, 2π y, 2π z] <= 3*10^-3, {{x, 0, 1}, {y, 0, 1}, {z, 0, 1}}];

reg = BoundaryDiscretizeRegion[R2, {{0, 1}, {0, 1}, {0, 1}}, MaxCellMeasure -> .0001]

enter image description here

0.015 * (2π)^2 * Volume[reg]

0.121605

If we clean up fn a bit, we get a noticeable speedup in BoundaryDiscretizeRegion:

f2[x_, y_, z_] = ComplexExpand[Rationalize[fn[x, y, z]]]^2;

R3 = ImplicitRegion[f2[3 x/200, 2π y, 2π z] <= (3*10^-3)^2, {{x, 0, 1}, {y, 0, 1}, {z, 0, 1}}];

Compare:

timeDisc[r_] := First[AbsoluteTiming[BoundaryDiscretizeRegion[r, {{0, 1}, {0, 1}, {0, 1}}, MaxCellMeasure -> .0001]]];

{timeDisc[R2], timeDisc[R3]}

{1.92265, 1.18681}

Note that this approach works on the original unscaled region, but there are strange fringes and it's very slow:

BoundaryDiscretizeRegion[R, {{0, 0.015}, {0, 2π}, {0, 2π}}, 
  BoxRatios -> {1, 1, 1}, MaxCellMeasure -> .0001]; // AbsoluteTiming

enter image description here

|improve this answer
$\endgroup$
  • $\begingroup$ Yes, that is what I got without rescaling. I am wondering why those defects disappear doing that... +1 $\endgroup$ – José Antonio Díaz Navas Feb 10 '18 at 22:03
  • $\begingroup$ My guess is since the y and z bounds are so much larger, they were sample more. The fact that the fringes are parallel to the x axis seem (to me) to support this. $\endgroup$ – Chip Hurst Feb 10 '18 at 22:09
  • $\begingroup$ @ José Antonio Díaz Navas @ Chip Hurst Thanks for your suggestions. It is true that If we set the MaxCellMeasure smaller, we can get a result of greater accuracy. Can I know the total error estimate of the volume when using the method of BoundaryDiscretizeRegion and Volume? $\endgroup$ – K.D Feb 11 '18 at 3:44
  • $\begingroup$ @Chip Hurst I tried the timeDisc some times with MMA 10.3 and found the result varied from about 10s to 50s even though I used ClearSystemCache before timeDisc, so I failed to compare timeDisc[R2] and timeDisc[R3]. Does it happen in your computer? $\endgroup$ – K.D Feb 11 '18 at 7:09
  • $\begingroup$ I am curious about what kind of algorithm is used by the Volume when the input of Volume is the output of BoundaryDiscretizeRegion as in the above case. $\endgroup$ – K.D Feb 11 '18 at 9:53
2
$\begingroup$

I have tried this with MMA 11.2 (macOs 10.13.3) in my old Intel I3. First, I have simplified a bit your fn as it involves Abs, and by setting the definition of ImplicitFucntion in a standard form:

fn[x_, y_, z_] := 8.4*10^-6 + 0.5804 x^2 + 0.00375659 x Sqrt[1 + 12500 x^2] Cos[y] + 
0.000885438 Sqrt[3 + 1000 x^2] (1. x Cos[y - z] + 0.00383326 Sqrt[1 + 12500 x^2] Cos[z])

R = ImplicitRegion[fn[x, y, z] <= 9/1000000 && 
                   0 <= x <= 0.015 && 0 <= y <= 2 Pi && 0 <= z <= 2 Pi, {x, y, z}];

NIntegrate[1, {x, y, z} \[Element] R,Method -> {"MonteCarlo", "MaxPoints" -> 10^6}, 
PrecisionGoal -> 5]

The high number of points used is to provide the required PrecisionGoal, and after 11.5 sec:

NIntegrate::maxp: The integral failed to converge after 1000100 integrand evaluations. 
NIntegrate obtained 0.12154680553578552` and 0.00023916058934647544` for the integral 
and error estimates.

(* 0.121547 *)

which is slow, and not enough.

However, by using BoundaryDiscretizeRegion, we can use Volume. This is quite slow, as we need to represent the region with enough accuracy (see also this), but provides a reference value for the volume calculated by using NIntegrate: 

R2 = BoundaryDiscretizeRegion[R, BoxRatios -> {1, 1, 1}, MaxCellMeasure -> 0.000075];
Volume[R2]

(* 0.121502 *)

enter image description here

For MaxCellMeasure->0.00006 the volume is 0.121509.

By means of:

NIntegrate[Boole[fn[x, y, z] <= 9/1000000], {x, 0, 0.005, 0.01, 0.015}, 
{y, 0, 2/3 Pi, 4/3 Pi, 2*Pi}, {z, 0, 2/3 Pi, 4/3 Pi, 2*Pi}, 
Method -> {"MonteCarlo", "MaxPoints" -> 10^7, "RandomSeed" -> 9}, 
PrecisionGoal -> 5] // AbsoluteTiming

I get (again, very slow):

NIntegrate::maxp: The integral failed to converge after 10000100 integrand evaluations. 
NIntegrate obtained 0.12149430710916392` and 0.000041400109821011765` for the integral 
and error estimates.

(* {20.1802, 0.121494} *)

Well, It seems that your problem is not so easy to make it more efficient. I always used in this case the MonteCarlo method as the faster one with NIntegrate. However, an expert on numerical integration can explore the multiple options to obtain a proper combination of them to speed up the calculation.

The use of Volume seems to be the faster as long as the process in obtaining the BoundaryDiscretizeRegion will be also fast. It was so to achieve the required PrecisionGoal in my computer, compared to NIntegrate .

Hope this helps.

|improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.