001/*
002 * $Id$
003 */
004
005package edu.jas.gbufd;
006
007
008import java.util.ArrayList;
009import java.util.List;
010
011import org.apache.logging.log4j.Logger;
012import org.apache.logging.log4j.LogManager; 
013
014import edu.jas.gb.GroebnerBaseAbstract;
015import edu.jas.gb.PairList;
016import edu.jas.gb.Reduction;
017import edu.jas.gb.ReductionSeq;
018import edu.jas.poly.ExpVector;
019import edu.jas.poly.GenPolynomial;
020import edu.jas.poly.GenPolynomialRing;
021import edu.jas.poly.PolyUtil;
022import edu.jas.poly.TermOrder;
023import edu.jas.structure.GcdRingElem;
024import edu.jas.structure.RingFactory;
025
026
027/**
028 * Groebner Base sequential FGLM algorithm. Implements Groebner base computation
029 * via FGLM algorithm.
030 * @param <C> coefficient type
031 * @author Jan Suess
032 *
033 * @see edu.jas.application.GBAlgorithmBuilder
034 * @see edu.jas.gbufd.GBFactory
035 */
036public class GroebnerBaseFGLM<C extends GcdRingElem<C>> extends GroebnerBaseAbstract<C> {
037
038
039    private static final Logger logger = LogManager.getLogger(GroebnerBaseFGLM.class);
040
041
042    //private static final boolean debug = logger.isDebugEnabled();
043
044
045    /**
046     * The backing GB algorithm implementation.
047     */
048    private GroebnerBaseAbstract<C> sgb;
049
050
051    /**
052     * Constructor.
053     */
054    public GroebnerBaseFGLM() {
055        super();
056        sgb = null;
057    }
058
059
060    /**
061     * Constructor.
062     * @param red Reduction engine
063     */
064    public GroebnerBaseFGLM(Reduction<C> red) {
065        super(red);
066        sgb = null;
067    }
068
069
070    /**
071     * Constructor.
072     * @param red Reduction engine
073     * @param pl pair selection strategy
074     */
075    public GroebnerBaseFGLM(Reduction<C> red, PairList<C> pl) {
076        super(red, pl);
077        sgb = null;
078    }
079
080
081    /**
082     * Constructor.
083     * @param red Reduction engine
084     * @param pl pair selection strategy
085     * @param gb backing GB algorithm.
086     */
087    public GroebnerBaseFGLM(Reduction<C> red, PairList<C> pl, GroebnerBaseAbstract<C> gb) {
088        super(red, pl);
089        sgb = gb;
090    }
091
092
093    /**
094     * Constructor.
095     * @param gb backing GB algorithm.
096     */
097    public GroebnerBaseFGLM(GroebnerBaseAbstract<C> gb) {
098        super();
099        sgb = gb;
100    }
101
102
103    /**
104     * Get the String representation with GB engine.
105     * @see java.lang.Object#toString()
106     */
107    @Override
108    public String toString() {
109        if (sgb == null) {
110            return "GroebnerBaseFGLM()";
111        }
112        return "GroebnerBaseFGLM( " + sgb.toString() + " )";
113    }
114
115
116    /**
117     * Groebner base using FGLM algorithm.
118     * @param modv module variable number.
119     * @param F polynomial list.
120     * @return GB(F) a inv lex term order Groebner base of F.
121     */
122    public List<GenPolynomial<C>> GB(int modv, List<GenPolynomial<C>> F) {
123        if (modv != 0) {
124            throw new UnsupportedOperationException("case modv != 0 not yet implemented");
125        }
126        if (F == null || F.size() == 0) {
127            return F;
128        }
129        List<GenPolynomial<C>> G = new ArrayList<GenPolynomial<C>>();
130        if (F.size() <= 1) {
131            GenPolynomial<C> p = F.get(0).monic();
132            G.add(p);
133            return G;
134        }
135        // convert to graded term order
136        List<GenPolynomial<C>> Fp = new ArrayList<GenPolynomial<C>>(F.size());
137        GenPolynomialRing<C> pfac = F.get(0).ring;
138        if (!pfac.coFac.isField()) {
139            throw new IllegalArgumentException("coefficients not from a field: " + pfac.coFac);
140        }
141        TermOrder tord = new TermOrder(TermOrder.IGRLEX);
142        GenPolynomialRing<C> gfac = new GenPolynomialRing<C>(pfac.coFac, pfac.nvar, tord, pfac.getVars());
143        for (GenPolynomial<C> p : F) {
144            GenPolynomial<C> g = gfac.copy(p); // change term order
145            Fp.add(g);
146        }
147        // compute graded term order Groebner base
148        if (sgb == null) {
149            sgb = GBFactory.<C> getImplementation(pfac.coFac, strategy);
150        }
151        List<GenPolynomial<C>> Gp = sgb.GB(modv, Fp);
152        logger.info("graded GB = {}", Gp);
153        if (tord.equals(pfac.tord)) {
154            return Gp;
155        }
156        if (Gp.size() == 0) {
157            return Gp;
158        }
159        if (Gp.size() == 1) { // also dimension -1
160            GenPolynomial<C> p = pfac.copy(Gp.get(0)); // change term order
161            G.add(p);
162            return G;
163        }
164        // check dimension zero
165        int z = commonZeroTest(Gp);
166        if (z > 0) {
167            logger.error("use Groebner Walk algorithm");
168            throw new IllegalArgumentException("ideal(G) not zero dimensional, dim =  " + z);
169        }
170        // compute invlex Groebner base via FGLM
171        G = convGroebnerToLex(Gp);
172        return G;
173    }
174
175
176    /**
177     * Algorithm CONVGROEBNER: Converts Groebner bases w.r.t. total degree
178     * termorder into Groebner base w.r.t to inverse lexicographical term order
179     * @return Groebner base w.r.t to inverse lexicographical term order
180     */
181    public List<GenPolynomial<C>> convGroebnerToLex(List<GenPolynomial<C>> groebnerBasis) {
182        if (groebnerBasis == null || groebnerBasis.size() == 0) {
183            throw new IllegalArgumentException("G may not be null or empty");
184        }
185        //Polynomial ring of input Groebnerbasis G
186        GenPolynomialRing<C> ring = groebnerBasis.get(0).ring;
187        int numberOfVariables = ring.nvar; //Number of Variables of the given Polynomial Ring
188        String[] ArrayOfVariables = ring.getVars(); //Variables of given polynomial ring w.r.t. to input G
189        RingFactory<C> cfac = ring.coFac;
190
191        //Main Algorithm
192        //Initialization
193
194        TermOrder invlex = new TermOrder(TermOrder.INVLEX);
195        //Polynomial ring of newGB with invlex order
196        GenPolynomialRing<C> ufac = new GenPolynomialRing<C>(cfac, numberOfVariables, invlex,
197                        ArrayOfVariables);
198
199        //Local Lists
200        List<GenPolynomial<C>> newGB = new ArrayList<GenPolynomial<C>>(); //Instantiate the return list of polynomials
201        List<GenPolynomial<C>> H = new ArrayList<GenPolynomial<C>>(); //Instantiate a help list of polynomials
202        List<GenPolynomial<C>> redTerms = new ArrayList<GenPolynomial<C>>();//Instantiate the return list of reduced terms
203
204        //Local Polynomials
205        GenPolynomial<C> t = ring.ONE; //Create ONE polynom of original polynomial ring
206        GenPolynomial<C> h; //Create help polynomial
207        GenPolynomial<GenPolynomial<C>> hh; //h as polynomial in rfac
208        GenPolynomial<GenPolynomial<C>> p; //Create another help polynomial
209        redTerms.add(t); //Add ONE to list of reduced terms
210
211        //create new indeterminate Y1
212        int indeterminates = 1; //Number of indeterminates, starting with Y1
213        GenPolynomialRing<C> cpfac = createRingOfIndeterminates(ring, indeterminates);
214        GenPolynomialRing<GenPolynomial<C>> rfac = new GenPolynomialRing<GenPolynomial<C>>(cpfac, ring);
215        GenPolynomial<GenPolynomial<C>> q = rfac.getZERO().sum(cpfac.univariate(0));
216
217        //Main while loop
218        int z = -1;
219        t = lMinterm(H, t);
220        while (t != null) {
221            //System.out.println("t = " + t);
222            h = red.normalform(groebnerBasis, t);
223            //System.out.println("Zwischennormalform h = " + h.toString());
224            hh = PolyUtil.<C> toRecursive(rfac, h);
225            p = hh.sum(q);
226            List<GenPolynomial<C>> Cf = new ArrayList<GenPolynomial<C>>(p.getMap().values());
227            Cf = red.irreducibleSet(Cf);
228            //System.out.println("Cf = " + Cf);
229            //System.out.println("Current Polynomial ring in Y_n: " + rfac.toString());
230
231            z = commonZeroTest(Cf);
232            //System.out.println("z = " + z);
233            if (z != 0) { //z=1 OR z=-1 --> Infinite number of solutions OR No solution
234                indeterminates++; //then, increase number of indeterminates by one
235                redTerms.add(t); //add current t to list of reduced terms
236                cpfac = addIndeterminate(cpfac);
237                rfac = new GenPolynomialRing<GenPolynomial<C>>(cpfac, ring);
238                hh = PolyUtil.<C> toRecursive(rfac, h);
239                GenPolynomial<GenPolynomial<C>> Yt = rfac.getZERO().sum(cpfac.univariate(0));
240                GenPolynomial<GenPolynomial<C>> Yth = hh.multiply(Yt);
241                q = PolyUtil.<C> extendCoefficients(rfac, q, 0, 0L);
242                q = Yth.sum(q);
243            } else { // z=0 --> one solution
244                GenPolynomial<C> g = ufac.getZERO();
245                for (GenPolynomial<C> pc : Cf) {
246                    ExpVector e = pc.leadingExpVector();
247                    //System.out.println("e = " + e);
248                    if (e == null) {
249                        continue;
250                    }
251                    int[] v = e.dependencyOnVariables();
252                    if (v == null || v.length == 0) {
253                        continue;
254                    }
255                    int vi = v[0];
256                    vi = indeterminates - vi;
257                    C tc = pc.trailingBaseCoefficient();
258                    if (!tc.isZERO()) {
259                        tc = tc.negate();
260                        GenPolynomial<C> csRedterm = redTerms.get(vi - 1).multiply(tc);
261                        //System.out.println("csRedterm = " + csRedterm);
262                        g = g.sum(csRedterm);
263                    }
264                }
265                g = g.sum(t);
266                g = ufac.copy(g);
267                H.add(t);
268                if (!g.isZERO()) {
269                    newGB.add(g);
270                    logger.info("new element for GB = {}", g.leadingExpVector());
271                }
272            }
273            t = lMinterm(H, t); // compute lMINTERM of current t (lexMinterm)
274        }
275        //logger.info("GB = {}", newGB);
276        return newGB;
277    }
278
279
280    /**
281     * Algorithm lMinterm: MINTERM algorithm for inverse lexicographical term
282     * order.
283     * @param t Term
284     * @param G Groebner basis
285     * @return Term that specifies condition (D) or null (Condition (D) in
286     *         "A computational approach to commutative algebra", Becker,
287     *         Weispfenning, Kredel 1993, p. 427)
288     */
289    public GenPolynomial<C> lMinterm(List<GenPolynomial<C>> G, GenPolynomial<C> t) {
290        //not ok: if ( G == null || G.size() == 0 ) ...
291        GenPolynomialRing<C> ring = t.ring;
292        int numberOfVariables = ring.nvar;
293        GenPolynomial<C> u = new GenPolynomial<C>(ring, t.leadingBaseCoefficient(), t.leadingExpVector()); //HeadTerm of of input polynomial
294        ReductionSeq<C> redHelp = new ReductionSeq<C>(); // Create instance of ReductionSeq to use method isReducible
295        //not ok: if ( redHelp.isTopReducible(G,u) ) ...
296        for (int i = numberOfVariables - 1; i >= 0; i--) { // Walk through all variables, starting with least w.r.t to lex-order
297            GenPolynomial<C> x = ring.univariate(i); // Create Linear Polynomial X_i
298            u = u.multiply(x); // Multiply current u with x
299            if (!redHelp.isTopReducible(G, u)) { // Check if any term in HT(G) divides current u
300                return u;
301            }
302            GenPolynomial<C> s = ring.univariate(i, u.degree(numberOfVariables - (i + 1))); //if not, eliminate variable x_i
303            u = u.divide(s);
304        }
305        return null;
306    }
307
308
309    /**
310     * Compute the residues to given polynomial list.
311     * @return List of reduced terms
312     */
313    public List<GenPolynomial<C>> redTerms(List<GenPolynomial<C>> groebnerBasis) {
314        if (groebnerBasis == null || groebnerBasis.size() == 0) {
315            throw new IllegalArgumentException("groebnerBasis may not be null or empty");
316        }
317        GenPolynomialRing<C> ring = groebnerBasis.get(0).ring;
318        int numberOfVariables = ring.nvar; //Number of Variables of the given Polynomial Ring
319        long[] degrees = new long[numberOfVariables]; //Array for the degree-limits for the reduced terms
320
321        List<GenPolynomial<C>> terms = new ArrayList<GenPolynomial<C>>(); //Instantiate the return object
322        for (GenPolynomial<C> g : groebnerBasis) { //For each polynomial of G
323            if (g.isONE()) {
324                terms.clear();
325                return terms; //If 1 e G, return empty list terms
326            }
327            ExpVector e = g.leadingExpVector(); //Take the exponent of the leading monomial             
328            if (e.totalDeg() == e.maxDeg()) { //and check, whether a variable x_i is isolated
329                for (int i = 0; i < numberOfVariables; i++) {
330                    long exp = e.getVal(i);
331                    if (exp > 0) {
332                        degrees[i] = exp; //if true, add the degree of univariate x_i to array degrees
333                    }
334                }
335            }
336        }
337        long max = maxArray(degrees); //Find maximum in Array degrees
338        for (int i = 0; i < degrees.length; i++) { //Set all zero grades to maximum of array "degrees"
339            if (degrees[i] == 0) {
340                logger.info("dimension not zero, setting degree to {}", max);
341                degrees[i] = max; //--> to "make" the ideal zero-dimensional
342            }
343        }
344        terms.add(ring.ONE); //Add the one-polynomial of the ring to the list of reduced terms
345        ReductionSeq<C> s = new ReductionSeq<C>(); //Create instance of ReductionSeq to use method isReducible
346
347        //Main Algorithm
348        for (int i = 0; i < numberOfVariables; i++) {
349            GenPolynomial<C> x = ring.univariate(i); //Create  Linear Polynomial X_i
350            List<GenPolynomial<C>> S = new ArrayList<GenPolynomial<C>>(terms); //Copy all entries of return list "terms" into list "S"
351            for (GenPolynomial<C> t : S) {
352                for (int l = 1; l <= degrees[i]; l++) {
353                    t = t.multiply(x); //Multiply current element t with Linear Polynomial X_i
354                    if (!s.isReducible(groebnerBasis, t)) { //Check, if t is irreducible mod groebnerbase
355                        terms.add(t); //Add t to return list terms
356                    }
357                }
358            }
359        }
360        return terms;
361    }
362
363
364    /**
365     * Internal method to create a polynomial ring in i indeterminates. Create
366     * new ring over coefficients of ring with i variables Y1,...,Yi
367     * (indeterminate)
368     * @return polynomial ring with variables Y1...Yi and coefficient of ring.
369     */
370    GenPolynomialRing<C> createRingOfIndeterminates(GenPolynomialRing<C> ring, int i) {
371        RingFactory<C> cfac = ring.coFac;
372        int indeterminates = i;
373        String[] stringIndeterminates = new String[indeterminates];
374        for (int j = 1; j <= indeterminates; j++) {
375            stringIndeterminates[j - 1] = ("Y" + j);
376        }
377        TermOrder invlex = new TermOrder(TermOrder.INVLEX);
378        GenPolynomialRing<C> cpfac = new GenPolynomialRing<C>(cfac, indeterminates, invlex,
379                        stringIndeterminates);
380        return cpfac;
381    }
382
383
384    /**
385     * Internal method to add new indeterminates. Add another variable
386     * (indeterminate) Y_{i+1} to existing ring
387     * @return polynomial ring with variables Y1,..,Yi,Yi+1 and coefficients of
388     *         ring.
389     */
390    GenPolynomialRing<C> addIndeterminate(GenPolynomialRing<C> ring) {
391        String[] stringIndeterminates = new String[1];
392        int number = ring.nvar + 1;
393        stringIndeterminates[0] = ("Y" + number);
394        ring = ring.extend(stringIndeterminates);
395        return ring;
396    }
397
398
399    /**
400     * Maximum of an array.
401     * @return maximum of an array
402     */
403    long maxArray(long[] t) {
404        if (t.length == 0) {
405            return 0L;
406        }
407        long maximum = t[0];
408        for (int i = 1; i < t.length; i++) {
409            if (t[i] > maximum) {
410                maximum = t[i];
411            }
412        }
413        return maximum;
414    }
415
416
417    /**
418     * Cleanup and terminate ThreadPool.
419     */
420    @Override
421    public void terminate() {
422        if (sgb == null) {
423            return;
424        }
425        sgb.terminate();
426    }
427
428
429    /**
430     * Cancel ThreadPool.
431     */
432    @Override
433    public int cancel() {
434        if (sgb == null) {
435            return 0;
436        }
437        return sgb.cancel();
438    }
439
440}