Extremely memory consuming Expand

Question

Expand floods all my 64GB RAM in MMA 12.1 (Windows) just by sorting the powers of 16 variables. Somebody with >64GB RAM could run it. A memory saving alternative would even be more interesting.

a=(x[1]-x[11])^2+(x[2]-x[12])^2+(x[3]-x[13])^2+(x[4]-x[14])^2;
b=(x[1]-x[21])^2+(x[2]-x[22])^2+(x[3]-x[23])^2+(x[4]-x[24])^2;
c=(x[1]-x[31])^2+(x[2]-x[32])^2+(x[3]-x[33])^2+(x[4]-x[34])^2;
d=(x[21]-x[31])^2+(x[22]-x[32])^2+(x[23]-x[33])^2+(x[24]-x[34])^2;
e=(x[11]-x[31])^2+(x[12]-x[32])^2+(x[13]-x[33])^2+(x[14]-x[34])^2;
f=(x[11]-x[21])^2+(x[12]-x[22])^2+(x[13]-x[23])^2+(x[14]-x[24])^2;
g=1/12*Sqrt[4*a*b*c-a*(b+c-d)^2-b*(a+c-e)^2-c*(a+b-f)^2+(b+c-d)*(a+c-e)*(a+b-f)];
g2=g /. x[i_] -> (z[i] - μ[i]) t + μ[i];
taylor= (Series[g2, {t,0,2}] // Normal) /. t -> 1;
(the next line exceeds 64 GB RAM)
taylor=Expand[taylor];
mean= taylor //. z[i_]^2 -> σ^2 + μ[i]^2;
mean= mean //. z[i_] -> μ[i];
FullSimplify[mean==g+[Sigma]^2 (4(ab+ac+bc+ed)-(a+b-f)^2-(a+c-e)^2-(b+c-d)^2-(e+d-f)^2)/(288g) /. x[i_] -> [Mu][i]]

Memory overflow happens after about 15min on a modern CPU. The code is adapted from here and here.

If the code outputs True then your code fix is right or you had enough memory.

In case you want to test the code and you need a running example without high memory demands then set in the 7th line g=Sqrt[a];.

The result will be different with and without Expand, meaning that this command cannot be omitted.

Answer 1

To avoid confusion I post another answer. The OP is right: the Expand is needed for the correct replacement of the square z terms. You can do the given task with 64 GB of RAM. One way is to do Expand not on the whole expression but in 'chunks', that is on the parts appearing within square roots, nominators and denominators separately. Begin with above code:

a = (x[1] - x[11])^2 + (x[2] - x[12])^2 + (x[3] - x[13])^2 + (x[4] - x[14])^2;
b = (x[1] - x[21])^2 + (x[2] - x[22])^2 + (x[3] - x[23])^2 + (x[4] - x[24])^2;
c = (x[1] - x[31])^2 + (x[2] - x[32])^2 + (x[3] - x[33])^2 + (x[4] - x[34])^2;
d = (x[21] - x[31])^2 + (x[22] - x[32])^2 + (x[23] - x[33])^2 + (x[24] - x[34])^2;
e = (x[11] - x[31])^2 + (x[12] - x[32])^2 + (x[13] - x[33])^2 + (x[14] - x[34])^2;
f = (x[11] - x[21])^2 + (x[12] - x[22])^2 + (x[13] - x[23])^2 + (x[14] - x[24])^2;
g = 1/12*Sqrt[4*a*b*c - a*(b + c - d)^2 - b*(a + c - e)^2 - 
 c*(a + b - f)^2 + (b + c - d)*(a + c - e)*(a + b - f)];
g = g /. x[i_] -> (z[i] - \[Mu][i]) t + \[Mu][i];
taylor = (Series[g, {t, 0, 2}] // Normal) /. t -> 1

the resulting expression appears as short output, -> Show All. Then identify nine terms, copy them 'by hand' into variables t1,...,t9. I worked from the back. Denominators, nominators and brackets within square roots. The result has the following structure:

(1/6)*Sqrt[t9] + (Sqrt[t8]*t7)/t6 + (1/12)*Sqrt[t5]*(-(t4/t3) + t2/t1).

Treat the variables all like

t1 = (((t1 // Expand) //. z[i_]^2 -> \[Sigma]^2 + \[Mu][i]^2) //.z[i_] -> \[Mu][i]);

The total result can be simplified as t9=t8=t5 to

Sqrt[t9]*(1/6 + t7/t6 + (1/12)*(t2/t1 - t4/t3))

and evaluate above expression. It will be VERY large, takes a long time to even display (the nb file will be about half a GB). I would not try Simplification on the whole expression, but maybe on the parts t1...t9 after the replacement of z^2 and z. So it can be done, but will it be useful?

Answer 2

Try this:

a = (x[1] - x[11])^2 + (x[2] - x[12])^2 + (x[3] - x[13])^2 + (x[4] - x[14])^2; 
b = (x[1] - x[21])^2 + (x[2] - x[22])^2 + (x[3] - x[23])^2 + (x[4] - x[24])^2; 
c = (x[1] - x[31])^2 + (x[2] - x[32])^2 + (x[3] - x[33])^2 + (x[4] - x[34])^2; 
d = (x[21] - x[31])^2 + (x[22] - x[32])^2 + (x[23] - x[33])^2 + (x[24] - x[34])^2; 
e = (x[11] - x[31])^2 + (x[12] - x[32])^2 + (x[13] - x[33])^2 + (x[14] - x[34])^2; 
f = (x[11] - x[21])^2 + (x[12] - x[22])^2 + (x[13] - x[23])^2 + (x[14] - x[24])^2; 
g = (1/12)*Sqrt[4*a*b*c - a*(b + c - d)^2 -  b*(a + c - e)^2 + (b + c - d)*
(a + c - e)*(a + b - f) - c*(a + b - f)^2]; 
g = g /. x[i_] -> (z[i] - \[Mu][i])*t + \[Mu][i]; 
taylor = Normal[Series[g, {t, 0, 2}]] /. t -> 1;
mean = taylor //. z[i_]^2 -> \[Sigma]^2 + \[Mu][i]^2; 
mean = Simplify[mean //. z[i_] -> \[Mu][i]]

just omit the expansion