Solving using Groebner basis computation

We load the needed packages, define the ring of polynomials we will need, and the polynomial system we will solve:

using AbstractAlgebra, Groebner, AlgebraicSolvers

R, (x,y,z) = QQ["x","y","z"]


P = [
    8*x^2*y^2+5*x*y^3+3*x^3*z+x^2*y*z,
    x^5+2*y^2*z^2+13*y^2*z^3+5*y*z^4,
    8*x^3+12*y^3+x*z^2+3,
    7*x^2*y^4+18*x*y^3*z^2+y^3*z^3
    ]

4-element Vector{AbstractAlgebra.Generic.MPoly{Rational{BigInt}}}:
 3*x^3*z + 8*x^2*y^2 + x^2*y*z + 5*x*y^3
 x^5 + 13*y^2*z^3 + 2*y^2*z^2 + 5*y*z^4
 8*x^3 + x*z^2 + 12*y^3 + 3
 7*x^2*y^4 + 18*x*y^3*z^2 + y^3*z^3

We solve the system P using the Groebner solver:

Xi, G, B = solve_groebner(P);

and we get 6 solutions (the columns of Xi):

Xi

3×6 Matrix{ComplexF64}:
       0.0+0.0im          …           0.0+0.0im
 -0.629961-1.664e-31im            0.31498+0.545562im
       0.0+4.22682e-10im     -1.18437e-10-1.9944e-10im

How does it work ?

First we compute the Groebner basis of P (for the degree reverse lexicographic ordering):

G = groebner(P)

3-element Vector{AbstractAlgebra.Generic.MPoly{Rational{BigInt}}}:
 z^2
 y^3 + 1//4
 x

Then we deduce the basis of quotient by the ideal $(P)$:

B, BIdx = quotient_basis(G); B

6-element Vector{AbstractAlgebra.Generic.MPoly{Rational{BigInt}}}:
 1
 z
 y
 y*z
 y^2
 y^2*z

(Here Bidx is a dictionary of monomials giving their index in the basis B).

Then we compute the matrices of multiplication by the variables in the basis Bof the quotient:

M = [mult_matrix(v, G, B) for v in [x,y,z]]

3-element Vector{Matrix{Rational{BigInt}}}:
 [0//1 0//1 … 0//1 0//1; 0//1 0//1 … 0//1 0//1; … ; 0//1 0//1 … 0//1 0//1; 0//1 0//1 … 0//1 0//1]
 [0//1 0//1 … -1//4 0//1; 0//1 0//1 … 0//1 -1//4; … ; 0//1 0//1 … 0//1 0//1; 0//1 0//1 … 0//1 0//1]
 [0//1 0//1 … 0//1 0//1; 1//1 0//1 … 0//1 0//1; … ; 0//1 0//1 … 0//1 0//1; 0//1 0//1 … 1//1 0//1]

Then we triangularise them jointly in the same basis and deduce the points Xifrom the values on the diagonal of triangularised $M_i$.

For that purpose, we compute the Schur factorization of a random combination of the matrices.

This gives us the points (possibly repeated with their multiplicity):

Xi

3×6 Matrix{ComplexF64}:
       0.0+0.0im          …           0.0+0.0im
 -0.629961-1.664e-31im            0.31498+0.545562im
       0.0+4.22682e-10im     -1.18437e-10-1.9944e-10im