0
$\begingroup$

I am currently working with the following integrals

\begin{equation} \int_{0}^{\infty} dk\thinspace \frac{k^{3}e^{-2kd}}{\omega^{2}+k^{4}} = \frac{1}{\omega^{2}d^{4}}\int_{0}^{\infty}d\epsilon\thinspace\frac{\epsilon^{3}e^{-2\epsilon}}{1+\frac{\epsilon^{4}}{\omega^{2}d^{4}}} \end{equation} which can be verified by the use of $\epsilon = k d$.

The problem is: depending on which side of the above equation I use, I obtain different results for its value as a function of $d$.

Here is how I define the integrals

SEND[\[Omega]_, d_] := NIntegrate[k*((k^2*Exp[-2*k*d])/((k^2)^2 + \[Omega]^2)), {k,0,Infinity}, AccuracyGoal -> 20, PrecisionGoal -> 20, MaxRecursion -> 20,   WorkingPrecision -> 50]

SEdimN[\[Omega]_, d_] := (1/(\[Omega]^2*d^4))*NIntegrate[(\[Epsilon]^3*Exp[-2*\[Epsilon]])/(1 + (1/(\[Omega]^2*d^4))*\[Epsilon]^4), {\[Epsilon], 0, Infinity},AccuracyGoal -> 20, PrecisionGoal -> 20,MaxRecursion -> 20, WorkingPrecision -> 50]

and here is how I plot them

LogLogPlot[{SEND[1, z], SEdimN[1, z]}, {z, 10^-10, 10^5}, LabelStyle -> Black, PlotLegends -> "Expressions", PlotRange -> {10^-6, 10^3}, PlotPoints -> 1000,WorkingPrecision -> 50]

The result I obtain is something like

enter image description here

from where we can clearly see that for most of the range both integrals yield the same value. So my question is: Whey they are yielding different results? How do I know which one is the correct one!?

Improve this question
$\endgroup$
2
  • 2
    $\begingroup$ This integral can be calculated symbolically Integrate[(k^3 Exp[-2 d k])/(\[Omega]^2 + k^4) , {k, 0, \[Infinity]}, Assumptions -> d > 0 && \[Omega] > 0]. The shape of the blue curve suggests that the numerical integral is calculated with insufficient accuracy. For related issues see e.g. Numerical underflow for a scaled error function. $\endgroup$
    – Artes
    Sep 28, 2021 at 2:41
  • 1
    $\begingroup$ The jump around $10^{-4}$ strongly suggests a numerical issue. $\endgroup$
    – Michael E2
    Sep 28, 2021 at 3:14

2 Answers 2

Reset to default
4
$\begingroup$

Better precision control is required.

Clear["Global`*"]

SEND[ωv_?NumericQ, dv_?NumericQ] :=
 Module[
  {ω = SetPrecision[ωv, 60], d = SetPrecision[dv, 60]},
  NIntegrate[
   k*((k^2*Exp[-2*k*d])/((k^2)^2 + ω^2)),
   {k, 0, Infinity},
   MaxRecursion -> 100,
   WorkingPrecision -> 50]]

SEdimN[ωv_?NumericQ, dv_?NumericQ] :=
 Module[
  {ω = SetPrecision[ωv, 60], d = SetPrecision[dv, 60]},
  (1/(ω^2*d^4))*
   NIntegrate[
    (ϵ^3*
       Exp[-2*ϵ])/(1 + (1/(ω^2*d^4))*ϵ^4),
    {ϵ, 0, Infinity},
    MaxRecursion -> 100,
    WorkingPrecision -> 50]]

{SEND[1, z], SEdimN[1, z]} /. z -> 10^-7

(* {14.847733027640972818058721332928476020439557891260, 
    14.847733027640972818058721332928476020439557891260} *)

LogLogPlot[
 {SEND[1, z], SEdimN[1, z]},
 {z, 10^-10, 100},
 LabelStyle -> Black,
 PlotStyle -> {Automatic, Dashed},
 PlotLegends -> Placed["Expressions", {0.35, 0.4}],
 PlotRange -> {10^-6, 10^2},
 PlotPoints -> 75,
 MaxRecursion -> 5,
 WorkingPrecision -> 50]

enter image description here

Improve this answer
$\endgroup$
3
$\begingroup$

Remove the AccuracyGoal -> 20 from SEdimN[] or from both (from both shown):

enter image description here

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.