# Integration analitycally and numerically

I've been trying to solve the next integral :

$$\int_{0}^{\infty}{\frac{652rSin[kr]}{72248.9k+ke^{1.85185r}}dr}$$

But if I do it analytically , I obtain a function which I can not work with:

ConditionalExpression[ 1/k (((0. + 330.575 I) HypergeometricPFQ[{1., 1. - (0. + 0.54 I) k,
1. - (0. + 0.54 I) k}, {2. - (0. + 0.54 I) k,
2. - (0. + 0.54 I) k}, -72248.9])/((0. + 1.85185 I) +
1. k)^2 - ((0. + 330.575 I) HypergeometricPFQ[{1.,
1. + (0. + 0.54 I) k,
1. + (0. + 0.54 I) k}, {2. + (0. + 0.54 I) k,
2. + (0. + 0.54 I) k}, -72248.9])/((0. + 1.85185 I) -
1. k)^2), Abs[Im[k]] < 1.85185]


And if I try to solve it numerically I got this :

"The relative error 7.769066430156647 is larger than expected for \the integrand (6611.51\(0. +r)\Sin[0.1\r])/(72248.9 +E^(1.85185 (0. \+r))) over {0,\[Infinity]} with DoubleExponentialOscillatory method \and automatic tuning parameters, TuningParameters -> {10,5}. The \integration will proceed with TuningParameters -> {1,5}. "


Several times!!! and I cannot fit a good function. The ideal thing here ( I guess is the goal) is to find a function because I have to put that function inside another integral but I don't know how to accomplish that. I really need help on this :D. Thank you in advance.

• The ConditionalExpression depends on the size of the imaginary part of k. Numerical integration also needs the value of k. Is k real? Do you have a value for k? What form do you absolutely have to have the result in? How precisely do you need to know the answer?
– Bill
Commented Aug 22, 2017 at 5:59

## 2 Answers

For this integral, you can use a trick where you replace a parameter with an exact numeric quantity, and then replace the exact numeric quantity with a variable. For example, lets replace your 72247.9 with Pi and 1.85185 with EulerGamma:

sint = Integrate[
(r Sin[k r])/(Pi k + k Exp[EulerGamma r]),
{r,0,Infinity},
Assumptions -> k>0
]


(2 EulerGamma)/(EulerGamma^2 + k^2)^2 + ( I π (LerchPhi[-π, 2, 2 - (I k)/EulerGamma] - LerchPhi[-π, 2, 2 + (I k)/EulerGamma]))/(2 EulerGamma^2 k)

(it is possible to evaluate the integral using s and t directly, but it takes a much longer time). Now, we need to replace Pi with s and EulerGamma with t:

sint /. {EulerGamma -> t, Pi -> s}


(2 t)/(k^2 + t^2)^2 + ( I s (LerchPhi[-s, 2, 2 - (I k)/t] - LerchPhi[-s, 2, 2 + (I k)/t]))/(2 k t^2)

Based on the above, it's reasonable to assume that the integral is:

int[s_, t_, k_] = sint /. {EulerGamma -> t, Pi -> s};


Or, in TeX:

$$\int _0^{\infty }\frac{r \sin (k r)}{k e^{r t}+k s}dr=\frac{2 t}{\left(k^2+t^2\right)^2}+\frac{i s \left(\Phi \left(-s,2,2-\frac{i k}{t}\right)-\Phi \left(-s,2,\frac{i k}{t}+2\right)\right)}{2 k t^2}$$

The NIntegrate version is:

Options[nint] = {WorkingPrecision->50};
nint[s_, t_, k_, OptionsPattern[]] := NIntegrate[
(r Sin[k r])/(s k + k Exp[t r]),
{r,0,Infinity},
WorkingPrecision->OptionValue[WorkingPrecision]
]


Here are some plots comparing int with nint:

Plot[int[1000, t, 2], {t, 0, 1}, WorkingPrecision->50,PlotRange->All]


Plot[nint[1000, t, 2], {t, 0, 1}, PlotRange->All, WorkingPrecision->50]


NIntegrate::precw: The precision of the argument function ((r Sin[2 r])/(2000+2 E^(0.0000204286 r))) is less than WorkingPrecision (50.).

NIntegrate::deodiv: DoubleExponentialOscillatory returns a finite integral estimate, but the integral might be divergent.

These plots take a long time to evaluate, so I will only show one more. Here is a plot of the function at the point of interest for the OP:

Plot[int[72248.950, 1.8518550, k], {k, 0, 2}, WorkingPrecision->40]


Plot[nint[72249.950, 1.8518550, k], {k, 0, 2}, WorkingPrecision->40]


NIntegrate::precw: The precision of the argument function ((r Sin[0.0000408571 r])/(2.95192 +0.0000408571 E^(1.8518500000000000000000000000000000000000000000000 r))) is less than WorkingPrecision (50.).

NIntegrate::deorela: The relative error 1.627388142932627660797142387729047612568727184720250. is larger than expected for the integrand (r Sin[0.000040857142857142848336223645011600069665291812270880 r])/(2.9519244857142847848763267393223941326141357421875+0.000040857142857142848336223645011600069665291812270880 2.7182818284590452353602874713526624977572470937000^(1.8518500000000000000000000000000000000000000000000 r)) over {0,[Infinity]} with DoubleExponentialOscillatory method and automatic tuning parameters, TuningParameters -> {10,5}. The integration will proceed with TuningParameters -> {1,5}.

• Thank you very much Carl! I will try to use this process to evaluate integral , I hope it works in a double integral that I'm pretty sure I have to calculate numercally. Have a nice day :D Commented Aug 23, 2017 at 2:45

This is your expression:

expr1 = (652*r*Sin[k r])/(72248.9*k + k*Exp[1.85*r]);


Let us replace k r->R so that there is Sin[R] instead of Sin[k r]:

expr2 = expr1 /. r -> R/k // Simplify

(*   (652. R Sin[R])/((7.22*10^4 + E^((1.85 R)/k)) k^2)  *)


Now one can integrate over R, not forget to divide once more by k. Try to use the Method->"LevinRule" option:

lst = Table[{k,
1/k^3 NIntegrate[(
652. R Sin[R])/(72248.9 + E^((1.85 R)/k)) , {R,
0, \[Infinity]}, Method -> "LevinRule"]}, {k,
Join[Range[0.01, 0.99, 0.01], Range[1, 5, 0.05]]}];


yielding this:

ListPlot[lst]
`

Have fun!

• Thank you very much for your help , Alexei ! I will try it right now and I'm pretty sure it will be very fun ! :D Commented Aug 23, 2017 at 2:43