1
$\begingroup$

I am working on a problem that has as part of its result an integral of the form

$\text{Integrate}\left[c,\left\{x_1,0,m_1\right\},\left\{x_2,x_1,m_2\right\},\text{...},\left\{x_n,x_{n-1},m_n\right\}\right]$

where

$m_1<=m_2<=... m_{n-1}<=m_n$ and $n>50$, all $m_x$ are non-negative rationals, with $c$ a rational non-negative constant.

So each subsequent integration limit has a lower limit equal to the prior limit's value.

I can evaluate this acceptably quickly using NIntegrate with AdaptiveQuasiMonteCarlo as the method, but I'd prefer if possible to get exact results.

Evaluation of the exact result takes a torturous amount of time as written above. Is there a technique to speed this up in Mathematica?

A small runnable (in reasonable time) example:

Integrate[1, {x1, 0, 2}, {x2, x1, 4}, {x3, x1, 6}, {x4, x1, 8}, {x5, 
  x1, 10}, {x6, x1, 12}, {x7, x1, 14}, {x8, x1, 16}, {x9, x1, 
  18}, {x10, x1, 20}, {x11, x1, 22}, {x12, x1, 24}, {x13, x1, 
  26}, {x14, x1, 28}, {x15, x1, 30}, {x16, x1, 32}, {x17, x1, 
  34}, {x18, x1, 36}, {x19, x1, 38}, {x20, x1, 40}, {x21, x1, 
  42}, {x22, x1, 44}, {x23, x1, 46}, {x24, x1, 48}, {x25, x1, 
  50}, {x26, x1, 52}, {x27, x1, 54}, {x28, x1, 56}, {x29, x1, 
  58}, {x30, x1, 60}, {x31, x1, 62}, {x32, x1, 64}, {x33, x1, 
  66}, {x34, x1, 68}, {x35, x1, 70}, {x36, x1, 72}, {x37, x1, 
  74}, {x38, x1, 76}, {x39, x1, 78}, {x40, x1, 80}, {x41, x1, 
  82}, {x42, x1, 84}, {x43, x1, 86}, {x44, x1, 88}, {x45, x1, 
  90}, {x46, x1, 92}, {x47, x1, 94}, {x48, x1, 96}, {x49, x1, 
  98}, {x50, x1, 100}]

This takes a couple of seconds on my machine, the NIntegrate version takes a few hundredths of a second. I'd be pleased getting the exact result in a few tenths of a second if possible.

Ideas?

Improve this question
$\endgroup$

2 Answers 2

Reset to default
3
$\begingroup$

One idea is to use ImplicitRegion and RegionMeasure. First, convert your variable ranges to an ImplicitRegion:

intervals = {
    {x1, 0, 2}, {x2, x1, 4}, {x3, x1, 6}, {x4, x1, 8}, {x5, x1, 10},
    {x6, x1, 12}, {x7, x1, 14}, {x8, x1, 16}, {x9, x1, 18},
    {x10, x1, 20}, {x11, x1, 22}, {x12, x1, 24}, {x13, x1, 26},
    {x14, x1, 28}, {x15, x1, 30}, {x16, x1, 32}, {x17, x1, 34},
    {x18, x1, 36}, {x19, x1, 38}, {x20, x1, 40}, {x21, x1, 42},
    {x22, x1, 44}, {x23, x1, 46}, {x24, x1, 48}, {x25, x1, 50},
    {x26, x1, 52}, {x27, x1, 54}, {x28, x1, 56}, {x29, x1, 58},
    {x30, x1, 60}, {x31, x1, 62}, {x32, x1, 64}, {x33, x1, 66},
    {x34, x1, 68}, {x35, x1, 70}, {x36, x1, 72}, {x37, x1, 74},
    {x38, x1, 76}, {x39, x1, 78}, {x40, x1, 80}, {x41, x1, 82},
    {x42, x1, 84}, {x43, x1, 86}, {x44, x1, 88}, {x45, x1, 90},
    {x46, x1, 92}, {x47, x1, 94}, {x48, x1, 96}, {x49, x1, 98},
    {x50, x1, 100}
};

reg = ImplicitRegion[
    Evaluate[And @@ Replace[intervals, {a_, l_, r_} -> l < a < r, {1}]],
    Evaluate[Cases[intervals, _Symbol, Infinity] //Union]
]

ImplicitRegion[ 0 < x1 < 2 && x1 < x2 < 4 && x1 < x3 < 6 && x1 < x4 < 8 && x1 < x5 < 10 && x1 < x6 < 12 && x1 < x7 < 14 && x1 < x8 < 16 && x1 < x9 < 18 && x1 < x10 < 20 && x1 < x11 < 22 && x1 < x12 < 24 && x1 < x13 < 26 && x1 < x14 < 28 && x1 < x15 < 30 && x1 < x16 < 32 && x1 < x17 < 34 && x1 < x18 < 36 && x1 < x19 < 38 && x1 < x20 < 40 && x1 < x21 < 42 && x1 < x22 < 44 && x1 < x23 < 46 && x1 < x24 < 48 && x1 < x25 < 50 && x1 < x26 < 52 && x1 < x27 < 54 && x1 < x28 < 56 && x1 < x29 < 58 && x1 < x30 < 60 && x1 < x31 < 62 && x1 < x32 < 64 && x1 < x33 < 66 && x1 < x34 < 68 && x1 < x35 < 70 && x1 < x36 < 72 && x1 < x37 < 74 && x1 < x38 < 76 && x1 < x39 < 78 && x1 < x40 < 80 && x1 < x41 < 82 && x1 < x42 < 84 && x1 < x43 < 86 && x1 < x44 < 88 && x1 < x45 < 90 && x1 < x46 < 92 && x1 < x47 < 94 && x1 < x48 < 96 && x1 < x49 < 98 && x1 < x50 < 100, {x1, x10, x11, x12, x13, x14, x15, x16, x17, x18, x19, x2, x20, x21, x22, x23, x24, x25, x26, x27, x28, x29, x3, x30, x31, x32, x33, x34, x35, x36, x37, x38, x39, x4, x40, x41, x42, x43, x44, x45, x46, x47, x48, x49, x5, x50, x6, x7, x8, x9}]

Now, use RegionMeasure:

res = RegionMeasure @ reg; //AbsoluteTiming
res
N @ res

{2.09289, Null}

2965552832547973527300539562780813113506017776171989352778554799616984916548386816/325

9.12478*10^78

Compare to the numerical answer:

NIntegrate[1, {x1, 0, 2}, {x2, x1, 4}, {x3, x1, 6}, {x4, x1, 8}, {x5, 
x1, 10}, {x6, x1, 12}, {x7, x1, 14}, {x8, x1, 16}, {x9, x1, 
18}, {x10, x1, 20}, {x11, x1, 22}, {x12, x1, 24}, {x13, x1, 
26}, {x14, x1, 28}, {x15, x1, 30}, {x16, x1, 32}, {x17, x1, 
34}, {x18, x1, 36}, {x19, x1, 38}, {x20, x1, 40}, {x21, x1, 
42}, {x22, x1, 44}, {x23, x1, 46}, {x24, x1, 48}, {x25, x1, 
50}, {x26, x1, 52}, {x27, x1, 54}, {x28, x1, 56}, {x29, x1, 
58}, {x30, x1, 60}, {x31, x1, 62}, {x32, x1, 64}, {x33, x1, 
66}, {x34, x1, 68}, {x35, x1, 70}, {x36, x1, 72}, {x37, x1, 
74}, {x38, x1, 76}, {x39, x1, 78}, {x40, x1, 80}, {x41, x1, 
82}, {x42, x1, 84}, {x43, x1, 86}, {x44, x1, 88}, {x45, x1, 
90}, {x46, x1, 92}, {x47, x1, 94}, {x48, x1, 96}, {x49, x1, 
98}, {x50, x1, 100}]

9.15451*10^78

Since the numerical answer uses MonteCarlo methods, it fluctuates every time it is evaluated.

Addendum

For this particular integral, all of the integration variables depend on the first integration, so you can just do:

Integrate[Times @@ (Range[4, 100, 2] - x), {x, 0, 2}] // AbsoluteTiming

{0.254912, 2965552832547973527300539562780813113506017776171989352778554799616984916548386816/325}

Improve this answer
$\endgroup$
6
  • $\begingroup$ +1 for interesting idea. This however takes about the same time (within a few tenths) of 2+ seconds. I'm curious what the timing is on my example in your environment: my laptop is no slouch, I find it hard to fathom your environment is nearly 1000x faster. $\endgroup$
    – ciao
    Nov 1, 2019 at 22:24
  • $\begingroup$ @ciao It takes about the same amount of time as yours. I think it is possible to come up with a symbolic result, and avoid integration entirely. Would that be better? $\endgroup$
    – Carl Woll
    Nov 1, 2019 at 22:27
  • $\begingroup$ A symbolic result could be interesting - but I'm still curious about the timings, unless I'm misreading your post & comments: On my machine, my integral and the regionmeasure results take ~2 seconds, but on your environment, the integral takes ~2 seconds but the regionmeasure takes ~0.03 seconds? Is that correct? $\endgroup$
    – ciao
    Nov 1, 2019 at 22:31
  • $\begingroup$ @ciao I see your question. It gets cached, so I must have put in the cached result by mistake. The first time it executes it takes ~2 seconds, but the second time it is much faster. $\endgroup$
    – Carl Woll
    Nov 1, 2019 at 22:36
  • $\begingroup$ Your observation in the addendum does the trick. Accepted, and thanks! $\endgroup$
    – ciao
    Nov 2, 2019 at 1:39
2
$\begingroup$

To make comparable timings, the definition of the ImplicitRegion should be included in the timing of the RegionMeasure approach.

Clear["Global`*"];

$HistoryLength = 0;

res1 = N[Integrate[1, 
    Sequence @@ 
     Table[{x[k], x[1] (1 - KroneckerDelta[k, 1]), 2 k}, {k, 50}]]] // 
  AbsoluteTiming

(* {2.89516, 9.12478*10^78} *)

res2 = Module[{reg = 
     ImplicitRegion[
      And @@ Table[x[1] (1 - KroneckerDelta[k, 1]) < x[k] < 2 k, {k, 50}], 
      Evaluate@Array[x, 50]]},
   RegionMeasure[reg] // N] // AbsoluteTiming

(* {2.47009, 9.12478*10^78} *)

Integrate is about 17% slower than RegionMeasure

res1[[1]]/res2[[1]]

1.17208

The results are identical

res1[[-1]] === res2[[-1]]

(* True *)
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.