Hypercube Line Picking

Author

Eric W. Weisstein
March 10, 2006

This notebook downloaded from http://mathworld.wolfram.com/notebooks/ComputationalGeometry/HypercubeLinePicking.nb.

For more information, see Eric's MathWorld entry http://mathworld.wolfram.com/HypercubeLinePicking.html.

©2006 Wolfram Research, Inc. except for portions noted otherwise

Additional Authors

1-D integral for n-D by
Michael Trott
February 23, 2005

Simplified form for Δ(4) and Δ(5) from
Jonathan Borwein
Ferbrary 13–15, 2006

Code

Delta[n_,opts___]:=Module[{v=Array[x,{2,n}]},
Integrate[Sqrt[Plus@@(Subtract@@v^2)],Sequence@@({#,0,1}&/@Flatten[v]),opts]
]
NDelta[n_,opts___]:=Module[{v=Array[x,{2,n}]},
NIntegrate[Evaluate[Sqrt[Plus@@(Subtract@@v^2)]],Evaluate[Sequence@@({#,0,1}&/@Flatten[v])],
opts]
]

By Michael Trott

cubeDistanceIntegrand[n_Integer?Positive,x_]:=Module[{i},i[p_Plus]:=i/@p;i[a_?NumericQ b_]:=a i[b];
(#-Normal[Series[#,{x,0,-1}]])&[FixedPoint[ExpandAll[#/.(i[a_ x^m_?(#<-1&)]:>-i[(1/(m+1) x^(m+1)) D[a,x]])]&,i[-((x (-(1/x^3)+x^(-3-2*n) (-1+E^(-x^2)+Sqrt[Pi] x Erf[x])^n))/Sqrt[Pi])]]/.i->Identity]]

Note that the following integral is equivalent to (Sqrt[s])

Integrate[-1/Sqrt[Pi]/2 p^(-3/2) (Exp[-p s]-1),{p,0,Infinity},Assumptions->Re[s]>0]

so the square root (which is what causes most of the difficulty in symbolic evaluation of the multiple integrals) can be eliminated at the
expense of introducing one additional integral.

A second trick involves rewriting the integrands using partial integration to place them in a form that Mathematica can do out of the box.

Numerical

int[n_,opts___]:=Module[{x,l},
l=Array[x,{2,n}];
NIntegrate[Evaluate[Sqrt[#.#]&[Subtract@@l]],Evaluate[Sequence@@({#,0,1}&/@Flatten[l])],opts]
]

QuasiMonteCarlo

Table[int[n,Method->QuasiMonteCarlo],{n,8}]//Timing

Table[int[n,Method->QuasiMonteCarlo,MaxPoints->10^5,PrecisionGoal->10],{n,8}]//Timing

Table[int[n,Method->QuasiMonteCarlo,MaxPoints->10^6,PrecisionGoal->10],{n,8}]//Timing

Table[int[n,Method->QuasiMonteCarlo,MaxPoints->10^7,PrecisionGoal->10],{n,8}]//Timing

int[5,Method->QuasiMonteCarlo,MaxPoints->5 10^6,PrecisionGoal->10]//Timing

int[5,Method->QuasiMonteCarlo,MaxPoints->7 10^6,PrecisionGoal->10]//Timing

int[5,Method->QuasiMonteCarlo,MaxPoints->10^7,PrecisionGoal->10]//Timing

MonteCarlo

Table[int[n,Method->MonteCarlo],{n,8}]//Timing

Table[int[n,Method->MonteCarlo,MaxPoints->10^5,PrecisionGoal->10],{n,8}]//Timing

Table[int[n,Method->MonteCarlo,MaxPoints->10^6,PrecisionGoal->10],{n,8}]//Timing

Table[int[n,Method->MonteCarlo,MaxPoints->10^7,PrecisionGoal->10],{n,8}]//Timing

Plots

<<Graphics`Colors`

Plot[{Sqrt[n]/3,Sqrt[n/6]Sqrt[(1+2Sqrt[1-3/(5n)])/3]},{n,1,8},
PlotStyle->{Blue,Red},
Epilog->{PointSize[.03],MapIndexed[Point[{#2[[1]],#1}]&,{1/3,.521,.6617,.77766,.87852,.96895,1.05159,1.12817}]},
AxesLabel->TraditionalForm/@{n,Δ[n]},
AxesOrigin->{1,1}
]

Δ(n)

Bailey et al.

Integrate[Exp[-(u^2)(((x-y))^2)],{x,0,1},{y,0,1}]

v5.2

d[u_]:=((Exp[-u^2]-1+(Sqrt[π])u Erf[u])/(u^2))

Table[With[{s=1},(s/(Gamma[1-(s/2)]))Integrate[((1-((d[u])^n))/(u^(s+1))),{u,0,∞}]],{n,3}]

N[%]

(((∫_0^∞)Apart[((1-((((-1+(E^(-(u^2)))+(Sqrt[π]) u Erf[u]))^2)/(u^4)))/(u^2))]\[DifferentialD]u)/(Sqrt[π]))

With[{s=1},(s/(Gamma[1-(s/2)]))Integrate[((1-((d[u])^n))/(u^(s+1))),{u,0,∞}]]

v6

d[u_]:=((Exp[-u^2]-1+(Sqrt[π])u Erf[u])/(u^2))

Table[With[{s=1},(s/(Gamma[1-(s/2)]))Integrate[((1-((d[u])^n))/(u^(s+1))),{u,0,∞}]],{n,3}]

N[%]

(((∫_0^∞)Apart[((1-((((-1+(E^(-(u^2)))+(Sqrt[π]) u Erf[u]))^2)/(u^4)))/(u^2))]\[DifferentialD]u)/(Sqrt[π]))

With[{s=1},(s/(Gamma[1-(s/2)]))Integrate[((1-((d[u])^n))/(u^(s+1))),{u,0,∞}]]

Trott

MakeBoxes[label,TraditionalForm]:=RowBox[{SubscriptBox["I","n"],"(","x",")"}]

Plot[Evaluate[cubeDistanceIntegrand[#,x]&/@Range[5]],{x,0,2},PlotStyle->{Red,Yellow,Green,Blue,Purple},
AxesLabel->TraditionalForm/@{x,label},PlotRange->{-1,4}]

Plot[Evaluate[Log[cubeDistanceIntegrand[#,x]&/@Range[5]]],{x,2,10},PlotStyle->{Red,Yellow,Green,Blue,Purple},
AxesLabel->TraditionalForm/@{x,label}]

(zero=(Sqrt[π])Limit[cubeDistanceIntegrand[#,x]&/@Range[50],x->0])//Timing

(zero=(Sqrt[π])Limit[cubeDistanceIntegrand[#,x]&/@Range[20],x->0])//Timing

cubeDistanceIntegrand[#,x]&/@Range[5]//ColumnForm//TraditionalForm

Integrate[cubeDistanceIntegrand[1,x],{x,0,∞}]

Integrate[cubeDistanceIntegrand[2,x],{x,0,∞}]

TrigToExp[Integrate[cubeDistanceIntegrand[3,x],{x,0,∞}]]//FullSimplify

Table[NIntegrate[cubeDistanceIntegrand[n,x],{x,0,∞},PrecisionGoal->11,WorkingPrecision->30],{n,8}]//Timing

Δ(3)

v4.2

v5

NIntegrate[(Sqrt[(x1-x2)^2+(y1-y2)^2+(z1-z2)^2+(w1-w2)^2]),{x1,0,1},{x2,0,1},{y1,0,1},{y2,0,1},{z1,0,1},{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo,PrecisionGoal->3,MaxPoints->10^6]

v6

1st Integral

i1=Simplify[Integrate[(Sqrt[(x1-x2)^2+(y1-y2)^2+(z1-z2)^2+(w1-w2)^2]),x1]]

int1=Simplify[(i1/.x1->1)-(i1/.x1->0)]

NIntegrate[int1,{x2,0,1},{y1,0,1},{y2,0,1},{z1,0,1},{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo,PrecisionGoal->3,MaxPoints->10^6]

2nd Integral

i2=Simplify[Integrate[int1,x2]]

int2=Simplify[(i2/.x2->1)-(i2/.x2->0)]

NIntegrate[int2,{y1,0,1},{y2,0,1},{z1,0,1},{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo,PrecisionGoal->3,MaxPoints->10^6]

3rd Integral

int2s=FullSimplify[int2]

int2s=(1/6) (1/(-(1/(4 (Sqrt[(((w1-w2))^2)+(((y1-y2))^2)+(((z1-z2))^2)])))-((Sqrt[1+(((w1-w2))^2)+(((y1-y2))^2)+(((z1-z2))^2)])/(4 ((((w1-w2))^2)+(((y1-y2))^2)+(((z1-z2))^2)))))+2 (Sqrt[1+(((w1-w2))^2)+(((y1-y2))^2)+(((z1-z2))^2)])+3 ((((w1-w2))^2)+(((y1-y2))^2)+(((z1-z2))^2)) (-2 Log[-1+(Sqrt[1+(((w1-w2))^2)+(((y1-y2))^2)+(((z1-z2))^2)])]+Log[(((w1-w2))^2)+(((y1-y2))^2)+(((z1-z2))^2)]));

NIntegrate[int2s,{y1,0,1},{y2,0,1},{z1,0,1},{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo,PrecisionGoal->3,MaxPoints->10^6]

i3=Simplify[Integrate[int2s,y1],0<y2<1&&0<z1<1&&0<z2<1&&0<w1<1&&0<w2<1]

int3=Simplify[(i3/.y1->1)-(i3/.y1->0),0<y2<1&&0<z1<1&&0<z2<1&&0<w1<1&&0<w2<1]

NIntegrate[int3,{y2,0,1},{z1,0,1},{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo,PrecisionGoal->3,MaxPoints->10^6]//Timing

4th Integral (in z1)

i4z1=Simplify[Integrate[int3,z1],0<y2<1&&0<z2<1&&0<w1<1&&0<w2<1]

...hangs...

4th Integral (in y2)

i4y2=Simplify[Integrate[int3,y2],0<z1<1&&0<z2<1&&0<w1<1&&0<w2<1]

int4y2=Simplify[(i4y2/.y2->1)-(i4y2/.y2->0),0<z1<1&&0<z2<1&&0<w1<1&&0<w2<1]

NIntegrate[int4y2,{z1,0,1},{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo]

ByteCount[int4y2]

l=ExpandAll[int4y2];

{ByteCount[l],Length[l]}

ints={};
Do[Print[{n,Timing[v=Integrate[FullSimplify[l[[n]]],z1]][[1]],If[FreeQ[v,Integrate],"OK","****"],Short[v]}];
AppendTo[ints,v],
{n,Length[l]}]

Gets stuck because we can't do these:

Timing[notdone=Select[ints,!FreeQ[#,Integrate]&]]

Type 1

notdone[[1,3]]

e=notdone[[1,3]]/.Integrate[(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)])],z1]:>(Integrate[(a Cosh[u])^2 ArcTan[a Cosh[u]],u]/.{a->w1-w2,u->ArcSinh[(z1-z2)/(w1-w2)]})

myint=(e/.z1->1)-(e/.z1->0);

ByteCount[myint]

NIntegrate[myint,{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo]

NIntegrate[(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)])],{z1,0,1},{z2,0,1},{w1,0,1},{w2,0,1}]

ByteCount[myint2=Simplify[myint,0<z2<1&&0<w1<1&&0<w2<1]]//Timing

Type 2

notdone[[4]]

e=notdone[[4]]/.(Integrate[z1^n_(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)])],z1]:>(Integrate[(a Cosh[u])^2(z2+a Sinh[u])^n ArcTan[a Cosh[u]],u]/.{a->w1-w2,u->ArcSinh[(z1-z2)/(w1-w2)]}))

ByteCount[myint=(e/.z1->1)-(e/.z1->0)]//Timing

NIntegrate[myint,{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo]//Timing

NIntegrate[-(2/3) (z1^2) (Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)])],{z1,0,1},{z2,0,1},{w1,0,1},{w2,0,1}]

Type 4

FullSimplify[(z1^2) (Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcCot[√(((((w1-w2))^2)+(((z1-z2))^2)) (2+(((w1-w2))^2)+(((z1-z2))^2)))]/.z1->z2+(w1-w2) Sinh[u]]

Simplify[(a^2) ((Cosh[u])^2) (2+(a^2)+(a^2) ((Sinh[u])^2)),u>0&&a>0]

Integrate[z^2((Sqrt[a^2+z^2]))ArcCot[(Sqrt[(b^2-1+z^2)(b^2+1+z^2)])],z]

Integrate[((Sqrt[a^2+z^2]))ArcCot[(Sqrt[(b^2-1+z^2)(b^2+1+z^2)])],z]

i5=Plus@@done+Plus@@nowdone;

ByteCount[int5=(i5/.z1->1)-(i5/.z1->0)];

,0<z2<1&&0<w1<1&&0<w2<1]

Type 3

NIntegrate[((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(1+(((w1-w2))^2)-(Sqrt[2+(((w1-w2))^2)+(((z1-z2))^2)])+(((z1-z2))^2))/((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]))]),{z1,0,1},{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo]//Timing

Simplify[((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(1+(((w1-w2))^2)-(Sqrt[2+(((w1-w2))^2)+(((z1-z2))^2)])+(((z1-z2))^2))/((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]))])/.{z1->z2+a Sinh[u],w1->w2+a},a>0&&u∈Reals]

i=Integrate[((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(1+(((w1-w2))^2)-(Sqrt[2+(((w1-w2))^2)+(((z1-z2))^2)])+(((z1-z2))^2))/((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]))]),z1]/.Integrate[((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(1+(((w1-w2))^2)-(Sqrt[2+(((w1-w2))^2)+(((z1-z2))^2)])+(((z1-z2))^2))/((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]))]),z1]:>Integrate[(a Cosh[u])^2ArcTan[(1/a)(Sech[u] (1+(a^2)+(a^2) ((Sinh[u])^2)-(Sqrt[2+(a^2)+(a^2) ((Sinh[u])^2)])))] ,u]/.{a->w1-w2,u->ArcSinh[(z1-z2)/(w1-w2)]}

i=(Integrate[((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(1+(((w1-w2))^2)-(Sqrt[2+(((w1-w2))^2)+(((z1-z2))^2)])+(((z1-z2))^2))/((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]))]),z1]/.Integrate[((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(1+(((w1-w2))^2)-(Sqrt[2+(((w1-w2))^2)+(((z1-z2))^2)])+(((z1-z2))^2))/((Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]))]),z1]:>Integrate[(a Cosh[u])^2ArcTan[(1/a)(Sech[u] (1+(a^2)+(a^2) ((Sinh[u])^2)-(Sqrt[2+(a^2)+(a^2) ((Sinh[u])^2)])))] ,u])

i1=(i/.z1->1)-(i/.z1->0);

NIntegrate[i1,{z2,0,1},{w1,0,1},{w2,0,1},Method->QuasiMonteCarlo,MaxPoints->10^5]//Timing

Together

nowdone=notdone/.{Integrate[(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)])],z1]:>(Integrate[(a Cosh[u])^(0+2) ArcTan[a Cosh[u]],u]/.{a->w1-w2,u->ArcSinh[(z1-z2)/(w1-w2)]}),
Integrate[z1^n_.(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)]) ArcTan[(Sqrt[(((w1-w2))^2)+(((z1-z2))^2)])],z1]:>(Integrate[(a Cosh[u])^2(z2+a Sinh[u])^n ArcTan[a Cosh[u]],u]/.{a->w1-w2,u->ArcSinh[(z1-z2)/(w1-w2)]})
};//Timing

done=Select[ints,FreeQ[#,Integrate]&];

Select[nowdone,!FreeQ[#,Integrate]&]

6th Integral

Δ(4)

Numeric

NIntegrate[cubeDistanceIntegrand[4,x],{x,0,Infinity},WorkingPrecision->30,PrecisionGoal->10]//Timing

Analytic

delta[4]=Integrate[#,{x,0,∞}]&/@cubeDistanceIntegrand[4,x]

{
NIntegrate[Exp[-2x^2]Erf[x]^2,{x,0,Infinity},WorkingPrecision->30],
N[1/2 Sqrt[2]ArcTan[Sqrt[2]/4]/Sqrt[Pi],20]
}

Integrate[(E^(-p (x^2)))  Erf[a x]Erf[b x],{x,0,∞}]==ArcTan[b c/Sqrt[(b^2+c^2+p)p]]/Sqrt[Pi p]

Integrate[Exp[-x^2]Erf[x]^3/x,{x,0,∞}]
==(-64*Catalan+4*(8-9*Sqrt[2])*Pi^2+72*Im[PolyLog[2,3*I-(2*I)*Sqrt[2]]]-72*Im[PolyLog[2,(3*I+(2*I)*Sqrt[2])^(-1)]]-288*Im[PolyLog[2,(1/2+I/2)*(-2+Sqrt[2])]]+288*Im[PolyLog[2,(-I)*(-1+Sqrt[2])]]-288*Im[PolyLog[2,I*(-1+Sqrt[2])]]+288*Im[PolyLog[2,(-1+I)/(2+Sqrt[2])]]-4*PolyGamma[1,1/12]+9*Sqrt[2]*PolyGamma[1,1/8]+9*Sqrt[2]*PolyGamma[1,3/8]-4*PolyGamma[1,5/12])/(24*Pi)
==(1/24Pi)*(-64*Catalan+4*(8-9*Sqrt[2])*Pi^2+144*Im[PolyLog[2,3*I-2*I*Sqrt[2]]]-4*PolyGamma[1,1/12]+9*Sqrt[2]*PolyGamma[1,1/8]+9*Sqrt[2]*PolyGamma[1,3/8]-4*(PolyGamma[1,5/12]+72*I*(PolyLog[2,(-1+I)-(-1)^(3/4)]-PolyLog[2,(1/2+I/2)*(-2+Sqrt[2])]+PolyLog[2,(-I)*(-1+Sqrt[2])]-PolyLog[2,I*(-1+Sqrt[2])])))

{NIntegrate[Exp[-x^2]Erf[x]^3/x,{x,0,∞}],
(-64*Catalan+4*(8-9*Sqrt[2])*Pi^2+72*Im[PolyLog[2,3*I-(2*I)*Sqrt[2]]]-72*Im[PolyLog[2,(3*I+(2*I)*Sqrt[2])^(-1)]]-288*Im[PolyLog[2,(1/2+I/2)*(-2+Sqrt[2])]]+288*Im[PolyLog[2,(-I)*(-1+Sqrt[2])]]-288*Im[PolyLog[2,I*(-1+Sqrt[2])]]+288*Im[PolyLog[2,(-1+I)/(2+Sqrt[2])]]-4*PolyGamma[1,1/12]+9*Sqrt[2]*PolyGamma[1,1/8]+9*Sqrt[2]*PolyGamma[1,3/8]-4*PolyGamma[1,5/12])/(24*Pi),
(1/24/Pi)*(-64*Catalan+4*(8-9*Sqrt[2])*Pi^2+144*Im[PolyLog[2,3*I-2*I*Sqrt[2]]]-4*PolyGamma[1,1/12]+9*Sqrt[2]*PolyGamma[1,1/8]+9*Sqrt[2]*PolyGamma[1,3/8]-4*(PolyGamma[1,5/12]+72*I*(PolyLog[2,(-1+I)-(-1)^(3/4)]-PolyLog[2,(1/2+I/2)*(-2+Sqrt[2])]+PolyLog[2,(-I)*(-1+Sqrt[2])]-PolyLog[2,I*(-1+Sqrt[2])])))}//N//Chop

delta4=delta[4]/.{
Integrate[c_ (E^(-2 (x^2)))  ((Erf[x])^2),{x,0,∞}]->c/2 Sqrt[2]ArcTan[Sqrt[2]/4]/Sqrt[Pi],
Integrate[c_ Exp[-x^2]Erf[x]^3/x,{x,0,∞}]->
c(1/24/Pi)*(-64*Catalan+4*(8-9*Sqrt[2])*Pi^2+144*Im[PolyLog[2,3*I-2*I*Sqrt[2]]]-4*PolyGamma[1,1/12]+9*Sqrt[2]*PolyGamma[1,1/8]+9*Sqrt[2]*PolyGamma[1,3/8]-4*(PolyGamma[1,5/12]+72*I*(PolyLog[2,(-1+I)-(-1)^(3/4)]-PolyLog[2,(1/2+I/2)*(-2+Sqrt[2])]+PolyLog[2,(-I)*(-1+Sqrt[2])]-PolyLog[2,I*(-1+Sqrt[2])])))
}

TraditionalForm[delta4=-(23/135)+(73/(315 (Sqrt[2])))+(8/(35 (Sqrt[3])))-((16 π)/315)+((ArcSinh[1])/14)+(104/105) ArcSinh[(1/(Sqrt[2]))]-(68/105) (Sqrt[2]) ArcTan[(1/(2 (Sqrt[2])))]+((197 Log[3])/420)+(1/180) (64 Catalan-4 (8-9 (Sqrt[2])) (π^2)-144 Im[PolyLog[2,3 I-2 I (Sqrt[2])]]+4 PolyGamma[1,(1/12)]-9 (Sqrt[2]) PolyGamma[1,(1/8)]-9 (Sqrt[2]) PolyGamma[1,(3/8)]+4 (PolyGamma[1,(5/12)]+72 I (PolyLog[2,(-1+I)-(((-1))^(3/4))]-PolyLog[2,((1/2)+(I/2)) (-2+(Sqrt[2]))]+PolyLog[2,-I (-1+(Sqrt[2]))]-PolyLog[2,I (-1+(Sqrt[2]))])))]

N[delta4,20]

Borwein's Simplified Version

delta4b=(14/15)Catalan+(1/10)(π^2)(Sqrt[2])-(34/105)π(Sqrt[2])-(16/315)π+(197/420)Log[3]+(52/105)Log[2+(Sqrt[3])]+(1/14)Log[1+(Sqrt[2])]+(8/105)(Sqrt[3])+(73/630)(Sqrt[2])-(23/135)+(136/105)(Sqrt[2])ArcCot[(Sqrt[2])]-((Sqrt[2])/40)(PolyGamma[1,(1/8)]+PolyGamma[1,(3/8)])-(8/5)Im[PolyLog[2,I(1-(Sqrt[2]))]+PolyLog[2,(1-I)((1/(Sqrt[2]))-1)]];

TraditionalForm[HoldForm[(14/15)Catalan+(1/10)(π^2)(Sqrt[2])-(34/105)π(Sqrt[2])-(16/315)π+(197/420)Log[3]+(52/105)Log[2+(Sqrt[3])]+(1/14)Log[1+(Sqrt[2])]+(8/105)(Sqrt[3])+(73/630)(Sqrt[2])-(23/135)+(136/105)(Sqrt[2])ArcCot[(Sqrt[2])]-((Sqrt[2])/40)(PolyGamma[1,(1/8)]+PolyGamma[1,(3/8)])-(8/5)Im[PolyLog[2,I(1-(Sqrt[2]))]+PolyLog[2,(1-I)((1/(Sqrt[2]))-1)]]]]

Timing[diff=FullSimplify[delta4-delta4b]]

Timing[diff=Simplify[delta4-delta4b]]

Block[{$MaxExtraPrecision=1000},N[diff,1000]]

Developer`ZeroQ[diff]

N[delta4b,20]

Borwein Again in Terms of (Subscript[L, -8])(2)

(14/15)Catalan-(34/105)π(Sqrt[2])-(16/315)π+(197/420)Log[3]+(52/105)Log[2+(Sqrt[3])]+(1/14)Log[1+(Sqrt[2])]+(8/105)(Sqrt[3])+(73/630)(Sqrt[2])-(23/135)+(136/105)(Sqrt[2])ArcTan[((Sqrt[2])/2)]-(4/5)(Sqrt[2])LSeries[-8,2]+(8/5)Im[PolyLog[2,I((2^(1/2))-1)]+PolyLog[2,(I+1)((2^(-1/2))-1)]]/.LSeries[-8,2]:>(1/64)(PolyGamma[1,1/8]+PolyGamma[1,3/8]-PolyGamma[1,5/8]-PolyGamma[1,7/8])//N

TraditionalForm[(14/15)Catalan-(34/105)π(Sqrt[2])-(16/315)π+(197/420)Log[3]+(52/105)Log[2+(Sqrt[3])]+(1/14)Log[1+(Sqrt[2])]+(8/105)(Sqrt[3])+(73/630)(Sqrt[2])-(23/135)+(136/105)(Sqrt[2])ArcTan[((Sqrt[2])/2)]-(4/5)(Sqrt[2])LSeries[-8,2]+(8/5)Im[PolyLog[2,I((2^(1/2))-1)]+PolyLog[2,(I+1)((2^(-1/2))-1)]]]

Borwein Take 3

N[26/15*Catalan-34/105*Pi*2^(1/2)-16/315*Pi+197/420*Log[3]+52/105*Log[2+3^(1/2)]+1/14*Log[1+2^(1/2)]+8/105*3^(1/2)+73/630*2^(1/2)-23/135+136/105*2^(1/2)*ArcTan[1/2*2^(1/2)]-1/5*Pi*Log[1+2^(1/2)]+4/5*Log[1+2^(1/2)]*ArcSin[2/3-1/6*2^(1/2)]+4/5*ClausenCl[2,-ArcSin[2/3-1/6*2^(1/2)]]-4/5*ClausenCl[2,ArcSin[2/3+1/6*2^(1/2)]],10]//Chop

HoldForm[26/15*Catalan-34/105*Pi*2^(1/2)-16/315*Pi+197/420*Log[3]+52/105*Log[2+3^(1/2)]+1/14*Log[1+2^(1/2)]+8/105*3^(1/2)+73/630*2^(1/2)-23/135+136/105*2^(1/2)*ArcTan[1/2*2^(1/2)]-1/5*Pi*Log[1+2^(1/2)]+4/5*Log[1+2^(1/2)]*ArcSin[2/3-1/6*2^(1/2)]+4/5*ClausenCl[2,-ArcSin[2/3-1/6*2^(1/2)]]-4/5*ClausenCl[2,ArcSin[2/3+1/6*2^(1/2)]]]//TraditionalForm

Borwein Take 4

N[expr=With[{α=-ArcSin[2/3-1/6*2^(1/2)]},26/15*Catalan-34/105*Pi*2^(1/2)-16/315*Pi+197/420*Log[3]+52/105*Log[2+3^(1/2)]+1/14*Log[1+2^(1/2)]+8/105*3^(1/2)+73/630*2^(1/2)-23/135+136/105*2^(1/2)*ArcTan[1/2*2^(1/2)]-1/5*Pi*Log[1+2^(1/2)]-4/5α*Log[1+2^(1/2)]+4/5*ClausenCl[2,α]-4/5*ClausenCl[2,α+π/2]],10]//Chop

N[expr=With[{α=ArcSin[1/6*2^(1/2)-2/3]},26/15*Catalan-34/105*Pi*2^(1/2)-16/315*Pi+197/420*Log[3]+52/105*Log[2+3^(1/2)]+1/14*Log[1+2^(1/2)]+8/105*3^(1/2)+73/630*2^(1/2)-23/135+136/105*2^(1/2)*ArcTan[1/2*2^(1/2)]-1/5*Pi*Log[1+2^(1/2)]-4/5α*Log[1+2^(1/2)]+4/5*ClausenCl[2,α]-4/5*ClausenCl[2,α+π/2]],10]//Chop

HoldForm[26/15*Catalan-34/105*Pi*2^(1/2)-16/315*Pi+197/420*Log[3]+52/105*Log[2+3^(1/2)]+1/14*Log[1+2^(1/2)]+8/105*3^(1/2)+73/630*2^(1/2)-23/135+136/105*2^(1/2)*ArcTan[1/2*2^(1/2)]-1/5*Pi*Log[1+2^(1/2)]-4/5α*Log[1+2^(1/2)]+4/5*ClausenCl[2,α]-4/5*ClausenCl[2,α+π/2]]//TraditionalForm

Δ(5)

Numeric

NIntegrate[cubeDistanceIntegrand[5,x],{x,0,Infinity},WorkingPrecision->30,PrecisionGoal->10]//Timing

Analytic

delta[5]=Integrate[#,{x,0,∞}]&/@cubeDistanceIntegrand[5,x]

delta5=delta[5]/.{
Integrate[c_. Exp[p_. x^2]Erf[x]^2,{x,0,∞}]:>c(1/(Sqrt[-π p]))ArcTan[(1/(Sqrt[-p(2-p)]))],
Integrate[c_ Exp[-x^2]Erf[x]^3/x,{x,0,∞}]->
c(1/24/Pi)*(-64*Catalan+4*(8-9*Sqrt[2])*Pi^2+144*Im[PolyLog[2,3*I-2*I*Sqrt[2]]]-4*PolyGamma[1,1/12]+9*Sqrt[2]*PolyGamma[1,1/8]+9*Sqrt[2]*PolyGamma[1,3/8]-4*(PolyGamma[1,5/12]+72*I*(PolyLog[2,(-1+I)-(-1)^(3/4)]-PolyLog[2,(1/2+I/2)*(-2+Sqrt[2])]+PolyLog[2,(-I)*(-1+Sqrt[2])]-PolyLog[2,I*(-1+Sqrt[2])])))
}

-(449/3465)+(3239/(31185 (Sqrt[2])))+(568/(1155 (Sqrt[3])))-((380 (Sqrt[5]))/6237)-((4 π)/189)+((4 (π^2))/315)+(64/189) ArcSinh[(1/2)]+((ArcSinh[1])/54)+(40/63) ArcSinh[(1/(Sqrt[2]))]-(73/63) (Sqrt[2]) ArcTan[(1/(2 (Sqrt[2])))]-((8 ArcTan[(1/(Sqrt[15]))])/(7 (Sqrt[3])))+(∫_0^∞)-((52 (E^(-2 (x^2))) π ((Erf[x])^3))/(63 x))\[DifferentialD]x+((295 Log[3])/252)-(1/1008)(5 (-64 Catalan+4 (8-9 (Sqrt[2])) (π^2)+144 Im[PolyLog[2,3 I-2 I (Sqrt[2])]]-4 PolyGamma[1,(1/12)]+9 (Sqrt[2]) PolyGamma[1,(1/8)]+9 (Sqrt[2]) PolyGamma[1,(3/8)]-4 (PolyGamma[1,(5/12)]+72 I (PolyLog[2,(-1+I)-(((-1))^(3/4))]-PolyLog[2,((1/2)+(I/2)) (-2+(Sqrt[2]))]+PolyLog[2,-I (-1+(Sqrt[2]))]-PolyLog[2,I (-1+(Sqrt[2]))]))))//N

TraditionalForm[delta5=-(449/3465)+(3239/(31185 (Sqrt[2])))+(568/(1155 (Sqrt[3])))-((380 (Sqrt[5]))/6237)-((4 π)/189)+((4 (π^2))/315)+(64/189) ArcSinh[(1/2)]+((ArcSinh[1])/54)+(40/63) ArcSinh[(1/(Sqrt[2]))]-(73/63) (Sqrt[2]) ArcTan[(1/(2 (Sqrt[2])))]-((8 ArcTan[(1/(Sqrt[15]))])/(7 (Sqrt[3])))+(∫_0^∞)-((52 (E^(-2 (x^2))) π ((Erf[x])^3))/(63 x))\[DifferentialD]x+((295 Log[3])/252)-(1/1008)(5 (-64 Catalan+4 (8-9 (Sqrt[2])) (π^2)+144 Im[PolyLog[2,3 I-2 I (Sqrt[2])]]-4 PolyGamma[1,(1/12)]+9 (Sqrt[2]) PolyGamma[1,(1/8)]+9 (Sqrt[2]) PolyGamma[1,(3/8)]-4 (PolyGamma[1,(5/12)]+72 I (PolyLog[2,(-1+I)-(((-1))^(3/4))]-PolyLog[2,((1/2)+(I/2)) (-2+(Sqrt[2]))]+PolyLog[2,-I (-1+(Sqrt[2]))]-PolyLog[2,I (-1+(Sqrt[2]))]))))]

Borwein Simplified Version

With[{α=ArcSin[((Sqrt[2])/6)-(2/3)]},(65/42)Catalan-(380/6237)(Sqrt[5])+(568/3465)(Sqrt[3])-(4/189)π-(449/3465)-(73/63)(Sqrt[2])ArcTan[(1/4)(Sqrt[2])]-(184/189)Log[2]+(64/189)Log[(Sqrt[5])+1]+(1/54)Log[1+(Sqrt[2])]+(40/63)Log[(Sqrt[2])+(Sqrt[6])]-(5/28)π Log[1+(Sqrt[2])]+(52/63)π Log[2]+(295/252)Log[3]+(4/315)(π^2)+(3239/62370)(Sqrt[2])-(8/21)(Sqrt[3])ArcCot[(Sqrt[15])]-(52/63)π Log[(Sqrt[2])+(Sqrt[6])]-(5/7)α Log[1+(Sqrt[2])]+(5/7)ClausenCl[2,α]-(5/7)ClausenCl[2,α+(π/2)]+(52/63)Integrate[((ArcSec[x])/(Sqrt[(x^2)-4x+3])),{x,3,4}]
]//N[#,10]&

Missing Piece

Plot[(( (E^(-2 (x^2))) ((Erf[x])^3))/x),{x,0,10},PlotStyle->Red,PlotRange->All]

NIntegrate[(((E^(-2 (x^2))) ((Erf[x])^3))/x),{x,0,∞}]

Integrate[(((E^(-2 (x^2))) ((Erf[x])^3))/x),x]

Integrate[(((E^(-2 (x^2))) ((Erf[x])^3))/x),{x,a,∞}]

(8/(π^(3/2)))NIntegrate[(((E^(-2 (x^2))) Exp[-(t1^2+t2^2+t3^2)])/x),{x,0,∞},{t1,0,x},{t2,0,x},{t3,0,x}]

Transform it

(((∫_3^4)(((sec^(-1))(x))/(Sqrt[(x-3) (x-1)]))\[DifferentialD]x)/π)

Plot[ArcSec[x]/Sqrt[(x-3) (x-1)],{x,3,4},PlotStyle->Red]

ArcCsch[Sqrt[2]]-1/Pi Integrate[ArcSec[x]/Sqrt[(x-3) (x-1)],{x,3,4}]//N

ArcCsch[Sqrt[2]]-1/Pi Integrate[ArcSec[x]/Sqrt[(x-3) (x-1)],{x,3,4}]//N

NIntegrate[ArcSec[x]/Sqrt[(x-3) (x-1)],{x,3,4}]

ArcCsch[Sqrt[2]]-1/Pi NIntegrate[ArcSec[x]/Sqrt[(x-3) (x-1)],{x,3,4},WorkingPrecision->220,MaxRecursion->15]//Timing

NIntegrate[ArcSec[x]/Sqrt[(x-3) (x-1)],{x,3,4}]

NIntegrate[ArcSec[2+Cosh[θ]],{θ,0,Im[ArcCos[2]]}]

ArcCsch[Sqrt[2]]-1/Pi NIntegrate[ArcSec[x]/Sqrt[(x-3) (x-1)],{x,3,4},WorkingPrecision->220,MaxRecursion->20]//Timing

NIntegrate[ArcSec[x]/Sqrt[(((x-2))^2)-1],{x,3,4}]

(∫_((sec^(-1))(3))^((sec^(-1))(4)))θ sec(θ) (Sqrt[1+(4/(sec(θ)-3))])\[DifferentialD]θ

Table[ArcSec[Cos[I θ]+2],{θ,0.,Im@ArcCos[2],.2}]

NIntegrate[θ Sec[θ] (Sqrt[1+(4/(-3+Sec[θ]))]),{θ,ArcSec[3],ArcSec[4]}]

-I NIntegrate[ArcSec[Cos[θ]+2],{θ,0,ArcCos[2]}]

Plot[ArcSec[Cos[I θ]+2],{θ,0,Im[ArcCos[2]]},PlotStyle->Red]

Plot[ArcSec[Cosh[θ]+2],{θ,0,Im[ArcCos[2]]},PlotStyle->Red]

NIntegrate[ArcSec[2+Cosh[θ]],{θ,0,Im[ArcCos[2]]}]

ArcCsch[Sqrt[2]]-1/Pi NIntegrate[ArcSec[2+Cosh[θ]],{θ,0,Im[ArcCos[2]]}]

ArcCsch[Sqrt[2]]-1/Pi NIntegrate[ArcSec[2+Cosh[θ]],{θ,0,Im[ArcCos[2]]},WorkingPrecision->220]//Timing

(∫_0^((cos^(-1))(2)))(1/2) (2 i log(((2 i)/(4+(e^(-i θ))+(e^(i θ))))+(Sqrt[1-(4/(((4+(e^(-i θ))+(e^(i θ))))^2))]))+π)dθ

Integrate[(1/2) (π+2 I Log[((2 I)/(4+(E^(-I θ))+(E^(I θ))))+(Sqrt[1-(4/(((4+(E^(-I θ))+(E^(I θ))))^2))])]),{θ,0,ArcCos[2]}]//Timing

Rest

N[delta5,20]

FullSimplify[DeleteCases[delta[5],_Integrate]]

TrigToExp[FullSimplify[TrigToExp[%]]]

Borwein Version

TraditionalForm[delta5b=(-25/84*Sqrt[2]+4/315)*Pi^2+(-52/63*Log[Sqrt[2]+Sqrt[6]]-4/189+52/63*Log[2])*Pi-449/3465+3239/62370*Sqrt[2]+568/3465*Sqrt[3]-380/6237*Sqrt[5]+1/54*Log[1+Sqrt[2]]+40/63*Log[Sqrt[2]+Sqrt[6]]-184/189*Log[2]+64/189*Log[Sqrt[5]+1]+295/252*Log[3]-80/21*Catalan-73/63*Sqrt[2]*ArcTan[1/4*Sqrt[2]]-8/21*ArcTan[1/15*Sqrt[15]]*Sqrt[3]+25/336*Sqrt[2]*PolyGamma[1,3/8]+100/21*Im[PolyLog[2,-1+1/2*Sqrt[2]+I*(1-1/2*Sqrt[2])]]+100/21*Im[PolyLog[2,-I*Sqrt[2]+I]]+25/336*Sqrt[2]*PolyGamma[1,1/8]+65/21*Im[PolyLog[2,3*I-2*I*Sqrt[2]]]+52/63*Integrate[ArcSec[x]/Sqrt[x^2-4*x+3],{x,3,4}]]

N[delta5b]

TraditionalForm[delta5=(∫_0^∞)-((52 (E^(-2 (x^2))) π ((erf(x))^3))/(63 x))\[DifferentialD]x-(10/7) I (Subscript[Li, 2])(I (-1+(Sqrt[2])))+(10/7) I (Subscript[Li, 2])(-I (-1+(Sqrt[2])))-(10/7) I (Subscript[Li, 2])(((1/2)+(I/2)) (-2+(Sqrt[2])))+(10/7) I (Subscript[Li, 2])((-1+I)-(((-1))^(3/4)))+(5/252) (ψ^((1)))((5/12))-((5 (ψ^((1)))((3/8)))/(56 (Sqrt[2])))-((5 (ψ^((1)))((1/8)))/(56 (Sqrt[2])))+(5/252) (ψ^((1)))((1/12))+((295 log(3))/252)-(5/7) Im((Subscript[Li, 2])(3 I-2 I (Sqrt[2])))-((8 (tan^(-1))((1/(Sqrt[15]))))/(7 (Sqrt[3])))-(73/63) (Sqrt[2]) (tan^(-1))((1/(2 (Sqrt[2]))))+(40/63) (sinh^(-1))((1/(Sqrt[2])))+(1/54) (sinh^(-1))(1)+(64/189) (sinh^(-1))((1/2))+((5 (π^2))/(14 (Sqrt[2])))-((46 (π^2))/315)-((4 π)/189)+((20 C)/63)-((380 (Sqrt[5]))/6237)+(568/(1155 (Sqrt[3])))+(3239/(31185 (Sqrt[2])))-(449/3465)]

N[delta5]//Chop

ByteCount/@{delta5,delta5b}

Subtract@@%

s=Simplify/@{delta5,delta5b};//Timing

ByteCount[delta5]

Δ(6)

Numeric

NIntegrate[cubeDistanceIntegrand[6,x],{x,0,Infinity},WorkingPrecision->30,PrecisionGoal->10]//Timing

Analytic

delta[6]=Integrate[#,{x,0,∞}]&/@cubeDistanceIntegrand[6,x]

N[%,20]

delta6=delta[6]/.{
Integrate[c_. Exp[p_. x^2]Erf[x]^2,{x,0,∞}]:>c(1/(Sqrt[-π p]))ArcTan[(1/(Sqrt[-p(2-p)]))],
Integrate[c_ Exp[-x^2]Erf[x]^3/x,{x,0,∞}]->
c(1/24/Pi)*(-64*Catalan+4*(8-9*Sqrt[2])*Pi^2+144*Im[PolyLog[2,3*I-2*I*Sqrt[2]]]-4*PolyGamma[1,1/12]+9*Sqrt[2]*PolyGamma[1,1/8]+9*Sqrt[2]*PolyGamma[1,3/8]-4*(PolyGamma[1,5/12]+72*I*(PolyLog[2,(-1+I)-(-1)^(3/4)]-PolyLog[2,(1/2+I/2)*(-2+Sqrt[2])]+PolyLog[2,(-I)*(-1+Sqrt[2])]-PolyLog[2,I*(-1+Sqrt[2])])))
}

Plus@@Cases[delta[6],_Integrate]

FullSimplify[DeleteCases[delta[6],_Integrate]]

TrigToExp[FullSimplify[TrigToExp[%]]]