Source code

001/*
002 * $Id: SolvableSyzygySeq.java 5870 2018-07-20 15:56:23Z kredel $
003 */
004
005package edu.jas.gbufd;
006
007
008import java.util.ArrayList;
009import java.util.Arrays;
010import java.util.Iterator;
011import java.util.List;
012
013import org.apache.logging.log4j.Logger;
014import org.apache.logging.log4j.LogManager; 
015
016import edu.jas.gb.SolvableExtendedGB;
017import edu.jas.gb.SolvableGroebnerBase;
018import edu.jas.poly.GenPolynomial;
019import edu.jas.poly.GenSolvablePolynomial;
020import edu.jas.poly.GenSolvablePolynomialRing;
021import edu.jas.poly.ModuleList;
022import edu.jas.poly.PolynomialList;
023import edu.jas.structure.GcdRingElem;
024import edu.jas.structure.RingFactory;
025
026
027/**
028 * Syzygy sequential class for solvable polynomials. Implements Syzygy
029 * computations and tests with Groebner bases.
030 * @param <C> coefficient type
031 * @author Heinz Kredel
032 */
033
034public class SolvableSyzygySeq<C extends GcdRingElem<C>> extends SolvableSyzygyAbstract<C> {
035
036
037    private static final Logger logger = LogManager.getLogger(SolvableSyzygySeq.class);
038
039
040    private static final boolean debug = logger.isDebugEnabled();
041
042
043    private static boolean assertEnabled = false;
044
045
046    static {
047        assert assertEnabled = true; // official hack to check assertions enabled
048    }
049
050
051    /**
052     * Groebner basis engine.
053     */
054    protected SolvableGroebnerBase<C> sbb;
055
056
057    /*
058     * Module Groebner basis engine.
059    //protected ModSolvableGroebnerBase<C> msbb;
060     */
061
062
063    /**
064     * Constructor.
065     * @param cf coefficient ring.
066     */
067    public SolvableSyzygySeq(RingFactory<C> cf) {
068        sbb = SGBFactory.getImplementation(cf);
069        //msbb = new ModSolvableGroebnerBaseSeq<C>(cf);
070    }
071
072
073    /**
074     * Resolution of a module. Only with direct GBs.
075     * @param M a module list of a Groebner basis.
076     * @return a resolution of M.
077     */
078    public List<SolvResPart<C>> resolution(ModuleList<C> M) {
079        List<SolvResPart<C>> R = new ArrayList<SolvResPart<C>>();
080        ModuleList<C> MM = M;
081        ModuleList<C> GM;
082        ModuleList<C> Z;
083        while (true) {
084            GM = sbb.leftGB(MM);
085            Z = leftZeroRelations(GM);
086            R.add(new SolvResPart<C>(MM, GM, Z));
087            if (Z == null || Z.list == null || Z.list.size() == 0) {
088                break;
089            }
090            MM = Z;
091        }
092        return R;
093    }
094
095
096    /**
097     * Resolution of a polynomial list. Only with direct GBs.
098     * @param F a polynomial list of a Groebner basis.
099     * @return a resolution of F.
100     */
101    @SuppressWarnings("unchecked")
102    public List/*<SolvResPart<C>|SolvResPolPart<C>>*/resolution(PolynomialList<C> F) {
103        List<List<GenSolvablePolynomial<C>>> Z;
104        ModuleList<C> Zm;
105        List<GenSolvablePolynomial<C>> G;
106        PolynomialList<C> Gl;
107
108        G = sbb.leftGB(F.castToSolvableList());
109        Z = leftZeroRelations(G);
110        Gl = new PolynomialList<C>((GenSolvablePolynomialRing<C>) F.ring, G);
111        Zm = new ModuleList<C>((GenSolvablePolynomialRing<C>) F.ring, Z);
112
113        List R = resolution(Zm);
114        R.add(0, new SolvResPolPart<C>(F, Gl, Zm));
115        return R;
116    }
117
118
119    /**
120     * Resolution of a module.
121     * @param M a module list of an arbitrary basis.
122     * @return a resolution of M.
123     */
124    public List<SolvResPart<C>> resolutionArbitrary(ModuleList<C> M) {
125        List<SolvResPart<C>> R = new ArrayList<SolvResPart<C>>();
126        ModuleList<C> MM = M;
127        ModuleList<C> GM = null;
128        ModuleList<C> Z;
129        while (true) {
130            //GM = sbb.leftGB(MM);
131            Z = leftZeroRelationsArbitrary(MM);
132            R.add(new SolvResPart<C>(MM, GM, Z));
133            if (Z == null || Z.list == null || Z.list.size() == 0) {
134                break;
135            }
136            MM = Z;
137        }
138        return R;
139    }
140
141
142    /**
143     * Resolution of a polynomial list.
144     * @param F a polynomial list of an arbitrary basis.
145     * @return a resolution of F.
146     */
147    @SuppressWarnings("unchecked")
148    public List/*<SolvResPart<C>|SolvResPolPart<C>>*/resolutionArbitrary(PolynomialList<C> F) {
149        List<List<GenSolvablePolynomial<C>>> Z;
150        ModuleList<C> Zm;
151        PolynomialList<C> Gl = null;
152        //List<GenSolvablePolynomial<C>> G = sbb.leftGB( F.castToSolvableList() );
153        //Gl = new PolynomialList<C>((GenSolvablePolynomialRing<C>)F.ring, G);
154        Z = leftZeroRelationsArbitrary(F.castToSolvableList());
155        Zm = new ModuleList<C>((GenSolvablePolynomialRing<C>) F.ring, Z);
156
157        List R = resolutionArbitrary(Zm);
158        R.add(0, new SolvResPolPart<C>(F, Gl, Zm));
159        return R;
160    }
161
162
163    /**
164     * Left syzygy module from arbitrary base.
165     * @param modv number of module variables.
166     * @param F a solvable polynomial list.
167     * @return syz(F), a basis for the module of left syzygies for F.
168     */
169    @SuppressWarnings("unchecked")
170    public List<List<GenSolvablePolynomial<C>>> leftZeroRelationsArbitrary(int modv,
171                    List<GenSolvablePolynomial<C>> F) {
172        if (F == null) {
173            return null; //leftZeroRelations( modv, F );
174        }
175        if (F.size() <= 1) {
176            return leftZeroRelations(modv, F);
177        }
178        final int lenf = F.size();
179        SolvableExtendedGB<C> exgb = sbb.extLeftGB(modv, F);
180        if (debug) {
181            logger.info("exgb = " + exgb);
182        }
183        if (assertEnabled) {
184            logger.info("check1 exgb start");
185            if (!sbb.isLeftReductionMatrix(exgb)) {
186                logger.error("is reduction matrix ? false");
187            }
188            logger.info("check1 exgb end");
189        }
190        List<GenSolvablePolynomial<C>> G = exgb.G;
191        List<List<GenSolvablePolynomial<C>>> G2F = exgb.G2F;
192        List<List<GenSolvablePolynomial<C>>> F2G = exgb.F2G;
193
194        List<List<GenSolvablePolynomial<C>>> sg = leftZeroRelations(modv, G);
195        GenSolvablePolynomialRing<C> ring = G.get(0).ring;
196        ModuleList<C> S = new ModuleList<C>(ring, sg);
197        if (debug) {
198            logger.info("syz = " + S);
199        }
200        if (assertEnabled) {
201            logger.info("check2 left syz start");
202            if (!isLeftZeroRelation(sg, G)) {
203                logger.error("is syzygy ? false");
204            }
205            logger.info("check2 left syz end");
206        }
207
208        List<List<GenSolvablePolynomial<C>>> sf;
209        sf = new ArrayList<List<GenSolvablePolynomial<C>>>(sg.size());
210        //List<GenPolynomial<C>> row;
211
212        for (List<GenSolvablePolynomial<C>> r : sg) {
213            Iterator<GenSolvablePolynomial<C>> it = r.iterator();
214            Iterator<List<GenSolvablePolynomial<C>>> jt = G2F.iterator();
215
216            List<GenSolvablePolynomial<C>> rf;
217            rf = new ArrayList<GenSolvablePolynomial<C>>(lenf);
218            for (int m = 0; m < lenf; m++) {
219                rf.add(ring.getZERO());
220            }
221            while (it.hasNext() && jt.hasNext()) {
222                GenSolvablePolynomial<C> si = it.next();
223                List<GenSolvablePolynomial<C>> ai = jt.next();
224                //System.out.println("si = " + si);
225                //System.out.println("ai = " + ai);
226                if (si == null || ai == null) {
227                    continue;
228                }
229                // pi has wrong type:
230                List<GenPolynomial<C>> pi = blas.scalarProduct(si, PolynomialList.<C> castToList(ai));
231                //System.out.println("pi = " + pi);
232                rf = PolynomialList.<C> castToSolvableList(blas.vectorAdd(PolynomialList.<C> castToList(rf),
233                                pi));
234            }
235            if (it.hasNext() || jt.hasNext()) {
236                logger.error("leftZeroRelationsArbitrary wrong sizes");
237            }
238            //System.out.println("\nrf = " + rf + "\n");
239            sf.add(rf);
240        }
241        if (assertEnabled) {
242            logger.info("check3 left syz start");
243            if (!isLeftZeroRelation(sf, F)) {
244                logger.error("is partial syz sf ? false");
245            }
246            logger.info("check3 left syz end");
247        }
248
249        List<List<GenSolvablePolynomial<C>>> M;
250        M = new ArrayList<List<GenSolvablePolynomial<C>>>(lenf);
251        for (List<GenSolvablePolynomial<C>> r : F2G) {
252            Iterator<GenSolvablePolynomial<C>> it = r.iterator();
253            Iterator<List<GenSolvablePolynomial<C>>> jt = G2F.iterator();
254
255            List<GenSolvablePolynomial<C>> rf;
256            rf = new ArrayList<GenSolvablePolynomial<C>>(lenf);
257            for (int m = 0; m < lenf; m++) {
258                rf.add(ring.getZERO());
259            }
260            while (it.hasNext() && jt.hasNext()) {
261                GenSolvablePolynomial<C> si = it.next();
262                List<GenSolvablePolynomial<C>> ai = jt.next();
263                //System.out.println("si = " + si);
264                //System.out.println("ai = " + ai);
265                if (si == null || ai == null) {
266                    continue;
267                }
268                //pi has wrong type, should be: List<GenSolvablePolynomial<C>>
269                List<GenPolynomial<C>> pi = blas.scalarProduct(si, PolynomialList.<C> castToList(ai));
270                //System.out.println("pi = " + pi);
271                rf = PolynomialList.<C> castToSolvableList(blas.vectorAdd(PolynomialList.<C> castToList(rf),
272                                pi));
273            }
274            if (it.hasNext() || jt.hasNext()) {
275                logger.error("zeroRelationsArbitrary wrong sizes");
276            }
277            //System.out.println("\nMg Mf = " + rf + "\n");
278            M.add(rf);
279        }
280        //ModuleList<C> ML = new ModuleList<C>( ring, M );
281        //System.out.println("syz ML = " + ML);
282        // debug only:
283        //List<GenSolvablePolynomial<C>> F2 = new ArrayList<GenSolvablePolynomial<C>>( F.size() );
284        /* not true in general
285        List<GenPolynomial<C>> Fp = PolynomialList.<C>castToList(F);
286        for ( List<GenSolvablePolynomial<C>> rr: M ) {
287            GenSolvablePolynomial<C> rrg = PolynomialList.<C>castToSolvableList(blas.scalarProduct(Fp,PolynomialList.<C>castToList(rr)));
288            F2.add( rrg );
289        }
290        PolynomialList<C> pF = new PolynomialList<C>( ring, F );
291        PolynomialList<C> pF2 = new PolynomialList<C>( ring, F2 );
292        if ( ! pF.equals( pF2 ) ) {
293           logger.error("is FAB = F ? false");
294           //System.out.println("pF  = " + pF.list.size());
295           //System.out.println("pF2 = " + pF2.list.size());
296        }
297        */
298        int sflen = sf.size();
299        List<List<GenSolvablePolynomial<C>>> M2;
300        M2 = new ArrayList<List<GenSolvablePolynomial<C>>>(lenf);
301        int i = 0;
302        for (List<GenSolvablePolynomial<C>> ri : M) {
303            List<GenSolvablePolynomial<C>> r2i;
304            r2i = new ArrayList<GenSolvablePolynomial<C>>(ri.size());
305            int j = 0;
306            for (GenSolvablePolynomial<C> rij : ri) {
307                GenSolvablePolynomial<C> p = null;
308                if (i == j) {
309                    p = (GenSolvablePolynomial<C>) ring.getONE().subtract(rij);
310                } else {
311                    if (rij != null) {
312                        p = (GenSolvablePolynomial<C>) rij.negate();
313                    }
314                }
315                r2i.add(p);
316                j++;
317            }
318            M2.add(r2i);
319            if (!blas.isZero(PolynomialList.<C> castToList(r2i))) {
320                sf.add(r2i);
321            }
322            i++;
323        }
324        ModuleList<C> M2L = new ModuleList<C>(ring, M2);
325        if (debug) {
326            logger.debug("syz M2L = " + M2L);
327        }
328
329        if (debug) {
330            ModuleList<C> SF = new ModuleList<C>(ring, sf);
331            logger.debug("syz sf = " + SF);
332            logger.debug("#syz " + sflen + ", " + sf.size());
333        }
334        if (assertEnabled) {
335            logger.info("check4 left syz start");
336            if (!isLeftZeroRelation(sf, F)) {
337                logger.error("is syz sf ? false");
338            }
339            logger.info("check4 left syz end");
340        }
341        return sf;
342    }
343
344
345    /**
346     * Left Ore condition. Generators for the left Ore condition of two solvable
347     * polynomials.
348     * @param a solvable polynomial
349     * @param b solvable polynomial
350     * @return [p,q] with p*a = q*b
351     */
352    @SuppressWarnings({ "cast", "unchecked" })
353    public GenSolvablePolynomial<C>[] leftOreCond(GenSolvablePolynomial<C> a, GenSolvablePolynomial<C> b) {
354        if (a == null || a.isZERO() || b == null || b.isZERO()) {
355            throw new IllegalArgumentException("a and b must be non zero");
356        }
357        GenSolvablePolynomialRing<C> pfac = a.ring;
358        GenSolvablePolynomial<C>[] oc = (GenSolvablePolynomial<C>[]) new GenSolvablePolynomial[2];
359        if (a.equals(b)) {
360            oc[0] = pfac.getONE();
361            oc[1] = pfac.getONE();
362            return oc;
363        }
364        if (a.isConstant()) {
365            if (pfac.coFac.isCommutative()) { // ??
366                oc[0] = b;
367                oc[1] = a;
368                return oc;
369            }
370            oc[1] = pfac.getONE();
371            C c = a.leadingBaseCoefficient().inverse();
372            oc[0] = b.multiply(c);
373            return oc;
374        }
375        if (b.isConstant()) {
376            if (pfac.coFac.isCommutative()) { // ??
377                oc[0] = b;
378                oc[1] = a;
379                return oc;
380            }
381            oc[0] = pfac.getONE();
382            C c = b.leadingBaseCoefficient().inverse();
383            oc[1] = a.multiply(c);
384            return oc;
385        }
386        logger.info("computing left Ore condition: " + a + ", " + b);
387        List<GenSolvablePolynomial<C>> F = new ArrayList<GenSolvablePolynomial<C>>(2);
388        F.add(a);
389        F.add(b);
390        List<List<GenSolvablePolynomial<C>>> Gz = leftZeroRelationsArbitrary(F);
391        /*
392        if (Gz.size() < 0) { // always false
393            //System.out.println("Gz = " + Gz);
394            ModuleList<C> M = new ModuleList<C>(pfac, Gz);
395            ModuleList<C> GM = sbb.leftGB(M);
396            //System.out.println("GM = " + GM);
397            Gz = GM.castToSolvableList();
398        }
399        */
400        List<GenSolvablePolynomial<C>> G1 = null;
401        GenSolvablePolynomial<C> g1 = null;
402        for (List<GenSolvablePolynomial<C>> Gi : Gz) {
403            //System.out.println("Gi = " + Gi);
404            if (Gi.get(0).isZERO()) {
405                continue;
406            }
407            if (G1 == null) {
408                G1 = Gi;
409            }
410            if (g1 == null) {
411                g1 = G1.get(0);
412            } else if (g1.compareTo(Gi.get(0)) > 0) { //g1.degree() > Gi.get(0).degree() 
413                G1 = Gi;
414                g1 = G1.get(0);
415            }
416        }
417        oc[0] = g1; //G1.get(0);
418        oc[1] = (GenSolvablePolynomial<C>) G1.get(1).negate();
419        //logger.info("Ore multiple: " + oc[0].multiply(a) + ", " + Arrays.toString(oc));
420        return oc;
421    }
422
423
424    /**
425     * Right Ore condition. Generators for the right Ore condition of two
426     * solvable polynomials.
427     * @param a solvable polynomial
428     * @param b solvable polynomial
429     * @return [p,q] with a*p = b*q
430     */
431    @SuppressWarnings({ "cast", "unchecked" })
432    public GenSolvablePolynomial<C>[] rightOreCond(GenSolvablePolynomial<C> a, GenSolvablePolynomial<C> b) {
433        if (a == null || a.isZERO() || b == null || b.isZERO()) {
434            throw new IllegalArgumentException("a and b must be non zero");
435        }
436        GenSolvablePolynomialRing<C> pfac = a.ring;
437        GenSolvablePolynomial<C>[] oc = (GenSolvablePolynomial<C>[]) new GenSolvablePolynomial[2];
438        if (a.equals(b)) {
439            oc[0] = pfac.getONE();
440            oc[1] = pfac.getONE();
441            return oc;
442        }
443        if (a.isConstant()) {
444            if (pfac.coFac.isCommutative()) { // ??
445                oc[0] = b;
446                oc[1] = a;
447                return oc;
448            }
449            oc[1] = pfac.getONE();
450            C c = a.leadingBaseCoefficient().inverse();
451            oc[0] = b.multiply(c);
452            return oc;
453        }
454        if (b.isConstant()) {
455            if (pfac.coFac.isCommutative()) { // ??
456                oc[0] = b;
457                oc[1] = a;
458                return oc;
459            }
460            oc[0] = pfac.getONE();
461            C c = b.leadingBaseCoefficient().inverse();
462            oc[1] = a.multiply(c);
463            return oc;
464        }
465        logger.info("computing right Ore condition: " + a + ", " + b);
466        List<GenSolvablePolynomial<C>> F = new ArrayList<GenSolvablePolynomial<C>>(2);
467        F.add(a);
468        F.add(b);
469        List<List<GenSolvablePolynomial<C>>> Gz = rightZeroRelationsArbitrary(F);
470        List<GenSolvablePolynomial<C>> G1 = null;
471        GenSolvablePolynomial<C> g1 = null;
472        for (List<GenSolvablePolynomial<C>> Gi : Gz) {
473            if (Gi.get(0).isZERO()) {
474                continue;
475            }
476            if (G1 == null) {
477                G1 = Gi;
478            }
479            if (g1 == null) {
480                g1 = G1.get(0);
481            } else if (g1.compareTo(Gi.get(0)) > 0) {
482                G1 = Gi;
483                g1 = G1.get(0);
484            }
485        }
486        oc[0] = G1.get(0);
487        oc[1] = (GenSolvablePolynomial<C>) G1.get(1).negate();
488        //logger.info("Ore multiple: " + a.multiply(oc[0]) + ", " + Arrays.toString(oc));
489        return oc;
490    }
491
492
493    /**
494     * Left simplifier. Method of Apel &amp; Lassner (1987).
495     * @param a solvable polynomial
496     * @param b solvable polynomial
497     * @return [p,q] with a/b = p/q and q is minimal and monic
498     */
499    @Override
500    @SuppressWarnings({ "unchecked" })
501    public GenSolvablePolynomial<C>[] leftSimplifier(GenSolvablePolynomial<C> a, GenSolvablePolynomial<C> b) {
502        if (a == null || a.isZERO() || b == null || b.isZERO()) {
503            throw new IllegalArgumentException("a and b must be non zero");
504        }
505        GenSolvablePolynomial<C>[] oc = null;
506        if (a.isConstant() || b.isConstant()) {
507            oc = new GenSolvablePolynomial[] { a, b };
508            return oc;
509        }
510        if (a.totalDegree() > 3 || b.totalDegree() > 3) { // how avoid too long running GBs ?
511            //if (a.totalDegree() + b.totalDegree() > 6) { 
512            // && a.length() < 10 && b.length() < 10
513            logger.info("skipping simplifier GB computation: degs = " + a.totalDegree() + ", " + b.totalDegree());
514            oc = new GenSolvablePolynomial[] { a, b };
515            return oc;
516        }
517        //GenSolvablePolynomialRing<C> pfac = a.ring;
518        oc = rightOreCond(a, b);
519        logger.info("oc = " + Arrays.toString(oc)); // + ", a = " + a + ", b = " + b);
520        List<GenSolvablePolynomial<C>> F = new ArrayList<GenSolvablePolynomial<C>>(oc.length);
521        // opposite order and undo negation
522        F.add((GenSolvablePolynomial<C>) oc[1].negate());
523        F.add(oc[0]);
524        //logger.info("F = " + F);
525        List<List<GenSolvablePolynomial<C>>> Gz = leftZeroRelationsArbitrary(F);
526        //logger.info("Gz: " + Gz);
527        List<GenSolvablePolynomial<C>> G1 = new ArrayList<GenSolvablePolynomial<C>>(Gz.size());
528        List<GenSolvablePolynomial<C>> G2 = new ArrayList<GenSolvablePolynomial<C>>(Gz.size());
529        for (List<GenSolvablePolynomial<C>> ll : Gz) {
530            if (!ll.get(0).isZERO()) { // && !ll.get(1).isZERO()
531                G1.add(ll.get(0)); // denominators
532                G2.add(ll.get(1)); // numerators
533            }
534        }
535        logger.info("G1(den): " + G1 + ", G2(num): " + G2);
536        SolvableExtendedGB<C> exgb = sbb.extLeftGB(G1);
537        logger.info("exgb.F: " + exgb.F + ", exgb.G: " + exgb.G);
538        List<GenSolvablePolynomial<C>> G = exgb.G;
539        int m = 0;
540        GenSolvablePolynomial<C> min = null;
541        for (int i = 0; i < G.size(); i++) {
542            if (min == null) {
543                min = G.get(i);
544                m = i;
545            } else if (min.compareTo(G.get(i)) > 0) {
546                min = G.get(i);
547                m = i;
548            }
549        }
550        //wrong: blas.scalarProduct(G2,exgb.G2F.get(m));
551        GenSolvablePolynomial<C> min2 = (GenSolvablePolynomial<C>) blas.scalarProduct(
552                        PolynomialList.<C> castToList(exgb.G2F.get(m)), PolynomialList.<C> castToList(G2));
553        logger.info("min(den): " + min + ", min(num): " + min2 + ", m = " + m + ", " + exgb.G2F.get(m));
554        // opposite order
555        GenSolvablePolynomial<C> n = min2; // nominator
556        GenSolvablePolynomial<C> d = min; // denominator
557        // normalize
558        if (d.signum() < 0) {
559            n = (GenSolvablePolynomial<C>) n.negate();
560            d = (GenSolvablePolynomial<C>) d.negate();
561        }
562        C lc = d.leadingBaseCoefficient();
563        if (!lc.isONE() && lc.isUnit()) {
564            lc = lc.inverse();
565            n = n.multiplyLeft(lc);
566            d = d.multiplyLeft(lc);
567        }
568        if (debug) {
569            int t = compare(a, b, n, d);
570            if (t != 0) {
571                oc[0] = a; // nominator
572                oc[1] = b; // denominator
573                throw new RuntimeException("simp wrong, giving up: t = " + t);
574                //logger.error("simp wrong, giving up: t = " + t);
575                //return oc;
576            }
577        }
578        oc[0] = n; // nominator
579        oc[1] = d; // denominator
580        return oc;
581    }
582
583}