EDIT #3 complete revision
This text is intended to exhibit the current status of the problem comprising the discussions.
Main issues are
- transformation of p and q to cartesian coordinates
- convergence of the integral in the origin (obvious after 1.)
- introduction of explicit cut-off for large p and q
- model case to judge the accuracy of NIntegrate
For the ease of reading I start again from the beginning.
Let us write down the twofold two-dimensional integral step by step.
Define length of vector x
abs[x_] := Sqrt[x.x]
Vectors pv, qv and kv (x and y are angles between 0 and 2 pi)
pv = p {Cos[x], Sin[x]};
qv = q {Cos[y], Sin[y]};
kv = {k, 0};
The volume elements are
d^2p = p dp dx
d^2 q = q dq dy
Hence the integrand of F(k) is
int = Simplify[
p q 1/(abs[pv] abs[pv - qv] ) (1/abs[kv - qv] - 1/q), {p > 0, q > 0}]
(* Out[1020]= (q - Sqrt[k^2 + q^2 - 2 k q Cos[y]])/(Sqrt[p^2 + q^2 - 2 p q Cos[x - y]] Sqrt[k^2 + q^2 - 2 k q Cos[y]]) *)
In Latex
$$int = \frac{q-\sqrt{k^2-2 k q \cos (y)+q^2}}{\sqrt{k^2-2 k q \cos (y)+q^2} \sqrt{p^2-2 p q \cos (x-y)+q^2}}$$
Typical behaviour of the integrand
r := RandomReal[];
qr = -Log[r]; pr = -Log[r]; xr = 2 \[Pi] r; yr = 2 \[Pi] r;
Plot[int /. {p -> pr, q -> qr, x -> xr, y -> yr}, {k, 0, 15},
PlotRange -> All,
PlotLabel -> "Monte-Carlo check of integrand\nas a function of k",
AxesLabel -> {"k", "int"}]
Polar coordinates in 2-dimensional p-q-space
It is convenient to consider p and q formally as cartesian coordinates and make a coordinate transformation to polar coordinates (notice the factor c from the volume element dpdq = dc du)
intc = Simplify[
c int /. {p -> c Cos[u], q -> c Sin[u]}, {c > 0, k > 0, 0 < x < 2 \[Pi],
0 < y < 2 \[Pi], 0 < u < \[Pi]/2}]
(* Out[1086]= (c Sin[u] - Sqrt[
k^2 - 2 c k Cos[y] Sin[u] +
c^2 Sin[u]^2])/Sqrt[(k^2 - 2 c k Cos[y] Sin[u] + c^2 Sin[u]^2) (Cos[u]^2 +
Sin[u]^2 - Cos[x - y] Sin[2 u])] *)
in Latex
$$intc = \frac{c \sin (u)-\sqrt{c^2 \sin ^2(u)-2 c k \sin (u) \cos (y)+k^2}}{\sqrt{(1-\sin (2 u) \cos (x-y)) \left(c^2 \sin ^2(u)-2 c k \sin (u) \cos (y)+k^2\right)}}$$
Monte-Carlo check of integrand
cr = -Log[r]; xr = 2 \[Pi] r; yr = 2 \[Pi] r; ur = r \[Pi]/2;
Plot[intc /. {c -> cr, x -> xr, y -> yr, u -> ur}, {k, 0, 30},
PlotRange -> All,
PlotLabel ->
"Monte-Carlo check of integrand\ntransformed to Cartesian coordinates\nas a \
function of k", AxesLabel -> {"k", "intc"}]
Sumarizing the (numerical) integral of the OP is given by the function
f[kk_, R_] :=
Integrate[intc /. k -> kk, {c, 0, R}, {u, 0, \[Pi]/2}, {x, 0, 2 \[Pi]}, {y,
0, 2 \[Pi]}]
fN[kk_, R_] :=
NIntegrate[
intc /. k -> kk, {c, 0, R}, {u, 0, \[Pi]/2}, {x, 0, 2 \[Pi]}, {y, 0,
2 \[Pi]}]
Here we have introduced a cutoff radis R in order to avoid the divergence at c=[Infinity].
A typical numerical result is
fN[30, 10]
During evaluation of In[1158]:= 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. >>
During evaluation of In[1158]:= NIntegrate::eincr: The global error of
the strategy GlobalAdaptive has increased more than 2000 times. The
global error is expected to decrease monotonically after a number of
integrand evaluations. Suspect one of the following: the working
precision is insufficient for the specified precision goal; the
integrand is highly oscillatory or it is not a (piecewise) smooth
function; or the true value of the integral is 0. Increasing the value
of the GlobalAdaptive option MaxErrorIncreases might lead to a
convergent numerical integration. NIntegrate obtained -745.356 and
27.317349155628584` for the integral and error estimates. >>
(* Out[1158]= -745.356 *)
Reliability
In order to get a feeling of the reliability of the numeral results we study the integral dxdydu over the integrand at c = 0
intc0 = -(1/Sqrt[1 - Cos[x - y] Sin[2 u]]);
The exact value of the integral over intc0 is (see proof below)
ic0 = -((2 \[Pi]^4)/(Gamma[5/8]^2 Gamma[7/8]^2));
% // N
(* Out[1139]= -79.7335 *)
This value can be used to ckeck the prcision of various numerical approaches
ic0N1 = NIntegrate[intc0, {u, 0, \[Pi]/2}, {x, 0, 2 \[Pi]}, {y, 0, 2 \[Pi]}]
During evaluation of In[1125]:= 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. >>
During evaluation of In[1125]:= NIntegrate::eincr: The global error of
the strategy GlobalAdaptive has increased more than 2000 times. The
global error is expected to decrease monotonically after a number of
integrand evaluations. Suspect one of the following: the working
precision is insufficient for the specified precision goal; the
integrand is highly oscillatory or it is not a (piecewise) smooth
function; or the true value of the integral is 0. Increasing the value
of the GlobalAdaptive option MaxErrorIncreases might lead to a
convergent numerical integration. NIntegrate obtained -79.4062 and
0.14460520038353433` for the integral and error estimates. >>
(* Out[1125]= -79.4062 *)
ic0N2 = NIntegrate[intc0, {u, 0, \[Pi]/2}, {x, 0, 2 \[Pi]}, {y, 0, 2 \[Pi]},
Method -> "MonteCarlo", PrecisionGoal -> 2, WorkingPrecision -> 10]
(* Out[1126]= -78.43086734 *)
We observe very good agreement using the first method and confirm the error estimate of Mathematica.
Calculation of the exact value of the integral over intc0
Brute force Mathematica gives
Integrate[intc0, {u, 0, \[Pi]/2}, {x, 0, 2 \[Pi]}, {y, 0,
2 \[Pi]}] (* Mathematica error *)
(* Out[1120]= 0 *)
Which is oviously wrong!
We can find the correct result expanding the square root into a series and integrating term by term
Sum[Binomial[-(1/2), n] (-1)^n Cos[x - y]^n Sin[2 u]^n, {n,
0, \[Infinity]}] // FullSimplify
(* Out[1133]= 1/Sqrt[1 - Cos[x - y] Sin[2 u]] *)
Integrate[Binomial[-(1/2), n] (-1)^n Cos[x - y]^n Sin[2 u]^n, {u, 0, \[Pi]/2},
Assumptions -> n > -1]
(* Out[1134]= ((-1)^n Sqrt[\[Pi]]
Binomial[-(1/2), n] Cos[x - y]^n Gamma[(1 + n)/2])/(2 Gamma[1 + n/2])
The x and y integral is
Simplify[2 \[Pi] Integrate[ Cos[x]^n , {x, 0, 2 \[Pi]}],
n \[Element] Integers]
(* Out[1113]= -(((-2)^(1 + n) (1 + (-1)^n) \[Pi]^3)/(Gamma[1/2 - n/2]^2 Gamma[1 + n])) *)
And the sum becomes
Sum[-(((-2)^(1 + n) (1 + (-1)^n) \[Pi]^3)/(
Gamma[1/2 - n/2]^2 Gamma[1 + n])) ((-1)^n Sqrt[\[Pi]]
Binomial[-(1/2), n] Gamma[(1 + n)/2])/(2 Gamma[1 + n/2]), {n,
0, \[Infinity]}]
(* Out[1136]= (2 \[Pi]^4)/(Gamma[5/8]^2 Gamma[7/8]^2) *)
% // N
(* Out[1137]= 79.7335 *)
Q.E.D.
My original post
This is not a solution but a clarification of my statement of the comment. The next step to be done is to show that the integral is divergent so that a numerical treatment is ruled out.
There is an inconsistency in the OP: the integrand in the definition of F(k) differs from the integrand used in NIntegrate. Here we calculate the correct Expression of the Integrand.
Let us write down the twofold two-dimensional integral step by step.
Define the length of vector x
abs[x_] := Sqrt[x.x]
Vectors pv, qv and kv (x and y are angles between 0 and 2 pi)
pv = p {Cos[x], Sin[x]};
qv = q {Cos[y], Sin[y]};
kv = {k, 0};
The volume elements of the integrals are
d^2 p = p dp dx
d^2 q = q dq dy
Hence the integrand of F(k) is
int = Simplify[
p q 1/(abs[pv] abs[pv - qv] ) (1/abs[kv - qv] - 1/q), {p > 0, q > 0}]
(* Out[133]= (q - Sqrt[k^2 + q^2 - 2 k q Cos[y]])/(Sqrt[p^2 + q^2 - 2 p q Cos[x - y]] Sqrt[ k^2 + q^2 - 2 k q Cos[y]]) *)
In Latex
$$int = \frac{q-\sqrt{k^2-2 k q \cos (y)+q^2}}{\sqrt{k^2-2 k q \cos (y)+q^2} \sqrt{p^2-2 p q \cos (x-y)+q^2}}$$
p q/p
reduces toq
-- is that a typo? $\endgroup$