001/*
002 * $Id: ComprehensiveGroebnerBaseSeq.java 4125 2012-08-19 19:05:22Z kredel $
003 */
004
005package edu.jas.application;
006
007
008import java.util.ArrayList;
009import java.util.Collections;
010import java.util.List;
011
012import org.apache.log4j.Logger;
013
014import edu.jas.gb.GroebnerBase;
015import edu.jas.gbufd.GBFactory;
016import edu.jas.poly.ExpVector;
017import edu.jas.poly.GenPolynomial;
018import edu.jas.poly.GenPolynomialRing;
019import edu.jas.structure.GcdRingElem;
020import edu.jas.structure.RingFactory;
021import edu.jas.ufd.SquarefreeAbstract;
022import edu.jas.ufd.SquarefreeFactory;
023
024
025/**
026 * Comprehensive Groebner Base sequential algorithm. Implements faithful
027 * comprehensive Groebner bases via Groebner systems and CGB test.
028 * @param <C> coefficient type
029 * @author Heinz Kredel
030 */
031
032public class ComprehensiveGroebnerBaseSeq<C extends GcdRingElem<C>>
033/* extends GroebnerBaseAbstract<GenPolynomial<C>> */{
034
035
036    private static final Logger logger = Logger.getLogger(ComprehensiveGroebnerBaseSeq.class);
037
038
039    private static final boolean debug = logger.isDebugEnabled();
040
041
042    /**
043     * Squarefree for coefficient content and primitive parts.
044     */
045    protected final SquarefreeAbstract<C> engine;
046
047
048    /*
049     * Flag if gcd engine should be used.
050     */
051    //private final boolean notFaithfull = false;
052
053
054    /**
055     * Comprehensive reduction engine.
056     */
057    protected final CReductionSeq<C> cred;
058
059
060    /**
061     * Polynomial coefficient ring factory.
062     */
063    protected final RingFactory<C> cofac;
064
065
066    /**
067     * Constructor.
068     * @param rf base coefficient ring factory.
069     */
070    public ComprehensiveGroebnerBaseSeq(RingFactory<C> rf) {
071        this(new CReductionSeq<C>(rf), rf);
072    }
073
074
075    /**
076     * Constructor.
077     * @param red C-pseudo-Reduction engine
078     * @param rf base coefficient ring factory.
079     */
080    @SuppressWarnings("unchecked")
081    public ComprehensiveGroebnerBaseSeq(CReductionSeq<C> red, RingFactory<C> rf) {
082        // super(null); // red not possible since type of type
083        cred = red;
084        cofac = rf;
085        // selection for C but used for R:
086        //engine e = GCDFactory.<C> getImplementation(cofac);
087        engine = SquarefreeFactory.<C> getImplementation(rf);
088    }
089
090
091    /**
092     * Comprehensive-Groebner base test.
093     * @param F polynomial list.
094     * @return true, if F is a Comprehensive-Groebner base, else false.
095     */
096    // @Override
097    public boolean isGB(List<GenPolynomial<GenPolynomial<C>>> F) {
098        return isGB(0, F);
099    }
100
101
102    /**
103     * Comprehensive-Groebner base test.
104     * @param modv module variable number.
105     * @param F polynomial list.
106     * @return true, if F is a Comprehensive-Groebner base, else false.
107     */
108    // @Override
109    public boolean isGB(int modv, List<GenPolynomial<GenPolynomial<C>>> F) {
110        // return isGBcol( modv, F );
111        return isGBsubst(modv, F);
112    }
113
114
115    /**
116     * Comprehensive-Groebner base test using colored systems.
117     * @param F polynomial list.
118     * @return true, if F is a Comprehensive-Groebner base, else false.
119     */
120    // @Override
121    public boolean isGBcol(List<GenPolynomial<GenPolynomial<C>>> F) {
122        return isGBcol(0, F);
123    }
124
125
126    /**
127     * Comprehensive-Groebner base test using colored systems.
128     * @param modv module variable number.
129     * @param F polynomial list.
130     * @return true, if F is a Comprehensive-Groebner base, else false.
131     */
132    // @Override
133    public boolean isGBcol(int modv, List<GenPolynomial<GenPolynomial<C>>> F) {
134        if (F == null || F.size() == 0) {
135            return true;
136        }
137        List<ColoredSystem<C>> CS = cred.determine(F);
138        return isGBsys(modv, CS);
139    }
140
141
142    /**
143     * Comprehensive-Groebner system test.
144     * @param CS list of colored systems.
145     * @return true, if CS is a Comprehensive-Groebner system, else false.
146     */
147    // @Override
148    public boolean isGBsys(List<ColoredSystem<C>> CS) {
149        return isGBsys(0, CS);
150    }
151
152
153    /**
154     * Comprehensive-Groebner system test.
155     * @param modv module variable number.
156     * @param CS list of colored systems.
157     * @return true, if CS is a Comprehensive-Groebner system, else false.
158     */
159    // @Override
160    public boolean isGBsys(int modv, List<ColoredSystem<C>> CS) {
161        if (CS == null || CS.size() == 0) {
162            return true;
163        }
164        ColorPolynomial<C> p, q, h, hp;
165        for (ColoredSystem<C> cs : CS) {
166            if (debug) {
167                if (!cs.isDetermined()) {
168                    System.out.println("not determined, cs = " + cs);
169                    return false;
170                }
171                if (!cs.checkInvariant()) {
172                    System.out.println("not invariant, cs = " + cs);
173                    return false;
174                }
175            }
176            Condition<C> cond = cs.condition;
177            List<ColorPolynomial<C>> S = cs.list;
178            int k = S.size();
179            for (int j = 0; j < k; j++) {
180                p = S.get(j);
181                for (int l = j + 1; l < k; l++) {
182                    q = S.get(l);
183                    h = cred.SPolynomial(p, q);
184                    // System.out.println("spol(a,b) = " + h);
185                    h = cred.normalform(cond, S, h);
186                    // System.out.println("NF(spol(a,b)) = " + h);
187                    if (debug) {
188                        if (!cred.isNormalform(S, h)) {
189                            System.out.println("not normalform, h = " + h);
190                            System.out.println("cs = " + cs);
191                            return false;
192                        }
193                    }
194                    if (!h.isZERO()) {
195                        hp = cond.reDetermine(h);
196                        if (!hp.isZERO()) {
197                            System.out.println("p = " + p);
198                            System.out.println("q = " + q);
199                            System.out.println("not zero:   NF(spol(p,q))  = " + h);
200                            System.out.println("redetermine(NF(spol(p,q))) = " + hp);
201                            System.out.println("cs = " + cs);
202                            return false;
203                        }
204                    }
205                }
206            }
207        }
208        return true;
209    }
210
211
212    /**
213     * Comprehensive-Groebner base test using substitution.
214     * @param F polynomial list.
215     * @return true, if F is a Comprehensive-Groebner base, else false.
216     */
217    // @Override
218    public boolean isGBsubst(List<GenPolynomial<GenPolynomial<C>>> F) {
219        return isGBsubst(0, F);
220    }
221
222
223    /**
224     * Comprehensive-Groebner base test using substitution.
225     * @param modv module variable number.
226     * @param F polynomial list.
227     * @return true, if F is a Comprehensive-Groebner base, else false.
228     */
229    // @Override
230    public boolean isGBsubst(int modv, List<GenPolynomial<GenPolynomial<C>>> F) {
231        if (F == null || F.size() == 0) {
232            return true;
233        }
234        GenPolynomial<GenPolynomial<C>> f = F.get(0); // assert non Zero
235        GenPolynomialRing<GenPolynomial<C>> cf = f.ring;
236
237        List<ColoredSystem<C>> CS = cred.determine(F);
238        if (logger.isDebugEnabled()) {
239            logger.info("determined polynomials =\n" + CS);
240        }
241        // substitute zero conditions into parameter coefficients and test
242        for (ColoredSystem<C> cs : CS) {
243            Ideal<C> id = cs.condition.zero;
244            ResidueRing<C> r = new ResidueRing<C>(id);
245            GenPolynomialRing<Residue<C>> rf = new GenPolynomialRing<Residue<C>>(r, cf);
246            List<GenPolynomial<Residue<C>>> list = PolyUtilApp.<C> toResidue(rf, F);
247            //GroebnerBase<Residue<C>> bb = new GroebnerBasePseudoSeq<Residue<C>>(r);
248            GroebnerBase<Residue<C>> bb = GBFactory.getImplementation(r);
249            boolean t = bb.isGB(list);
250            if (!t) {
251                System.out.println("test condition = " + cs.condition);
252                System.out.println("no GB for residue coefficients = " + list);
253                return false;
254            }
255        }
256
257        // substitute random ideal into parameter coefficients and test
258        GenPolynomialRing<C> ccf = (GenPolynomialRing<C>) cf.coFac;
259        int nv = ccf.nvar - 2;
260        if (nv < 1) {
261            nv = 1;
262        }
263        List<GenPolynomial<C>> il = new ArrayList<GenPolynomial<C>>();
264        int i = 0;
265        int j = 1;
266        while (i < nv) {
267            j++;
268            GenPolynomial<C> p = ccf.random(j + 1);
269            // System.out.println("p = " + p);
270            if (p.isConstant()) {
271                continue;
272            }
273            if (p.isZERO()) {
274                continue;
275            }
276            p = engine.squarefreePart(p);
277            il.add(p);
278            i++;
279        }
280        logger.info("random ideal = " + il);
281        Ideal<C> id = new Ideal<C>(ccf, il);
282        ResidueRing<C> r = new ResidueRing<C>(id);
283        GenPolynomialRing<Residue<C>> rf = new GenPolynomialRing<Residue<C>>(r, cf);
284        List<GenPolynomial<Residue<C>>> list = PolyUtilApp.<C> toResidue(rf, F);
285        logger.info("random residue = " + r.ideal.getList());
286        //GroebnerBase<Residue<C>> bb = new GroebnerBasePseudoSeq<Residue<C>>(r);
287        GroebnerBase<Residue<C>> bb = GBFactory.getImplementation(r);
288        boolean t = bb.isGB(list);
289        if (!t) {
290            System.out.println("no GB for residue coefficients = " + list);
291            return false;
292        }
293        return true;
294    }
295
296
297    /**
298     * Comprehensive-Groebner system test.
299     * @param F Groebner system.
300     * @return true, if F is a Comprehensive-Groebner system, else false.
301     */
302    // @Override
303    public boolean isGBsys(GroebnerSystem<C> F) {
304        return isGBsys(0, F.list);
305    }
306
307
308    /**
309     * Comprehensive-Groebner base test.
310     * @param F Groebner system.
311     * @return true, if F is a Comprehensive-Groebner base, else false.
312     */
313    // @Override
314    public boolean isCGB(GroebnerSystem<C> F) {
315        return isGB(F.getCGB());
316    }
317
318
319    /**
320     * Comprehensive-Groebner system and base test.
321     * @param F Groebner system.
322     * @return true, if F is a Comprehensive-Groebner system and base, else
323     *         false.
324     */
325    // @Override
326    public boolean isGB(GroebnerSystem<C> F) {
327        return isGBsys(0, F.list) && isGB(F.getCGB());
328    }
329
330
331    /**
332     * Comprehensive Groebner base system using pairlist class.
333     * @param F polynomial list.
334     * @return GBsys(F) a Comprehensive Groebner system of F.
335     */
336    // @Override
337    // @SuppressWarnings("unchecked")
338    public GroebnerSystem<C> GBsys(List<GenPolynomial<GenPolynomial<C>>> F) {
339        if (F == null) {
340            return null;
341        }
342        List<ColoredSystem<C>> CSp = new ArrayList<ColoredSystem<C>>();
343        if (F.size() == 0) {
344            return new GroebnerSystem<C>(CSp);
345        }
346        // extract coefficient factory
347        GenPolynomial<GenPolynomial<C>> f = F.get(0);
348        GenPolynomialRing<GenPolynomial<C>> fac = f.ring;
349        // determine polynomials
350        List<ColoredSystem<C>> CS = cred.determine(F);
351        // System.out.println("CS = " + CS);
352        // CS.remove(0); // empty colored system
353        if (logger.isInfoEnabled()) {
354            logger.info("determined polynomials =\n" + CS);
355        }
356
357        // setup pair lists
358        List<ColoredSystem<C>> CSs = new ArrayList<ColoredSystem<C>>();
359        ColoredSystem<C> css;
360        for (ColoredSystem<C> cs : CS) {
361            OrderedCPairlist<C> pairlist = new OrderedCPairlist<C>(fac);
362            for (ColorPolynomial<C> p : cs.list) {
363                // System.out.println("p = " + p);
364                pairlist.put(p);
365            }
366            css = new ColoredSystem<C>(cs.condition, cs.list, pairlist);
367            CSs.add(css);
368        }
369
370        // main loop
371        List<ColoredSystem<C>> CSb = new ArrayList<ColoredSystem<C>>();
372        List<ColoredSystem<C>> ncs;
373        List<ColoredSystem<C>> CSh; //, CSbh;
374        ColoredSystem<C> cs;
375        List<ColorPolynomial<C>> G;
376        OrderedCPairlist<C> pairlist;
377        Condition<C> cond;
378        int si = 0;
379        while (CSs.size() > 0) {
380            cs = CSs.get(0); // remove(0);
381            si++;
382            logger.info("poped GBsys number    " + si + " with condition = " + cs.condition);
383            logger.info("poped GBsys (remaining " + (CSs.size() - 1) + ") with pairlist  = " + cs.pairlist);
384            if (!cs.isDetermined()) {
385                cs = cs.reDetermine();
386            }
387            pairlist = cs.pairlist;
388            G = cs.list;
389            cond = cs.condition;
390            // logger.info( pairlist.toString() );
391
392            CPair<C> pair;
393            ColorPolynomial<C> pi;
394            ColorPolynomial<C> pj;
395            ColorPolynomial<C> S;
396            // GenPolynomial<GenPolynomial<C>> H;
397            ColorPolynomial<C> H;
398            while (pairlist.hasNext()) {
399                pair = pairlist.removeNext();
400                if (pair == null)
401                    continue;
402
403                pi = pair.pi;
404                pj = pair.pj;
405                if (debug) {
406                    logger.info("pi    = " + pi);
407                    logger.info("pj    = " + pj);
408                }
409
410                S = cred.SPolynomial(pi, pj);
411                if (S.isZERO()) {
412                    pair.setZero();
413                    continue;
414                }
415                if (debug) {
416                    // logger.info("ht(S) = " + S.leadingExpVector() );
417                    logger.info("S = " + S);
418                }
419
420                H = cred.normalform(cond, G, S);
421                if (H.isZERO()) {
422                    pair.setZero();
423                    continue;
424                }
425                if (debug) {
426                    logger.info("ht(H) = " + H.leadingExpVector());
427                }
428
429                H = H.abs();
430                if (debug) {
431                    logger.debug("H = " + H);
432                }
433                logger.info("H = " + H);
434                if (!H.isZERO()) {
435                    //CSh = new ArrayList<ColoredSystem<C>>();
436                    ncs = determineAddPairs(cs, H);
437                    if (ncs.size() == 0) {
438                        continue;
439                    }
440                    cs = ncs.remove(0); // remove other?
441                    pairlist = cs.pairlist;
442                    G = cs.list;
443                    cond = cs.condition;
444                    logger.info("replaced main branch = " + cond);
445                    logger.info("#new systems       = " + ncs.size());
446                    int yi = CSs.size();
447                    for (ColoredSystem<C> x : ncs) {
448                        if (!x.isDetermined()) {
449                            x = x.reDetermine();
450                        }
451                        CSs = x.addToList(CSs);
452                    }
453                    logger.info("#new systems added = " + (CSs.size() - yi));
454                }
455            }
456            // all s-pols reduce to zero in this branch
457            if (!cs.isDetermined()) {
458                cs = cs.reDetermine();
459            }
460            CSb.add(cs);
461            CSs.remove(0);
462            logger.info("done with = " + cs.condition);
463        }
464        // all branches done
465        CSh = new ArrayList<ColoredSystem<C>>();
466        for (ColoredSystem<C> x : CSb) {
467            // System.out.println("G = " + x.list );
468            if (!x.isDetermined()) {
469                x = x.reDetermine();
470            }
471            cs = minimalGB(x);
472            // System.out.println("min(G) = " + cs.list );
473            if (!cs.isDetermined()) {
474                cs = cs.reDetermine();
475            }
476            // cs = new ColoredSystem<C>( x.condition, G, x.pairlist );
477            CSh.add(cs);
478            logger.info("#sequential done = " + x.condition);
479            logger.info(x.pairlist.toString());
480        }
481        CSb = new ArrayList<ColoredSystem<C>>(CSh);
482        return new GroebnerSystem<C>(CSb);
483    }
484
485
486    /**
487     * Determine polynomial relative to a condition of a colored system and add
488     * pairs.
489     * @param cs a colored system.
490     * @param A color polynomial.
491     * @return list of colored systems, the conditions extending the condition
492     *         of cs.
493     */
494    public List<ColoredSystem<C>> determineAddPairs(ColoredSystem<C> cs, ColorPolynomial<C> A) {
495        List<ColoredSystem<C>> NCS = new ArrayList<ColoredSystem<C>>();
496        if (A == null || A.isZERO()) {
497            // NCS.add( cs );
498            return NCS;
499        }
500        List<ColorPolynomial<C>> S = cs.list;
501        Condition<C> cond = cs.condition; // .clone(); done in Condition
502        // itself
503        OrderedCPairlist<C> pl = cs.pairlist;
504
505        List<ColorPolynomial<C>> Sp;
506        ColorPolynomial<C> nz;
507        ColoredSystem<C> NS;
508        // if ( A.isDetermined() ) { ... } // dont use this
509        // System.out.println("to determine = " + A);
510        GenPolynomial<GenPolynomial<C>> Ap = A.getPolynomial();
511        List<Condition<C>> cd = cred.caseDistinction(cond, Ap);
512        logger.info("# cases = " + cd.size());
513        for (Condition<C> cnd : cd) {
514            //nz = cnd.determine(Ap);
515            nz = cnd.reDetermine(A);
516            if (nz.isZERO()) {
517                logger.info("zero determined nz = " + nz);
518                Sp = new ArrayList<ColorPolynomial<C>>(S);
519                OrderedCPairlist<C> PL = pl.copy();
520                NS = new ColoredSystem<C>(cnd, Sp, PL);
521                try {
522                    if (!NS.isDetermined()) {
523                        NS = NS.reDetermine();
524                    }
525                } catch (RuntimeException e) {
526                    System.out.println("Contradiction in NS_0 = " + NS);
527                    //e.printStackTrace();
528                    continue;
529                }
530                NCS = NS.addToList(NCS);
531                continue;
532            }
533            if (S.contains(nz)) {
534                System.out.println("*** S.contains(nz) ***");
535                continue;
536            }
537            logger.info("new determined nz = " + nz);
538            Sp = new ArrayList<ColorPolynomial<C>>(S);
539            Sp.add(nz);
540            OrderedCPairlist<C> PL = pl.copy();
541            PL.put(nz);
542            NS = new ColoredSystem<C>(cnd, Sp, PL);
543            try {
544                if (!NS.isDetermined()) {
545                    NS = NS.reDetermine();
546                }
547            } catch (RuntimeException e) {
548                System.out.println("Contradiction in NS = " + NS);
549                //e.printStackTrace();
550                continue;
551            }
552            NCS = NS.addToList(NCS);
553        }
554        // System.out.println("new determination = " + NCS);
555        return NCS;
556    }
557
558
559    /**
560     * Comprehensive Groebner base via Groebner system.
561     * @param F polynomial list.
562     * @return GB(F) a Comprehensive Groebner base of F.
563     */
564    // @Override
565    // @SuppressWarnings("unchecked")
566    public List<GenPolynomial<GenPolynomial<C>>> GB(List<GenPolynomial<GenPolynomial<C>>> F) {
567        if (F == null) {
568            return F;
569        }
570        // compute Groebner system
571        GroebnerSystem<C> gs = GBsys(F);
572        // System.out.println("\n\nGBsys = " + gs);
573        return gs.getCGB();
574    }
575
576
577    /**
578     * Minimal ordered Groebner basis.
579     * @param cs colored system.
580     * @return a reduced Groebner base of Gp.
581     */
582    // @Override
583    public ColoredSystem<C> minimalGB(ColoredSystem<C> cs) {
584        // List<ColorPolynomial<C>> Gp ) {
585        if (cs == null || cs.list == null || cs.list.size() <= 1) {
586            return cs;
587        }
588        // remove zero polynomials
589        List<ColorPolynomial<C>> G = new ArrayList<ColorPolynomial<C>>(cs.list.size());
590        for (ColorPolynomial<C> a : cs.list) {
591            if (a != null && !a.isZERO()) { // always true in GB()
592                // already positive a = a.abs();
593                G.add(a);
594            }
595        }
596        if (G.size() <= 1) {
597            return new ColoredSystem<C>(cs.condition, G, cs.pairlist);
598        }
599        // System.out.println("G check " + G);
600        // remove top reducible polynomials
601        Condition<C> cond = cs.condition;
602        ColorPolynomial<C> a, b;
603        List<ColorPolynomial<C>> F;
604        F = new ArrayList<ColorPolynomial<C>>(G.size());
605        while (G.size() > 0) {
606            a = G.remove(0);
607            b = a;
608            // System.out.println("check " + b);
609            //if (false) {
610            //    if (a.red.leadingBaseCoefficient().isConstant()) { // dont drop
611            //        // these
612            //        F.add(a);
613            //        continue;
614            //    }
615            //}
616            if (cred.isTopReducible(G, a) || cred.isTopReducible(F, a)) {
617                // drop polynomial
618                if (debug) {
619                    // System.out.println("trying to drop " + a);
620                    List<ColorPolynomial<C>> ff;
621                    ff = new ArrayList<ColorPolynomial<C>>(G);
622                    ff.addAll(F);
623                    a = cred.normalform(cond, ff, a);
624                    try {
625                        a = cond.reDetermine(a);
626                    } catch (RuntimeException ignored) {
627                    }
628                    if (!a.isZERO()) {
629                        logger.error("nf(a) != 0 " + b + ", " + a);
630                        F.add(b);
631                    }
632                }
633            } else {
634                F.add(a);
635            }
636        }
637        G = F;
638        if (G.size() <= 1) {
639            return new ColoredSystem<C>(cs.condition, G, cs.pairlist);
640        }
641        Collections.reverse(G); // important for lex GB
642        // reduce remaining polynomials
643        int len = G.size();
644        int i = 0;
645        while (i < len) {
646            a = G.remove(0);
647            b = a;
648            ExpVector e = a.red.leadingExpVector();
649            // System.out.println("reducing " + a);
650            a = cred.normalform(cond, G, a); // unchanged by top reduction
651            // System.out.println("reduced " + a);
652            try {
653                a = cond.reDetermine(a);
654            } catch (RuntimeException ignored) {
655            }
656            ExpVector f = a.red.leadingExpVector();
657            // a = ##engine.basePrimitivePart(a); //a.monic(); was not required
658            // a = a.abs();
659            // a = red.normalform( F, a );
660            if (e.equals(f)) {
661                G.add(a); // adds as last
662            } else { // should not happen
663                if (debug) {
664                    logger.error("nf(a) not determined " + b + ", " + a);
665                }
666                G.add(b); // adds as last
667            }
668            i++;
669        }
670        return new ColoredSystem<C>(cs.condition, G, cs.pairlist);
671    }
672
673}