Can this type of function be plotted faster using NDSolve instead of NIntegrate?

In this post How to speed up the plot of NIntegrate? it was suggested that we can speed up the plotting of NIntegrate by using NDSolve instead. Can a similar way be done to plot this type of function

f[x_] := NIntegrate[Exp[-x y z - y^(1/2) - z^(1/2)], {y, 0, ∞}, {z, 0, ∞}];

Plot[f[x], {x, 0, 10}] // Timing


How can I apply NDSolve to this function?

How can I apply NDSolve to this function?

Do you want to apply NDSolve or do you want to plot the integral faster?

The plot can be produced ≈ 20 faster using appropriate option settings. Here is an example:

Clear[f2];
f2[x_?NumericQ] :=
NIntegrate[Exp[-x y z - y^(1/2) - z^(1/2)], {y, 0, Infinity}, {z, 0, Infinity}, PrecisionGoal -> 2, AccuracyGoal -> 3, Method -> {Automatic, "SymbolicProcessing" -> 0}];

Plot[f2[x], {x, 0, 10}, PlotRange -> All]; // AbsoluteTiming
(* 0.694531 *)


• To answer your question I just want to plot the integral faster because some functions I have are more complicated than this and they take much longer. Commented Sep 27, 2022 at 23:27
• Ok. Do the options above help? For complicated integrals using less precision and accuracy might make NIntegrate "skip" certain singularities. If the integrals have more or less the same shape then "en bloc" numerical integration can be used; see mathematica.stackexchange.com/a/126041 . Commented Sep 28, 2022 at 1:21

Note: This is an extended comment.

This integral can be expressed in terms of the MeijerG special function:

16/π*MeijerG[{{1},{}},{{1,1,3/2,3/2},{}},1/(16*x)]


Plotting using this expression takes about 6 seconds, compared to 12 seconds using OPs code. To obtain this expression, I first made a change of variables, replacing y by u/z where u is a new variable, and then evaluated the resulting integral

Integrate[Exp[-x*u-(u/z)^(1/2)-z^(1/2)]/z,{u,0,∞},{z,0,∞},
Assumptions->{x>0}]


which takes about 10 seconds and gives the MeijerG result.

Having a symbolic result has advantages besides plotting, for example if OP wanted to know about $$x \to \infty$$ asymptotics, they could use

Series[16/π*MeijerG[{{1},{}},{{1,1,3/2,3/2},{}},1/(16*x)],{x,∞,1}]
(* (-EulerGamma+Log[16]+Log[x]+2 PolyGamma[0,1/2])/x+O[1/x]^(3/2) *)