# Certain integral over torus

Let $F(t_1,s_1,t_2,s_2)=$ $$\big((2+\cos t_1)\cos s_1 - (2+\cos t_2)\cos s_2\big)^2 + \big((2+\cos t_1)\sin s_1-(2+\cos t_2)\sin s_2\big)^2 + (\sin t_1 -\sin t_2 )^2.$$ I am interested in computing the following integral:

$$\int_0^{2\pi}\int_0^{2\pi}\int_0^{2\pi}\int_0^{2\pi}\frac{(2+\cos t_1)(2+\cos t_2)\,dt_1\,ds_1\,dt_2\,ds_2}{F(t_1,s_1,t_2,s_2)}$$ with Mathematica using the following

NIntegrate
[((2+Cos[t])(2+Cos[u]))/((2+Cos[t])Cos[s] - (2+Cos[u])Cos[v])^2
+ ((2+Cos[t])Sin[s]-(2+Cos[u])Sin[v])^2+(Sin[t]-Sin[u])^2,
{t,0, 2Pi}, {s, 0, 2Pi}, {u, 0, 2Pi}, {v,0, 2Pi}]


, but so far no success.

• To me it gives the result 1422.46, together with the warning that this result may not be reliable. Do you have a way to verify it?
– glS
Commented Nov 23, 2015 at 5:13
• Actually its the other way around.i was going to check with Mathematica to see if it converges or not
– BigM
Commented Nov 23, 2015 at 5:17

I am by no means expert in numerical integration techniques, but using the simplicistic rule of thumb high dimensional integrals -> Monte Carlo techniques work better I tried using a Monte Carlo method for the integration (you can find here an explanation of all numerical techniques employed by Mathematica). Using this and incrementing MaxPoint to around 107 seems to compute without warnings:

f[t1_, s1_, t2_,
s2_] := ((2 + Cos[t1]) Cos[s1] - (2 + Cos[t2] Cos[s2]))^2 +
((2 + Cos[t1]) Sin[s1] - (2 + Cos[t2] Sin[s2]))^2 +
(Sin[t1] - Sin[t2])^2
NIntegrate[
(2 + Cos[t1]) (2 + Cos[t2])/f[t1, s1, t2, s2],
{t1, 0, 2 Pi}, {s1, 0, 2 Pi}, {t2, 0, 2 Pi}, {s2, 0, 2 Pi},
Method -> {"MonteCarlo", "MaxPoints" -> 10^7}
]


Which gives a value $\approx 1400$, though varying a bit in each run.

• Very interesting.this integral came up while I was studying certain potentials over tori. its counterpart diverges over sphere ,so I thought this one also must diverge but it seems i doesn't.
– BigM
Commented Nov 23, 2015 at 5:32
• @BigM please take this result with a grain of salt though. While it seems to converge, it would probably be best to find an analytical confirmation. You could try that on math.SE (but I don't know if there is a simple technique to assess convergence in such a case).
– glS
Commented Nov 23, 2015 at 5:51
• That sounds like a great idea. I'll post it on math.SE.
– BigM
Commented Nov 23, 2015 at 5:54
• I wouldn't call 4 D "high dimensional", but MC is definitely a safe option in any number of dimensions. If the result is stable across runs and for different values of points, it is likely reliable. You should be careful for $t_1 = t_2$ and $s_1=s_2$, because the integrand diverges, so I'd sample that region more thoroughly to be safer. Commented Nov 23, 2015 at 9:22

The integrand has singularities occasionally when the coordinates are a multiple of Pi/4. If we subdivide the domain at multiples of Pi/4, we seem to get divergence.

integrand = ((2 + Cos[t]) (2 + Cos[u]))/(
((2 + Cos[t]) Cos[s] - (2 + Cos[u]) Cos[v])^2 +
((2 + Cos[t]) Sin[s] - (2 + Cos[u]) Sin[v])^2 +
(Sin[t] - Sin[u])^2);

NIntegrate[integrand,
{t, 0, π/4, π/2, (3 π)/4, π, (5 π)/4, (3 π)/2, (7 π)/4, 2 Pi},
{s, 0, π/4, π/2, (3 π)/4, π, (5 π)/4, (3 π)/2, (7 π)/4, 2 Pi},
{u, 0, π/4, π/2, (3 π)/4, π, (5 π)/4, (3 π)/2, (7 π)/4, 2 Pi},
{v, 0, π/4, π/2, (3 π)/4, π, (5 π)/4, (3 π)/2, (7 π)/4, 2 Pi}]


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. >>

NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 36 recursive bisections in t near {t,s,u,v} = {0.392699,3.01917,0.662973,3.01917}. NIntegrate obtained 5.391691817606044*^12 and 5.934913000757257*^13 for the integral and error estimates. >>

(*  5.39169*10^12  *)


Over one of the subregions, it seems to diverge, too.

NIntegrate[integrand,
{t, 0, π/4}, {s, 0, π/4}, {u, 0, π/4}, {v, 0, π/4}]


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. >>

NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 36 recursive bisections in t near {t,s,u,v} = {0.392699,0.662973,0.662973,0.662973}. NIntegrate obtained 5.390318262293057*^12 and 5.933648796990452*^13 for the integral and error estimates. >>

(*  5.39032*10^12  *)


Local adaptive method gives 1423.97 without errors

NIntegrate[(2 + Cos[t1]) (2 + Cos[t2])/f[t1, s1, t2, s2], {t1, 0,
2 Pi}, {s1, 0, 2 Pi}, {t2, 0, 2 Pi}, {s2, 0, 2 Pi}, Method -> "LocalAdaptive"]
(* 1423.97 *)


Global adaptive method also converges to this value with option Method -> {"GlobalAdaptive", "MaxErrorIncreases" -> bigNumber} but it is much slower.

This value has expected order of magnitude because your function is of the order of 1 and the integration volume is

(2 π)^4 // N
(* 1558.55 *)