2
$\begingroup$

I'm trying to evaluate an integral numerically, but I get the following warning message

NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small. >>

and get no results. I tried to adjust the PrecisionGoal and MaxRecursion parameters for several values, but it still doesn't work. How can I fix this?

The integration mathematically is

\begin{equation} \varepsilon_s=\frac{1}{2}-\frac{A^2\alpha_Q}{2(K-1)!\alpha_p}\int_0^{\infty}\frac{1}{x}e^{-x}I(x)\,dx, \end{equation} where \begin{equation} I(x)=(-1)^{n-1}\left[\frac{A}{\alpha_p}+\frac{A^2\alpha_Q}{x\alpha_p}\right]^ne^{\frac{A}{\alpha_p}+\frac{A^2\alpha_Q}{x\alpha_p}}\text{Ei}\left(-\left[\frac{A}{\alpha_p}+\frac{A^2\alpha_Q}{x\alpha_p}\right]\right)+\sum_{k=1}^n(k-1)!\left(-\left[\frac{A}{\alpha_p}+\frac{A^2\alpha_Q}{x\gamma_p}\right]\right)^{n-k}, \end{equation} where $n=K-1$. The mathematical code is

alphap = 10^(0/10); 
ACont = 1; 
KPU = 2; 
NPU = KPU - 1; 
For[alphaQdB = -10, alphaQdB <= 5, alphaQdB++; alphaQ = 10^(alphaQdB/10); 
   Result = 0.5 - 0.5*((ACont^2*alphaQ)/((KPU - 1)!*alphap))*(-1)^(NPU - 1)*
      NIntegrate[(1/x)*(ACont/alphap + (ACont^2*alphaQ)/(x*alphap))^NPU*
         Exp[ACont/alphap + (ACont^2*alphaQ)/(x*alphap) - x]*
         ExpIntegralE[1, ACont/alphap + (ACont^2*alphaQ)/(x*alphap)] + 
        Sum[(kk - 1)!*(-(ACont/alphap + (ACont^2*alphaQ)/(x*alphap)))^
           (NPU - kk), {kk, 1, NPU}], {x, 0, Infinity}, 
       PrecisionGoal -> 20, MaxRecursion -> 50]; Print[Result]; ]

I hope the code format is OK. I copies and pasted it directly from Mathematica editor. Also, sorry, the integration expression isn't that pleasant, and that's why I hesitated to put it.

PS: $$-Ei(-x)=E_1(x)$$

Thanks

$\endgroup$
1
$\begingroup$

too long for a comment.. lets look at your integrand:

alphap = 10^(0/10);
ACont = 1;
KPU = 2;
NPU = KPU - 1;
alphaQdB = -10;
alphaQ = 10^(alphaQdB/10);
(1/x)*(ACont/alphap + (ACont^2*alphaQ)/(x*alphap))^NPU*
  Exp[ACont/alphap + (ACont^2*alphaQ)/(x*alphap) - x]*
  ExpIntegralE[1, ACont/alphap + (ACont^2*alphaQ)/(x*alphap)] + 
 Sum[(kk - 1)!*(-(ACont/alphap + (ACont^2*alphaQ)/(x*alphap)))^(NPU - 
      kk), {kk, 1, NPU}]

enter image description here

That constant 1 causes the infinite integral to diverge and as noted in comments arises because have failed to multiply the the Exp[-x]/x through the entire I expression.

fixing that,

enter image description here

you get a convergent integral on the infinite domain, so long as you avoid zero.

Without knowing the application its not clear what the point would be of setting the lower bound to some arbitrary finite value since the result will depend strongly on the selected value.

$\endgroup$
  • $\begingroup$ OMG, you are right, I just noticed that. The second part that includes the summation in I isn't multiplied by Exp[-x]/x. I don't have access to Mathematica right now, but you're saying that the integral still doesn't give results in the range between 0 and infinity because of the zero lower bound, is that correct? The idea of arbitrary small value for the lower bound is to evaluate the integral as close as possible to the zero lower bound such that it's integrable to get results. $\endgroup$ – EngDavid Aug 5 '16 at 19:04
  • $\begingroup$ right you can do NIntegrate[ .. , {x,10^-6,Infinity}] and it converges (~25 for this example) $\endgroup$ – george2079 Aug 5 '16 at 19:32
  • $\begingroup$ Thanks. I will check this when I have access to my computer. $\endgroup$ – EngDavid Aug 5 '16 at 19:51
  • $\begingroup$ Thank you. I just checked it, and the integral is evaluated to a finite value now, even with the zero lower bound. Thanks again for notifying me on the the missing term. $\endgroup$ – EngDavid Aug 8 '16 at 15:06
  • $\begingroup$ OK, now when I change KPU to 3 I get warning message: overflow, indeterminate, or infinity!! Why with KPU=2 it works while with KPU=3 it doesn't?! $\endgroup$ – EngDavid Aug 8 '16 at 15:33

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.