Optimal way of performing a high-dimensional numerical integral of a sharped function

Question

I am trying to solve the following numerical integral

fun=E^(1/48 (-113 k[1]^2+102 Cos[x1-x3] k[1] k[3]-113 k[3]^2-2 Cos[x1-x5] k[1] k[5]+102 Cos[x3-x5] k[3] k[5]-113 k[5]^2+102 Cos[x1-x7] k[1] k[7]-2 Cos[x3-x7] k[3] k[7]+102 Cos[x5-x7] k[5] k[7]-113 k[7]^2))/(36864 \[Pi]^4);
Timing[NIntegrate[k[1]k[3]k[5]k[7]fun Cos[2(x1+x3-x5-x7)],{k[1],0,Infinity},{k[3],0,Infinity},{k[5],0,Infinity},{k[7],0,Infinity},{x1,0,2Pi},{x3,0,2Pi},{x5,0,2Pi},{x7,0,2Pi}]]

which is basically a 8-dimensional gaussian function times a cosine.

The result that I'm obtaining is

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::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 5.956319478409388`*^-7 and 8.378969793277401`*^-7 for the integral and error estimates.
{43.3281,5.95632*10^-7}

My problem is that the error estimates of the integral is bigger than its value which makes me doubt about the accuracy of the result and also I would like to minimize the timing of computation.

Is there a method of computing this integral in less than 40 seconds and obtaining a more or less accurate result?

I don't mind if the result is $6 \cdot 10^{-7}$ or $5 \cdot 10^{-7}$, I just want an approximate value in the less time possible.

Answer 1

... a more or less accurate result ? Have you tried Monte-Carlo integration? It's well suited for higher dimensional integrals and it cuts the time down a lot with about the same result:

Timing[NIntegrate[
  k[1] k[3] k[5] k[7] fun Cos[2 (x1 + x3 - x5 - x7)], {k[1], 0, 
   Infinity}, {k[3], 0, Infinity}, {k[5], 0, Infinity}, {k[7], 0, 
   Infinity}, {x1, 0, 2 Pi}, {x3, 0, 2 Pi}, {x5, 0, 2 Pi}, {x7, 0, 2 Pi},
   Method -> {"MonteCarlo", "MaxPoints" -> 10^6, "RandomSeed" -> 12345}]]

(* result: {2.09375, 5.51509*10^-7} *)