1

I've been familiarising myself with the different options available in the AdaptiveMonteCarlo method for NIntegrate but one of the things that I can't really understand is what the numbers correspond to in the RandomSeed option? Naively I thought that they correspond to which set of pseudo-random numbers that Mathematica generates, but when I change the number it can change the result quite a bit (at least for the integration that I'm carrying out). The documentation doesn't go very deep into how RandomSeed other than that it enables reproducibility of results, but how can I trust which result is correct when chaning the RandomSeed number seems to change the result so much? 

g = 5.*10.^-7;

λ = 5.*10.^-1;

m = 0.1;

f = 1.;

φ = 10.^5;

ε = 10.^-4;

ϕ[t_] = φ*Cos[f*t];


p[P_, α_, β_] := {P*Sin[α]*Cos[β], P*Sin[α]*Sin[β], P*Cos[α]};

q[Q_, a_] := {Q*Sin[a], 0, Q*Cos[a]};

k[X_] := {0, 0, X};

ω[x_] := Sqrt[x.x + m^2];
ωbar[x_] := Sqrt[x.x + m^2 + (g*φ^2)/2] // FullSimplify;

omega[x_, t_] := Sqrt[x.x + m^2 + g/2*ϕ[t]^2] // FullSimplify;

omegadot[x_, t_] = D[omega[x, t], t] // FullSimplify;


ϰ[k_, p_, q_, s_] := λ^2/(96*ωbar[k]*ωbar[p]*ωbar[q]*ωbar[s]);

ζ[x_, t_] = I*omegadot[x, t]/(4*omega[x, t]) // FullSimplify;

delta[ε_, x_] := 1/(2*Sqrt[Pi*ε])*Exp[-(x^2/(4*ε))];


solution[X_] := NDSolve[{y1'[j] == omegadot[k[X], j]/(2*omega[k[X], j])*(y2[j] + y3[j]), y2'[j] == -2*I*omega[k[X], j]*y2[j] + omegadot[k[X], j]/(2*omega[k[X], j])*(2*y1[j] + 1), y3'[j] == 2*I*omega[k[X], j]*y3[j] + omegadot[k[X], j]/(2*omega[k[X], j])*(2*y1[j] + 1), y1[0] == 0, y2[0] == 0, y3[0] == 0}, {y1, y2, y3}, {j, 0, 4*Pi}];

nfreetable = Table[{X, (y1 /. solution[X][[1]])[Pi]}, {X, 0., 12., 10.^-1}];

nfree = Interpolation[nfreetable];


A5[k_, p_, q_, s_] := (1 - (g*φ^2)/(8*ωbar[k]^2))*ωbar[k] + (1 - (g*φ^2)/(8*ωbar[q]^2))*ωbar[q] - (1 - (g*φ^2)/(8*ωbar[p]^2))*ωbar[p] - (1 - (g*φ^2)/(8*ωbar[s]^2))*ωbar[s];

FA5[nk_, np_, nq_, nl_] := nk*(1 + np)*nq*(1 + nl) - (1 + nk)*np*(1 + np)*nl;

FA7[nk_, np_, nq_, nl_] := nk*(1 + nl)*(1 + nq)*np - (1 + nk)*nl*nq*(1 + np);

ncomp5[k_, p_, q_, l_, α_?NumericQ, a_?NumericQ, nk_, np_, nq_, nl_] := ϰ[k, p, q, l]*Sin[α]*Sin[a]*(p.p)*(q.q)*FA5[nk, np, nq, nl]*delta[ε, A5[k, p, q, l]] // Simplify;

ncomp7[k_, p_, q_, l_, α_?NumericQ, a_?NumericQ, nk_, np_, nq_, nl_] := ϰ[k, p, q, l]*Sin[α]*Sin[a]*(p.p)*(q.q)*FA7[nk, np, nq, nl]*delta[ε, A5[k, p, q, l]] // Simplify;


norder2[k_, nk_, np_, nq_, nl_] := Re[NIntegrate[If[Sqrt[(k - p[P, α, β] - q[Q, a]).(k - p[P, α, β] - q[Q, a])] <= 12, 
                                               (2*ncomp5[k, p[P, α, β], q[Q, a], (k - p[P, α, β] - q[Q, a]), α, a, nk, np, nq, nl] + ncomp7[k, p[P, α, β], q[Q, a], (k - p[P, α, β] - q[Q, a]), α, a, nk, np, nq, nl]),
                                               (2*ncomp5[k, p[P, α, β], q[Q, a], (k - p[P, α, β] - q[Q, a]), α, a, nk, np, nq, 0] + ncomp7[k, p[P, α, β], q[Q, a], (k - p[P, α, β] - q[Q, a]), α, a, nk, np, nq, 0])],
                                               {P, 0., 12.}, {Q, 0., 12.}, {α, 0., N[Pi, MachinePrecision]}, {β, 0., N[2*Pi, MachinePrecision]}, {a, 0., N[Pi, MachinePrecision]},
                                               Method -> {"AdaptiveMonteCarlo", "RandomSeed" -> 10, "BisectionDithering" -> 1/100, "MaxPoints" -> 10^6, "SymbolicProcessing" -> 0}, AccuracyGoal -> 3, PrecisionGoal -> 3]];

norder2[k[5], nfree[Sqrt[k[5].k[5]]], nfree[Sqrt[p[P, α, β].p[P, α, β]]], nfree[Sqrt[q[Q, a].q[Q, a]]], nfree[Sqrt[(k[5]-p[P, α, β]-q[Q,a]).((k[5]-p[P, α, β]-q[Q,a]))]]]
user35305
  • 291
  • 1
  • 4
  • 1
    I think this answer can help you understand how AdaptiveMonteCarlo works with randomly generated points. – Anton Antonov Sep 07 '17 at 12:03
  • @AntonAntonov Thanks for the link. Your answer is helpful, but I'm still unsure as to exactly what RandomSeed does?! – user35305 Sep 07 '17 at 12:45
  • The number (or expression in general) fed to "RandomSeed" is used by the internal PRNG as a starting point to generate a sequence of variates; of course, different seeds will yield different sequences, and thus slightly different Monte Carlo results. – J. M.'s missing motivation Sep 08 '17 at 03:16
  • @J.M. Thanks for the explanation. So is there any benefit of choosing a bigger number, e.g. 10 over a smaller one e.g. 2? – user35305 Sep 08 '17 at 10:21
  • The accuracy of the Monte Carlo results is not really correlated to the size of the number used as seed. – J. M.'s missing motivation Sep 08 '17 at 10:23
  • @J.M. Ah ok, so what is the point really of having different numbers to choose as a seed? – user35305 Sep 08 '17 at 11:15
  • It allows you to compare different Monte Carlo evaluations. Theoretically, you should be getting consistent results independent of the seed. In practice, you check consistency by performing different runs on the same seed, or different seeds. – J. M.'s missing motivation Sep 08 '17 at 12:08
  • @J.M. I see, that makes more sense. I don't seem to be able to get consistent results for different seeds, but I do if I use the same seed. – user35305 Sep 08 '17 at 12:28

1 Answers1

1

"RandomSeed"->seed sets the seed used to generate random number. It is a bit like surroudning the call to NIntegrate in BlockRandom and putting an explict call to SeedRandom[seed].

Many numeric and maching learning functions in 11.2 take a new RandomSeeding option which is quite similar in spirit.

Itai Seggev
  • 14,113
  • 60
  • 84