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