# Nintegrate and non numerical values

I would like to calculate the next integral for different values to the parameter $w$

$\int_{0}^{\infty}{du}\int_{0}^{\infty}{dv} \frac{u}{u^2+v^2}\left(\int_{0}^{0.8}\int_{0}^{\infty} xe^{-2w\sqrt{x^2+y^2}} J_{0}(ux)[cos(\pi y/0.8)]^2 cos(yv) {dy}{dx}\right)^2$

I am trying first to define the double integral Inside of the square:

    f[u_, v_, w_?NumericQ] :=NIntegrate[x*Exp[-2*w*Sqrt[x^2 + y^2]]*BesselJ[0, u*x]*Cos[y*v]*(Cos[Pi*y/0.30040290])^2, {x, 0, Infinity}, {y, 0, 0.30040290/2},Method -> "QuasiMonteCarlo"]


and after to integrate this expression with respect to $u$ and $v$ for a given value of $w$, for example 4

    NIntegrate[ (u/(u^2 + v^2))*((f[u_, v_,4])^2), {u, 0, Infinity}, {v, 0,  Infinity}, Method -> "QuasiMonteCarlo"]


Of course something is wrong in my reasoning because I get this warning message:

    NIntegrate::inumr: The integrand E^(-4 Sqrt[x^2+y^2]) x BesselJ[0,x u_] Cos[10.4579 y]^2 Cos[y v_] has evaluated to non-numerical values for all sampling points in the region with boundaries {{0,1},{0,1}}.


Thank you for all your help

• Do Clear[f] . You probably have an old def for f without the NumericQ Mar 18, 2018 at 13:58
• Just definef to have the signature f[u_?NumericQ, v_?NumericQ, w_?NumericQ]. Mar 18, 2018 at 14:34
• And for the second integration? What is the sintaxis to integrate the square of this f?
– rafa
Mar 18, 2018 at 14:38
• please fix the noted errors and update the question with whatever new problem you have Mar 18, 2018 at 15:21
• "example 4" suggests there's a reference for this problem. If so, what is it? Mar 18, 2018 at 15:23

The following definitions produce no messages, but the computations are slow and I did not wait for NIntegrate to produce a result.

Clear[f]
f[u_?NumericQ, v_?NumericQ, w_?NumericQ] :=
NIntegrate[
x*Exp[-2*w*Sqrt[x^2 + y^2]]*BesselJ[0, u*x]*
Cos[y*v]*(Cos[Pi*y/0.30040290])^2, {x, 0, Infinity}, {y, 0,
0.30040290/2}, Method -> "GlobalAdaptive", MaxRecursion -> 2,
AccuracyGoal -> 2, PrecisionGoal -> 2];

NIntegrate[(u/(u^2 + v^2))*((f[u, v, 4])^2), {u, 0, Infinity}, {v, 0,
Infinity}, Method -> "GlobalAdaptive", MaxRecursion -> 2,
PrecisionGoal -> 4]


The computations are slow because of BesselJ[0,x_]. If we use approximation we get very quickly very high values as a result together with non-convergence messages.

Clear[f]
f[u_?NumericQ, v_?NumericQ, w_?NumericQ] :=
Block[{BesselJ0},
BesselJ0[x_] :=
1 - x^2/4 + x^4/64 - x^6/2304 + x^8/147456 - x^10/14745600 + x^12/
2123366400;
NIntegrate[
x*Exp[-2*w*Sqrt[x^2 + y^2]]*BesselJ0[u*x]*
Cos[y*v]*(Cos[Pi*y/0.30040290])^2, {x, 0, Infinity}, {y, 0,
0.30040290/2}, Method -> "GlobalAdaptive", MaxRecursion -> 4,
PrecisionGoal -> 2]
];

NIntegrate[(u/(u^2 + v^2))*((f[u, v, 4])^2), {u, 0, Infinity}, {v, 0,
Infinity}, Method -> "GlobalAdaptive", MaxRecursion -> 5,
PrecisionGoal -> 4]

(* 1.674809585047211*10^4159650 *)

• Since the integral goes to infinity, I am not very sure replacing the Bessel function with a series would capture the oscillatory behavior for large arguments. It might be profitable to split the integral and use that polynomial expression you have in the small x region, and the asymptotic expression for large x. Mar 18, 2018 at 14:59
• @J.M. Thank you for your comment! I am not saying the approximation is valid, I just wanted to demonstrate where the slowness comes from. Mar 18, 2018 at 15:00
– rafa
Mar 18, 2018 at 15:14
• @rafa You can experiment with different option values. My answer shows definitions that do not have the NIntegrate::inumr: message in your question. Mar 18, 2018 at 15:20

No need to define a intermediate function by using NIntegrate, as it assumes that u and v are known numbers. Further, you have a syntax error when calculating NIntegrate: your should write f[u,v,4] NOT f[u_,v_,4].

I am also assuming that your code is correct although the integral limits are not equal for y:

NIntegrate[(u/(u^2 + v^2))*((x*Exp[-2*4*Sqrt[x^2 + y^2]]*
BesselJ[0, u*x]*Cos[y*v]*(Cos[Pi*y/0.30040290])^2)^2), {u, 0,
Infinity}, {v, 0, Infinity}, {x, 0, Infinity}, {y, 0, 0.30040290/2},
Method -> {"AdaptiveMonteCarlo", "BisectionDithering" -> 1/10}]

(* 0.000433598 *)


I suggest to use AdaptiveMonteCarlo as you integrand is quite peaked, together with BisectionDithering option to sample it properly. You can play with other integration strategies and/or change the value for this option to see whether your result is reliable.

Oh, I am using macOS 10.13.3 and MMA 11.3 on an old Intel i3.

## EDIT

I missed previously the fact that the integration f[u,v,w] is squared. When correcting it, the integration is too slow as noted in other answers.

So I propose another workaround to be checked. Being defined f, we can sample it in u and v :

pts = ParallelTable[{{u, v}, (f[u, v, 4])^2}, {u, 0, 20}, {v, 0, 20}];


The restricted range has to do with the low values of the integral. The, we can interpolate the points, and, finally, integrate:

fun = Interpolation[Flatten[pts, 1]];
NIntegrate[(u/(u^2 + v^2))*fun[u, v], {u, 0, 20}, {v, 0, 20},
Method -> {"AdaptiveMonteCarlo", "BisectionDithering" -> 1/10}]

(* 0.000039325 *)


After almost 12 hours of computation:

NIntegrate[(u/(u^2 + v^2))*((f[u, v, 4.])^2), {u, 0, Infinity}, {v, 0,
Infinity}, Method -> {"AdaptiveMonteCarlo", "BisectionDithering" -> 1/10}]

(* 0.0000401483 *)


which shows that my previous approach is not so bad and quite faster.

• Thak you for your reply. I think it is necessary to define this intermediate function with Nintegrate, because is not the same: the square of the integral, that the integral of the square....isnt it?
– rafa
Mar 18, 2018 at 13:59
• + for noticing the syntax thing, but he's right you cant pull the integral out of the square. Mar 18, 2018 at 14:03
• Oh, I misunderstood that! Mar 18, 2018 at 14:22
• @rafa take a look at my last edit Mar 19, 2018 at 7:49