Source code

001/*
002 * $Id: GroebnerBaseSigSeqIter.java 5998 2020-03-17 15:40:05Z kredel $
003 */
004
005package edu.jas.gb;
006
007
008import java.util.ArrayList;
009import java.util.List;
010
011import org.apache.logging.log4j.LogManager;
012import org.apache.logging.log4j.Logger;
013
014import edu.jas.poly.ExpVector;
015import edu.jas.poly.GenPolynomial;
016import edu.jas.poly.GenPolynomialRing;
017import edu.jas.poly.OrderedPolynomialList;
018import edu.jas.poly.PolyUtil;
019import edu.jas.structure.RingElem;
020
021
022/**
023 * Groebner Base signature based sequential iterative algorithm. Implements
024 * Groebner bases after the paper "Signature-based Algorithms to Compute Gröbner
025 * Bases" by Christian Eder and John Perry, ISSAC 2011. Compare the jython+JAS
026 * code in examples/basic_sigbased_gb.py. Originally the Python+Sage code is
027 * from http://www.math.usm.edu/perry/Research/basic_sigbased_gb.py
028 * 
029 * @param <C> coefficient type
030 * @author Heinz Kredel
031 * 
032 * @see edu.jas.application.GBAlgorithmBuilder
033 * @see edu.jas.gbufd.GBFactory
034 * @see edu.jas.gb.GroebnerBaseGGVSigSeqIter
035 * @see edu.jas.gb.GroebnerBaseArriSigSeqIter
036 * @see edu.jas.gb.GroebnerBaseF5zSigSeqIter
037 */
038
039public class GroebnerBaseSigSeqIter<C extends RingElem<C>> extends GroebnerBaseAbstract<C> {
040
041
042    private static final Logger logger = LogManager.getLogger(GroebnerBaseSigSeqIter.class);
043
044
045    private static final boolean debug = logger.isDebugEnabled();
046
047
048    final SigReductionSeq<C> sred;
049
050
051    /**
052     * Constructor.
053     */
054    public GroebnerBaseSigSeqIter() {
055        this(new SigReductionSeq<C>());
056    }
057
058
059    /**
060     * Constructor.
061     * @param red Reduction engine
062     */
063    public GroebnerBaseSigSeqIter(SigReductionSeq<C> red) {
064        super();
065        sred = red;
066    }
067
068
069    /**
070     * Groebner base signature iterative algorithm.
071     * @param modv module variable number.
072     * @param F polynomial list.
073     * @return GB(F) a Groebner base of F.
074     */
075    public List<GenPolynomial<C>> GB(int modv, List<GenPolynomial<C>> F) {
076        List<GenPolynomial<C>> G = normalizeZerosOnes(F);
077        G = PolyUtil.<C> monic(G);
078        if (G.size() <= 1) {
079            return G;
080        }
081        // sort, no reverse
082        //  G = OrderedPolynomialList.<C> sort(G);
083        G = OrderedPolynomialList.<C> sortDegree(G);
084        //no: Collections.reverse(G);
085        logger.info("G-sort = " + G);
086        List<GenPolynomial<C>> Gp = new ArrayList<GenPolynomial<C>>();
087        for (GenPolynomial<C> p : G) {
088            if (logger.isInfoEnabled()) {
089                logger.info("p = " + p);
090            }
091            GenPolynomial<C> pp = red.normalform(Gp, p);
092            if (pp.isZERO()) {
093                continue;
094            }
095            Gp = GB(modv, Gp, p);
096            if (Gp.size() > 0) {
097                if (Gp.get(0).isONE()) {
098                    return Gp;
099                }
100            }
101        }
102        return Gp;
103    }
104
105
106    /**
107     * Groebner base iterated.
108     * @param modv module variable number.
109     * @param G polynomial list of a Groebner base.
110     * @param f polynomial.
111     * @return GB(G,f) a Groebner base of G+(f).
112     */
113    public List<GenPolynomial<C>> GB(int modv, List<GenPolynomial<C>> G, GenPolynomial<C> f) {
114        List<GenPolynomial<C>> F = new ArrayList<GenPolynomial<C>>(G);
115        GenPolynomial<C> g = f.monic();
116        if (F.isEmpty()) {
117            F.add(g);
118            return F; // commutative
119        }
120        if (g.isZERO()) {
121            return F;
122        }
123        if (g.isONE()) {
124            F.clear();
125            F.add(g);
126            return F;
127        }
128        GenPolynomialRing<C> ring = F.get(0).ring;
129        if (!ring.coFac.isField()) {
130            throw new IllegalArgumentException("coefficients not from a field");
131        }
132        if (modv != 0) {
133            throw new UnsupportedOperationException("motv != 0 not implemented");
134        }
135        // add signatures
136        List<SigPoly<C>> Gs = new ArrayList<SigPoly<C>>();
137        for (GenPolynomial<C> p : F) {
138            Gs.add(new SigPoly<C>(ring.getZERO(), p));
139        }
140        SigPoly<C> gs = new SigPoly<C>(ring.getONE(), g);
141        Gs.add(gs);
142        //logger.info("Gs = " + Gs);
143        // construct critical pair list
144        List<SigPair<C>> pairlist = new ArrayList<SigPair<C>>();
145        for (SigPoly<C> p : Gs) { // F via continue
146            if (p.poly.equals(g)) {
147                continue;
148            }
149            pairlist.add(newPair(gs, p, Gs));
150        }
151        //logger.info("start " + pairlist.size());
152        logger.info("start " + Gs);
153
154        List<ExpVector> syz = initializeSyz(F, Gs);
155        List<SigPoly<C>> done = new ArrayList<SigPoly<C>>();
156
157        SigPair<C> pair;
158        //SigPoly<C> pi, pj;
159        GenPolynomial<C> S, H, sigma;
160        while (!pairlist.isEmpty()) {
161            pairlist = pruneP(pairlist, syz);
162            if (pairlist.isEmpty()) {
163                continue;
164            }
165            List<SigPair<C>>[] spl = sred.minDegSubset(pairlist);
166            List<SigPair<C>> Sl = spl[0];
167            long mdeg = sred.minimalSigDegree(Sl);
168            pairlist = spl[1];
169            logger.info("treating " + Sl.size() + " signatures of degree " + mdeg);
170            //logger.info("Sl(" + mdeg + ") = " + Sl);
171            while (!Sl.isEmpty()) {
172                //logger.info("Sl_full = " + sred.sigmas(Sl));
173                Sl = pruneS(Sl, syz, done, Gs);
174                if (Sl.isEmpty()) {
175                    continue;
176                }
177                Sl = sred.sortSigma(Sl);
178                //logger.info("Sl_sort = " + Sl);
179                pair = Sl.remove(0);
180                if (pair == null) {
181                    continue;
182                }
183                //logger.info("sigma = " + pair.sigma);
184                S = SPolynomial(pair);
185                SigPoly<C> Ss = new SigPoly<C>(pair.sigma, S);
186                if (S.isZERO()) {
187                    updateSyz(syz, Ss);
188                    done.add(Ss);
189                    continue;
190                }
191                if (debug) {
192                    logger.debug("ht(S) = " + S.leadingExpVector());
193                }
194
195                SigPoly<C> Hs = sigNormalform(F, Gs, Ss);
196                H = Hs.poly;
197                sigma = Hs.sigma;
198                if (debug) {
199                    logger.info("new polynomial = " + Hs); //.leadingExpVector() );
200                }
201                if (H.isZERO()) {
202                    updateSyz(syz, Hs);
203                    done.add(Hs);
204                    continue;
205                }
206                H = H.monic();
207                if (debug) {
208                    logger.info("ht(H) = " + H.leadingExpVector());
209                }
210
211                if (H.isONE()) {
212                    G.clear();
213                    G.add(H);
214                    logger.info("end " + pairlist);
215                    return G; // since no threads are activated
216                }
217                if (sred.isSigRedundant(Gs, Hs)) {
218                    continue;
219                }
220                if (logger.isInfoEnabled()) {
221                    //logger.info("sigma::h = " + sigma + " :: " + ring.toScript(H.leadingExpVector()));
222                    logger.info("sigma::h = " + sigma + " :: " + H.leadingExpVector());
223                }
224                if (H.length() > 0) {
225                    for (SigPoly<C> p : Gs) {
226                        if (p.poly.isZERO()) {
227                            continue;
228                        }
229                        GenPolynomial<C> tau = p.sigma;
230                        GenPolynomial<C>[] mult = SPolynomialFactors(Hs, p);
231                        //System.out.print("sigma = " + sigma + ", tau = " + tau);
232                        //System.out.println(", mult  = " + Arrays.toString(mult));
233                        ExpVector se = sigma.leadingExpVector();
234                        ExpVector te = tau.leadingExpVector();
235                        if (mult[0].multiply(se).equals(mult[1].multiply(te))) {
236                            //logger.info("skip by sigma");
237                            continue;
238                        }
239                        SigPair<C> pp;
240                        //boolean xy = mult[0].multiply(se).compareTo(mult[1].multiply(te)) > 0;
241                        if (mult[0].multiply(se).compareTo(mult[1].multiply(te)) > 0) {
242                            pp = newPair(sigma.multiply(mult[0]), Hs, p, Gs);
243                        } else {
244                            pp = newPair(tau.multiply(mult[1]), p, Hs, Gs);
245                        }
246                        //System.out.println("new_pair " + pp.sigma + ", xy = " + xy + ", sigma = " + sigma + ", tau = " + tau + ", mult  = " + Arrays.toString(mult) + ", m0*se = " + mult[0].multiply(se) + ", m1*te = " + mult[1].multiply(te) );
247                        if (pp.sigma.degree() == mdeg) { // mdeg is sigma.degree()
248                            Sl.add(pp); // do not check contains
249                        } else {
250                            pairlist.add(pp); // do not check contains
251                        }
252                    }
253                    Gs.add(Hs);
254                    done.add(Hs);
255                }
256            }
257        }
258        logger.info("#sequential list before reduction = " + Gs.size());
259        List<GenPolynomial<C>> Gp = sred.polys(Gs);
260        //logger.info("G_full = " + Gp);
261        G = minimalGB(Gp);
262        //G = red.irreducibleSet(Gp);
263        //G = OrderedPolynomialList.<C> sortDegree(G);
264        //logger.info("G_reduced = " + G);
265        logger.info("end " + pairlist);
266        return G;
267    }
268
269
270    /**
271     * S-Polynomial.
272     * @param A monic polynomial.
273     * @param B monic polynomial.
274     * @return spol(A,B) the S-polynomial of the A and B.
275     */
276    GenPolynomial<C> SPolynomial(SigPoly<C> A, SigPoly<C> B) {
277        return sred.SPolynomial(A, B);
278    }
279
280
281    /**
282     * S-Polynomial.
283     * @param P pair.
284     * @return spol(A,B) the S-polynomial of the pair (A,B).
285     */
286    GenPolynomial<C> SPolynomial(SigPair<C> P) {
287        return sred.SPolynomial(P.pi, P.pj);
288    }
289
290
291    /**
292     * S-Polynomial polynomial factors.
293     * @param A monic polynomial.
294     * @param B monic polynomial.
295     * @return polynomials [e,f] such that spol(A,B) = e*a - f*B.
296     */
297    GenPolynomial<C>[] SPolynomialFactors(SigPoly<C> A, SigPoly<C> B) {
298        return sred.SPolynomialFactors(A, B);
299    }
300
301
302    /**
303     * Pair with signature.
304     * @param A polynomial with signature.
305     * @param B polynomial with signature.
306     * @param G polynomial ith signature list.
307     * @return signature pair according to algorithm.
308     */
309    SigPair<C> newPair(SigPoly<C> A, SigPoly<C> B, List<SigPoly<C>> G) {
310        ExpVector e = A.poly.leadingExpVector().lcm(B.poly.leadingExpVector())
311                        .subtract(A.poly.leadingExpVector());
312        return new SigPair<C>(e, A, B, G);
313    }
314
315
316    /**
317     * Pair with signature.
318     * @param s signature for pair.
319     * @param A polynomial with signature.
320     * @param B polynomial with signature.
321     * @param G polynomial ith signature list.
322     * @return signature pair according to algorithm.
323     */
324    SigPair<C> newPair(GenPolynomial<C> s, SigPoly<C> A, SigPoly<C> B, List<SigPoly<C>> G) {
325        return new SigPair<C>(s, A, B, G);
326    }
327
328
329    /**
330     * Top normalform.
331     * @param A polynomial.
332     * @param F polynomial list.
333     * @param G polynomial list.
334     * @return nf(A) with respect to F and G.
335     */
336    SigPoly<C> sigNormalform(List<GenPolynomial<C>> F, List<SigPoly<C>> G, SigPoly<C> A) {
337        return sred.sigNormalform(F, G, A);
338    }
339
340
341    /**
342     * Prune total pair list P.
343     * @param P pair list.
344     * @param syz list of exponent vectors representing syzygies.
345     * @return updated pair list.
346     */
347    List<SigPair<C>> pruneP(List<SigPair<C>> P, List<ExpVector> syz) {
348        if (debug) {
349            logger.debug("unused " + syz);
350        }
351        return P;
352    }
353
354
355    /**
356     * Prune pair list of degree d.
357     * @param S pair list.
358     * @param syz list of exponent vectors representing syzygies.
359     * @param done list of treated polynomials.
360     * @param G polynomial with signature list.
361     * @return updated pair list.
362     */
363    List<SigPair<C>> pruneS(List<SigPair<C>> S, List<ExpVector> syz, List<SigPoly<C>> done,
364                    List<SigPoly<C>> G) {
365        if (debug) {
366            logger.debug("unused " + syz + " " + done + " " + G);
367        }
368        return S;
369    }
370
371
372    /**
373     * Initializes syzygy list.
374     * @param F polynomial list.
375     * @param G polynomial with signature list.
376     * @return list of exponent vectors representing syzygies.
377     */
378    List<ExpVector> initializeSyz(List<GenPolynomial<C>> F, List<SigPoly<C>> G) {
379        if (debug) {
380            logger.debug("unused " + G + " " + F);
381        }
382        List<ExpVector> P = new ArrayList<ExpVector>();
383        return P;
384    }
385
386
387    /**
388     * Update syzygy list.
389     * @param syz list of exponent vectors representing syzygies.
390     * @param r polynomial. <b>Note:</b> szy is modified to represent updated
391     *            list of exponent vectors.
392     */
393    void updateSyz(List<ExpVector> syz, SigPoly<C> r) {
394        if (debug) {
395            logger.debug("unused " + syz + " " + r);
396        }
397        return;
398    }
399
400}