How to calculate the double integral Numerically in Mathematica?

Question

How to calculate the double integral Numerically in Mathematica ?

Integrate[Exp[-0.099308 s]
       * Integrate[Exp[0.041657423 u] Exp[-3.1413 s + 3.12 u]
             * ((u/(s - u))^(1/2) BesselI[1,2 (u (s - u))^(1/2)]
             + 0.293 BesselI[0,2 (u (s - u))^(1/2)]),{u,0,s}],{s,0,10}]

enter image description here

possible (near) duplicate of [How to calculate the double integral Numerically in MATLAB?](http://stackoverflow.com/questions/9554617/how-to-calculate-the-double-integral-numerically-in-matlab) — High Performance Mark, Mar 06 '12 at 14:28
How could a Mathematica question be a dupe of a Matlab question? That would be like a Java question and a C# question. — Verbeia, Mar 10 '12 at 21:51
@Verbeia you have a point, yet the OP asked how to approach the same integral in both languages. — rcollyer, Mar 11 '12 at 02:50

score 4 · Accepted Answer · edited May 23 '17 at 10:30

Two things. First of all, Integrate accepts multiple "iterators", i.e. {x, x1, x2}, so you can specify a multiple integral without nesting them, as follows

Integrate[x y, {x, 0, 1}, {y, 0, x}]

integrates x y over the triangle bounded by y == x, x == 0, and x == 1. Note, the order of the limits, they go from outer to inner, so the integration is performed from right to left. Then, your integral becomes

Integrate[Exp[-0.099308 s] Exp[0.041657423 u] Exp[-3.1413 s + 3.12 u]
   * ((u/(s - u))^(1/2) BesselI[1,2 (u (s - u))^(1/2)]
   + 0.293 BesselI[0,2 (u (s - u))^(1/2)]),
 {s,0,10}, {u,0,s}]

Second, Mathematica has a number of numerical equivalents to its standard algorithms, like NSolve, NDSolve, NSum, and NIntegrate. They are all identifiable by the leading N, which is a function in its own right, too. The nice thing about these functions is that they have the same signature as their analytical equivalent. So, to numerically integrate your integral, you simply change Integrate to NIntegrate, as follows

NIntegrate[Exp[-0.099308 s] Exp[0.041657423 u] Exp[-3.1413 s + 3.12 u]
   * ((u/(s - u))^(1/2) BesselI[1,2 (u (s - u))^(1/2)]
   + 0.293 BesselI[0,2 (u (s - u))^(1/2)]),
 {s,0,10}, {u,0,s}]

which gives 27.4182, as noted by tkott, but without any warnings generated.

score 1 · Answer 2 · edited Apr 13 '17 at 12:55

If you are looking to numerically integrate something in Mathematica, you should use NIntegrate instead. This will give you a numerical approximation rather than the definite integral.

NIntegrate[ Exp[-0.099308*s]*
  NIntegrate[Exp[0.041657423*u]*(Exp[-3.1413*s + 3.12*u])*((u/(s - u))^(1/2)*
       BesselI[1, 2*(u*(s - u))^(1/2)] + 0.293*BesselI[0, 2*(u*(s - u))^(1/2)])
  , {u, 0, s}]
, {s, 0, 10}]

It complains because the inner integrate is evaluated first, and cannot be evaluated without specifying s, but gives

27.4182

No idea if this is correct or not!

BTW: you might get a better response over at https://mathematica.stackexchange.com/

How to calculate the double integral Numerically in Mathematica?

2 Answers2