Mfcttrf

I want to calculate
$$
int_mathbbR^n e^-t_1^4-t_2^4- ...-t_n^4-(1-t_1-t_2-...-t_n)^4 textdt_1 textdt_2 ... textdt_n
$$
with the highest possible n.

I tried with

Module[n = 3, intVars, poly,
 intVars = Table[Subscript[t, i], -∞, ∞, i, 1, n - 1];
 poly = Sum[Subscript[t, i]^4, i, 1, n - 1] + (1 - Sum[Subscript[t, i], i, 1, n - 1])^4;

 Print[E^-poly];

 NIntegrate[E^-poly, ##, Method -> Automatic, "SymbolicProcessing" -> 0] & @@ intVars 
]

which seems to work well only for n<6.
Is there any trick I can use?

One idea I like to use for these kinds of integrals is to add an auxiliary variable and a Dirac delta function, convert the Dirac delta function to it's integral formulation, and then do a bunch of simple 1D integrals. For your case, this would proceed as follows, starting from:

$$undersettin mathbbR^nint e^-t_1^4-t_2^4-ldots
-t_n^4-left(1-t_1-t_2-ldots -t_nright)^4$$

Introduce the auxiliary variable $s$ and add a Dirac delta function:

$$undersettin mathbbR^nint undersetsin
mathbbRint e^-t_1^4-t_2^4-ldots -t_n^4-s^4 delta(1-t_1-t_2-ldots -t_n-s)$$

Next, introduce the integral formulation of the Dirac delta function:

$$frac12 pi undersettin
mathbbR^nint undersetsin mathbbRint undersetuin
mathbbRint e^-t_1^4-t_2^4-ldots -t_n^4-s^4 e^i u-i t_1 u-i t_2
u-ldots -i t_n u-i s u$$

Finally, we can do all of the $t$ and $s$ integrals to obtain:

$$frac12 pi int_-infty
^infty e^i u g(u)^n+1 , du$$

where:

$$g(u)=int_-infty ^infty e^-i t u-t^4 , dt$$

Now, let's have Mathematica do these integrals:

g[u_] = Sqrt[2 Pi] FourierTransform[Exp[-t^4], t, u]

2 Gamma[5/4] HypergeometricPFQ[, 1/2, 3/4, u^4/256] -
1/4 u^2 Gamma[3/4] HypergeometricPFQ[, 5/4, 3/2, u^4/256]

So, the desired integral has become:

int[n_, opts:OptionsPattern[NIntegrate]] := (1/(2 Pi)) NIntegrate[
 Cos[u] g[u]^(n+1), 
 u, -Infinity, Infinity,
 opts,
 Method->"GlobalAdaptive", Method->"GaussKronrodRule"
] 

Now, g[u] is real:

Refine[g[u] ∈ Reals, u ∈ Reals]

True

so I use Cos[u] instead of Exp[I u] since the integral is real. I also customize the integration method. Let's check:

int[2, WorkingPrecision->20]
int[2, WorkingPrecision->40]
int[2, WorkingPrecision->60]

1.4733172914977911077

1.473317291497785926905017339845596712841

1.47331729149778592690501733984559670949096610342311667206502

The result seems correct, and improves with higher working precision. Now, for higher orders:

int[4, WorkingPrecision -> 20]
int[5, WorkingPrecision -> 20]
int[6, WorkingPrecision -> 20]
int[40, WorkingPrecision -> 20]

4.4732305211180348293

7.7543594355221995796

13.461688085347942892

4.0351905913672630176*10^9

We can increase n, reducing the accuracy of calculations, for example

f[m_, p_] := 
 Block[n = m, $MinPrecision = p, $MaxPrecision = p, 
 intVars = 
 Table[Subscript[t, i], -[Infinity], [Infinity], i, 1, n];
 poly = Sum[
 Subscript[t, i]^4, i, 1, 
 n] + (1 - Sum[Subscript[t, i], i, 1, n])^4;
 NIntegrate[E^-poly, ##, PrecisionGoal -> p/2, 
 AccuracyGoal -> p/2] & @@ intVars]
f[4, 6]

(* 4.47148*)

f[5, 4]

(* 7.75615*)

 f[6, 2]

(* 15.8289*)