001/*
002 * $Id$
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.LogManager;
014import org.apache.logging.log4j.Logger;
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        //System.out.println("isLeftZeroRelation leftGB: " + sbb.isLeftGB(exgb.G));
181        if (debug) {
182            logger.info("exgb = {}", exgb);
183        }
184        if (assertEnabled) {
185            logger.info("check1 exgb start");
186            if (!sbb.isLeftReductionMatrix(exgb)) {
187                logger.error("is reduction matrix ? false");
188            }
189            logger.info("check1 exgb end");
190        }
191        List<GenSolvablePolynomial<C>> G = exgb.G;
192        List<List<GenSolvablePolynomial<C>>> G2F = exgb.G2F;
193        List<List<GenSolvablePolynomial<C>>> F2G = exgb.F2G;
194
195        List<List<GenSolvablePolynomial<C>>> sg = leftZeroRelations(modv, G);
196        GenSolvablePolynomialRing<C> ring = G.get(0).ring;
197        ModuleList<C> S = new ModuleList<C>(ring, sg);
198        if (debug) {
199            logger.info("syz = {}", S);
200        }
201        if (assertEnabled) {
202            logger.info("check2 left syz start");
203            if (!isLeftZeroRelation(sg, G)) {
204                logger.error("is syzygy ? false");
205            }
206            logger.info("check2 left syz end");
207        }
208        //System.out.println("isLeftZeroRelation extGB: " + isLeftZeroRelation(sg, G));
209
210        List<List<GenSolvablePolynomial<C>>> sf;
211        sf = new ArrayList<List<GenSolvablePolynomial<C>>>(sg.size());
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(
233                                blas.vectorAdd(PolynomialList.<C> castToList(rf), 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        //System.out.println("isLeftZeroRelation backtrans: " + isLeftZeroRelation(sf, F));
249
250        List<List<GenSolvablePolynomial<C>>> M;
251        M = new ArrayList<List<GenSolvablePolynomial<C>>>(lenf);
252        for (List<GenSolvablePolynomial<C>> r : F2G) {
253            Iterator<GenSolvablePolynomial<C>> it = r.iterator();
254            Iterator<List<GenSolvablePolynomial<C>>> jt = G2F.iterator();
255
256            List<GenSolvablePolynomial<C>> rf;
257            rf = new ArrayList<GenSolvablePolynomial<C>>(lenf);
258            for (int m = 0; m < lenf; m++) {
259                rf.add(ring.getZERO());
260            }
261            while (it.hasNext() && jt.hasNext()) {
262                GenSolvablePolynomial<C> si = it.next();
263                List<GenSolvablePolynomial<C>> ai = jt.next();
264                //System.out.println("si = " + si);
265                //System.out.println("ai = " + ai);
266                if (si == null || ai == null) {
267                    continue;
268                }
269                //pi has wrong type, should be: List<GenSolvablePolynomial<C>>
270                List<GenPolynomial<C>> pi = blas.scalarProduct(si, PolynomialList.<C> castToList(ai));
271                //--List<GenPolynomial<C>> pi = blas.scalarProduct(PolynomialList.<C> castToList(ai), si);
272                //System.out.println("pi = " + pi);
273                rf = PolynomialList.<C> castToSolvableList(
274                                blas.vectorAdd(PolynomialList.<C> castToList(rf), pi));
275            }
276            if (it.hasNext() || jt.hasNext()) {
277                logger.error("zeroRelationsArbitrary wrong sizes");
278            }
279            //System.out.println("\nMg Mf = " + rf + "\n");
280            M.add(rf);
281        }
282        /* not true in general, debug only:
283        List<GenSolvablePolynomial<C>> F2 = new ArrayList<GenSolvablePolynomial<C>>( F.size() );
284        List<GenPolynomial<C>> Fp = PolynomialList.<C>castToList(F);
285        for ( List<GenSolvablePolynomial<C>> rr: M ) {
286            GenSolvablePolynomial<C> rrg = (GenSolvablePolynomial<C>) blas.scalarProduct(PolynomialList.<C>castToList(rr), Fp);
287            F2.add( rrg );
288        }
289        if ( ! F.equals( F2 ) ) {
290           logger.error("is FAB = F ? false");
291           System.out.println("F  = " + F);
292           System.out.println("F2 = " + F2);
293        }
294        System.out.println("isLeftZeroRelation F == F2: " + F.equals(F2));
295        */
296
297        int sflen = sf.size();
298        int i = 0;
299        for (List<GenSolvablePolynomial<C>> ri : M) {
300            List<GenSolvablePolynomial<C>> r2i;
301            r2i = new ArrayList<GenSolvablePolynomial<C>>(ri.size());
302            int j = 0;
303            for (GenSolvablePolynomial<C> rij : ri) {
304                GenSolvablePolynomial<C> p = null;
305                if (i == j) {
306                    p = (GenSolvablePolynomial<C>) ring.getONE().subtract(rij);
307                } else {
308                    if (rij != null) {
309                        p = (GenSolvablePolynomial<C>) rij.negate();
310                    }
311                }
312                r2i.add(p);
313                j++;
314            }
315            if (!blas.isZero( PolynomialList.<C> castToList(r2i)) ) {
316                sf.add(r2i);
317            }
318            i++;
319        }
320        if (assertEnabled) {
321            logger.info("check4 left syz start");
322            if (!isLeftZeroRelation(sf, F)) {
323                logger.error("is syz sf ? false");
324            }
325            logger.info("check4 left syz end");
326        }
327        //System.out.println("isLeftZeroRelation diagonal sf: " + isLeftZeroRelation(sf, F));
328        return sf;
329    }
330
331
332    /**
333     * Left Ore condition. Generators for the left Ore condition of two solvable
334     * polynomials.
335     * @param a solvable polynomial
336     * @param b solvable polynomial
337     * @return [p,q] with p*a = q*b
338     */
339    @SuppressWarnings({ "cast", "unchecked" })
340    public GenSolvablePolynomial<C>[] leftOreCond(GenSolvablePolynomial<C> a, GenSolvablePolynomial<C> b) {
341        if (a == null || a.isZERO() || b == null || b.isZERO()) {
342            throw new IllegalArgumentException("a and b must be non zero");
343        }
344        GenSolvablePolynomialRing<C> pfac = a.ring;
345        GenSolvablePolynomial<C>[] oc = (GenSolvablePolynomial<C>[]) new GenSolvablePolynomial[2];
346        if (a.equals(b)) {
347            oc[0] = pfac.getONE();
348            oc[1] = pfac.getONE();
349            return oc;
350        }
351        if (a.equals(b.negate())) {
352            oc[0] = pfac.getONE();
353            oc[1] = (GenSolvablePolynomial<C>) pfac.getONE().negate();
354            return oc;
355        }
356        if (a.isConstant()) {
357            if (pfac.coFac.isCommutative()) { // ??
358                oc[0] = b;
359                oc[1] = a;
360                return oc;
361            }
362            C c = a.leadingBaseCoefficient().inverse();
363            oc[0] = b.multiply(c); // was Left
364            oc[1] = pfac.getONE();
365            return oc;
366        }
367        if (b.isConstant()) {
368            if (pfac.coFac.isCommutative()) { // ??
369                oc[0] = b;
370                oc[1] = a;
371                return oc;
372            }
373            oc[0] = pfac.getONE();
374            C c = b.leadingBaseCoefficient().inverse();
375            oc[1] = a.multiply(c); // was Left
376            return oc;
377        }
378        logger.info("computing left Ore condition: {}, {}", a, b);
379        List<GenSolvablePolynomial<C>> F = new ArrayList<GenSolvablePolynomial<C>>(2);
380        F.add(a);
381        F.add(b);
382        List<List<GenSolvablePolynomial<C>>> Gz = leftZeroRelationsArbitrary(F);
383        boolean zr = isLeftZeroRelation(Gz, F);
384        //System.out.println("isLeftZeroRelation: " + zr); //isLeftZeroRelation(Gz, F));
385        // if (false) {
386        //     //System.out.println("Gz = " + Gz);
387        //     ModuleList<C> M = new ModuleList<C>(pfac, Gz);
388        //     System.out.println("M  = " + M.list);
389        //     ModuleList<C> GM = sbb.leftGB(M);
390        //     System.out.println("GM = " + GM.list);
391        //     Gz = GM.castToSolvableList();
392        // }
393        List<GenSolvablePolynomial<C>> G1 = null;
394        GenSolvablePolynomial<C> g1 = null;
395        for (List<GenSolvablePolynomial<C>> Gi : Gz) {
396            //System.out.println("Gil = " + Gi);
397            if (Gi.get(0).isZERO()) {
398                continue;
399            }
400            if (G1 == null) {
401                G1 = Gi;
402            }
403            if (g1 == null) {
404                g1 = G1.get(0);
405            } else if (g1.compareTo(Gi.get(0)) > 0) { //g1.degree() > Gi.get(0).degree() 
406                G1 = Gi;
407                g1 = G1.get(0);
408            }
409        }
410        oc[0] = g1; //G1.get(0);
411        oc[1] = (GenSolvablePolynomial<C>) G1.get(1).negate();
412        //logger.info("Ore multiple: {}, {}", oc[0].multiply(a), Arrays.toString(oc));
413        return oc;
414    }
415
416
417    /**
418     * Right Ore condition. Generators for the right Ore condition of two
419     * solvable polynomials.
420     * @param a solvable polynomial
421     * @param b solvable polynomial
422     * @return [p,q] with a*p = b*q
423     */
424    @SuppressWarnings({ "cast", "unchecked" })
425    public GenSolvablePolynomial<C>[] rightOreCond(GenSolvablePolynomial<C> a, GenSolvablePolynomial<C> b) {
426        if (a == null || a.isZERO() || b == null || b.isZERO()) {
427            throw new IllegalArgumentException("a and b must be non zero");
428        }
429        GenSolvablePolynomialRing<C> pfac = a.ring;
430        GenSolvablePolynomial<C>[] oc = (GenSolvablePolynomial<C>[]) new GenSolvablePolynomial[2];
431        if (a.equals(b)) {
432            oc[0] = pfac.getONE();
433            oc[1] = pfac.getONE();
434            return oc;
435        }
436        if (a.equals(b.negate())) {
437            oc[0] = pfac.getONE();
438            oc[1] = (GenSolvablePolynomial<C>) pfac.getONE().negate();
439            return oc;
440        }
441        if (a.isConstant()) {
442            if (pfac.coFac.isCommutative()) { // ??
443                oc[0] = b;
444                oc[1] = a;
445                return oc;
446            }
447            C c = a.leadingBaseCoefficient().inverse();
448            oc[0] = b.multiplyLeft(c);
449            oc[1] = pfac.getONE();
450            return oc;
451        }
452        if (b.isConstant()) {
453            if (pfac.coFac.isCommutative()) { // ??
454                oc[0] = b;
455                oc[1] = a;
456                return oc;
457            }
458            oc[0] = pfac.getONE();
459            C c = b.leadingBaseCoefficient().inverse();
460            oc[1] = a.multiplyLeft(c);
461            return oc;
462        }
463        logger.info("computing right Ore condition: {}, {}", a, b);
464        List<GenSolvablePolynomial<C>> F = new ArrayList<GenSolvablePolynomial<C>>(2);
465        F.add(a);
466        F.add(b);
467        List<List<GenSolvablePolynomial<C>>> Gz = rightZeroRelationsArbitrary(F);
468        //System.out.println("isRightZeroRelation: " + isRightZeroRelation(Gz, F));
469        List<GenSolvablePolynomial<C>> G1 = null;
470        GenSolvablePolynomial<C> g1 = null;
471        for (List<GenSolvablePolynomial<C>> Gi : Gz) {
472            //System.out.println("Gir = " + Gi);
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; //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(),
514                            b.totalDegree());
515            oc = new GenSolvablePolynomial[] { a, b };
516            return oc;
517        }
518        //GenSolvablePolynomialRing<C> pfac = a.ring;
519        oc = rightOreCond(a, b);
520        logger.info("oc = {}", Arrays.toString(oc)); // + ", a = {}, b = {}", a, b);
521        List<GenSolvablePolynomial<C>> F = new ArrayList<GenSolvablePolynomial<C>>(oc.length);
522        // opposite order and undo negation
523        F.add((GenSolvablePolynomial<C>) oc[1].negate());
524        F.add(oc[0]);
525        //logger.info("F = {}", F);
526        List<List<GenSolvablePolynomial<C>>> Gz = leftZeroRelationsArbitrary(F);
527        //logger.info("Gz: {}", Gz);
528        List<GenSolvablePolynomial<C>> G1 = new ArrayList<GenSolvablePolynomial<C>>(Gz.size());
529        List<GenSolvablePolynomial<C>> G2 = new ArrayList<GenSolvablePolynomial<C>>(Gz.size());
530        for (List<GenSolvablePolynomial<C>> ll : Gz) {
531            if (!ll.get(0).isZERO()) { // && !ll.get(1).isZERO()
532                G1.add(ll.get(0)); // denominators
533                G2.add(ll.get(1)); // numerators
534            }
535        }
536        logger.info("G1(den): {}, G2(num): {}", G1, G2);
537        SolvableExtendedGB<C> exgb = sbb.extLeftGB(G1);
538        logger.info("exgb.F: {}, exgb.G: {}", exgb.F, exgb.G);
539        List<GenSolvablePolynomial<C>> G = exgb.G;
540        int m = 0;
541        GenSolvablePolynomial<C> min = null;
542        for (int i = 0; i < G.size(); i++) {
543            if (min == null) {
544                min = G.get(i);
545                m = i;
546            } else if (min.compareTo(G.get(i)) > 0) {
547                min = G.get(i);
548                m = i;
549            }
550        }
551        //wrong: blas.scalarProduct(G2,exgb.G2F.get(m));
552        GenSolvablePolynomial<C> min2 = (GenSolvablePolynomial<C>) blas.scalarProduct(
553                        PolynomialList.<C> castToList(exgb.G2F.get(m)), PolynomialList.<C> castToList(G2));
554        logger.info("min(den): {}, min(num): {}, m = {}, {}", min, min2, m, exgb.G2F.get(m));
555        // opposite order
556        GenSolvablePolynomial<C> n = min2; // nominator
557        GenSolvablePolynomial<C> d = min; // denominator
558        // normalize
559        if (d.signum() < 0) {
560            n = (GenSolvablePolynomial<C>) n.negate();
561            d = (GenSolvablePolynomial<C>) d.negate();
562        }
563        C lc = d.leadingBaseCoefficient();
564        if (!lc.isONE() && lc.isUnit()) {
565            lc = lc.inverse();
566            n = n.multiplyLeft(lc);
567            d = d.multiplyLeft(lc);
568        }
569        if (debug) {
570            int t = compare(a, b, n, d);
571            if (t != 0) {
572                oc[0] = a; // nominator
573                oc[1] = b; // denominator
574                throw new RuntimeException("simp wrong, giving up: t = " + t);
575                //logger.error("simp wrong, giving up: t = {}", t);
576                //return oc;
577            }
578        }
579        oc[0] = n; // nominator
580        oc[1] = d; // denominator
581        return oc;
582    }
583
584}