Ask Your Question
0

Solving polynomial equations using Grobner basis technique

asked 2024-06-18 22:25:46 +0100

hamouda gravatar image

updated 2024-06-19 02:26:29 +0100

Given a list of 24 polynomials over $\mathbb{Z}[x,y,z]$ namely Lf, computed by sage using the code

R.<x,y,z>=QQ['x, y, z']
I=(x*y^3-z)*R
Q= R.quo(I)
N=15
p=3
q=5
L=[(0, 0, 0),(0, 0, 1),(0, 0, 2),(0, 0, 3),(0, 0, 4),(0, 0, 5),(0, 0, 6),(0, 0, 7),(0, 1, 0),(0, 1, 1),(0, 1, 2),(0, 2, 0),(0,2,1), 
(0, 2, 2),(1, 0, 0),(1, 0, 1),(1, 0, 2),(1, 0, 3),(1, 0, 4),(1, 1, 0),(1, 1, 1),(1, 1, 2),(2, 0, 0),(2, 0, 1)]
X=floor(N^2)
Y=floor(N^2)
Z=X*Y^3
phi=((p^4-1)*(q^4-1))/((p-1)*(q-1))
e=79
d=79
f=x*(y^3+(N+1)*y^2+(N^2-2*N+1)*y+N^3-N^2-N+1)+1
M=Matrix(0,len(L))
for (k,i,j) in L:
             g=Q(x^i*y^j*f^k*e^(2-k)).lift()
             h=g(x*X,y*Y,z*Z)
             L1=vector(h[x^i*y^j*z^k] for (k,i,j) in L)
             M = M.stack(L1)
N=M.LLL()  

L2=[x^i*y^j*z^k for (k,i,j) in L]
r=0
Lf=[]
Li=[] 
for i in srange(24):
         D=N[i,:]*transpose(matrix([L2]))
         r=D[0,0]
         Li=[r]
         Lf.extend(Li)
print(Lf)

In the theory, I know that there are three polynomials from the list Lf say $p(x,y,z)$, $q(x,y,z)$, and $r(x,y,z)$ having the same root (x,y,z). To extract such root, it suffices to compute the grobner basis of the three polynomials which yields three polynomilas and solving them to obtain the common root. I have computed the grobner basis by sage, here it is denoted by $E$.

for i in srange(24):
     for j in srange(i+1,24):
          for k in srange(j+1,24):
               J = R.ideal([Lf[i],Lf[j],Lf[k]])
               E=J.groebner_basis()
               if len(E)==3:
                   print(E)

So far, E consists of the three polynomials we want. My probelm is what is the command to solve them.?

edit retag flag offensive close merge delete

Comments

What do you want to do if J.dimension() > 0 (i.e. there are not only common points but there is a common curve or surface)? And where do you seek the point? Over the complex numbers CC?

rburing gravatar imagerburing ( 2024-06-19 12:03:18 +0100 )edit

2 Answers

Sort by » oldest newest most voted
2

answered 2024-06-19 17:34:02 +0100

Max Alekseyev gravatar image

updated 2024-06-19 20:23:08 +0100

The set of solutions, when it's finite, can be computed with .variety() method. In your problem, there are many cases when the ideal does not depend on one of the variables, say z, essentially making it free. To take care of those cases, variety can be computed over the subring depending only on the rest of the variables (say, x and y). The corresponding code is given below.

For example, when it prints {x: 0, y: 0} it means that solutions are given by x=0, y=0 and any z.

In cases when variety still fails, it reports an error like

ERROR: Dimension: 1 Basis: [x - 3796875/8*z - 1/5400, y]
ERROR: Dimension: 1 Basis: [x^2 - 21515625/143*x*z + 14416259765625/143*z^2 - 62/7425*x + 11250/143*z + 1/65154375, y]

In some cases, one variable can be expressed in terms of another (like x in terms of z in the first error above), while others (like the second error above) may lead to quadratic Pell-type, elliptic, or higher-degree equations. One may try to address them with an appropriate functionality in Sage if needed.

from sage.misc.verbose import set_verbose
set_verbose(-1)

for fff in Combinations(Lf,3):
    B = R.ideal(fff).groebner_basis()
    V = set.union(*(set(f.variables()) for f in B))
    R1 = PolynomialRing(QQ,len(V),', '.join(map(str,V)))
    J = R1.ideal(B)
    try:
        S = J.variety()
    except:
        print('ERROR: Dimension:',J.dimension(),'Basis:',B)
        continue
    if S:
        print(S)

ADDED. Here is a reduced list of solutions (without repeats) resulted from cases not involving quadratic or higher-degree equations in each of the variables:

{y: 0}
{y: 0, x: 0}
{y: 0, x: -16/2475}
{y: 0, x: -1/525}
{y: 0, x: -23/225}
{y: 0, x: -19/4050}
{y: 0, x: 3796875/8*z + 1/5400}
{z: -21499/568107421875, y: 0, x: 0}
{z: -79/7688671875, y: 0, x: 0}
{z: -1/2562890625, y: 0, x: 0}
{z: 24964/1606932421875, y: 0, x: -16/2475}
{z: 6241/414333984375, y: 0, x: 2401/327375}
{z: (-6733125*x^2 + 472050*x + 32)/(20182763671875*x - 82012500000), y: 0}
{z: 118579/2145139453125, y: 0, x: -23/225}
{z: 68651/2493692578125, y: 0, x: 967/72975}
{z: -79/5638359375, y: 0, x: -16/2475}
{z: 755161/57249850781250, y: 0, x: -19/4050}
{z: -79/7688671875, y: 0, x: -19/4050}
{z: -553/2562890625, y: 0, x: -23/225}
{z: 6241/843191015625, y: 0, x: -1/525}
{z: -79/17940234375, y: 0, x: -1/525}
{z: -1/2562890625, y: -19/3600, x: 0}
{z: 493/48694921875, y: 19/3600, x: 0}
{z: -2209/1355769140625, y: -23/2700, x: 0}
{z: -79/7688671875, y: -4/1425, x: 0}
{z: 1381/2050312500, y: 2/675, x: 3/25}
{z: -12637/5207793750000, y: -4/675, x: 3/3175}
edit flag offensive delete link more
1

answered 2024-06-20 14:43:33 +0100

Emmanuel Charpentier gravatar image

updated 2024-06-21 13:54:01 +0100

A lazy alternative to @Max Alekseyev's excellent answer is to use one of Sage's symbolic solvers, which already take care of non-zero-dimension solutions, and to filter out the undeterminates returned by these solvers.

After running :

SymUnk=list(map(SR, (x, y, z)))
def foo():
    Sols=[]
    for fff in Combinations(Lf, 3):
        if len((B:=RR.ideal(fff).groebner_basis()))==3:
            SSol=solve([SR(u)==0 for u in B],
                           SymUnk,
                           solution_dict=True)
            Sol=[]
            for S in SSol:
                SRR={}
                for u in S.keys():
                    try:
                        SRR|={RR(u):RR(S[u])}
                    except:
                        continue
                Sol+=[SRR]
            Sols+=[Sol]
    return Sols

Running :

sage: %time Sols=foo()
CPU times: user 2min 1s, sys: 105 ms, total: 2min 1s
Wall time: 1min 48s
sage: len(Sols)
245

returns 245 solutions. However, many of them are duplicates, as can be shown by indexing them on the groebner_basis() solved :

def bar():
    Sols={}
    for fff in Combinations(Lf, 3):
        if len((B:=RR.ideal(fff).groebner_basis()))==3:
            SSol=solve([SR(u)==0 for u in B],
                           SymUnk,
                           solution_dict=True)
            Sol=[]
            for S in SSol:
                SRR={}
                for u in S.keys():
                    try:
                        SRR|={RR(u):RR(S[u])}
                    except:
                        continue
                Sol+=[SRR]
            Sols|={B:Sol}
    return Sols

sage: %time DSols=bar()
CPU times: user 2min 2s, sys: 84.4 ms, total: 2min 2s
Wall time: 1min 49s
sage: len(DSols)
156

These solutions to 156 distinct Groebner bases give 260 solutions :

sage: Set(flatten(list(DSols.values()))).cardinality()
260

In the present case, the laziness isn't too costly (2 minutes is way less than the time necessary to rig up a general solution with polynomial rings tools...).

UPDATE :

It turns out that Sympy's solver does a better job than Sage's for dimension>0 ideals. An alternative way to write the solutions is :

 # List of pairs (triplets, base) with Groebnar base of length 3
Triplets = [[u, b]
            for u in Combinations(Lf, 3)
            if len(b:=Rr.ideal(u).groebner_basis())==3]
# List of *distinct* bases
Bases = list(uniq([u[1] for u in Triplets]))
# Solve them *IN SR* using sympy :
# Precompute the unknowns
Unk = list(map(lambda u:SR(u), Rr.gens()))
# Solve
SSols = [solve([SR(u) for u in v], Unk, algorithm="sympy") for v in Bases]

FWIW :

sage: %time SSols = [solve([SR(u) for u in v], Unk, algorithm="sympy") for v in Bases]
CPU times: user 3.27 s, sys: 0 ns, total: 3.27 s
Wall time: 3.27 s

Sympy's solver seems also faster than Sage's...

Since neither the triplets nor the bases can be used as dictionary keys ("unhashable"), we need a bit of record-keeping :

# Record the base-solutions link
BS=list(map(list, zip(Bases, SSols)))
# Link triplets to their solutions
TBS=[[u[0]]+BS[Bases.index(u[1])] for u in Triplets]

Expressing the SR solutions in the original ring is lazily left to the reader...

HTH

edit flag offensive delete link more

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 2024-06-18 22:25:46 +0100

Seen: 253 times

Last updated: Jun 21