Source code

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