Finding singularities for Math SE 779275

mx779275.txt

sage: # http://math.stackexchange.com/a/779318/35416
sage: P1.<x,y,z,w> = QQ[]
sage: p = x^3+y^3+z^3+w*x*y*z
sage: px, py, pz = [P1(p.polynomial(v).differentiate()) for v in [x, y, z]]
sage: px
y*z*w + 3*x^2
sage: pxy = px.resultant(py, x)
sage: pxz = px.resultant(pz, x)
sage: pxyz = pxy.resultant(pxz, y)
sage: f = pxyz.factor()
sage: f
(-27) * (w - 3)^3 * (w + 3)^3 * z^16 * (w^2 - 3*w + 9)^3 * (w^2 + 3*w + 9)^3
sage: ws = flatten([f[i][0].univariate_polynomial().roots(QQbar, False) for i in [0, 1, 3, 4]])
sage: for i in ws:
...       print(i.radical_expression())
3
-3
-3/2*I*sqrt(3) + 3/2
3/2*I*sqrt(3) + 3/2
-3/2*I*sqrt(3) - 3/2
3/2*I*sqrt(3) - 3/2
sage: w6 = [3*QQbar.zeta(6)^k for k in range(6)]
sage: set(w6) == set(ws)
True
sage: pxy
y*z^3*w^3 + 27*y^4
sage: zi = 1
sage: for wi in w6:
...       assert pxyz(z=zi, w=wi).is_zero()
...       for yi in pxy(z=zi, w=wi).univariate_polynomial().roots(QQbar, False):
...           for xi in py(y=yi, z=zi, w=wi).univariate_polynomial().roots(QQbar, False):
...               if p(x=xi, y=yi, z=zi, w=wi).is_zero():
...                   print((wi, xi, yi, zi))
(1.500000000000000? + 2.598076211353316?*I, -0.500000000000000? + 0.866025403784439?*I, 1, 1)
(1.500000000000000? + 2.598076211353316?*I, -0.500000000000000? - 0.866025403784439?*I, -0.500000000000000? - 0.866025403784439?*I, 1)
(1.500000000000000? + 2.598076211353316?*I, 1, -0.500000000000000? + 0.866025403784439?*I, 1)
(-3, 1, 1, 1)
(-3, -0.500000000000000? + 0.866025403784439?*I, -0.500000000000000? - 0.866025403784439?*I, 1)
(-3, -0.500000000000000? - 0.866025403784439?*I, -0.500000000000000? + 0.866025403784439?*I, 1)
(1.500000000000000? - 2.598076211353316?*I, -0.500000000000000? - 0.866025403784439?*I, 1, 1)
(1.500000000000000? - 2.598076211353316?*I, 1, -0.500000000000000? - 0.866025403784439?*I, 1)
(1.500000000000000? - 2.598076211353316?*I, -0.500000000000000? + 0.866025403784439?*I, -0.500000000000000? + 0.866025403784439?*I, 1)
sage: pxy(z=1, w=3)
27*y^4 + 27*y
sage: pxz(z=1, w=3)
27*y^3 + 27
sage: pxy(z=1,w=3).univariate_polynomial().roots(QQbar, False)
[-1, 0, 0.500000000000000? - 0.866025403784439?*I, 0.500000000000000? + 0.866025403784439?*I]
sage: px(z=1, w=3, y=-1)
3*x^2 - 3
sage: py(z=1, w=3, y=-1)
3*x + 3
sage: pz(z=1, w=3, y=-1)
-3*x + 3
sage: (w-3)*(w+3)*(w^2-3*w+9)*(w^2+3*w+9)
w^6 - 729
sage: (w^6-3^6)
w^6 - 729