# required package
using Oscar;

# polynomial ring in variables X,Y,Z 
F = QQ; # F = GF(3); # F = GF(2);
R,(X,Y,Z) = polynomial_ring(F,["X","Y","Z"]);

f = X + Y + Z - 1;
g = X^2 + Y^2 + Z^2 - 1;
h = X^3 + Y^3 + Z^3 - 1;

##### solve polynomial equations by Gröbner bases #####

I = ideal(R, [f,g,h]);
gb = groebner_basis(I, ordering = lex(R), complete_reduction = true)
# Gröbner basis with elements
#   1: Z^3 - Z^2                # = Z^2(Z-1)
#   2: Y^2 + Y*Z - Y + Z^2 - Z  # = Y^2 + (Z-1)Y + (Z-1)Z
#   3: X + Y + Z - 1
# with respect to the ordering lex([X, Y, Z]), that is X > Y > Z > 1

a = h - X*g           # = -X*Y^2 - X*Z^2 + X + Y^3 + Z^3 - 1
a = a + (Y^2+Z^2-1)*f # = 2*Y^3 + Y^2*Z - Y^2 + Y*Z^2 - Y + 2*Z^3 - Z^2 - Z

b = h - X^2*f         # = -X^2*Y - X^2*Z + X^2 + Y^3 + Z^3 - 1 
b = b + (X*Y+X*Z-X -Y^2-2*Y*Z-Z^2 +2*Y+2*Z -1)*f
b = b // 3            # = -Y^2*Z + Y^2 - Y*Z^2 + 2*Y*Z - Y + Z^2 - Z

c = g - X*f           # = -X*Y - X*Z + X + Y^2 + Z^2 - 1
c = c + (Y+Z-1)*f
c = c // 2            # = Y^2 + Y*Z - Y + Z^2 - Z
# hence I = (c,b,f) 

d = b + (Z-1)*c       # = Z^3 - Z^2
e = a + (-2*Y+Z-1)*c 
e = e // 3            # = Z^3 - Z^2
# hence I = (d,c,f) 

##### apply the Newton identities ##### 
e2 = (f^2-g +2*f)//2 # = X*Y + X*Z + Y*Z
e3 = (h -f^3+3*f*e2 -3*g-3*e2 +3*f)//3 # = X*Y*Z