Source code

001/*
002 * $Id: GroebnerBaseSeq.java 5908 2018-08-25 11:49:29Z kredel $
003 */
004
005package edu.jas.gb;
006
007import java.util.ArrayList;
008import java.util.List;
009import java.util.ListIterator;
010
011import org.apache.logging.log4j.Logger;
012import org.apache.logging.log4j.LogManager; 
013
014import edu.jas.structure.RingElem;
015import edu.jas.gb.OrderedPairlist;
016import edu.jas.poly.GenPolynomial;
017import edu.jas.poly.GenPolynomialRing;
018import edu.jas.poly.PolyUtil;
019
020
021/**
022 * Groebner Base sequential algorithm.
023 * Implements Groebner bases and GB test.
024 * @param <C> coefficient type
025 * @author Heinz Kredel
026 *
027 * @see edu.jas.application.GBAlgorithmBuilder
028 * @see edu.jas.gbufd.GBFactory
029 */
030
031public class GroebnerBaseSeq<C extends RingElem<C>> 
032    extends GroebnerBaseAbstract<C>  {
033
034    private static final Logger logger = LogManager.getLogger(GroebnerBaseSeq.class);
035    private static final boolean debug = logger.isDebugEnabled();
036
037
038    /**
039     * Constructor.
040     */
041    public GroebnerBaseSeq() {
042        super();
043    }
044
045
046    /**
047     * Constructor.
048     * @param red Reduction engine
049     */
050    public GroebnerBaseSeq(Reduction<C> red) {
051        super(red);
052    }
053
054
055    /**
056     * Constructor.
057     * @param pl pair selection strategy
058     */
059    public GroebnerBaseSeq(PairList<C> pl) {
060        super(pl);
061    }
062
063
064    /**
065     * Constructor.
066     * @param red Reduction engine
067     * @param pl pair selection strategy
068     */
069    public GroebnerBaseSeq(Reduction<C> red, PairList<C> pl) {
070        super(red,pl);
071    }
072
073
074    /**
075     * Groebner base using pairlist class.
076     * @param modv module variable number.
077     * @param F polynomial list.
078     * @return GB(F) a Groebner base of F.
079     */
080    public List<GenPolynomial<C>> GB( int modv, List<GenPolynomial<C>> F ) {  
081        List<GenPolynomial<C>> G = normalizeZerosOnes(F);
082        G = PolyUtil.<C> monic(G);
083        if ( G.size() <= 1 ) {
084            return G;
085        }
086        GenPolynomialRing<C> ring = G.get(0).ring;
087        if ( ! ring.coFac.isField() ) {
088            throw new IllegalArgumentException("coefficients not from a field");
089        }
090        PairList<C> pairlist = strategy.create( modv, ring ); 
091        pairlist.put(G);
092        logger.info("start " + pairlist); 
093
094        Pair<C> pair;
095        GenPolynomial<C> pi, pj, S, H;
096        while ( pairlist.hasNext() ) {
097            pair = pairlist.removeNext();
098            //logger.debug("pair = " + pair);
099            if ( pair == null ) {
100                continue; 
101            }
102            pi = pair.pi; 
103            pj = pair.pj; 
104            if ( /*false &&*/ debug ) {
105                logger.debug("pi    = " + pi );
106                logger.debug("pj    = " + pj );
107            }
108
109            S = red.SPolynomial( pi, pj );
110            if ( S.isZERO() ) {
111                pair.setZero();
112                continue;
113            }
114            if ( debug ) {
115                logger.debug("ht(S) = " + S.leadingExpVector() );
116            }
117
118            H = red.normalform( G, S );
119            if ( debug ) {
120                //logger.info("pair = " + pair); 
121                //logger.info("ht(S) = " + S.monic()); //.leadingExpVector() );
122                logger.info("ht(H) = " + H.monic()); //.leadingExpVector() );
123            }
124            if ( H.isZERO() ) {
125                pair.setZero();
126                continue;
127            }
128            H = H.monic();
129            if ( debug ) {
130                logger.info("ht(H) = " + H.leadingExpVector() );
131            }
132
133            H = H.monic();
134            if ( H.isONE() ) {
135                G.clear(); G.add( H );
136                pairlist.putOne();
137                logger.info("end " + pairlist); 
138                return G; // since no threads are activated
139            }
140            if ( debug ) {
141                logger.info("H = " + H );
142            }
143            if ( H.length() > 0 ) {
144                //l++;
145                G.add( H );
146                pairlist.put( H );
147            }
148        }
149        logger.debug("#sequential list = " + G.size());
150        G = minimalGB(G);
151        logger.info("end " + pairlist); 
152        return G;
153    }
154
155
156    /**
157     * Extended Groebner base using critical pair class.
158     * @param modv module variable number.
159     * @param F polynomial list.
160     * @return a container for an extended Groebner base of F.
161     */
162    @Override
163    public ExtendedGB<C> 
164        extGB( int modv, List<GenPolynomial<C>> F ) {  
165        if ( F == null || F.isEmpty() ) {
166            throw new IllegalArgumentException("null or empty F not allowed");
167        }
168        List<GenPolynomial<C>> G = new ArrayList<GenPolynomial<C>>();
169        List<List<GenPolynomial<C>>> F2G = new ArrayList<List<GenPolynomial<C>>>();
170        List<List<GenPolynomial<C>>> G2F = new ArrayList<List<GenPolynomial<C>>>();
171        PairList<C> pairlist = null; 
172        boolean oneInGB = false;
173        int len = F.size();
174
175        List<GenPolynomial<C>> row = null;
176        List<GenPolynomial<C>> rows = null;
177        List<GenPolynomial<C>> rowh = null;
178        GenPolynomialRing<C> ring = null;
179        GenPolynomial<C> H;
180        GenPolynomial<C> p;
181
182        int nzlen = 0;
183        for ( GenPolynomial<C> f : F ) { 
184            if ( f.length() > 0 ) {
185                nzlen++;
186            }
187            if ( ring == null ) {
188                ring = f.ring;
189            }
190        }
191        GenPolynomial<C> mone = ring.getONE(); //.negate();
192        int k = 0;
193        ListIterator<GenPolynomial<C>> it = F.listIterator();
194        while ( it.hasNext() ) { 
195            p = it.next();
196            if ( p.length() > 0 ) {
197                row = new ArrayList<GenPolynomial<C>>( nzlen );
198                for ( int j = 0; j < nzlen; j++ ) {
199                    row.add(null);
200                }
201                //C c = p.leadingBaseCoefficient();
202                //c = c.inverse();
203                //p = p.multiply( c );
204                row.set( k, mone ); //.multiply(c) );
205                k++;
206                if ( p.isUnit() ) {
207                    G.clear(); G.add( p );
208                    G2F.clear(); G2F.add( row );
209                    oneInGB = true;
210                    break;
211                }
212                G.add( p );
213                G2F.add( row );
214                if ( pairlist == null ) {
215                    //pairlist = new OrderedPairlist<C>( modv, p.ring );
216                    pairlist = strategy.create( modv, p.ring );
217                    if ( ! p.ring.coFac.isField() ) {
218                        throw new RuntimeException("coefficients not from a field");
219                    }
220                }
221                // putOne not required
222                pairlist.put( p );
223            } else { 
224                len--;
225            }
226        }
227        ExtendedGB<C> exgb;
228        if ( len <= 1 || oneInGB ) {
229            // adjust F2G
230            for ( GenPolynomial<C> f : F ) {
231                row = new ArrayList<GenPolynomial<C>>( G.size() );
232                for ( int j = 0; j < G.size(); j++ ) {
233                    row.add(null);
234                }
235                H = red.normalform( row, G, f );
236                if ( ! H.isZERO() ) {
237                    logger.error("nonzero H = " + H );
238                }
239                F2G.add( row );
240            }
241            exgb = new ExtendedGB<C>(F,G,F2G,G2F);
242            //System.out.println("exgb 1 = " + exgb);
243            return exgb;
244        }
245
246        Pair<C> pair;
247        int i, j;
248        GenPolynomial<C> pi;
249        GenPolynomial<C> pj;
250        GenPolynomial<C> S;
251        GenPolynomial<C> x;
252        GenPolynomial<C> y;
253        //GenPolynomial<C> z;
254        while ( pairlist.hasNext() && ! oneInGB ) {
255            pair = pairlist.removeNext();
256            if ( pair == null ) { 
257                continue; 
258            }
259            i = pair.i; 
260            j = pair.j; 
261            pi = pair.pi; 
262            pj = pair.pj; 
263            if ( debug ) {
264                logger.info("i, pi    = " + i + ", " + pi );
265                logger.info("j, pj    = " + j + ", " + pj );
266            }
267
268            rows = new ArrayList<GenPolynomial<C>>( G.size() );
269            for ( int m = 0; m < G.size(); m++ ) {
270                rows.add(null);
271            }
272            S = red.SPolynomial( rows, i, pi, j, pj );
273            if ( debug ) {
274                logger.debug("is reduction S = " 
275                             + red.isReductionNF( rows, G, ring.getZERO(), S ) );
276            }
277            if ( S.isZERO() ) {
278                // do not add to G2F
279                continue;
280            }
281            if ( debug ) {
282                logger.debug("ht(S) = " + S.leadingExpVector() );
283            }
284
285            rowh = new ArrayList<GenPolynomial<C>>( G.size() );
286            for ( int m = 0; m < G.size(); m++ ) {
287                rowh.add(null);
288            }
289            H = red.normalform( rowh, G, S );
290            if ( debug ) {
291                logger.debug("is reduction H = " 
292                             + red.isReductionNF( rowh, G, S, H ) );
293            }
294            if ( H.isZERO() ) {
295                // do not add to G2F
296                continue;
297            }
298            if ( debug ) {
299                logger.debug("ht(H) = " + H.leadingExpVector() );
300            }
301
302            row = new ArrayList<GenPolynomial<C>>( G.size()+1 );
303            for ( int m = 0; m < G.size(); m++ ) {
304                x = rows.get(m);
305                if ( x != null ) {
306                    //System.out.println("ms = " + m + " " + x);
307                    x = x.negate();
308                }
309                y = rowh.get(m);
310                if ( y != null ) {
311                    y = y.negate();
312                    //System.out.println("mh = " + m + " " + y);
313                }
314                if ( x == null ) {
315                    x = y;
316                } else {
317                    x = x.sum( y );
318                }
319                //System.out.println("mx = " + m + " " + x);
320                row.add( x );
321            }
322            if ( debug ) {
323                logger.debug("is reduction 0+sum(row,G) == H : " 
324                             + red.isReductionNF( row, G, H, ring.getZERO() ) );
325            }
326            row.add( null );
327
328
329            //  H = H.monic();
330            C c = H.leadingBaseCoefficient();
331            c = c.inverse();
332            H = H.multiply( c );
333            row = blas.scalarProduct( mone.multiply(c), row );
334            row.set( G.size(), mone );
335            if ( H.isONE() ) {
336                // G.clear(); 
337                G.add( H );
338                G2F.add( row );
339                oneInGB = true;
340                break; 
341            }
342            if ( debug ) {
343                logger.debug("H = " + H );
344            }
345            G.add( H );
346            pairlist.put( H );
347            G2F.add( row );
348        }
349        if ( debug ) {
350            exgb = new ExtendedGB<C>(F,G,F2G,G2F);
351            logger.info("exgb unnorm = " + exgb);
352        }
353        G2F = normalizeMatrix( F.size(), G2F );
354        if ( debug ) {
355            exgb = new ExtendedGB<C>(F,G,F2G,G2F);
356            logger.info("exgb nonmin = " + exgb);
357            boolean t2 = isReductionMatrix( exgb );
358            logger.info("exgb t2 = " + t2);
359        }
360        exgb = minimalExtendedGB(F.size(),G,G2F);
361        G = exgb.G;
362        G2F = exgb.G2F;
363        logger.debug("#sequential list = " + G.size());
364        logger.info("" + pairlist); 
365        // setup matrices F and F2G
366        for ( GenPolynomial<C> f : F ) {
367            row = new ArrayList<GenPolynomial<C>>( G.size() );
368            for ( int m = 0; m < G.size(); m++ ) {
369                row.add(null);
370            }
371            H = red.normalform( row, G, f );
372            if ( ! H.isZERO() ) {
373                logger.error("nonzero H = " + H );
374            }
375            F2G.add( row );
376        }
377        exgb = new ExtendedGB<C>(F,G,F2G,G2F);
378        if ( debug ) {
379            logger.info("exgb nonmin = " + exgb);
380            boolean t2 = isReductionMatrix( exgb );
381            logger.info("exgb t2 = " + t2);
382        }
383        return exgb;
384    }
385
386}