an integral containing BesselJ
dimm...@yahoo.com
In[1]:=
\!\(Integrate[\(Log[x]\/\@\(1 + x\^2\)\) BesselJ[0, x], {x, 0, 8}] //
Timing\)
Out[1]=
{46.171 Second,0}
which is a special case of:
In[2]:=
\!\(Integrate[\(Log[a\ x]\/\@\(1 + x\^2\)\) BesselJ[0,
x], {x, 0, 8}] // Timing\)
Out[2]=
\!\({32.18800000000001`\ Second, BesselI[0, 1\/2]\ BesselK[0, 1\/2]\
Log[a]}\)
However performing the numerical integration gave the result:
In[4]:=
\!\(NIntegrate[\(Log[x]\/\@\(1 +
x\^2\)\) BesselJ[0, x], {x, 0, 8}, Method \[Rule]
Oscillatory] // Timing\)
Out[4]=
{0.75 Second,-0.997939}
I also perfrorm the symbolic integration in Mathematica 4.0 and took a
lengthy result, which agrees numerically with the obtained result from
the numerical integration:
In[1]:=
\!\(Integrate[\(Log[a\ x]\/\@\(1 + x\^2\)\) BesselJ[0, x], {x,
0, \[Infinity]}] // Timing\)
Out[1]=
\!\(\*
RowBox[{"{",
RowBox[{\(18.594`\ Second\), ",",
RowBox[{\(1\/12\ \((\(-1\) +
HypergeometricPFQ[{1\/2}, {1, 1},
1\/4])\)\ \((6\ EulerGamma\^2 - \[Pi]\^2 -
6\ EulerGamma\ Log[4])\)\),
"+", \(1\/12\ \((6\ EulerGamma\^2 - \[Pi]\^2 -
6\ EulerGamma\ Log[4] +
12\ BesselI[0, 1\/2]\ BesselK[0, 1\/2]\ Log[a])\)\), "+",
RowBox[{\(1\/\(2\ \@\[Pi]\)\),
RowBox[{"(",
RowBox[{\(Log[4]\), " ",
RowBox[{"(",
RowBox[{\(1\/2\ \@\[Pi]\ \((\(-1\) +
HypergeometricPFQ[{1\/2}, {1, 1}, 1\/4])\)\
Log[4]\),
"+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14174 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(\(2\^\(\(-1\) - 2\ K$14174\)\ Gamma[
1\/2 + K$14174]\ PolyGamma[0,
1\/2 + K$14174]\)\/Gamma[1 +
K$14174]\^3\)}], "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14174 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(-\(\(2\^\(\(-1\) - 2\ K$14174\)\ Gamma[
1\/2 + K$14174]\ PolyGamma[0,
1 + K$14174]\)\/Gamma[1 +
K$14174]\^3\)\)}]}],
")"}]}], ")"}]}], "+",
FractionBox[
RowBox[{\((\(-2\)\ EulerGamma + Log[4])\), " ",
RowBox[{"(",
RowBox[{\(EulerGamma\ \@\[Pi]\ \((\(-1\) +
HypergeometricPFQ[{1\/2}, {1, 1}, 1\/4])\)\),
"+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14447 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(\(4\^\(-K$14447\)\ Gamma[
1\/2 + K$14447]\ PolyGamma[0,
1 + K$14447]\)\/Gamma[1 + K$14447]\^3\)}]}],
")"}]}], \(2\ \@\[Pi]\)], "+",
RowBox[{\(1\/\(3\ \@\[Pi]\)\),
RowBox[{"(",
RowBox[{\(1\/8\ \[Pi]\^\(5/2\)\ \((\(-1\) +
HypergeometricPFQ[{1\/2}, {1, 1}, 1\/4])\)\),
"-", \(3\/8\ \@\[Pi]\ \((\(-1\) +
HypergeometricPFQ[{1\/2}, {1, 1}, 1\/4])\)\
Log[4]\^2\),
"+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14601 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(-\(\(3\ 2\^\(\(-3\) - 2\ K$14601\)\ Gamma[
1\/2 + K$14601]\ Log[16]\ PolyGamma[0,
1\/2 + K$14601]\)\/Gamma[1 +
K$14601]\^3\)\)}], "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14601 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(-\(\(3\ 2\^\(\(-3\) - 2\ K$14601\)\ Gamma[
1\/2 +
K$14601]\ PolyGamma[0, 1\/2 + \
K$14601]\^2\)\/Gamma[1 + K$14601]\^3\)\)}], "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14601 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(\(3\ 2\^\(\(-3\) - 2\ K$14601\)\ Gamma[
1\/2 + K$14601]\ Log[16]\ PolyGamma[0,
1 + K$14601]\)\/Gamma[1 + K$14601]\^3\)}], "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14601 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(\(3\ 2\^\(\(-2\) - 2\ K$14601\)\ Gamma[
1\/2 + K$14601]\ PolyGamma[0,
1\/2 + K$14601]\ PolyGamma[0,
1 + K$14601]\)\/Gamma[1 + K$14601]\^3\)}], "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14601 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(-\(\(3\ 2\^\(\(-3\) - 2\ K$14601\)\ Gamma[
1\/2 +
K$14601]\ PolyGamma[0, 1 +
K$14601]\^2\)\/Gamma[1 \
+ K$14601]\^3\)\)}], "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14601 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(-\(\(3\ 2\^\(\(-3\) - 2\ K$14601\)\ Gamma[
1\/2 + K$14601]\ PolyGamma[1,
1\/2 + K$14601]\)\/Gamma[1 +
K$14601]\^3\)\)}], "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$14601 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(\(3\ 2\^\(\(-3\) - 2\ K$14601\)\ Gamma[
1\/2 + K$14601]\ PolyGamma[1,
1 + K$14601]\)\/Gamma[1 + K$14601]\^3\)}]}],
")"}]}],
"+",
RowBox[{\(1\/\(3\ \@\[Pi]\)\),
RowBox[{"(",
RowBox[{\(3\/2\ EulerGamma\^2\ \@\[Pi]\ \((\(-1\) +
HypergeometricPFQ[{1\/2}, {1, 1}, 1\/4])\)\),
"+", \(1\/8\ \[Pi]\^\(5/2\)\ \((\(-1\) +
HypergeometricPFQ[{1\/2}, {1, 1}, 1\/4])\)\), "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$15304 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(\(3\ 4\^\(-K$15304\)\ EulerGamma\ Gamma[
1\/2 + K$15304]\ PolyGamma[0,
1 + K$15304]\)\/Gamma[1 + K$15304]\^3\)}], "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$15304 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(\(3\ 2\^\(\(-1\) - 2\ K$15304\)\ Gamma[
1\/2 +
K$15304]\ PolyGamma[0, 1 +
K$15304]\^2\)\/Gamma[1 + \
K$15304]\^3\)}], "+",
RowBox[{
UnderoverscriptBox["\[Sum]", \(K$15304 = 1\),
InterpretationBox["\[Infinity]",
DirectedInfinity[
1]]], \(-\(\(3\ 2\^\(\(-2\) - 2\ K$15304\)\ Gamma[
1\/2 + K$15304]\ PolyGamma[1,
1 + K$15304]\)\/Gamma[1 +
K$15304]\^3\)\)}]}],
")"}]}]}]}], "}"}]\)
In[2]:=
N[%[[2]]]
Out[2]=
0.0580341\[InvisibleSpace] + 0.0833333 (-12.6717 + 11.7973 Log[a])
In[3]:=
% /. a -> 1
Out[3]=
-0.997939
Can anyone help me?
ab_...@prontomail.com
Here's a link to a paper that explains how such integrals can be
evaluated: http://www.cs.cmu.edu/~adamchik/articles/integr/mier.nb
(also have a look at issac90.pdf). First we find the transforms of the
factors:
In[1]:= Integrate[Log[1/x]/Sqrt[1 + 1/x^2]*x^(s - 1), {x, 0, Infinity}]
In[2]:= Integrate[BesselJ[0, x]*x*x^(s - 1), {x, 0, Infinity}]
In[3]:= fs = Simplify[%%*%, -1 < s < 0]
Out[3]= (2^(-2 + s)*Gamma[-(s/2)]*Gamma[(1 + s)/2]^2*(PolyGamma[0,
-(s/2)] - PolyGamma[0, (1 + s)/2]))/(Sqrt[Pi]*Gamma[1/2 - s/2])
The original integral is equal to 1/(2*Pi*I) times the integral of fs
along the line Re[s] == -1/2:
In[4]:= 1/(2*Pi*I)*NIntegrate[fs, {s, -I*Infinity, -1/2, I*Infinity}]
Out[4]= -0.9979393746397847 - 2.6504622345529306*^-17*I
We can represent fs as a sum of the derivatives of gamma function
ratios:
In[5]:= fs == D[2^(-2 + s)*Gamma[-s/2 + a]*Gamma[(1 + s)/2]^2/
(Sqrt[Pi]*Gamma[(1 - s)/2]) -
1/2*2^(-2 + s)*Gamma[-s/2]*Gamma[(1 + s)/2 + a]^2/
(Sqrt[Pi]*Gamma[(1 - s)/2]),
a] /. a -> 0 // Simplify
Out[5]= True
Since the integral of the ratio of gamma functions is a Meijer
function, the answer is expressible in terms of the derivatives of
MeijerG:
In[6]:= D[1/(4*Sqrt[Pi])*MeijerG[{{1/2 - a}, {}}, {{0, 0}, {0}}, 1/4] -
1/(8*Sqrt[Pi])*MeijerG[{{1/2}, {}}, {{a, a}, {0}}, 1/4],
a] /. a -> 0 // Simplify
Out[6]= (1/(8*Sqrt[Pi]))*(-Derivative[{{0}, {}}, {{0, 1}, {0}},
0][MeijerG][{{1/2}, {}}, {{0, 0}, {0}}, 1/4] - Derivative[{{0}, {}},
{{1, 0}, {0}}, 0][MeijerG][{{1/2}, {}}, {{0, 0}, {0}}, 1/4] -
2*Derivative[{{1}, {}}, {{0, 0}, {0}}, 0][MeijerG][{{1/2}, {}}, {{0,
0}, {0}}, 1/4])
This is probably not very useful, because Mathematica cannot evaluate
N[%6], and FunctionExpand[%6] gives Indeterminate.
Using the same technique of differentiation with respect to a
parameter, we can express the integral in terms of the derivatives of
HypergeometricPFQ:
In[7]:= int = (1/24)*(-2*BesselI[0, 1/2]^2*(Pi^2 +
6*EulerGamma*(-EulerGamma + Log[4])) + (9*EulerGamma -
3*Log[16])*Derivative[{0}, {0, 1}, 0][HypergeometricPFQ][{1/2}, {1, 1},
1/4] + 3*Log[16]*Derivative[{1}, {0, 0}, 0][HypergeometricPFQ][{1/2},
{1, 1}, 1/4] + 3*Derivative[{1}, {0, 1}, 0][HypergeometricPFQ][{1/2},
{1, 1}, 1/4] - 3*Derivative[{2}, {0, 0}, 0][HypergeometricPFQ][{1/2},
{1, 1}, 1/4] - 15*EulerGamma*Derivative[{0, 1}, {0, 0, 0},
0][HypergeometricPFQ][{1/2, 1}, {1, 1, 1}, 1/4] - 12*Derivative[{0, 1},
{0, 0, 1}, 0][HypergeometricPFQ][{1/2, 1}, {1, 1, 1}, 1/4] -
3*Derivative[{0, 2}, {0, 0, 0}, 0][HypergeometricPFQ][{1/2, 1}, {1, 1,
1}, 1/4] + 9*Derivative[{1, 1}, {0, 0, 0}, 0][HypergeometricPFQ][{1/2,
1}, {1, 1, 1}, 1/4]);
Mathematica seems to get the approximations for mixed partial
derivatives wrong:
In[8]:= N[Derivative[{1, 1}, {0, 0, 0}, 0][
HypergeometricPFQ][{1/2, 1}, {1, 1, 1}, 1/4],
20]
Out[8]= 4.523987519373309478408412473406`20.*^40
This looks like a bug, because even the simplest finite difference
formula can give a good approximation:
In[9]:= Block[{$MaxExtraPrecision = 100},
int /.
Derivative[{si: (0)..., 1, si2: (0)..., 1, si3: (0)...},
{(0)..}, 0][HypergeometricPFQ][a_, b_, z_] :>
((f[h, k] - f[h, -k] - f[-h, k] + f[-h, -k])/(4*h*k) /.
f -> Function @@ {HypergeometricPFQ[
a + {si, #, si2, #2, si3}, b, z]}) /.
{h -> 10^-25, k -> 10^-25} //
N[#, 50]&
]
Out[9]=
-0.9979393746360745614717131506926624475792306906555576417821764083894`50.
Generalizing to Log[a*x] is easy:
Integrate[Log[a*x]/Sqrt[1 + x^2]*BesselJ[0, x], {x, 0, Infinity}] ==
Integrate[(Log[x] + Log[a])/Sqrt[1 + x^2]*BesselJ[0, x],
{x, 0, Infinity}] ==
int + BesselI[0, 1/2]*BesselK[0, 1/2]*Log[a]
Maxim Rytin
m...@inbox.ru