I have a long polynomial in four variables x, y, z, w:
((x^2+y^2+z^2+w^2+145/3)^2-4*(9*z^2+16*w^2))^2*((x^2+y^2+z^2+w^2+145/3)^2+296*(x^2+y^2)-4*(9*z^2+16*w^2)) -16*(x^2+y^2)*(x^2+y^2+z^2+w^2+145/3)^2*(37*(x^2+y^2+z^2+w^2+145/3)^2-1369*(x^2+y^2)-7*(225*z^2+448*w^2)) -16*sqrt(3)/9*(x^3-3*x*y^2)*(110*(x^2+y^2+z^2+w^2+145/3)^3 -148*(x^2+y^2+z^2+w^2+145/3)*(110*x^2+110*y^2-297*z^2+480*w^2)) -64*(x^2+y^2)*(3*(729*z^4+4096*w^4)+168*(x^2+y^2)*(15*z^2-22*w^2)) +64*(12100/27*(x^3-3*x*y^2)^2 -7056*(3*x^2*y-y^3)^2) -592240896*z^2*w^2
I’m working with R. I want to use the caracas package (a wrapper of Sympy) to get this expression as a polynomial after doing a change of variables. Namely, I want to substitue x, y, z and w by
a*x - b*y - c*z - d*w, a*y + b*x + c*w - d*z, a*z - b*w + c*x + d*y, a*w + b*z - c*y + d*x
respectively. I tried subs with no luck. Here is the only working way I found:
library(caracas) def_sym(x, y, z, w, a, b, c, d) X <- a*x - b*y - c*z - d*w Y <- a*y + b*x + c*w - d*z Z <- a*z - b*w + c*x + d*y W <- a*w + b*z - c*y + d*x expr <- ((X^2+Y^2+Z^2+W^2+145/3)^2-4*(9*Z^2+16*W^2))^2*((X^2+Y^2+Z^2+W^2+145/3)^2+296*(X^2+Y^2)-4*(9*Z^2+16*W^2)) -16*(X^2+Y^2)*(X^2+Y^2+Z^2+W^2+145/3)^2*(37*(X^2+Y^2+Z^2+W^2+145/3)^2-1369*(X^2+Y^2)-7*(225*Z^2+448*W^2)) -16*sqrt(3)/9*(X^3-3*X*Y^2)*(110*(X^2+Y^2+Z^2+W^2+145/3)^3 -148*(X^2+Y^2+Z^2+W^2+145/3)*(110*X^2+110*Y^2-297*Z^2+480*W^2)) -64*(X^2+Y^2)*(3*(729*Z^4+4096*W^4)+168*(X^2+Y^2)*(15*Z^2-22*W^2)) +64*(12100/27*(X^3-3*X*Y^2)^2 -7056*(3*X^2*Y-Y^3)^2) -592240896*Z^2*W^2 poly <- sympy_func( expr, "Poly", domain = "QQ[a,b,c,d]" )
But after 30 minutes the computation of poly is not finished. Is there a more efficient way?
Advertisement
Answer
Rather than generating the full expression and then requesting coefficients of the expansion when you are done, you can take it term-by-term. I have done so twice with the full expression, but you can do so with the toy expression that is followed by your full expression as a comment:
from sympy import *
from sympy.parsing.sympy_parser import *
from sympy.abc import x,y,z,w,a,b,c,d
eq = x*y + 1#parse_expr('((x^2+y^2+z^2+w^2+145/3)^2-4*(9*z^2+16*w^2))^2*((x^2+y^2+z^2+w^2+145/3)^2+296*(x^2+y^2)-4*(9*z^2+16*w^2)) -16*(x^2+y^2)*(x^2+y^2+z^2+w^2+145/3)^2*(37*(x^2+y^2+z^2+w^2+145/3)^2-1369*(x^2+y^2)-7*(225*z^2+448*w^2)) -16*sqrt(3)/9*(x^3-3*x*y^2)*(110*(x^2+y^2+z^2+w^2+145/3)^3 -148*(x^2+y^2+z^2+w^2+145/3)*(110*x^2+110*y^2-297*z^2+480*w^2)) -64*(x^2+y^2)*(3*(729*z^4+4096*w^4)+168*(x^2+y^2)*(15*z^2-22*w^2)) +64*(12100/27*(x^3-3*x*y^2)^2 -7056*(3*x^2*y-y^3)^2) -592240896*z^2*w^2', transformations=T[:])
reps = {x: a*x - b*y - c*z - d*w,y:a*y + b*x + c*w - d*z,z:a*z - b*w + c*x + d*y,w:a*w + b*z - c*y + d*x}
eq = eq.xreplace(reps)
c = {}
for i in Add.make_args(eq):
f = i.xreplace(reps).expand()
for s in Add.make_args(f):
co, mo = s.as_coeff_mul(x,y,z,w)
c.setdefault(Mul._from_args(mo), []).append(co)
for k in c:
print(k,Add(*c[k])))