001/* 002 * $Id: PolyUtilApp.java 4125 2012-08-19 19:05:22Z kredel $ 003 */ 004 005package edu.jas.application; 006 007 008import java.util.ArrayList; 009import java.util.Arrays; 010import java.util.List; 011import java.util.Map; 012import java.util.Set; 013import java.util.SortedMap; 014import java.util.TreeMap; 015import java.util.TreeSet; 016 017import org.apache.log4j.Logger; 018 019import edu.jas.arith.BigDecimal; 020import edu.jas.arith.Product; 021import edu.jas.arith.ProductRing; 022import edu.jas.arith.Rational; 023import edu.jas.poly.AlgebraicNumber; 024import edu.jas.poly.AlgebraicNumberRing; 025import edu.jas.poly.Complex; 026import edu.jas.poly.ComplexRing; 027import edu.jas.poly.ExpVector; 028import edu.jas.poly.GenPolynomial; 029import edu.jas.poly.GenPolynomialRing; 030import edu.jas.poly.PolyUtil; 031import edu.jas.poly.PolynomialList; 032import edu.jas.poly.TermOrder; 033import edu.jas.root.ComplexRoots; 034import edu.jas.root.ComplexRootsAbstract; 035import edu.jas.root.ComplexRootsSturm; 036import edu.jas.root.Interval; 037import edu.jas.root.InvalidBoundaryException; 038import edu.jas.root.RealAlgebraicNumber; 039import edu.jas.root.RealAlgebraicRing; 040import edu.jas.root.RealRootTuple; 041import edu.jas.root.RealRootsAbstract; 042import edu.jas.root.RealRootsSturm; 043import edu.jas.root.Rectangle; 044import edu.jas.root.RootFactory; 045import edu.jas.structure.GcdRingElem; 046import edu.jas.structure.RingElem; 047import edu.jas.structure.RingFactory; 048import edu.jas.structure.UnaryFunctor; 049import edu.jas.util.ListUtil; 050 051 052/** 053 * Polynomial utilities for applications, for example conversion ExpVector to 054 * Product or zero dimensional ideal root computation. 055 * @param <C> coefficient type 056 * @author Heinz Kredel 057 */ 058public class PolyUtilApp<C extends RingElem<C>> { 059 060 061 private static final Logger logger = Logger.getLogger(PolyUtilApp.class); 062 063 064 private static boolean debug = logger.isDebugEnabled(); 065 066 067 /** 068 * Product representation. 069 * @param <C> coefficient type. 070 * @param pfac polynomial ring factory. 071 * @param L list of polynomials to be represented. 072 * @return Product represenation of L in the polynomial ring pfac. 073 */ 074 public static <C extends GcdRingElem<C>> List<GenPolynomial<Product<Residue<C>>>> toProductRes( 075 GenPolynomialRing<Product<Residue<C>>> pfac, List<GenPolynomial<GenPolynomial<C>>> L) { 076 077 List<GenPolynomial<Product<Residue<C>>>> list = new ArrayList<GenPolynomial<Product<Residue<C>>>>(); 078 if (L == null || L.size() == 0) { 079 return list; 080 } 081 GenPolynomial<Product<Residue<C>>> b; 082 for (GenPolynomial<GenPolynomial<C>> a : L) { 083 b = toProductRes(pfac, a); 084 list.add(b); 085 } 086 return list; 087 } 088 089 090 /** 091 * Product representation. 092 * @param <C> coefficient type. 093 * @param pfac polynomial ring factory. 094 * @param A polynomial to be represented. 095 * @return Product represenation of A in the polynomial ring pfac. 096 */ 097 public static <C extends GcdRingElem<C>> GenPolynomial<Product<Residue<C>>> toProductRes( 098 GenPolynomialRing<Product<Residue<C>>> pfac, GenPolynomial<GenPolynomial<C>> A) { 099 100 GenPolynomial<Product<Residue<C>>> P = pfac.getZERO().copy(); 101 if (A == null || A.isZERO()) { 102 return P; 103 } 104 RingFactory<Product<Residue<C>>> rpfac = pfac.coFac; 105 ProductRing<Residue<C>> fac = (ProductRing<Residue<C>>) rpfac; 106 Product<Residue<C>> p; 107 for (Map.Entry<ExpVector, GenPolynomial<C>> y : A.getMap().entrySet()) { 108 ExpVector e = y.getKey(); 109 GenPolynomial<C> a = y.getValue(); 110 p = toProductRes(fac, a); 111 if (!p.isZERO()) { 112 P.doPutToMap(e, p); 113 } 114 } 115 return P; 116 } 117 118 119 /** 120 * Product representation. 121 * @param <C> coefficient type. 122 * @param pfac product ring factory. 123 * @param c coefficient to be represented. 124 * @return Product represenation of c in the ring pfac. 125 */ 126 public static <C extends GcdRingElem<C>> Product<Residue<C>> toProductRes(ProductRing<Residue<C>> pfac, 127 GenPolynomial<C> c) { 128 129 SortedMap<Integer, Residue<C>> elem = new TreeMap<Integer, Residue<C>>(); 130 for (int i = 0; i < pfac.length(); i++) { 131 RingFactory<Residue<C>> rfac = pfac.getFactory(i); 132 ResidueRing<C> fac = (ResidueRing<C>) rfac; 133 Residue<C> u = new Residue<C>(fac, c); 134 //fac.fromInteger( c.getVal() ); 135 if (!u.isZERO()) { 136 elem.put(i, u); 137 } 138 } 139 return new Product<Residue<C>>(pfac, elem); 140 } 141 142 143 /** 144 * Product residue representation. 145 * @param <C> coefficient type. 146 * @param CS list of ColoredSystems from comprehensive GB system. 147 * @return Product residue represenation of CS. 148 */ 149 public static <C extends GcdRingElem<C>> List<GenPolynomial<Product<Residue<C>>>> toProductRes( 150 List<ColoredSystem<C>> CS) { 151 152 List<GenPolynomial<Product<Residue<C>>>> list = new ArrayList<GenPolynomial<Product<Residue<C>>>>(); 153 if (CS == null || CS.size() == 0) { 154 return list; 155 } 156 GenPolynomialRing<GenPolynomial<C>> pr = null; 157 List<RingFactory<Residue<C>>> rrl = new ArrayList<RingFactory<Residue<C>>>(CS.size()); 158 for (ColoredSystem<C> cs : CS) { 159 Ideal<C> id = cs.condition.zero; 160 ResidueRing<C> r = new ResidueRing<C>(id); 161 if (!rrl.contains(r)) { 162 rrl.add(r); 163 } 164 if (pr == null) { 165 if (cs.list.size() > 0) { 166 pr = cs.list.get(0).green.ring; 167 } 168 } 169 } 170 ProductRing<Residue<C>> pfac; 171 pfac = new ProductRing<Residue<C>>(rrl); 172 //System.out.println("pfac = " + pfac); 173 GenPolynomialRing<Product<Residue<C>>> rf = new GenPolynomialRing<Product<Residue<C>>>(pfac, pr.nvar, 174 pr.tord, pr.getVars()); 175 GroebnerSystem<C> gs = new GroebnerSystem<C>(CS); 176 List<GenPolynomial<GenPolynomial<C>>> F = gs.getCGB(); 177 list = PolyUtilApp.<C> toProductRes(rf, F); 178 return list; 179 } 180 181 182 /** 183 * Residue coefficient representation. 184 * @param pfac polynomial ring factory. 185 * @param L list of polynomials to be represented. 186 * @return Represenation of L in the polynomial ring pfac. 187 */ 188 public static <C extends GcdRingElem<C>> List<GenPolynomial<Residue<C>>> toResidue( 189 GenPolynomialRing<Residue<C>> pfac, List<GenPolynomial<GenPolynomial<C>>> L) { 190 List<GenPolynomial<Residue<C>>> list = new ArrayList<GenPolynomial<Residue<C>>>(); 191 if (L == null || L.size() == 0) { 192 return list; 193 } 194 GenPolynomial<Residue<C>> b; 195 for (GenPolynomial<GenPolynomial<C>> a : L) { 196 b = toResidue(pfac, a); 197 if (!b.isZERO()) { 198 list.add(b); 199 } 200 } 201 return list; 202 } 203 204 205 /** 206 * Residue coefficient representation. 207 * @param pfac polynomial ring factory. 208 * @param A polynomial to be represented. 209 * @return Represenation of A in the polynomial ring pfac. 210 */ 211 public static <C extends GcdRingElem<C>> GenPolynomial<Residue<C>> toResidue( 212 GenPolynomialRing<Residue<C>> pfac, GenPolynomial<GenPolynomial<C>> A) { 213 GenPolynomial<Residue<C>> P = pfac.getZERO().copy(); 214 if (A == null || A.isZERO()) { 215 return P; 216 } 217 RingFactory<Residue<C>> rpfac = pfac.coFac; 218 ResidueRing<C> fac = (ResidueRing<C>) rpfac; 219 Residue<C> p; 220 for (Map.Entry<ExpVector, GenPolynomial<C>> y : A.getMap().entrySet()) { 221 ExpVector e = y.getKey(); 222 GenPolynomial<C> a = y.getValue(); 223 p = new Residue<C>(fac, a); 224 if (!p.isZERO()) { 225 P.doPutToMap(e, p); 226 } 227 } 228 return P; 229 } 230 231 232 /** 233 * Product slice. 234 * @param <C> coefficient type. 235 * @param L list of polynomials with product coefficients. 236 * @return Slices represenation of L. 237 */ 238 public static <C extends GcdRingElem<C>> Map<Ideal<C>, PolynomialList<GenPolynomial<C>>> productSlice( 239 PolynomialList<Product<Residue<C>>> L) { 240 241 Map<Ideal<C>, PolynomialList<GenPolynomial<C>>> map; 242 RingFactory<Product<Residue<C>>> fpr = L.ring.coFac; 243 ProductRing<Residue<C>> pr = (ProductRing<Residue<C>>) fpr; 244 int s = pr.length(); 245 map = new TreeMap<Ideal<C>, PolynomialList<GenPolynomial<C>>>(); 246 List<GenPolynomial<GenPolynomial<C>>> slist; 247 248 List<GenPolynomial<Product<Residue<C>>>> plist = L.list; 249 PolynomialList<GenPolynomial<C>> spl; 250 251 for (int i = 0; i < s; i++) { 252 RingFactory<Residue<C>> r = pr.getFactory(i); 253 ResidueRing<C> rr = (ResidueRing<C>) r; 254 Ideal<C> id = rr.ideal; 255 GenPolynomialRing<C> cof = rr.ring; 256 GenPolynomialRing<GenPolynomial<C>> pfc; 257 pfc = new GenPolynomialRing<GenPolynomial<C>>(cof, L.ring); 258 slist = fromProduct(pfc, plist, i); 259 spl = new PolynomialList<GenPolynomial<C>>(pfc, slist); 260 PolynomialList<GenPolynomial<C>> d = map.get(id); 261 if (d != null) { 262 throw new RuntimeException("ideal exists twice " + id); 263 } 264 map.put(id, spl); 265 } 266 return map; 267 } 268 269 270 /** 271 * Product slice at i. 272 * @param <C> coefficient type. 273 * @param L list of polynomials with product coeffients. 274 * @param i index of slice. 275 * @return Slice of of L at i. 276 */ 277 public static <C extends GcdRingElem<C>> PolynomialList<GenPolynomial<C>> productSlice( 278 PolynomialList<Product<Residue<C>>> L, int i) { 279 280 RingFactory<Product<Residue<C>>> fpr = L.ring.coFac; 281 ProductRing<Residue<C>> pr = (ProductRing<Residue<C>>) fpr; 282 List<GenPolynomial<GenPolynomial<C>>> slist; 283 284 List<GenPolynomial<Product<Residue<C>>>> plist = L.list; 285 PolynomialList<GenPolynomial<C>> spl; 286 287 RingFactory<Residue<C>> r = pr.getFactory(i); 288 ResidueRing<C> rr = (ResidueRing<C>) r; 289 GenPolynomialRing<C> cof = rr.ring; 290 GenPolynomialRing<GenPolynomial<C>> pfc; 291 pfc = new GenPolynomialRing<GenPolynomial<C>>(cof, L.ring); 292 slist = fromProduct(pfc, plist, i); 293 spl = new PolynomialList<GenPolynomial<C>>(pfc, slist); 294 return spl; 295 } 296 297 298 /** 299 * From product representation. 300 * @param <C> coefficient type. 301 * @param pfac polynomial ring factory. 302 * @param L list of polynomials to be converted from product representation. 303 * @param i index of product representation to be taken. 304 * @return Represenation of i-slice of L in the polynomial ring pfac. 305 */ 306 public static <C extends GcdRingElem<C>> List<GenPolynomial<GenPolynomial<C>>> fromProduct( 307 GenPolynomialRing<GenPolynomial<C>> pfac, List<GenPolynomial<Product<Residue<C>>>> L, 308 int i) { 309 310 List<GenPolynomial<GenPolynomial<C>>> list = new ArrayList<GenPolynomial<GenPolynomial<C>>>(); 311 312 if (L == null || L.size() == 0) { 313 return list; 314 } 315 GenPolynomial<GenPolynomial<C>> b; 316 for (GenPolynomial<Product<Residue<C>>> a : L) { 317 b = fromProduct(pfac, a, i); 318 if (!b.isZERO()) { 319 b = b.abs(); 320 if (!list.contains(b)) { 321 list.add(b); 322 } 323 } 324 } 325 return list; 326 } 327 328 329 /** 330 * From product representation. 331 * @param <C> coefficient type. 332 * @param pfac polynomial ring factory. 333 * @param P polynomial to be converted from product representation. 334 * @param i index of product representation to be taken. 335 * @return Represenation of i-slice of P in the polynomial ring pfac. 336 */ 337 public static <C extends GcdRingElem<C>> GenPolynomial<GenPolynomial<C>> fromProduct( 338 GenPolynomialRing<GenPolynomial<C>> pfac, GenPolynomial<Product<Residue<C>>> P, int i) { 339 340 GenPolynomial<GenPolynomial<C>> b = pfac.getZERO().copy(); 341 if (P == null || P.isZERO()) { 342 return b; 343 } 344 345 for (Map.Entry<ExpVector, Product<Residue<C>>> y : P.getMap().entrySet()) { 346 ExpVector e = y.getKey(); 347 Product<Residue<C>> a = y.getValue(); 348 Residue<C> r = a.get(i); 349 if (r != null && !r.isZERO()) { 350 GenPolynomial<C> p = r.val; 351 if (!p.isZERO()) { 352 b.doPutToMap(e, p); 353 } 354 } 355 } 356 return b; 357 } 358 359 360 /** 361 * Product slice to String. 362 * @param <C> coefficient type. 363 * @param L list of polynomials with to be represented. 364 * @return Product represenation of L in the polynomial ring pfac. 365 */ 366 public static <C extends GcdRingElem<C>> String productSliceToString( 367 Map<Ideal<C>, PolynomialList<GenPolynomial<C>>> L) { 368 Set<GenPolynomial<GenPolynomial<C>>> sl = new TreeSet<GenPolynomial<GenPolynomial<C>>>(); 369 PolynomialList<GenPolynomial<C>> pl = null; 370 StringBuffer sb = new StringBuffer(); //"\nproductSlice ----------------- begin"); 371 for (Map.Entry<Ideal<C>, PolynomialList<GenPolynomial<C>>> en : L.entrySet()) { 372 sb.append("\n\ncondition == 0:\n"); 373 sb.append(en.getKey().list.toScript()); 374 pl = en.getValue(); //L.get(id); 375 sl.addAll(pl.list); 376 sb.append("\ncorresponding ideal:\n"); 377 sb.append(pl.toScript()); 378 } 379 //List<GenPolynomial<GenPolynomial<C>>> sll 380 // = new ArrayList<GenPolynomial<GenPolynomial<C>>>( sl ); 381 //pl = new PolynomialList<GenPolynomial<C>>(pl.ring,sll); 382 // sb.append("\nunion = " + pl.toString()); 383 //sb.append("\nproductSlice ------------------------- end\n"); 384 return sb.toString(); 385 } 386 387 388 /** 389 * Product slice to String. 390 * @param <C> coefficient type. 391 * @param L list of polynomials with product coefficients. 392 * @return string represenation of slices of L. 393 */ 394 public static <C extends GcdRingElem<C>> String productToString(PolynomialList<Product<Residue<C>>> L) { 395 Map<Ideal<C>, PolynomialList<GenPolynomial<C>>> M; 396 M = productSlice(L); 397 String s = productSliceToString(M); 398 return s; 399 } 400 401 402 /** 403 * Construct superset of complex roots for zero dimensional ideal(G). 404 * @param I zero dimensional ideal. 405 * @param eps desired precision. 406 * @return list of coordinates of complex roots for ideal(G) 407 */ 408 public static <D extends GcdRingElem<D> & Rational> List<List<Complex<BigDecimal>>> complexRootTuples( 409 Ideal<D> I, D eps) { 410 List<GenPolynomial<D>> univs = I.constructUnivariate(); 411 if (logger.isInfoEnabled()) { 412 logger.info("univs = " + univs); 413 } 414 return complexRoots(I, univs, eps); 415 } 416 417 418 /** 419 * Construct superset of complex roots for zero dimensional ideal(G). 420 * @param I zero dimensional ideal. 421 * @param univs list of univariate polynomials. 422 * @param eps desired precision. 423 * @return list of coordinates of complex roots for ideal(G) 424 */ 425 public static <D extends GcdRingElem<D> & Rational> List<List<Complex<BigDecimal>>> complexRoots( 426 Ideal<D> I, List<GenPolynomial<D>> univs, D eps) { 427 List<List<Complex<BigDecimal>>> croots = new ArrayList<List<Complex<BigDecimal>>>(); 428 RingFactory<D> cf = (RingFactory<D>) I.list.ring.coFac; 429 ComplexRing<D> cr = new ComplexRing<D>(cf); 430 ComplexRootsAbstract<D> cra = new ComplexRootsSturm<D>(cr); 431 List<GenPolynomial<Complex<D>>> cunivs = new ArrayList<GenPolynomial<Complex<D>>>(); 432 for (GenPolynomial<D> p : univs) { 433 GenPolynomialRing<Complex<D>> pfac = new GenPolynomialRing<Complex<D>>(cr, p.ring); 434 //System.out.println("pfac = " + pfac.toScript()); 435 GenPolynomial<Complex<D>> cp = PolyUtil.<D> toComplex(pfac, p); 436 cunivs.add(cp); 437 //System.out.println("cp = " + cp); 438 } 439 for (int i = 0; i < I.list.ring.nvar; i++) { 440 List<Complex<BigDecimal>> cri = cra.approximateRoots(cunivs.get(i), eps); 441 //System.out.println("cri = " + cri); 442 croots.add(cri); 443 } 444 croots = ListUtil.<Complex<BigDecimal>> tupleFromList(croots); 445 return croots; 446 } 447 448 449 /** 450 * Construct superset of complex roots for zero dimensional ideal(G). 451 * @param Il list of zero dimensional ideals with univariate polynomials. 452 * @param eps desired precision. 453 * @return list of coordinates of complex roots for ideal(cap_i(G_i)) 454 */ 455 public static <D extends GcdRingElem<D> & Rational> List<List<Complex<BigDecimal>>> complexRootTuples( 456 List<IdealWithUniv<D>> Il, D eps) { 457 List<List<Complex<BigDecimal>>> croots = new ArrayList<List<Complex<BigDecimal>>>(); 458 for (IdealWithUniv<D> I : Il) { 459 List<List<Complex<BigDecimal>>> cr = complexRoots(I.ideal, I.upolys, eps); 460 croots.addAll(cr); 461 } 462 return croots; 463 } 464 465 466 /** 467 * Construct superset of complex roots for zero dimensional ideal(G). 468 * @param Il list of zero dimensional ideals with univariate polynomials. 469 * @param eps desired precision. 470 * @return list of ideals with coordinates of complex roots for 471 * ideal(cap_i(G_i)) 472 */ 473 public static <D extends GcdRingElem<D> & Rational> List<IdealWithComplexRoots<D>> complexRoots( 474 List<IdealWithUniv<D>> Il, D eps) { 475 List<IdealWithComplexRoots<D>> Ic = new ArrayList<IdealWithComplexRoots<D>>(Il.size()); 476 for (IdealWithUniv<D> I : Il) { 477 List<List<Complex<BigDecimal>>> cr = complexRoots(I.ideal, I.upolys, eps); 478 IdealWithComplexRoots<D> ic = new IdealWithComplexRoots<D>(I, cr); 479 Ic.add(ic); 480 } 481 return Ic; 482 } 483 484 485 /** 486 * Construct superset of complex roots for zero dimensional ideal(G). 487 * @param G list of polynomials of a of zero dimensional ideal. 488 * @param eps desired precision. 489 * @return list of ideals with coordinates of complex roots for ideal(G) 490 */ 491 public static <D extends GcdRingElem<D> & Rational> List<IdealWithComplexRoots<D>> complexRoots( 492 Ideal<D> G, D eps) { 493 List<IdealWithUniv<D>> Il = G.zeroDimDecomposition(); 494 return complexRoots(Il, eps); 495 } 496 497 498 /** 499 * Construct superset of real roots for zero dimensional ideal(G). 500 * @param I zero dimensional ideal. 501 * @param eps desired precision. 502 * @return list of coordinates of real roots for ideal(G) 503 */ 504 public static <D extends GcdRingElem<D> & Rational> List<List<BigDecimal>> realRootTuples(Ideal<D> I, 505 D eps) { 506 List<GenPolynomial<D>> univs = I.constructUnivariate(); 507 if (logger.isInfoEnabled()) { 508 logger.info("univs = " + univs); 509 } 510 return realRoots(I, univs, eps); 511 } 512 513 514 /** 515 * Construct superset of real roots for zero dimensional ideal(G). 516 * @param I zero dimensional ideal. 517 * @param univs list of univariate polynomials. 518 * @param eps desired precision. 519 * @return list of coordinates of real roots for ideal(G) 520 */ 521 public static <D extends GcdRingElem<D> & Rational> List<List<BigDecimal>> realRoots(Ideal<D> I, 522 List<GenPolynomial<D>> univs, D eps) { 523 List<List<BigDecimal>> roots = new ArrayList<List<BigDecimal>>(); 524 //RingFactory<D> cf = (RingFactory<D>) I.list.ring.coFac; 525 RealRootsAbstract<D> rra = new RealRootsSturm<D>(); 526 for (int i = 0; i < I.list.ring.nvar; i++) { 527 List<BigDecimal> rri = rra.approximateRoots(univs.get(i), eps); 528 //System.out.println("rri = " + rri); 529 roots.add(rri); 530 } 531 //System.out.println("roots-1 = " + roots); 532 roots = ListUtil.<BigDecimal> tupleFromList(roots); 533 //System.out.println("roots-2 = " + roots); 534 return roots; 535 } 536 537 538 /** 539 * Construct superset of real roots for zero dimensional ideal(G). 540 * @param Il list of zero dimensional ideals with univariate polynomials. 541 * @param eps desired precision. 542 * @return list of coordinates of real roots for ideal(cap_i(G_i)) 543 */ 544 public static <D extends GcdRingElem<D> & Rational> List<List<BigDecimal>> realRootTuples( 545 List<IdealWithUniv<D>> Il, D eps) { 546 List<List<BigDecimal>> rroots = new ArrayList<List<BigDecimal>>(); 547 for (IdealWithUniv<D> I : Il) { 548 List<List<BigDecimal>> rr = realRoots(I.ideal, I.upolys, eps); 549 rroots.addAll(rr); 550 } 551 return rroots; 552 } 553 554 555 /** 556 * Construct superset of real roots for zero dimensional ideal(G). 557 * @param Il list of zero dimensional ideals with univariate polynomials. 558 * @param eps desired precision. 559 * @return list of ideals with coordinates of real roots for 560 * ideal(cap_i(G_i)) 561 */ 562 public static <D extends GcdRingElem<D> & Rational> List<IdealWithRealRoots<D>> realRoots( 563 List<IdealWithUniv<D>> Il, D eps) { 564 List<IdealWithRealRoots<D>> Ir = new ArrayList<IdealWithRealRoots<D>>(Il.size()); 565 for (IdealWithUniv<D> I : Il) { 566 List<List<BigDecimal>> rr = realRoots(I.ideal, I.upolys, eps); 567 IdealWithRealRoots<D> ir = new IdealWithRealRoots<D>(I, rr); 568 Ir.add(ir); 569 } 570 return Ir; 571 } 572 573 574 /** 575 * Construct superset of real roots for zero dimensional ideal(G). 576 * @param G list of polynomials of a of zero dimensional ideal. 577 * @param eps desired precision. 578 * @return list of ideals with coordinates of real roots for ideal(G) 579 */ 580 public static <D extends GcdRingElem<D> & Rational> List<IdealWithRealRoots<D>> realRoots(Ideal<D> G, 581 D eps) { 582 List<IdealWithUniv<D>> Il = G.zeroDimDecomposition(); 583 return realRoots(Il, eps); 584 } 585 586 587 /** 588 * Test for real roots of zero dimensional ideal(L). 589 * @param L list of polynomials. 590 * @param roots list of real roots for ideal(G). 591 * @param eps desired precision. 592 * @return true if root is a list of coordinates of real roots for ideal(L) 593 */ 594 public static boolean isRealRoots(List<GenPolynomial<BigDecimal>> L, List<List<BigDecimal>> roots, 595 BigDecimal eps) { 596 if (L == null || L.size() == 0) { 597 return true; 598 } 599 // polynomials with decimal coefficients 600 BigDecimal dc = BigDecimal.ONE; 601 GenPolynomialRing<BigDecimal> dfac = L.get(0).ring; 602 //System.out.println("dfac = " + dfac); 603 for (GenPolynomial<BigDecimal> dp : L) { 604 //System.out.println("dp = " + dp); 605 for (List<BigDecimal> r : roots) { 606 //System.out.println("r = " + r); 607 BigDecimal ev = PolyUtil.<BigDecimal> evaluateAll(dc, dfac, dp, r); 608 if (ev.abs().compareTo(eps) > 0) { 609 System.out.println("ev = " + ev); 610 return false; 611 } 612 } 613 } 614 return true; 615 } 616 617 618 /** 619 * Test for complex roots of zero dimensional ideal(L). 620 * @param L list of polynomials. 621 * @param roots list of real roots for ideal(G). 622 * @param eps desired precision. 623 * @return true if root is a list of coordinates of complex roots for 624 * ideal(L) 625 */ 626 public static boolean isComplexRoots(List<GenPolynomial<Complex<BigDecimal>>> L, 627 List<List<Complex<BigDecimal>>> roots, BigDecimal eps) { 628 if (L == null || L.size() == 0) { 629 return true; 630 } 631 // polynomials with decimal coefficients 632 BigDecimal dc = BigDecimal.ONE; 633 ComplexRing<BigDecimal> dcc = new ComplexRing<BigDecimal>(dc); 634 GenPolynomialRing<Complex<BigDecimal>> dfac = L.get(0).ring; 635 //System.out.println("dfac = " + dfac); 636 for (GenPolynomial<Complex<BigDecimal>> dp : L) { 637 //System.out.println("dp = " + dp); 638 for (List<Complex<BigDecimal>> r : roots) { 639 //System.out.println("r = " + r); 640 Complex<BigDecimal> ev = PolyUtil.<Complex<BigDecimal>> evaluateAll(dcc, dfac, dp, r); 641 if (ev.norm().getRe().compareTo(eps) > 0) { 642 System.out.println("ev = " + ev); 643 return false; 644 } 645 } 646 } 647 return true; 648 } 649 650 651 /** 652 * Construct real roots for zero dimensional ideal(G). 653 * @param I zero dimensional ideal with univariate irreducible polynomials 654 * and bi-variate polynomials. 655 * @return real algebraic roots for ideal(G) 656 */ 657 public static <D extends GcdRingElem<D> & Rational> IdealWithRealAlgebraicRoots<D> realAlgebraicRoots( 658 IdealWithUniv<D> I) { 659 List<List<RealAlgebraicNumber<D>>> ran = new ArrayList<List<RealAlgebraicNumber<D>>>(); 660 if (I == null) { 661 throw new IllegalArgumentException("null ideal not permitted"); 662 } 663 if (I.ideal == null || I.upolys == null) { 664 throw new IllegalArgumentException("null ideal components not permitted " + I); 665 } 666 if (I.ideal.isZERO() || I.upolys.size() == 0) { 667 return new IdealWithRealAlgebraicRoots<D>(I, ran); 668 } 669 GenPolynomialRing<D> fac = I.ideal.list.ring; 670 // case i == 0: 671 GenPolynomial<D> p0 = I.upolys.get(0); 672 GenPolynomial<D> p0p = PolyUtil.<D> selectWithVariable(I.ideal.list.list, fac.nvar - 1); 673 if (p0p == null) { 674 throw new RuntimeException("no polynomial found in " + (fac.nvar - 1) + " of " + I.ideal); 675 } 676 //System.out.println("p0 = " + p0); 677 if (logger.isInfoEnabled()) { 678 logger.info("p0p = " + p0p); 679 } 680 int[] dep0 = p0p.degreeVector().dependencyOnVariables(); 681 //System.out.println("dep0 = " + Arrays.toString(dep0)); 682 if (dep0.length != 1) { 683 throw new RuntimeException("wrong number of variables " + Arrays.toString(dep0)); 684 } 685 List<RealAlgebraicNumber<D>> rra = RootFactory.<D> realAlgebraicNumbersIrred(p0); 686 if (logger.isInfoEnabled()) { 687 List<Interval<D>> il = new ArrayList<Interval<D>>(); 688 for (RealAlgebraicNumber<D> rr : rra) { 689 il.add(rr.ring.getRoot()); 690 } 691 logger.info("roots(p0) = " + il); 692 } 693 for (RealAlgebraicNumber<D> rr : rra) { 694 List<RealAlgebraicNumber<D>> rl = new ArrayList<RealAlgebraicNumber<D>>(); 695 rl.add(rr); 696 ran.add(rl); 697 } 698 // case i > 0: 699 for (int i = 1; i < I.upolys.size(); i++) { 700 List<List<RealAlgebraicNumber<D>>> rn = new ArrayList<List<RealAlgebraicNumber<D>>>(); 701 GenPolynomial<D> pi = I.upolys.get(i); 702 GenPolynomial<D> pip = PolyUtil.selectWithVariable(I.ideal.list.list, fac.nvar - 1 - i); 703 if (pip == null) { 704 throw new RuntimeException("no polynomial found in " + (fac.nvar - 1 - i) + " of " + I.ideal); 705 } 706 //System.out.println("i = " + i); 707 //System.out.println("pi = " + pi); 708 if (logger.isInfoEnabled()) { 709 logger.info("pi = " + pi); 710 logger.info("pip = " + pip); 711 } 712 int[] depi = pip.degreeVector().dependencyOnVariables(); 713 //System.out.println("depi = " + Arrays.toString(depi)); 714 if (depi.length < 1 || depi.length > 2) { 715 throw new RuntimeException("wrong number of variables " + Arrays.toString(depi)); 716 } 717 rra = RootFactory.<D> realAlgebraicNumbersIrred(pi); 718 if (logger.isInfoEnabled()) { 719 List<Interval<D>> il = new ArrayList<Interval<D>>(); 720 for (RealAlgebraicNumber<D> rr : rra) { 721 il.add(rr.ring.getRoot()); 722 } 723 logger.info("roots(pi) = " + il); 724 } 725 if (depi.length == 1) { 726 // all combinations are roots of the ideal I 727 for (RealAlgebraicNumber<D> rr : rra) { 728 //System.out.println("rr.ring = " + rr.ring); 729 for (List<RealAlgebraicNumber<D>> rx : ran) { 730 //System.out.println("rx = " + rx); 731 List<RealAlgebraicNumber<D>> ry = new ArrayList<RealAlgebraicNumber<D>>(); 732 ry.addAll(rx); 733 ry.add(rr); 734 rn.add(ry); 735 } 736 } 737 } else { // depi.length == 2 738 // select roots of the ideal I 739 GenPolynomial<D> pip2 = PolyUtil.<D> removeUnusedUpperVariables(pip); 740 //System.out.println("pip2 = " + pip2.ring); 741 GenPolynomialRing<D> ufac = pip2.ring.contract(1); 742 TermOrder to = new TermOrder(TermOrder.INVLEX); 743 GenPolynomialRing<GenPolynomial<D>> rfac = new GenPolynomialRing<GenPolynomial<D>>(ufac, 1, 744 to); 745 GenPolynomial<GenPolynomial<D>> pip2r = PolyUtil.<D> recursive(rfac, pip2); 746 int ix = fac.nvar - 1 - depi[depi.length - 1]; 747 //System.out.println("ix = " + ix); 748 for (RealAlgebraicNumber<D> rr : rra) { 749 //System.out.println("rr.ring = " + rr.ring); 750 Interval<D> rroot = rr.ring.getRoot(); 751 GenPolynomial<D> pip2el = PolyUtil.<D> evaluateMainRecursive(ufac, pip2r, rroot.left); 752 GenPolynomial<D> pip2er = PolyUtil.<D> evaluateMainRecursive(ufac, pip2r, rroot.right); 753 GenPolynomialRing<D> upfac = I.upolys.get(ix).ring; 754 GenPolynomial<D> pip2elc = convert(upfac, pip2el); 755 GenPolynomial<D> pip2erc = convert(upfac, pip2er); 756 //System.out.println("pip2elc = " + pip2elc); 757 //System.out.println("pip2erc = " + pip2erc); 758 for (List<RealAlgebraicNumber<D>> rx : ran) { 759 //System.out.println("rx = " + rx); 760 RealAlgebraicRing<D> rar = rx.get(ix).ring; 761 //System.out.println("rar = " + rar.toScript()); 762 RealAlgebraicNumber<D> rel = new RealAlgebraicNumber<D>(rar, pip2elc); 763 RealAlgebraicNumber<D> rer = new RealAlgebraicNumber<D>(rar, pip2erc); 764 int sl = rel.signum(); 765 int sr = rer.signum(); 766 //System.out.println("sl = " + sl + ", sr = " + sr + ", sl*sr = " + (sl*sr)); 767 if (sl * sr <= 0) { 768 //System.out.println("sl * sr <= 0: rar = " + rar.toScript()); 769 List<RealAlgebraicNumber<D>> ry = new ArrayList<RealAlgebraicNumber<D>>(); 770 ry.addAll(rx); 771 ry.add(rr); 772 rn.add(ry); 773 } 774 } 775 } 776 } 777 ran = rn; 778 } 779 if (logger.isInfoEnabled()) { 780 for (List<RealAlgebraicNumber<D>> rz : ran) { 781 List<Interval<D>> il = new ArrayList<Interval<D>>(); 782 for (RealAlgebraicNumber<D> rr : rz) { 783 il.add(rr.ring.getRoot()); 784 } 785 logger.info("root-tuple = " + il); 786 } 787 } 788 IdealWithRealAlgebraicRoots<D> Ir = new IdealWithRealAlgebraicRoots<D>(I, ran); 789 return Ir; 790 } 791 792 793 /** 794 * Construct real roots for zero dimensional ideal(G). 795 * @param I list of zero dimensional ideal with univariate irreducible 796 * polynomials and bi-variate polynomials. 797 * @return list of real algebraic roots for all ideal(I_i) 798 */ 799 public static <D extends GcdRingElem<D> & Rational> List<IdealWithRealAlgebraicRoots<D>> realAlgebraicRoots( 800 List<IdealWithUniv<D>> I) { 801 List<IdealWithRealAlgebraicRoots<D>> lir = new ArrayList<IdealWithRealAlgebraicRoots<D>>(I.size()); 802 for (IdealWithUniv<D> iu : I) { 803 IdealWithRealAlgebraicRoots<D> iur = PolyUtilApp.<D> realAlgebraicRoots(iu); 804 //System.out.println("iur = " + iur); 805 lir.add(iur); 806 } 807 return lir; 808 } 809 810 811 /** 812 * Construct complex roots for zero dimensional ideal(G). 813 * @param I zero dimensional ideal with univariate irreducible polynomials 814 * and bi-variate polynomials. 815 * @return complex algebraic roots for ideal(G) <b>Note:</b> implementation 816 * contains errors, do not use. 817 */ 818 public static <D extends GcdRingElem<D> & Rational> IdealWithComplexAlgebraicRoots<D> complexAlgebraicRootsWrong( // Wrong 819 IdealWithUniv<D> I) { 820 List<List<Complex<edu.jas.application.RealAlgebraicNumber<D>>>> can; 821 can = new ArrayList<List<Complex<edu.jas.application.RealAlgebraicNumber<D>>>>(); 822 if (I == null) { 823 throw new IllegalArgumentException("null ideal not permitted"); 824 } 825 if (I.ideal == null || I.upolys == null) { 826 throw new IllegalArgumentException("null ideal components not permitted " + I); 827 } 828 if (I.ideal.isZERO() || I.upolys.size() == 0) { 829 return new IdealWithComplexAlgebraicRoots<D>(I, can); 830 } 831 GenPolynomialRing<D> fac = I.ideal.list.ring; 832 if (fac.nvar == 0) { 833 return new IdealWithComplexAlgebraicRoots<D>(I, can); 834 } 835 if (fac.nvar != I.upolys.size()) { 836 throw new IllegalArgumentException("ideal not zero dimnsional: " + I); 837 } 838 // case i == 0: 839 GenPolynomial<D> p0 = I.upolys.get(0); 840 GenPolynomial<D> p0p = PolyUtil.<D> selectWithVariable(I.ideal.list.list, fac.nvar - 1); 841 if (p0p == null) { 842 throw new RuntimeException("no polynomial found in " + (fac.nvar - 1) + " of " + I.ideal); 843 } 844 if (logger.isInfoEnabled()) { 845 logger.info("p0 = " + p0); 846 logger.info("p0p = " + p0p); 847 } 848 int[] dep0 = p0p.degreeVector().dependencyOnVariables(); 849 //System.out.println("dep0 = " + Arrays.toString(dep0)); 850 if (dep0.length != 1) { 851 throw new RuntimeException("wrong number of variables " + Arrays.toString(dep0)); 852 } 853 RingFactory<D> cfac = p0.ring.coFac; 854 ComplexRing<D> ccfac = new ComplexRing<D>(cfac); 855 GenPolynomialRing<Complex<D>> facc = new GenPolynomialRing<Complex<D>>(ccfac, p0.ring); 856 GenPolynomial<Complex<D>> p0c = PolyUtil.<D> complexFromAny(facc, p0); 857 List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cra; 858 cra = edu.jas.application.RootFactory.<D> complexAlgebraicNumbersSquarefree(p0c); 859 logger.info("#roots(p0c) = " + cra.size()); 860 if (debug) { 861 boolean t = edu.jas.application.RootFactory.<D> isRoot(p0c, cra); 862 if (!t) { 863 throw new RuntimeException("no roots of " + p0c); 864 } 865 } 866 for (Complex<edu.jas.application.RealAlgebraicNumber<D>> cr : cra) { 867 List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cl; 868 cl = new ArrayList<Complex<edu.jas.application.RealAlgebraicNumber<D>>>(); 869 cl.add(cr); 870 can.add(cl); 871 } 872 if (fac.nvar == 1) { 873 return new IdealWithComplexAlgebraicRoots<D>(I, can); 874 } 875 // case i > 0: 876 for (int i = 1; i < I.upolys.size(); i++) { 877 List<List<Complex<edu.jas.application.RealAlgebraicNumber<D>>>> cn; 878 cn = new ArrayList<List<Complex<edu.jas.application.RealAlgebraicNumber<D>>>>(); 879 GenPolynomial<D> pi = I.upolys.get(i); 880 GenPolynomial<D> pip = PolyUtil.selectWithVariable(I.ideal.list.list, fac.nvar - 1 - i); 881 if (pip == null) { 882 throw new RuntimeException("no polynomial found in " + (fac.nvar - 1 - i) + " of " + I.ideal); 883 } 884 if (logger.isInfoEnabled()) { 885 logger.info("pi(" + i + ") = " + pi); 886 logger.info("pip = " + pip); 887 } 888 facc = new GenPolynomialRing<Complex<D>>(ccfac, pi.ring); 889 GenPolynomial<Complex<D>> pic = PolyUtil.<D> complexFromAny(facc, pi); 890 int[] depi = pip.degreeVector().dependencyOnVariables(); 891 //System.out.println("depi = " + Arrays.toString(depi)); 892 if (depi.length < 1 || depi.length > 2) { 893 throw new RuntimeException("wrong number of variables " + Arrays.toString(depi) + " for " 894 + pip); 895 } 896 cra = edu.jas.application.RootFactory.<D> complexAlgebraicNumbersSquarefree(pic); 897 logger.info("#roots(pic) = " + cra.size()); 898 if (debug) { 899 boolean t = edu.jas.application.RootFactory.<D> isRoot(pic, cra); 900 if (!t) { 901 throw new RuntimeException("no roots of " + pic); 902 } 903 } 904 if (depi.length == 1) { 905 // all combinations are roots of the ideal I 906 for (Complex<edu.jas.application.RealAlgebraicNumber<D>> cr : cra) { 907 //System.out.println("cr.ring = " + cr.ring); 908 for (List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cx : can) { 909 //System.out.println("cx = " + cx); 910 List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cy; 911 cy = new ArrayList<Complex<edu.jas.application.RealAlgebraicNumber<D>>>(); 912 cy.addAll(cx); 913 cy.add(cr); 914 cn.add(cy); 915 } 916 } 917 } else { // depi.length == 2 918 // select roots of the ideal I 919 GenPolynomial<D> pip2 = PolyUtil.<D> removeUnusedUpperVariables(pip); 920 GenPolynomialRing<GenPolynomial<D>> rfac = pip2.ring.recursive(1); 921 GenPolynomialRing<D> ufac = pip2.ring.contract(1); 922 GenPolynomialRing<Complex<D>> ucfac = new GenPolynomialRing<Complex<D>>(ccfac, ufac); 923 GenPolynomialRing<Complex<D>> c2fac = new GenPolynomialRing<Complex<D>>(ccfac, pip2.ring); 924 GenPolynomial<Complex<D>> pip2c = PolyUtil.<D> complexFromAny(c2fac, pip2); 925 //System.out.println("pip2c = " + pip2c); 926 GenPolynomialRing<GenPolynomial<Complex<D>>> rcfac; 927 rcfac = new GenPolynomialRing<GenPolynomial<Complex<D>>>(ucfac, rfac); 928 GenPolynomial<GenPolynomial<Complex<D>>> pip2cr = PolyUtil.<Complex<D>> recursive(rcfac, 929 pip2c); 930 //System.out.println("pip2cr = " + pip2cr); 931 932 int ix = fac.nvar - 1 - depi[depi.length - 1]; 933 //System.out.println("ix = " + ix); 934 for (Complex<edu.jas.application.RealAlgebraicNumber<D>> cr : cra) { 935 System.out.println("cr = " + toString(cr)); // <---------------------------------- 936 edu.jas.application.RealAlgebraicRing<D> cring = (edu.jas.application.RealAlgebraicRing<D>) cr.ring.ring; 937 RealRootTuple<D> rroot = cring.getRoot(); 938 List<RealAlgebraicNumber<D>> rlist = rroot.tuple; 939 Interval<D> vr = rlist.get(0).ring.getRoot(); 940 Interval<D> vi = rlist.get(1).ring.getRoot(); 941 logger.info("vr = " + vr + ", vi = " + vi); 942 if (vr.length().isZERO()) { 943 D e = vr.left.factory().parse("1/2"); 944 D m = vr.left; //middle(); 945 vr = new Interval<D>(m.subtract(e), m.sum(e)); 946 logger.info("|vr| == 0: " + vr); 947 } 948 if (vi.length().isZERO()) { 949 D e = vi.left.factory().parse("1/2"); 950 D m = vi.left; //middle(); 951 vi = new Interval<D>(m.subtract(e), m.sum(e)); 952 logger.info("|vi| == 0: " + vi); 953 } 954 Complex<D> sw = new Complex<D>(ccfac, vr.left, vi.left); 955 Complex<D> ne = new Complex<D>(ccfac, vr.right, vi.right); 956 logger.info("sw = " + toString1(sw) + ", ne = " + toString1(ne)); 957 GenPolynomial<Complex<D>> pip2cesw, pip2cene; 958 pip2cesw = PolyUtil.<Complex<D>> evaluateMainRecursive(ucfac, pip2cr, sw); 959 pip2cene = PolyUtil.<Complex<D>> evaluateMainRecursive(ucfac, pip2cr, ne); 960 GenPolynomialRing<D> upfac = I.upolys.get(ix).ring; 961 GenPolynomialRing<Complex<D>> upcfac = new GenPolynomialRing<Complex<D>>(ccfac, upfac); 962 //System.out.println("upfac = " + upfac); 963 //System.out.println("upcfac = " + upcfac); 964 GenPolynomial<Complex<D>> pip2eswc = convertComplexComplex(upcfac, pip2cesw); 965 GenPolynomial<Complex<D>> pip2enec = convertComplexComplex(upcfac, pip2cene); 966 //System.out.println("pip2eswc = " + pip2eswc); 967 //System.out.println("pip2enec = " + pip2enec); 968 for (List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cx : can) { 969 //System.out.println("cxi = " + toString(cx.get(ix))); 970 Complex<edu.jas.application.RealAlgebraicNumber<D>> cax = cx.get(ix); 971 ComplexRing<edu.jas.application.RealAlgebraicNumber<D>> car = cax.ring; 972 edu.jas.application.RealAlgebraicRing<D> rar = (edu.jas.application.RealAlgebraicRing<D>) car.ring; 973 //System.out.println("car = " + car); 974 //System.out.println("rar = " + rar); 975 TermOrder to = new TermOrder(TermOrder.INVLEX); 976 String vvr = rar.algebraic.ring.getVars()[0]; 977 String vvi = rar.algebraic.ring.getVars()[1]; 978 String[] vars = new String[] { vvr, vvi }; 979 GenPolynomialRing<Complex<D>> tfac = new GenPolynomialRing<Complex<D>>(ccfac, to, 980 vars); 981 GenPolynomial<Complex<D>> t = tfac.univariate(1, 1L).sum( 982 tfac.univariate(0, 1L).multiply(ccfac.getIMAG())); 983 //System.out.println("t = " + t); // t = x + i y 984 GenPolynomialRing<D> rtfac = new GenPolynomialRing<D>(cfac, tfac); 985 GenPolynomial<Complex<D>> su; 986 GenPolynomial<D> re, im; 987 su = PolyUtil.<Complex<D>> substituteUnivariate(pip2eswc, t); 988 //su = su.monic(); not here 989 re = PolyUtil.<D> realPartFromComplex(rtfac, su); 990 im = PolyUtil.<D> imaginaryPartFromComplex(rtfac, su); 991 //System.out.println("re = " + re); 992 //System.out.println("im = " + im); 993 edu.jas.application.RealAlgebraicNumber<D> resw, imsw, rene, imne; 994 resw = new edu.jas.application.RealAlgebraicNumber<D>(rar, re); 995 //System.out.println("resw = " + resw); 996 int sswr = resw.signum(); 997 imsw = new edu.jas.application.RealAlgebraicNumber<D>(rar, im); 998 //System.out.println("imsw = " + imsw); 999 int sswi = imsw.signum(); 1000 su = PolyUtil.<Complex<D>> substituteUnivariate(pip2enec, t); 1001 //su = su.monic(); not here 1002 re = PolyUtil.<D> realPartFromComplex(rtfac, su); 1003 im = PolyUtil.<D> imaginaryPartFromComplex(rtfac, su); 1004 //System.out.println("re = " + re); 1005 //System.out.println("im = " + im); 1006 rene = new edu.jas.application.RealAlgebraicNumber<D>(rar, re); 1007 //System.out.println("rene = " + rene); 1008 int sner = rene.signum(); 1009 imne = new edu.jas.application.RealAlgebraicNumber<D>(rar, im); 1010 //System.out.println("imne = " + imne); 1011 int snei = imne.signum(); 1012 //System.out.println("sswr = " + sswr + ", sswi = " + sswi); 1013 //System.out.println("sner = " + sner + ", snei = " + snei); 1014 if ((sswr * sner <= 0 && sswi * snei <= 0)) { // wrong ! 1015 logger.info(" hit, cxi = " + toString(cx.get(ix)) + ", cr = " + toString(cr)); 1016 List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cy; 1017 cy = new ArrayList<Complex<edu.jas.application.RealAlgebraicNumber<D>>>(); 1018 cy.addAll(cx); 1019 cy.add(cr); 1020 cn.add(cy); 1021 } else { 1022 logger.info("no hit, cxi = " + toString(cx.get(ix)) + ", cr = " + toString(cr)); 1023 } 1024 } 1025 } 1026 } 1027 can = cn; 1028 } 1029 IdealWithComplexAlgebraicRoots<D> Ic = new IdealWithComplexAlgebraicRoots<D>(I, can); 1030 return Ic; 1031 } 1032 1033 1034 /** 1035 * Construct complex roots for zero dimensional ideal(G). 1036 * @param I zero dimensional ideal with univariate irreducible polynomials 1037 * and bi-variate polynomials. 1038 * @return complex algebraic roots for ideal(G) <b>Note:</b> not jet 1039 * completed oin all cases. 1040 */ 1041 public static <D extends GcdRingElem<D> & Rational> IdealWithComplexAlgebraicRoots<D> complexAlgebraicRoots( 1042 IdealWithUniv<D> I) { 1043 List<List<Complex<edu.jas.application.RealAlgebraicNumber<D>>>> can; 1044 can = new ArrayList<List<Complex<edu.jas.application.RealAlgebraicNumber<D>>>>(); 1045 if (I == null) { 1046 throw new IllegalArgumentException("null ideal not permitted"); 1047 } 1048 if (I.ideal == null || I.upolys == null) { 1049 throw new IllegalArgumentException("null ideal components not permitted " + I); 1050 } 1051 if (I.ideal.isZERO() || I.upolys.size() == 0) { 1052 return new IdealWithComplexAlgebraicRoots<D>(I, can); 1053 } 1054 GenPolynomialRing<D> fac = I.ideal.list.ring; 1055 if (fac.nvar == 0) { 1056 return new IdealWithComplexAlgebraicRoots<D>(I, can); 1057 } 1058 if (fac.nvar != I.upolys.size()) { 1059 throw new IllegalArgumentException("ideal not zero dimnsional: " + I); 1060 } 1061 // case i == 0: 1062 GenPolynomial<D> p0 = I.upolys.get(0); 1063 GenPolynomial<D> p0p = PolyUtil.<D> selectWithVariable(I.ideal.list.list, fac.nvar - 1); 1064 if (p0p == null) { 1065 throw new RuntimeException("no polynomial found in " + (fac.nvar - 1) + " of " + I.ideal); 1066 } 1067 if (logger.isInfoEnabled()) { 1068 logger.info("p0 = " + p0); 1069 logger.info("p0p = " + p0p); 1070 } 1071 int[] dep0 = p0p.degreeVector().dependencyOnVariables(); 1072 //System.out.println("dep0 = " + Arrays.toString(dep0)); 1073 if (dep0.length != 1) { 1074 throw new RuntimeException("wrong number of variables " + Arrays.toString(dep0)); 1075 } 1076 RingFactory<D> cfac = p0.ring.coFac; 1077 ComplexRing<D> ccfac = new ComplexRing<D>(cfac); 1078 GenPolynomialRing<Complex<D>> facc = new GenPolynomialRing<Complex<D>>(ccfac, p0.ring); 1079 GenPolynomial<Complex<D>> p0c = PolyUtil.<D> complexFromAny(facc, p0); 1080 List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cra; 1081 cra = edu.jas.application.RootFactory.<D> complexAlgebraicNumbersSquarefree(p0c); 1082 logger.info("#roots(p0c) = " + cra.size()); 1083 for (Complex<edu.jas.application.RealAlgebraicNumber<D>> cr : cra) { 1084 List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cl; 1085 cl = new ArrayList<Complex<edu.jas.application.RealAlgebraicNumber<D>>>(); 1086 cl.add(cr); 1087 can.add(cl); 1088 } 1089 if (fac.nvar == 1) { 1090 return new IdealWithComplexAlgebraicRoots<D>(I, can); 1091 } 1092 // case i > 0: 1093 for (int i = 1; i < I.upolys.size(); i++) { 1094 List<List<Complex<edu.jas.application.RealAlgebraicNumber<D>>>> cn; 1095 cn = new ArrayList<List<Complex<edu.jas.application.RealAlgebraicNumber<D>>>>(); 1096 GenPolynomial<D> pi = I.upolys.get(i); 1097 GenPolynomial<D> pip = PolyUtil.selectWithVariable(I.ideal.list.list, fac.nvar - 1 - i); 1098 if (pip == null) { 1099 throw new RuntimeException("no polynomial found in " + (fac.nvar - 1 - i) + " of " + I.ideal); 1100 } 1101 if (logger.isInfoEnabled()) { 1102 logger.info("pi(" + i + ") = " + pi); 1103 logger.info("pip = " + pip); 1104 } 1105 facc = new GenPolynomialRing<Complex<D>>(ccfac, pi.ring); 1106 GenPolynomial<Complex<D>> pic = PolyUtil.<D> complexFromAny(facc, pi); 1107 int[] depi = pip.degreeVector().dependencyOnVariables(); 1108 //System.out.println("depi = " + Arrays.toString(depi)); 1109 if (depi.length < 1 || depi.length > 2) { 1110 throw new RuntimeException("wrong number of variables " + Arrays.toString(depi) + " for " 1111 + pip); 1112 } 1113 cra = edu.jas.application.RootFactory.<D> complexAlgebraicNumbersSquarefree(pic); 1114 logger.info("#roots(pic) = " + cra.size()); 1115 if (depi.length == 1) { 1116 // all combinations are roots of the ideal I 1117 for (Complex<edu.jas.application.RealAlgebraicNumber<D>> cr : cra) { 1118 //System.out.println("cr.ring = " + cr.ring); 1119 for (List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cx : can) { 1120 //System.out.println("cx = " + cx); 1121 List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cy; 1122 cy = new ArrayList<Complex<edu.jas.application.RealAlgebraicNumber<D>>>(); 1123 cy.addAll(cx); 1124 cy.add(cr); 1125 cn.add(cy); 1126 } 1127 } 1128 } else { // depi.length == 2 1129 // select roots of the ideal I 1130 GenPolynomial<D> pip2 = PolyUtil.<D> removeUnusedUpperVariables(pip); 1131 pip2 = PolyUtil.<D> removeUnusedLowerVariables(pip2); 1132 pip2 = PolyUtil.<D> removeUnusedMiddleVariables(pip2); 1133 GenPolynomialRing<GenPolynomial<D>> rfac = pip2.ring.recursive(1); 1134 GenPolynomialRing<D> ufac = pip2.ring.contract(1); 1135 GenPolynomialRing<Complex<D>> ucfac = new GenPolynomialRing<Complex<D>>(ccfac, ufac); 1136 GenPolynomialRing<Complex<D>> c2fac = new GenPolynomialRing<Complex<D>>(ccfac, pip2.ring); 1137 GenPolynomial<Complex<D>> pip2c = PolyUtil.<D> complexFromAny(c2fac, pip2); 1138 GenPolynomialRing<GenPolynomial<Complex<D>>> rcfac; 1139 rcfac = new GenPolynomialRing<GenPolynomial<Complex<D>>>(ucfac, rfac); 1140 GenPolynomial<GenPolynomial<Complex<D>>> pip2cr = PolyUtil.<Complex<D>> recursive(rcfac, 1141 pip2c); 1142 //System.out.println("pip2cr = " + pip2cr); 1143 int ix = fac.nvar - 1 - depi[depi.length - 1]; 1144 //System.out.println("ix = " + ix); 1145 for (Complex<edu.jas.application.RealAlgebraicNumber<D>> cr : cra) { 1146 //System.out.println("cr = " + toString(cr)); 1147 edu.jas.application.RealAlgebraicRing<D> cring = (edu.jas.application.RealAlgebraicRing<D>) cr.ring.ring; 1148 RealRootTuple<D> rroot = cring.getRoot(); 1149 List<RealAlgebraicNumber<D>> rlist = rroot.tuple; 1150 //System.out.println("rlist = " + rlist); 1151 Interval<D> vr = rlist.get(0).ring.getRoot(); 1152 Interval<D> vi = rlist.get(1).ring.getRoot(); 1153 //logger.info("vr = " + vr + ", vi = " + vi); 1154 edu.jas.application.RealAlgebraicNumber<D> vrl, vil, vrr, vir; 1155 vrl = new edu.jas.application.RealAlgebraicNumber<D>(cring, vr.left); 1156 vil = new edu.jas.application.RealAlgebraicNumber<D>(cring, vi.left); 1157 vrr = new edu.jas.application.RealAlgebraicNumber<D>(cring, vr.right); 1158 vir = new edu.jas.application.RealAlgebraicNumber<D>(cring, vi.right); 1159 ComplexRing<edu.jas.application.RealAlgebraicNumber<D>> crr; 1160 crr = new ComplexRing<edu.jas.application.RealAlgebraicNumber<D>>(cring); 1161 Complex<edu.jas.application.RealAlgebraicNumber<D>> csw, cne; 1162 csw = new Complex<edu.jas.application.RealAlgebraicNumber<D>>(crr, vrl, vil); 1163 cne = new Complex<edu.jas.application.RealAlgebraicNumber<D>>(crr, vrr, vir); 1164 //logger.info("csw = " + toString(csw) + ", cne = " + toString(cne)); 1165 Rectangle<edu.jas.application.RealAlgebraicNumber<D>> rec; 1166 rec = new Rectangle<edu.jas.application.RealAlgebraicNumber<D>>(csw, cne); 1167 //System.out.println("rec = " + rec); 1168 for (List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cx : can) { 1169 Complex<edu.jas.application.RealAlgebraicNumber<D>> cax = cx.get(ix); 1170 //System.out.println("cax = " + toString(cax)); 1171 ComplexRing<edu.jas.application.RealAlgebraicNumber<D>> car = cax.ring; 1172 //System.out.println("car = " + car); 1173 GenPolynomialRing<Complex<edu.jas.application.RealAlgebraicNumber<D>>> pcrfac; 1174 pcrfac = new GenPolynomialRing<Complex<edu.jas.application.RealAlgebraicNumber<D>>>( 1175 car, rcfac); 1176 GenPolynomial<Complex<edu.jas.application.RealAlgebraicNumber<D>>> pcr; 1177 pcr = evaluateToComplexRealCoefficients(pcrfac, pip2cr, cax); 1178 //System.out.println("pcr = " + pcr); 1179 ComplexRoots<edu.jas.application.RealAlgebraicNumber<D>> rengine; 1180 rengine = new ComplexRootsSturm<edu.jas.application.RealAlgebraicNumber<D>>(car); 1181 long nr = 0; 1182 try { 1183 nr = rengine.complexRootCount(rec, pcr); 1184 //logger.info("rootCount = " + nr); 1185 } catch (InvalidBoundaryException e) { 1186 e.printStackTrace(); 1187 } 1188 if (nr == 1) { // one root 1189 logger.info(" hit, cxi = " + toString(cx.get(ix)) + ", cr = " + toString(cr)); 1190 List<Complex<edu.jas.application.RealAlgebraicNumber<D>>> cy; 1191 cy = new ArrayList<Complex<edu.jas.application.RealAlgebraicNumber<D>>>(); 1192 cy.addAll(cx); 1193 cy.add(cr); 1194 cn.add(cy); 1195 } else if (nr > 1) { 1196 logger.error("to many roots, cxi = " + toString(cx.get(ix)) + ", cr = " 1197 + toString(cr)); 1198 } else { // no root 1199 logger.info("no hit, cxi = " + toString(cx.get(ix)) + ", cr = " + toString(cr)); 1200 } 1201 } 1202 } 1203 } 1204 can = cn; 1205 } 1206 IdealWithComplexAlgebraicRoots<D> Ic = new IdealWithComplexAlgebraicRoots<D>(I, can); 1207 return Ic; 1208 } 1209 1210 1211 /** 1212 * String representation of a deximal approximation of a complex number. 1213 * @param c compelx number. 1214 * @return String representation of c 1215 */ 1216 public static <D extends GcdRingElem<D> & Rational> String toString( 1217 Complex<edu.jas.application.RealAlgebraicNumber<D>> c) { 1218 edu.jas.application.RealAlgebraicNumber<D> re = c.getRe(); 1219 edu.jas.application.RealAlgebraicNumber<D> im = c.getIm(); 1220 String s = re.decimalMagnitude().toString(); 1221 if (!im.isZERO()) { 1222 s = s + "i" + im.decimalMagnitude(); 1223 } 1224 return s; 1225 } 1226 1227 1228 /** 1229 * String representation of a deximal approximation of a complex number. 1230 * @param c compelx number. 1231 * @return String representation of c 1232 */ 1233 public static <D extends GcdRingElem<D> & Rational> String toString1(Complex<D> c) { 1234 D re = c.getRe(); 1235 D im = c.getIm(); 1236 String s = new BigDecimal(re.getRational()).toString(); 1237 if (!im.isZERO()) { 1238 s = s + "i" + new BigDecimal(im.getRational()); 1239 } 1240 return s; 1241 } 1242 1243 1244 /** 1245 * Construct complex roots for zero dimensional ideal(G). 1246 * @param I list of zero dimensional ideal with univariate irreducible 1247 * polynomials and bi-variate polynomials. 1248 * @return list of complex algebraic roots for ideal(G) 1249 */ 1250 public static <D extends GcdRingElem<D> & Rational> List<IdealWithComplexAlgebraicRoots<D>> complexAlgebraicRoots( 1251 List<IdealWithUniv<D>> I) { 1252 List<IdealWithComplexAlgebraicRoots<D>> lic = new ArrayList<IdealWithComplexAlgebraicRoots<D>>(); 1253 for (IdealWithUniv<D> iu : I) { 1254 IdealWithComplexAlgebraicRoots<D> iuc = PolyUtilApp.<D> complexAlgebraicRoots(iu); 1255 //System.out.println("iuc = " + iuc); 1256 lic.add(iuc); 1257 } 1258 return lic; 1259 } 1260 1261 1262 /** 1263 * Construct exact set of complex roots for zero dimensional ideal(G). 1264 * @param I zero dimensional ideal. 1265 * @return list of coordinates of complex roots for ideal(G) 1266 */ 1267 public static <D extends GcdRingElem<D> & Rational> List<IdealWithComplexAlgebraicRoots<D>> complexAlgebraicRoots( 1268 Ideal<D> I) { 1269 List<IdealWithUniv<D>> Ir = I.zeroDimRootDecomposition(); 1270 //System.out.println("Ir = " + Ir); 1271 List<IdealWithComplexAlgebraicRoots<D>> roots = PolyUtilApp.<D> complexAlgebraicRoots(Ir); 1272 return roots; 1273 } 1274 1275 1276 /* 1277 * Convert to a polynomial in given ring. 1278 * @param fac result polynomial ring. 1279 * @param p polynomial. 1280 * @return polynomial in ring fac <b>Note: </b> if p can not be represented 1281 * in fac then the results are unpredictable. 1282 */ 1283 static <C extends RingElem<C>> GenPolynomial<C> convert(GenPolynomialRing<C> fac, GenPolynomial<C> p) { 1284 if (fac.equals(p.factory())) { 1285 return p; 1286 } 1287 GenPolynomial<C> q = fac.parse(p.toString()); 1288 if (!q.toString().equals(p.toString())) { 1289 throw new RuntimeException("convert(" + p + ") = " + q); 1290 } 1291 return q; 1292 } 1293 1294 1295 /* 1296 * Convert to a polynomial in given ring. 1297 * @param fac result polynomial ring. 1298 * @param p polynomial. 1299 * @return polynomial in ring fac <b>Note: </b> if p can not be represented 1300 * in fac then the results are unpredictable. 1301 */ 1302 static <C extends RingElem<C>> GenPolynomial<Complex<C>> convertComplex( 1303 GenPolynomialRing<Complex<C>> fac, GenPolynomial<C> p) { 1304 GenPolynomial<Complex<C>> q = fac.parse(p.toString()); 1305 if (!q.toString().equals(p.toString())) { 1306 throw new RuntimeException("convert(" + p + ") = " + q); 1307 } 1308 return q; 1309 } 1310 1311 1312 /* 1313 * Convert to a polynomial in given ring. 1314 * @param fac result polynomial ring. 1315 * @param p complex polynomial. 1316 * @return polynomial in ring fac <b>Note: </b> if p can not be represented 1317 * in fac then the results are unpredictable. 1318 */ 1319 static <C extends RingElem<C>> GenPolynomial<Complex<C>> convertComplexComplex( 1320 GenPolynomialRing<Complex<C>> fac, GenPolynomial<Complex<C>> p) { 1321 if (fac.equals(p.factory())) { 1322 return p; 1323 } 1324 GenPolynomial<Complex<C>> q = fac.parse(p.toString()); 1325 if (!q.toString().equals(p.toString())) { 1326 throw new RuntimeException("convert(" + p + ") = " + q); 1327 } 1328 return q; 1329 } 1330 1331 1332 /** 1333 * Construct exact set of real roots for zero dimensional ideal(G). 1334 * @param I zero dimensional ideal. 1335 * @return list of coordinates of real roots for ideal(G) 1336 */ 1337 public static <D extends GcdRingElem<D> & Rational> List<IdealWithRealAlgebraicRoots<D>> realAlgebraicRoots( 1338 Ideal<D> I) { 1339 List<IdealWithUniv<D>> Ir = I.zeroDimRootDecomposition(); 1340 //System.out.println("Ir = " + Ir); 1341 List<IdealWithRealAlgebraicRoots<D>> roots = PolyUtilApp.<D> realAlgebraicRoots(Ir); 1342 return roots; 1343 } 1344 1345 1346 /** 1347 * Construct primitive element for double field extension. 1348 * @param a algebraic number ring with squarefree monic minimal polynomial 1349 * @param b algebraic number ring with squarefree monic minimal polynomial 1350 * @return primitive element container with algebraic number ring c, with 1351 * Q(c) = Q(a,b) 1352 */ 1353 public static <C extends GcdRingElem<C>> PrimitiveElement<C> primitiveElement(AlgebraicNumberRing<C> a, 1354 AlgebraicNumberRing<C> b) { 1355 GenPolynomial<C> ap = a.modul; 1356 GenPolynomial<C> bp = b.modul; 1357 1358 // setup bivariate polynomial ring 1359 String[] cv = new String[2]; 1360 cv[0] = ap.ring.getVars()[0]; 1361 cv[1] = bp.ring.getVars()[0]; 1362 TermOrder to = new TermOrder(TermOrder.INVLEX); 1363 GenPolynomialRing<C> cfac = new GenPolynomialRing<C>(ap.ring.coFac, 2, to, cv); 1364 GenPolynomial<C> as = ap.extendUnivariate(cfac, 0); 1365 GenPolynomial<C> bs = bp.extendUnivariate(cfac, 1); 1366 List<GenPolynomial<C>> L = new ArrayList<GenPolynomial<C>>(2); 1367 L.add(as); 1368 L.add(bs); 1369 List<GenPolynomial<C>> Op = new ArrayList<GenPolynomial<C>>(); 1370 1371 Ideal<C> id = new Ideal<C>(cfac, L); 1372 //System.out.println("id = " + id); 1373 IdealWithUniv<C> iu = id.normalPositionFor(0, 1, Op); 1374 //System.out.println("iu = " + iu); 1375 1376 // extract result polynomials 1377 List<GenPolynomial<C>> Np = iu.ideal.getList(); 1378 //System.out.println("Np = " + Np); 1379 as = PolyUtil.<C> selectWithVariable(Np, 1); 1380 bs = PolyUtil.<C> selectWithVariable(Np, 0); 1381 GenPolynomial<C> cs = PolyUtil.<C> selectWithVariable(Np, 2); 1382 //System.out.println("as = " + as); 1383 //System.out.println("bs = " + bs); 1384 //System.out.println("cs = " + cs); 1385 String[] ev = new String[] { cs.ring.getVars()[0] }; 1386 GenPolynomialRing<C> efac = new GenPolynomialRing<C>(ap.ring.coFac, 1, to, ev); 1387 //System.out.println("efac = " + efac); 1388 cs = cs.contractCoeff(efac); 1389 //System.out.println("cs = " + cs); 1390 as = as.reductum().contractCoeff(efac); 1391 as = as.negate(); 1392 //System.out.println("as = " + as); 1393 bs = bs.reductum().contractCoeff(efac); 1394 bs = bs.negate(); 1395 //System.out.println("bs = " + bs); 1396 AlgebraicNumberRing<C> c = new AlgebraicNumberRing<C>(cs); 1397 AlgebraicNumber<C> ab = new AlgebraicNumber<C>(c, as); 1398 AlgebraicNumber<C> bb = new AlgebraicNumber<C>(c, bs); 1399 PrimitiveElement<C> pe = new PrimitiveElement<C>(c, ab, bb, a, b); 1400 if (logger.isInfoEnabled()) { 1401 logger.info("primitive element = " + c); 1402 } 1403 return pe; 1404 } 1405 1406 1407 /** 1408 * Convert to primitive element ring. 1409 * @param cfac primitive element ring. 1410 * @param A algebraic number representing the generating element of a in the 1411 * new ring. 1412 * @param a algebraic number to convert. 1413 * @return a converted to the primitive element ring 1414 */ 1415 public static <C extends GcdRingElem<C>> AlgebraicNumber<C> convertToPrimitiveElem( 1416 AlgebraicNumberRing<C> cfac, AlgebraicNumber<C> A, AlgebraicNumber<C> a) { 1417 GenPolynomialRing<C> aufac = a.ring.ring; 1418 GenPolynomialRing<AlgebraicNumber<C>> ar = new GenPolynomialRing<AlgebraicNumber<C>>(cfac, aufac); 1419 GenPolynomial<AlgebraicNumber<C>> aps = PolyUtil.<C> convertToAlgebraicCoefficients(ar, a.val); 1420 AlgebraicNumber<C> ac = PolyUtil.<AlgebraicNumber<C>> evaluateMain(cfac, aps, A); 1421 return ac; 1422 } 1423 1424 1425 /** 1426 * Convert coefficients to primitive element ring. 1427 * @param cfac primitive element ring. 1428 * @param A algebraic number representing the generating element of a in the 1429 * new ring. 1430 * @param a polynomial with coefficients algebraic number to convert. 1431 * @return a with coefficients converted to the primitive element ring 1432 */ 1433 public static <C extends GcdRingElem<C>> GenPolynomial<AlgebraicNumber<C>> convertToPrimitiveElem( 1434 AlgebraicNumberRing<C> cfac, AlgebraicNumber<C> A, GenPolynomial<AlgebraicNumber<C>> a) { 1435 GenPolynomialRing<AlgebraicNumber<C>> cr = new GenPolynomialRing<AlgebraicNumber<C>>(cfac, a.ring); 1436 return PolyUtil.<AlgebraicNumber<C>, AlgebraicNumber<C>> map(cr, a, new CoeffConvertAlg<C>(cfac, A)); 1437 } 1438 1439 1440 /** 1441 * Convert to primitive element ring. 1442 * @param cfac primitive element ring. 1443 * @param A algebraic number representing the generating element of a in the 1444 * new ring. 1445 * @param a recursive algebraic number to convert. 1446 * @return a converted to the primitive element ring 1447 */ 1448 public static <C extends GcdRingElem<C>> AlgebraicNumber<C> convertToPrimitiveElem( 1449 AlgebraicNumberRing<C> cfac, AlgebraicNumber<C> A, AlgebraicNumber<C> B, 1450 AlgebraicNumber<AlgebraicNumber<C>> a) { 1451 GenPolynomial<AlgebraicNumber<C>> aps = PolyUtilApp.<C> convertToPrimitiveElem(cfac, A, a.val); 1452 AlgebraicNumber<C> ac = PolyUtil.<AlgebraicNumber<C>> evaluateMain(cfac, aps, B); 1453 return ac; 1454 } 1455 1456 1457 /** 1458 * Construct primitive element for double field extension. 1459 * @param b algebraic number ring with squarefree monic minimal polynomial 1460 * over Q(a) 1461 * @return primitive element container with algebraic number ring c, with 1462 * Q(c) = Q(a)(b) 1463 */ 1464 public static <C extends GcdRingElem<C>> PrimitiveElement<C> primitiveElement( 1465 AlgebraicNumberRing<AlgebraicNumber<C>> b) { 1466 GenPolynomial<AlgebraicNumber<C>> bp = b.modul; 1467 AlgebraicNumberRing<C> a = (AlgebraicNumberRing<C>) b.ring.coFac; 1468 GenPolynomial<C> ap = a.modul; 1469 1470 // setup bivariate polynomial ring 1471 String[] cv = new String[2]; 1472 cv[0] = ap.ring.getVars()[0]; 1473 cv[1] = bp.ring.getVars()[0]; 1474 TermOrder to = new TermOrder(TermOrder.INVLEX); 1475 GenPolynomialRing<C> cfac = new GenPolynomialRing<C>(ap.ring.coFac, 2, to, cv); 1476 GenPolynomialRing<GenPolynomial<C>> rfac = new GenPolynomialRing<GenPolynomial<C>>(a.ring, 1, 1477 bp.ring.getVars()); 1478 GenPolynomial<C> as = ap.extendUnivariate(cfac, 0); 1479 GenPolynomial<GenPolynomial<C>> bss = PolyUtil.<C> fromAlgebraicCoefficients(rfac, bp); 1480 GenPolynomial<C> bs = PolyUtil.<C> distribute(cfac, bss); 1481 List<GenPolynomial<C>> L = new ArrayList<GenPolynomial<C>>(2); 1482 L.add(as); 1483 L.add(bs); 1484 List<GenPolynomial<C>> Op = new ArrayList<GenPolynomial<C>>(); 1485 1486 Ideal<C> id = new Ideal<C>(cfac, L); 1487 //System.out.println("id = " + id); 1488 IdealWithUniv<C> iu = id.normalPositionFor(0, 1, Op); 1489 //System.out.println("iu = " + iu); 1490 1491 // extract result polynomials 1492 List<GenPolynomial<C>> Np = iu.ideal.getList(); 1493 as = PolyUtil.<C> selectWithVariable(Np, 1); 1494 bs = PolyUtil.<C> selectWithVariable(Np, 0); 1495 GenPolynomial<C> cs = PolyUtil.<C> selectWithVariable(Np, 2); 1496 //System.out.println("as = " + as); 1497 //System.out.println("bs = " + bs); 1498 //System.out.println("cs = " + cs); 1499 String[] ev = new String[] { cs.ring.getVars()[0] }; 1500 GenPolynomialRing<C> efac = new GenPolynomialRing<C>(ap.ring.coFac, 1, to, ev); 1501 // System.out.println("efac = " + efac); 1502 cs = cs.contractCoeff(efac); 1503 // System.out.println("cs = " + cs); 1504 as = as.reductum().contractCoeff(efac); 1505 as = as.negate(); 1506 // System.out.println("as = " + as); 1507 bs = bs.reductum().contractCoeff(efac); 1508 bs = bs.negate(); 1509 //System.out.println("bs = " + bs); 1510 AlgebraicNumberRing<C> c = new AlgebraicNumberRing<C>(cs); 1511 AlgebraicNumber<C> ab = new AlgebraicNumber<C>(c, as); 1512 AlgebraicNumber<C> bb = new AlgebraicNumber<C>(c, bs); 1513 PrimitiveElement<C> pe = new PrimitiveElement<C>(c, ab, bb); // missing ,a,b); 1514 if (logger.isInfoEnabled()) { 1515 logger.info("primitive element = " + pe); 1516 } 1517 return pe; 1518 } 1519 1520 1521 /** 1522 * Convert to primitive element ring. 1523 * @param cfac primitive element ring. 1524 * @param A algebraic number representing the generating element of a in the 1525 * new ring. 1526 * @param a polynomial with recursive algebraic number coefficients to 1527 * convert. 1528 * @return a converted to the primitive element ring 1529 */ 1530 public static <C extends GcdRingElem<C>> GenPolynomial<AlgebraicNumber<C>> convertToPrimitiveElem( 1531 AlgebraicNumberRing<C> cfac, AlgebraicNumber<C> A, AlgebraicNumber<C> B, 1532 GenPolynomial<AlgebraicNumber<AlgebraicNumber<C>>> a) { 1533 GenPolynomialRing<AlgebraicNumber<C>> cr = new GenPolynomialRing<AlgebraicNumber<C>>(cfac, a.ring); 1534 return PolyUtil.<AlgebraicNumber<AlgebraicNumber<C>>, AlgebraicNumber<C>> map(cr, a, 1535 new CoeffRecConvertAlg<C>(cfac, A, B)); 1536 } 1537 1538 1539 /** 1540 * Convert to RealAlgebraicNumber coefficients. Represent as polynomial with 1541 * RealAlgebraicNumber<C> coefficients from package 1542 * 1543 * <pre> 1544 * edu.jas.root 1545 * </pre> 1546 * 1547 * . 1548 * @param afac result polynomial factory. 1549 * @param A polynomial with RealAlgebraicNumber<C> coefficients to be 1550 * converted. 1551 * @return polynomial with RealAlgebraicNumber<C> coefficients. 1552 */ 1553 public static <C extends GcdRingElem<C> & Rational> GenPolynomial<edu.jas.root.RealAlgebraicNumber<C>> realAlgFromRealCoefficients( 1554 GenPolynomialRing<edu.jas.root.RealAlgebraicNumber<C>> afac, 1555 GenPolynomial<edu.jas.application.RealAlgebraicNumber<C>> A) { 1556 edu.jas.root.RealAlgebraicRing<C> cfac = (edu.jas.root.RealAlgebraicRing<C>) afac.coFac; 1557 return PolyUtil.<edu.jas.application.RealAlgebraicNumber<C>, edu.jas.root.RealAlgebraicNumber<C>> map( 1558 afac, A, new ReAlgFromRealCoeff<C>(cfac)); 1559 } 1560 1561 1562 /** 1563 * Convert to RealAlgebraicNumber coefficients. Represent as polynomial with 1564 * RealAlgebraicNumber<C> coefficients from package 1565 * 1566 * <pre> 1567 * edu.jas.application 1568 * </pre> 1569 * 1570 * . 1571 * @param rfac result polynomial factory. 1572 * @param A polynomial with RealAlgebraicNumber<C> coefficients to be 1573 * converted. 1574 * @return polynomial with RealAlgebraicNumber<C> coefficients. 1575 */ 1576 public static <C extends GcdRingElem<C> & Rational> GenPolynomial<edu.jas.application.RealAlgebraicNumber<C>> realFromRealAlgCoefficients( 1577 GenPolynomialRing<edu.jas.application.RealAlgebraicNumber<C>> rfac, 1578 GenPolynomial<edu.jas.root.RealAlgebraicNumber<C>> A) { 1579 edu.jas.application.RealAlgebraicRing<C> cfac = (edu.jas.application.RealAlgebraicRing<C>) rfac.coFac; 1580 return PolyUtil.<edu.jas.root.RealAlgebraicNumber<C>, edu.jas.application.RealAlgebraicNumber<C>> map( 1581 rfac, A, new RealFromReAlgCoeff<C>(cfac)); 1582 } 1583 1584 1585 /** 1586 * Convert to Complex<RealAlgebraicNumber> coefficients. Represent as 1587 * polynomial with Complex<RealAlgebraicNumber> coefficients, C is 1588 * e.g. BigRational. 1589 * @param pfac result polynomial factory. 1590 * @param A polynomial with Complex coefficients to be converted. 1591 * @return polynomial with Complex<RealAlgebraicNumber> coefficients. 1592 */ 1593 public static <C extends GcdRingElem<C> & Rational> GenPolynomial<Complex<edu.jas.application.RealAlgebraicNumber<C>>> convertToComplexRealCoefficients( 1594 GenPolynomialRing<Complex<edu.jas.application.RealAlgebraicNumber<C>>> pfac, 1595 GenPolynomial<Complex<C>> A) { 1596 ComplexRing<edu.jas.application.RealAlgebraicNumber<C>> afac; 1597 afac = (ComplexRing<edu.jas.application.RealAlgebraicNumber<C>>) pfac.coFac; 1598 return PolyUtil.<Complex<C>, Complex<edu.jas.application.RealAlgebraicNumber<C>>> map(pfac, A, 1599 new CoeffToComplexReal<C>(afac)); 1600 } 1601 1602 1603 /** 1604 * Evaluate to Complex<RealAlgebraicNumber> coefficients. Represent as 1605 * polynomial with Complex<RealAlgebraicNumber> coefficients, C is 1606 * e.g. BigRational. 1607 * @param pfac result polynomial factory. 1608 * @param A = A(x,Y) a recursive polynomial with 1609 * GenPolynomial<Complex> coefficients to be converted. 1610 * @param r Complex<RealAlgebraicNumber> to be evaluated at. 1611 * @return A(r,Y), a polynomial with Complex<RealAlgebraicNumber> 1612 * coefficients. 1613 */ 1614 public static <C extends GcdRingElem<C> & Rational> GenPolynomial<Complex<edu.jas.application.RealAlgebraicNumber<C>>> evaluateToComplexRealCoefficients( 1615 GenPolynomialRing<Complex<edu.jas.application.RealAlgebraicNumber<C>>> pfac, 1616 GenPolynomial<GenPolynomial<Complex<C>>> A, 1617 Complex<edu.jas.application.RealAlgebraicNumber<C>> r) { 1618 return PolyUtil.<GenPolynomial<Complex<C>>, Complex<edu.jas.application.RealAlgebraicNumber<C>>> map( 1619 pfac, A, new EvaluateToComplexReal<C>(pfac, r)); 1620 } 1621 1622 1623} 1624 1625 1626/** 1627 * Coefficient to convert algebriac functor. 1628 */ 1629class CoeffConvertAlg<C extends GcdRingElem<C>> implements 1630 UnaryFunctor<AlgebraicNumber<C>, AlgebraicNumber<C>> { 1631 1632 1633 final protected AlgebraicNumberRing<C> afac; 1634 1635 1636 final protected AlgebraicNumber<C> A; 1637 1638 1639 public CoeffConvertAlg(AlgebraicNumberRing<C> fac, AlgebraicNumber<C> a) { 1640 if (fac == null || a == null) { 1641 throw new IllegalArgumentException("fac and a must not be null"); 1642 } 1643 afac = fac; 1644 A = a; 1645 } 1646 1647 1648 public AlgebraicNumber<C> eval(AlgebraicNumber<C> c) { 1649 if (c == null) { 1650 return afac.getZERO(); 1651 } 1652 return PolyUtilApp.<C> convertToPrimitiveElem(afac, A, c); 1653 } 1654} 1655 1656 1657/** 1658 * Coefficient recursive to convert algebriac functor. 1659 */ 1660class CoeffRecConvertAlg<C extends GcdRingElem<C>> implements 1661 UnaryFunctor<AlgebraicNumber<AlgebraicNumber<C>>, AlgebraicNumber<C>> { 1662 1663 1664 final protected AlgebraicNumberRing<C> afac; 1665 1666 1667 final protected AlgebraicNumber<C> A; 1668 1669 1670 final protected AlgebraicNumber<C> B; 1671 1672 1673 public CoeffRecConvertAlg(AlgebraicNumberRing<C> fac, AlgebraicNumber<C> a, AlgebraicNumber<C> b) { 1674 if (fac == null || a == null || b == null) { 1675 throw new IllegalArgumentException("fac, a and b must not be null"); 1676 } 1677 afac = fac; 1678 A = a; 1679 B = b; 1680 } 1681 1682 1683 public AlgebraicNumber<C> eval(AlgebraicNumber<AlgebraicNumber<C>> c) { 1684 if (c == null) { 1685 return afac.getZERO(); 1686 } 1687 return PolyUtilApp.<C> convertToPrimitiveElem(afac, A, B, c); 1688 } 1689} 1690 1691 1692/** 1693 * Coefficient to real algebriac from real algebraic functor. 1694 */ 1695class ReAlgFromRealCoeff<C extends GcdRingElem<C> & Rational> implements 1696 UnaryFunctor<edu.jas.application.RealAlgebraicNumber<C>, edu.jas.root.RealAlgebraicNumber<C>> { 1697 1698 1699 final protected edu.jas.root.RealAlgebraicRing<C> afac; 1700 1701 1702 public ReAlgFromRealCoeff(edu.jas.root.RealAlgebraicRing<C> fac) { 1703 if (fac == null) { 1704 throw new IllegalArgumentException("fac must not be null"); 1705 } 1706 afac = fac; 1707 } 1708 1709 1710 public edu.jas.root.RealAlgebraicNumber<C> eval(edu.jas.application.RealAlgebraicNumber<C> c) { 1711 if (c == null) { 1712 return afac.getZERO(); 1713 } 1714 return (edu.jas.root.RealAlgebraicNumber<C>) (Object) c.number; // force ignore recursion 1715 } 1716} 1717 1718 1719/** 1720 * Coefficient to real algebriac from algebraic functor. 1721 */ 1722class RealFromReAlgCoeff<C extends GcdRingElem<C> & Rational> implements 1723 UnaryFunctor<edu.jas.root.RealAlgebraicNumber<C>, edu.jas.application.RealAlgebraicNumber<C>> { 1724 1725 1726 final protected edu.jas.application.RealAlgebraicRing<C> rfac; 1727 1728 1729 public RealFromReAlgCoeff(edu.jas.application.RealAlgebraicRing<C> fac) { 1730 if (fac == null) { 1731 throw new IllegalArgumentException("fac must not be null"); 1732 } 1733 rfac = fac; 1734 } 1735 1736 1737 public edu.jas.application.RealAlgebraicNumber<C> eval(edu.jas.root.RealAlgebraicNumber<C> c) { 1738 if (c == null) { 1739 return rfac.getZERO(); 1740 } 1741 edu.jas.root.RealAlgebraicNumber<edu.jas.root.RealAlgebraicNumber<C>> rrc = (edu.jas.root.RealAlgebraicNumber<edu.jas.root.RealAlgebraicNumber<C>>) (Object) c; // force resurrect recursion 1742 return new edu.jas.application.RealAlgebraicNumber<C>(rfac, rrc); 1743 } 1744} 1745 1746 1747/** 1748 * Coefficient to complex real algebriac functor. 1749 */ 1750class CoeffToComplexReal<C extends GcdRingElem<C> & Rational> implements 1751 UnaryFunctor<Complex<C>, Complex<edu.jas.application.RealAlgebraicNumber<C>>> { 1752 1753 1754 final protected ComplexRing<edu.jas.application.RealAlgebraicNumber<C>> cfac; 1755 1756 1757 final edu.jas.application.RealAlgebraicRing<C> afac; 1758 1759 1760 final GenPolynomialRing<C> pfac; 1761 1762 1763 public CoeffToComplexReal(ComplexRing<edu.jas.application.RealAlgebraicNumber<C>> fac) { 1764 if (fac == null) { 1765 throw new IllegalArgumentException("fac must not be null"); 1766 } 1767 cfac = fac; 1768 afac = (edu.jas.application.RealAlgebraicRing<C>) cfac.ring; 1769 pfac = afac.univs.ideal.getRing(); 1770 } 1771 1772 1773 public Complex<edu.jas.application.RealAlgebraicNumber<C>> eval(Complex<C> c) { 1774 if (c == null) { 1775 return cfac.getZERO(); 1776 } 1777 GenPolynomial<C> pr, pi; 1778 pr = new GenPolynomial<C>(pfac, c.getRe()); 1779 pi = new GenPolynomial<C>(pfac, c.getIm()); 1780 //System.out.println("pr = " + pr); 1781 //System.out.println("pi = " + pi); 1782 edu.jas.application.RealAlgebraicNumber<C> re, im; 1783 re = new edu.jas.application.RealAlgebraicNumber<C>(afac, pr); 1784 im = new edu.jas.application.RealAlgebraicNumber<C>(afac, pi); 1785 //System.out.println("re = " + re); 1786 //System.out.println("im = " + im); 1787 return new Complex<edu.jas.application.RealAlgebraicNumber<C>>(cfac, re, im); 1788 } 1789} 1790 1791 1792/** 1793 * Polynomial coefficient to complex real algebriac evaluation functor. 1794 */ 1795class EvaluateToComplexReal<C extends GcdRingElem<C> & Rational> implements 1796 UnaryFunctor<GenPolynomial<Complex<C>>, Complex<edu.jas.application.RealAlgebraicNumber<C>>> { 1797 1798 1799 final protected GenPolynomialRing<Complex<edu.jas.application.RealAlgebraicNumber<C>>> pfac; 1800 1801 1802 final protected ComplexRing<edu.jas.application.RealAlgebraicNumber<C>> cfac; 1803 1804 1805 final protected Complex<edu.jas.application.RealAlgebraicNumber<C>> root; 1806 1807 1808 public EvaluateToComplexReal(GenPolynomialRing<Complex<edu.jas.application.RealAlgebraicNumber<C>>> fac, 1809 Complex<edu.jas.application.RealAlgebraicNumber<C>> r) { 1810 if (fac == null) { 1811 throw new IllegalArgumentException("fac must not be null"); 1812 } 1813 if (r == null) { 1814 throw new IllegalArgumentException("r must not be null"); 1815 } 1816 pfac = fac; 1817 cfac = (ComplexRing<edu.jas.application.RealAlgebraicNumber<C>>) fac.coFac; 1818 root = r; 1819 //System.out.println("cfac = " + cfac); 1820 //System.out.println("root = " + root); 1821 } 1822 1823 1824 public Complex<edu.jas.application.RealAlgebraicNumber<C>> eval(GenPolynomial<Complex<C>> c) { 1825 if (c == null) { 1826 return cfac.getZERO(); 1827 } 1828 //System.out.println("c = " + c); 1829 GenPolynomial<Complex<edu.jas.application.RealAlgebraicNumber<C>>> cp; 1830 cp = PolyUtilApp.<C> convertToComplexRealCoefficients(pfac, c); 1831 Complex<edu.jas.application.RealAlgebraicNumber<C>> cr; 1832 cr = PolyUtil.<Complex<edu.jas.application.RealAlgebraicNumber<C>>> evaluateMain(cfac, cp, root); 1833 return cr; 1834 } 1835}