001/* 002 * $Id: MultiVarPowerSeries.java 4125 2012-08-19 19:05:22Z kredel $ 003 */ 004 005package edu.jas.ps; 006 007 008import java.util.BitSet; 009import java.util.HashMap; 010import java.util.HashSet; 011import java.util.Iterator; 012import java.util.List; 013import java.util.Map; 014import java.util.Set; 015import java.util.TreeMap; 016 017import edu.jas.poly.AlgebraicNumber; 018import edu.jas.poly.ExpVector; 019import edu.jas.poly.GenPolynomial; 020import edu.jas.structure.BinaryFunctor; 021import edu.jas.structure.RingElem; 022import edu.jas.structure.Selector; 023import edu.jas.structure.UnaryFunctor; 024import edu.jas.util.MapEntry; 025 026 027/** 028 * Multivariate power series implementation. Uses inner classes and lazy 029 * evaluated generating function for coefficients. All ring element methods use 030 * lazy evaluation except where noted otherwise. Eager evaluated methods are 031 * <code>toString()</code>, <code>compareTo()</code>, <code>equals()</code>, 032 * <code>evaluate()</code>, or methods which use the <code>order()</code> or 033 * <code>orderExpVector()</code> methods, like <code>signum()</code>, 034 * <code>abs()</code>, <code>divide()</code>, <code>remainder()</code> and 035 * <code>gcd()</code>. <b>Note: </b> Currently the term order is fixed to the 036 * order defined by the iterator over exponent vectors in class 037 * <code>ExpVectorIterator</code>. 038 * @param <C> ring element type 039 * @author Heinz Kredel 040 */ 041 042public class MultiVarPowerSeries<C extends RingElem<C>> implements RingElem<MultiVarPowerSeries<C>> { 043 044 045 /** 046 * Power series ring factory. 047 */ 048 public final MultiVarPowerSeriesRing<C> ring; 049 050 051 /** 052 * Data structure / generating function for coeffcients. Cannot be final 053 * because of fixPoint, must be accessible in factory. 054 */ 055 /*package*/MultiVarCoefficients<C> lazyCoeffs; 056 057 058 /** 059 * Truncation of computations. 060 */ 061 private int truncate; 062 063 064 /** 065 * Order of power series. 066 */ 067 private int order = -1; // == unknown 068 069 070 /** 071 * ExpVector of order of power series. 072 */ 073 private ExpVector evorder = null; // == unknown 074 075 076 /** 077 * Private constructor. 078 */ 079 @SuppressWarnings("unused") 080 private MultiVarPowerSeries() { 081 throw new IllegalArgumentException("do not use no-argument constructor"); 082 } 083 084 085 /** 086 * Package constructor. Use in fixPoint only, must be accessible in factory. 087 * @param ring power series ring. 088 */ 089 /*package*/MultiVarPowerSeries(MultiVarPowerSeriesRing<C> ring) { 090 this.ring = ring; 091 this.lazyCoeffs = null; 092 this.truncate = ring.truncate; 093 } 094 095 096 /** 097 * Constructor. 098 * @param ring power series ring. 099 * @param lazyCoeffs generating function for coefficients. 100 */ 101 public MultiVarPowerSeries(MultiVarPowerSeriesRing<C> ring, MultiVarCoefficients<C> lazyCoeffs) { 102 this(ring, lazyCoeffs, ring.truncate); 103 } 104 105 106 /** 107 * Constructor. 108 * @param ring power series ring. 109 * @param lazyCoeffs generating function for coefficients. 110 * @param trunc truncate parameter for this power series. 111 */ 112 public MultiVarPowerSeries(MultiVarPowerSeriesRing<C> ring, MultiVarCoefficients<C> lazyCoeffs, int trunc) { 113 if (lazyCoeffs == null || ring == null) { 114 throw new IllegalArgumentException("null not allowed: ring = " + ring + ", lazyCoeffs = " 115 + lazyCoeffs); 116 } 117 this.ring = ring; 118 this.lazyCoeffs = lazyCoeffs; 119 this.truncate = Math.min(trunc, ring.truncate); 120 } 121 122 123 /** 124 * Get the corresponding element factory. 125 * @return factory for this Element. 126 * @see edu.jas.structure.Element#factory() 127 */ 128 public MultiVarPowerSeriesRing<C> factory() { 129 return ring; 130 } 131 132 133 /** 134 * Clone this power series. 135 * @see java.lang.Object#clone() 136 */ 137 @Override 138 public MultiVarPowerSeries<C> copy() { 139 return new MultiVarPowerSeries<C>(ring, lazyCoeffs); 140 } 141 142 143 /** 144 * String representation of power series. 145 * @see java.lang.Object#toString() 146 */ 147 @Override 148 public String toString() { 149 return toString(truncate); 150 } 151 152 153 /** 154 * To String with given truncate. 155 * @param trunc truncate parameter for this power series. 156 * @return string representation of this to given truncate. 157 */ 158 public String toString(int trunc) { 159 StringBuffer sb = new StringBuffer(); 160 MultiVarPowerSeries<C> s = this; 161 String[] vars = ring.vars; 162 //System.out.println("cache1 = " + s.lazyCoeffs.coeffCache); 163 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, trunc)) { 164 C c = s.coefficient(i); 165 //System.out.println("i = " + i + ", c = " +c); 166 int si = c.signum(); 167 if (si != 0) { 168 if (si > 0) { 169 if (sb.length() > 0) { 170 sb.append(" + "); 171 } 172 } else { 173 c = c.negate(); 174 sb.append(" - "); 175 } 176 if (!c.isONE() || i.isZERO()) { 177 if (c instanceof GenPolynomial || c instanceof AlgebraicNumber) { 178 sb.append("{ "); 179 } 180 sb.append(c.toString()); 181 if (c instanceof GenPolynomial || c instanceof AlgebraicNumber) { 182 sb.append(" }"); 183 } 184 if (!i.isZERO()) { 185 sb.append(" * "); 186 } 187 } 188 if (i.isZERO()) { 189 //skip; sb.append(" "); 190 } else { 191 sb.append(i.toString(vars)); 192 } 193 //sb.append(c.toString() + ", "); 194 } 195 //System.out.println("cache = " + s.coeffCache); 196 } 197 if (sb.length() == 0) { 198 sb.append("0"); 199 } 200 sb.append(" + BigO( (" + ring.varsToString() + ")^" + (trunc + 1) + "(" + (ring.truncate + 1) + ") )"); 201 //sb.append("..."); 202 //System.out.println("cache2 = " + s.lazyCoeffs.coeffCache); 203 return sb.toString(); 204 } 205 206 207 /** 208 * Get a scripting compatible string representation. 209 * @return script compatible representation for this Element. 210 * @see edu.jas.structure.Element#toScript() 211 */ 212 //JAVA6only: @Override 213 public String toScript() { 214 // Python case 215 StringBuffer sb = new StringBuffer(""); 216 MultiVarPowerSeries<C> s = this; 217 String[] vars = ring.vars; 218 //System.out.println("cache = " + s.coeffCache); 219 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, truncate)) { 220 C c = s.coefficient(i); 221 int si = c.signum(); 222 if (si != 0) { 223 if (si > 0) { 224 if (sb.length() > 0) { 225 sb.append(" + "); 226 } 227 } else { 228 c = c.negate(); 229 sb.append(" - "); 230 } 231 if (!c.isONE() || i.isZERO()) { 232 if (c instanceof GenPolynomial || c instanceof AlgebraicNumber) { 233 sb.append("( "); 234 } 235 sb.append(c.toScript()); 236 if (c instanceof GenPolynomial || c instanceof AlgebraicNumber) { 237 sb.append(" )"); 238 } 239 if (!i.isZERO()) { 240 sb.append(" * "); 241 } 242 } 243 if (i.isZERO()) { 244 //skip; sb.append(" "); 245 } else { 246 sb.append(i.toScript(vars)); 247 } 248 //sb.append(c.toString() + ", "); 249 } 250 //System.out.println("cache = " + s.coeffCache); 251 } 252 if (sb.length() == 0) { 253 sb.append("0"); 254 } 255 sb.append(" + BigO( (" + ring.varsToString() + ")**" + (truncate + 1) + " )"); 256 // sb.append("," + truncate + ""); 257 return sb.toString(); 258 } 259 260 261 /** 262 * Get a scripting compatible string representation of the factory. 263 * @return script compatible representation for this ElemFactory. 264 * @see edu.jas.structure.Element#toScriptFactory() 265 */ 266 //JAVA6only: @Override 267 public String toScriptFactory() { 268 // Python case 269 return factory().toScript(); 270 } 271 272 273 /** 274 * Get coefficient. 275 * @param index number of requested coefficient. 276 * @return coefficient at index. 277 */ 278 public C coefficient(ExpVector index) { 279 if (index == null) { 280 throw new IndexOutOfBoundsException("null index not allowed"); 281 } 282 return lazyCoeffs.get(index); 283 } 284 285 286 /** 287 * Homogeneous part. 288 * @param tdeg requested degree. 289 * @return polynomial part of given degree. 290 */ 291 public GenPolynomial<C> homogeneousPart(long tdeg) { 292 return lazyCoeffs.getHomPart(tdeg); 293 } 294 295 296 /** 297 * Get a GenPolynomial<C> from this. 298 * @return a GenPolynomial<C> from this up to truncate homogeneous 299 * parts. 300 */ 301 public GenPolynomial<C> asPolynomial() { 302 GenPolynomial<C> p = homogeneousPart(0L); 303 for (int i = 1; i <= truncate; i++) { 304 p = p.sum(homogeneousPart(i)); 305 } 306 return p; 307 } 308 309 310 /** 311 * Leading base coefficient. 312 * @return first coefficient. 313 */ 314 public C leadingCoefficient() { 315 return coefficient(ring.EVZERO); 316 } 317 318 319 /** 320 * Reductum. 321 * @param r variable for taking the reductum. 322 * @return this - leading monomial in the direction of r. 323 */ 324 public MultiVarPowerSeries<C> reductum(final int r) { 325 if (r < 0 || ring.nvar < r) { 326 throw new IllegalArgumentException("variable index out of bound"); 327 } 328 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 329 330 331 @Override 332 public C generate(ExpVector i) { 333 ExpVector e = i.subst(r, i.getVal(r) + 1); 334 return coefficient(e); 335 } 336 }); 337 } 338 339 340 /** 341 * Prepend a new leading coefficient. 342 * @param r variable for the direction. 343 * @param h new coefficient. 344 * @return new power series. 345 */ 346 public MultiVarPowerSeries<C> prepend(final C h, final int r) { 347 if (r < 0 || ring.nvar < r) { 348 throw new IllegalArgumentException("variable index out of bound"); 349 } 350 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 351 352 353 @Override 354 public C generate(ExpVector i) { 355 if (i.isZERO()) { 356 return h; 357 } 358 ExpVector e = i.subst(r, i.getVal(r) - 1); 359 if (e.signum() < 0) { 360 return pfac.coFac.getZERO(); 361 } 362 return coefficient(e); 363 } 364 }); 365 } 366 367 368 /** 369 * Shift coefficients. 370 * @param k shift index. 371 * @param r variable for the direction. 372 * @return new power series with coefficient(i) = old.coefficient(i-k). 373 */ 374 public MultiVarPowerSeries<C> shift(final int k, final int r) { 375 if (r < 0 || ring.nvar < r) { 376 throw new IllegalArgumentException("variable index out of bound"); 377 } 378 int nt = Math.min(truncate + k, ring.truncate); 379 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 380 381 382 @Override 383 public C generate(ExpVector i) { 384 long d = i.getVal(r); 385 if (d - k < 0) { 386 return ring.coFac.getZERO(); 387 } 388 ExpVector e = i.subst(r, i.getVal(r) - k); 389 return coefficient(e); 390 } 391 }, nt); 392 } 393 394 395 /** 396 * Reductum. 397 * @return this - leading monomial. 398 */ 399 public MultiVarPowerSeries<C> reductum() { 400 Map.Entry<ExpVector, C> m = orderMonomial(); 401 ExpVector e = m.getKey(); 402 long d = e.totalDeg(); 403 MultiVarCoefficients<C> mc = lazyCoeffs; 404 HashMap<Long, GenPolynomial<C>> cc = new HashMap<Long, GenPolynomial<C>>(mc.coeffCache); 405 GenPolynomial<C> p = cc.get(d); 406 if (p != null && !p.isZERO()) { 407 p = p.subtract(m.getValue(), e); // p contains this term after orderMonomial() 408 cc.put(d, p); 409 } 410 HashSet<ExpVector> z = new HashSet<ExpVector>(mc.zeroCache); 411 if (!mc.homCheck.get((int) d)) { 412 z.add(e); 413 //System.out.println("e = " + e); 414 } 415 416 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(mc.pfac, cc, z, mc.homCheck) { 417 418 419 @Override 420 public C generate(ExpVector i) { 421 return coefficient(i); 422 } 423 }); 424 } 425 426 427 /** 428 * Shift coefficients. Multiply by exponent vector. 429 * @param k shift ExpVector. 430 * @return new power series with coefficient(i) = old.coefficient(i-k). 431 */ 432 public MultiVarPowerSeries<C> shift(final ExpVector k) { 433 if (k == null) { 434 throw new IllegalArgumentException("null ExpVector not allowed"); 435 } 436 if (k.signum() == 0) { 437 return this; 438 } 439 int nt = Math.min(truncate + (int) k.totalDeg(), ring.truncate); 440 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 441 442 443 @Override 444 public C generate(ExpVector i) { 445 ExpVector d = i.subtract(k); 446 if (d.signum() < 0) { 447 return ring.coFac.getZERO(); 448 } 449 return coefficient(d); 450 } 451 }, nt); 452 } 453 454 455 /** 456 * Multiply by exponent vector and coefficient. 457 * @param k shift ExpVector. 458 * @param c coefficient multiplier. 459 * @return new power series with coefficient(i) = old.coefficient(i-k)*c. 460 */ 461 public MultiVarPowerSeries<C> multiply(final C c, final ExpVector k) { 462 if (k == null) { 463 throw new IllegalArgumentException("null ExpVector not allowed"); 464 } 465 if (k.signum() == 0) { 466 return this.multiply(c); 467 } 468 if (c.signum() == 0) { 469 return ring.getZERO(); 470 } 471 int nt = Math.min(ring.truncate, truncate + (int) k.totalDeg()); 472 473 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 474 475 476 @Override 477 public C generate(ExpVector i) { 478 ExpVector d = i.subtract(k); 479 if (d.signum() < 0) { 480 return ring.coFac.getZERO(); 481 } 482 long tdegd = d.totalDeg(); 483 if (lazyCoeffs.homCheck.get((int) tdegd)) { 484 GenPolynomial<C> p = homogeneousPart(tdegd).multiply(c, k); 485 long tdegi = i.totalDeg(); 486 coeffCache.put(tdegi, p); // overwrite 487 homCheck.set((int) tdegi); 488 C b = p.coefficient(i); 489 //System.out.println("b = " + b + ", i = " + i + ", tdegi = " + tdegi+ ", tdegd = " + tdegd); 490 //System.out.println("p = " + p + ", i = " + i); 491 return b; 492 } 493 return coefficient(d).multiply(c); 494 } 495 }, nt); 496 } 497 498 499 /** 500 * Sum monomial. 501 * @param m ExpVector , coeffcient pair 502 * @return this + ONE.multiply(m.coefficient,m.exponent). 503 */ 504 public MultiVarPowerSeries<C> sum(Map.Entry<ExpVector, C> m) { 505 if (m == null) { 506 throw new IllegalArgumentException("null Map.Entry not allowed"); 507 } 508 return sum(m.getValue(), m.getKey()); 509 } 510 511 512 /** 513 * Sum exponent vector and coefficient. 514 * @param k ExpVector. 515 * @param c coefficient. 516 * @return this + ONE.multiply(c,k). 517 */ 518 public MultiVarPowerSeries<C> sum(final C c, final ExpVector k) { 519 if (k == null) { 520 throw new IllegalArgumentException("null ExpVector not allowed"); 521 } 522 if (c.signum() == 0) { 523 return this; 524 } 525 long d = k.totalDeg(); 526 MultiVarCoefficients<C> mc = lazyCoeffs; 527 HashMap<Long, GenPolynomial<C>> cc = new HashMap<Long, GenPolynomial<C>>(mc.coeffCache); 528 GenPolynomial<C> p = cc.get(d); 529 if (p == null) { 530 p = mc.pfac.getZERO(); 531 } 532 p = p.sum(c, k); 533 //System.out.println("p = " + p); 534 cc.put(d, p); 535 HashSet<ExpVector> z = new HashSet<ExpVector>(mc.zeroCache); 536 //System.out.println("z = " + z); 537 if (p.coefficient(k).isZERO() && !mc.homCheck.get((int) d)) { 538 z.add(k); 539 } 540 541 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(mc.pfac, cc, z, mc.homCheck) { 542 543 544 @Override 545 public C generate(ExpVector i) { 546 return coefficient(i); 547 } 548 }); 549 } 550 551 552 /** 553 * Subtract exponent vector and coefficient. 554 * @param k ExpVector. 555 * @param c coefficient. 556 * @return this - ONE.multiply(c,k). 557 */ 558 public MultiVarPowerSeries<C> subtract(final C c, final ExpVector k) { 559 if (k == null) { 560 throw new IllegalArgumentException("null ExpVector not allowed"); 561 } 562 if (c.signum() == 0) { 563 return this; 564 } 565 long d = k.totalDeg(); 566 MultiVarCoefficients<C> mc = lazyCoeffs; 567 HashMap<Long, GenPolynomial<C>> cc = new HashMap<Long, GenPolynomial<C>>(mc.coeffCache); 568 GenPolynomial<C> p = cc.get(d); 569 if (p == null) { 570 p = mc.pfac.getZERO(); 571 } 572 p = p.subtract(c, k); 573 cc.put(d, p); 574 HashSet<ExpVector> z = new HashSet<ExpVector>(mc.zeroCache); 575 //System.out.println("z = " + z); 576 if (p.coefficient(k).isZERO() && !mc.homCheck.get((int) d)) { 577 z.add(k); 578 } 579 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(mc.pfac, cc, z, mc.homCheck) { 580 581 582 @Override 583 public C generate(ExpVector i) { 584 return coefficient(i); 585 } 586 }); 587 } 588 589 590 /** 591 * Sum exponent vector and coefficient. 592 * @param mvc cached coefficients. 593 * @return this + mvc. 594 */ 595 public MultiVarPowerSeries<C> sum(MultiVarCoefficients<C> mvc) { 596 MultiVarCoefficients<C> mc = lazyCoeffs; 597 TreeMap<Long, GenPolynomial<C>> cc = new TreeMap<Long, GenPolynomial<C>>(mc.coeffCache); 598 TreeMap<Long, GenPolynomial<C>> ccv = new TreeMap<Long, GenPolynomial<C>>(mvc.coeffCache); 599 long d1 = (cc.size() > 0 ? cc.lastKey() : 0); 600 long d2 = (ccv.size() > 0 ? ccv.lastKey() : 0); 601 HashSet<ExpVector> z = new HashSet<ExpVector>(mc.zeroCache); 602 z.addAll(mvc.zeroCache); 603 long d = Math.max(d1, d2); 604 BitSet hc = new BitSet((int) d); 605 for (long i = 0; i <= d; i++) { 606 GenPolynomial<C> p1 = cc.get(i); 607 GenPolynomial<C> p2 = mvc.coeffCache.get(i); 608 if (p1 == null) { 609 p1 = mc.pfac.getZERO(); 610 } 611 if (p2 == null) { 612 p2 = mc.pfac.getZERO(); 613 } 614 GenPolynomial<C> p = p1.sum(p2); 615 //System.out.println("p = " + p); 616 cc.put(i, p); 617 if (mc.homCheck.get((int) i) && mvc.homCheck.get((int) i)) { 618 hc.set((int) i); 619 } else { 620 Set<ExpVector> ev = new HashSet<ExpVector>(p1.getMap().keySet()); 621 ev.addAll(p2.getMap().keySet()); 622 ev.removeAll(p.getMap().keySet()); 623 z.addAll(ev); 624 } 625 } 626 //System.out.println("cc = " + cc); 627 628 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(mc.pfac, 629 new HashMap<Long, GenPolynomial<C>>(cc), z, hc) { 630 631 632 @Override 633 public C generate(ExpVector i) { 634 return coefficient(i); 635 } 636 }); 637 } 638 639 640 /** 641 * Select coefficients. 642 * @param sel selector functor. 643 * @return new power series with selected coefficients. 644 */ 645 public MultiVarPowerSeries<C> select(final Selector<? super C> sel) { 646 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 647 648 649 @Override 650 public C generate(ExpVector i) { 651 C c = coefficient(i); 652 if (sel.select(c)) { 653 return c; 654 } 655 return ring.coFac.getZERO(); 656 } 657 }); 658 } 659 660 661 /** 662 * Shift select coefficients. Not selected coefficients are removed from the 663 * result series. 664 * @param sel selector functor. 665 * @return new power series with shifted selected coefficients. 666 */ 667 public MultiVarPowerSeries<C> shiftSelect(final Selector<? super C> sel) { 668 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 669 670 671 ExpVectorIterable ib = new ExpVectorIterable(ring.nvar, true, truncate); 672 673 674 Iterator<ExpVector> pos = ib.iterator(); 675 676 677 @Override 678 public C generate(ExpVector i) { 679 C c; 680 if (i.signum() > 0) { 681 int[] deps = i.dependencyOnVariables(); 682 ExpVector x = i.subst(deps[0], i.getVal(deps[0]) - 1L); 683 c = get(x); 684 } 685 do { 686 c = null; 687 if (pos.hasNext()) { 688 ExpVector e = pos.next(); 689 c = coefficient(e); 690 } else { 691 break; 692 } 693 } while (!sel.select(c)); 694 if (c == null) { // not correct because not known 695 c = ring.coFac.getZERO(); 696 } 697 return c; 698 } 699 }); 700 } 701 702 703 /** 704 * Map a unary function to this power series. 705 * @param f evaluation functor. 706 * @return new power series with coefficients f(this(i)). 707 */ 708 public MultiVarPowerSeries<C> map(final UnaryFunctor<? super C, C> f) { 709 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 710 711 712 @Override 713 public C generate(ExpVector i) { 714 return f.eval(coefficient(i)); 715 } 716 }); 717 } 718 719 720 /** 721 * Map a binary function to this and another power series. 722 * @param f evaluation functor with coefficients f(this(i),other(i)). 723 * @param ps other power series. 724 * @return new power series. 725 */ 726 public MultiVarPowerSeries<C> zip(final BinaryFunctor<? super C, ? super C, C> f, 727 final MultiVarPowerSeries<C> ps) { 728 int m = Math.min(ring.truncate, Math.max(truncate, ps.truncate())); 729 730 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 731 732 733 @Override 734 public C generate(ExpVector i) { 735 return f.eval(coefficient(i), ps.coefficient(i)); 736 } 737 }, m); 738 } 739 740 741 /** 742 * Sum of two power series, using zip(). 743 * @param ps other power series. 744 * @return this + ps. 745 */ 746 public MultiVarPowerSeries<C> sumZip(MultiVarPowerSeries<C> ps) { 747 return zip(new BinaryFunctor<C, C, C>() { 748 749 750 //JAVA6only: @Override 751 public C eval(C c1, C c2) { 752 return c1.sum(c2); 753 } 754 }, ps); 755 } 756 757 758 /** 759 * Subtraction of two power series, using zip(). 760 * @param ps other power series. 761 * @return this - ps. 762 */ 763 public MultiVarPowerSeries<C> subtractZip(MultiVarPowerSeries<C> ps) { 764 return zip(new BinaryFunctor<C, C, C>() { 765 766 767 //JAVA6only: @Override 768 public C eval(C c1, C c2) { 769 return c1.subtract(c2); 770 } 771 }, ps); 772 } 773 774 775 /** 776 * Multiply by coefficient. 777 * @param a coefficient. 778 * @return this * a. 779 */ 780 public MultiVarPowerSeries<C> multiply(final C a) { 781 if (a.isZERO()) { 782 return ring.getZERO(); 783 } 784 if (a.isONE()) { 785 return this; 786 } 787 return map(new UnaryFunctor<C, C>() { 788 789 790 //JAVA6only: @Override 791 public C eval(C c) { 792 return c.multiply(a); 793 } 794 }); 795 } 796 797 798 /** 799 * Monic. 800 * @return 1/orderCoeff() * this. 801 */ 802 public MultiVarPowerSeries<C> monic() { 803 ExpVector e = orderExpVector(); 804 if (e == null) { 805 return this; 806 } 807 C a = coefficient(e); 808 if (a.isONE()) { 809 return this; 810 } 811 if (a.isZERO()) { // sic 812 return this; 813 } 814 final C b = a.inverse(); 815 return map(new UnaryFunctor<C, C>() { 816 817 818 //JAVA6only: @Override 819 public C eval(C c) { 820 return b.multiply(c); 821 } 822 }); 823 } 824 825 826 /** 827 * Negate. 828 * @return - this. 829 */ 830 public MultiVarPowerSeries<C> negate() { 831 return map(new UnaryFunctor<C, C>() { 832 833 834 //JAVA6only: @Override 835 public C eval(C c) { 836 return c.negate(); 837 } 838 }); 839 } 840 841 842 /** 843 * Absolute value. 844 * @return abs(this). 845 */ 846 public MultiVarPowerSeries<C> abs() { 847 if (signum() < 0) { 848 return negate(); 849 } 850 return this; 851 } 852 853 854 /** 855 * Evaluate at given point. 856 * @return ps(a). 857 */ 858 public C evaluate(List<C> a) { 859 C v = ring.coFac.getZERO(); 860 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, truncate)) { 861 C c = coefficient(i); 862 if (c.isZERO()) { 863 continue; 864 } 865 c = c.multiply(i.evaluate(ring.coFac, a)); 866 v = v.sum(c); 867 } 868 return v; 869 } 870 871 872 /** 873 * Order. 874 * @return index of first non zero coefficient. 875 */ 876 public int order() { 877 if (order >= 0) { 878 return order; 879 } 880 // must compute it 881 GenPolynomial<C> p = null; 882 int t = 0; 883 while (lazyCoeffs.homCheck.get(t)) { 884 p = lazyCoeffs.coeffCache.get(t); 885 if (p == null || p.isZERO()) { // ?? 886 t++; 887 continue; 888 } 889 order = t; 890 evorder = p.trailingExpVector(); 891 //System.out.println("order = " + t); 892 return order; 893 } 894 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, truncate)) { 895 if (!coefficient(i).isZERO()) { 896 order = (int) i.totalDeg(); //ord; 897 evorder = i; 898 //System.out.println("order = " + order + ", evorder = " + evorder); 899 return order; 900 } 901 } 902 order = truncate + 1; 903 // evorder is null 904 return order; 905 } 906 907 908 /** 909 * Order ExpVector. 910 * @return ExpVector of first non zero coefficient. 911 */ 912 @SuppressWarnings("unused") 913 public ExpVector orderExpVector() { 914 int x = order(); // ensure evorder is set 915 return evorder; 916 } 917 918 919 /** 920 * Order monomial. 921 * @return first map entry or null. 922 */ 923 public Map.Entry<ExpVector, C> orderMonomial() { 924 ExpVector e = orderExpVector(); 925 if (e == null) { 926 return null; 927 } 928 //JAVA6only: 929 //return new AbstractMap.SimpleImmutableEntry<ExpVector, C>(e, coefficient(e)); 930 return new MapEntry<ExpVector, C>(e, coefficient(e)); 931 } 932 933 934 /** 935 * Truncate. 936 * @return truncate index of power series. 937 */ 938 public int truncate() { 939 return truncate; 940 } 941 942 943 /** 944 * Set truncate. 945 * @param t new truncate index. 946 * @return old truncate index of power series. 947 */ 948 public int setTruncate(int t) { 949 if (t < 0) { 950 throw new IllegalArgumentException("negative truncate not allowed"); 951 } 952 int ot = truncate; 953 if (order >= 0) { 954 if (order > truncate) { 955 order = -1; // reset 956 evorder = null; 957 } 958 } 959 truncate = t; 960 return ot; 961 } 962 963 964 /** 965 * Ecart. 966 * @return ecart. 967 */ 968 public long ecart() { 969 ExpVector e = orderExpVector(); 970 if (e == null) { 971 return 0L; 972 } 973 long d = e.totalDeg(); 974 long hd = d; 975 for (long i = d + 1L; i <= truncate; i++) { 976 if (!homogeneousPart(i).isZERO()) { 977 hd = i; 978 } 979 } 980 //System.out.println("d = " + d + ", hd = " + hd + ", coeffCache = " + lazyCoeffs.coeffCache + ", this = " + this); 981 return hd - d; 982 } 983 984 985 /** 986 * Signum. 987 * @return sign of first non zero coefficient. 988 */ 989 public int signum() { 990 ExpVector ev = orderExpVector(); 991 if (ev != null) { 992 return coefficient(ev).signum(); 993 } 994 return 0; 995 } 996 997 998 /** 999 * Compare to. <b>Note: </b> compare only up to max(truncates). 1000 * @return sign of first non zero coefficient of this-ps. 1001 */ 1002 //JAVA6only: @Override 1003 public int compareTo(MultiVarPowerSeries<C> ps) { 1004 final int m = truncate(); 1005 final int n = ps.truncate(); 1006 final int pos = Math.min(ring.truncate, Math.min(m, n)); 1007 int s = 0; 1008 //System.out.println("coeffCache_c1 = " + lazyCoeffs.coeffCache); 1009 //System.out.println("coeffCache_c2 = " + ps.lazyCoeffs.coeffCache); 1010 // test homogeneous parts first is slower 1011 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, pos)) { 1012 s = coefficient(i).compareTo(ps.coefficient(i)); 1013 if (s != 0) { 1014 //System.out.println("i = " + i + ", coeff = " + coefficient(i) + ", ps.coeff = " + ps.coefficient(i)); 1015 return s; 1016 } 1017 } 1018 for (int j = pos + 1; j <= Math.min(ring.truncate, Math.max(m, n)); j++) { 1019 for (ExpVector i : new ExpVectorIterable(ring.nvar, j)) { 1020 s = coefficient(i).compareTo(ps.coefficient(i)); 1021 //System.out.println("i = " + i + ", coeff = " + coefficient(i) + ", ps.coeff = " + ps.coefficient(i)); 1022 if (s != 0) { 1023 //System.out.println("i = " + i + ", coeff = " + coefficient(i) + ", ps.coeff = " + ps.coefficient(i)); 1024 return s; 1025 } 1026 } 1027 } 1028 return s; 1029 } 1030 1031 1032 /** 1033 * Is power series zero. <b>Note: </b> compare only up to truncate. 1034 * @return If this is 0 then true is returned, else false. 1035 * @see edu.jas.structure.RingElem#isZERO() 1036 */ 1037 public boolean isZERO() { 1038 return (signum() == 0); 1039 } 1040 1041 1042 /** 1043 * Is power series one. <b>Note: </b> compare only up to truncate. 1044 * @return If this is 1 then true is returned, else false. 1045 * @see edu.jas.structure.RingElem#isONE() 1046 */ 1047 public boolean isONE() { 1048 if (!leadingCoefficient().isONE()) { 1049 return false; 1050 } 1051 return (compareTo(ring.ONE) == 0); 1052 //return reductum().isZERO(); 1053 } 1054 1055 1056 /** 1057 * Comparison with any other object. <b>Note: </b> compare only up to 1058 * truncate. 1059 * @see java.lang.Object#equals(java.lang.Object) 1060 */ 1061 @Override 1062 @SuppressWarnings("unchecked") 1063 public boolean equals(Object B) { 1064 if (!(B instanceof MultiVarPowerSeries)) { 1065 return false; 1066 } 1067 MultiVarPowerSeries<C> a = null; 1068 try { 1069 a = (MultiVarPowerSeries<C>) B; 1070 } catch (ClassCastException ignored) { 1071 } 1072 if (a == null) { 1073 return false; 1074 } 1075 return compareTo(a) == 0; 1076 } 1077 1078 1079 /** 1080 * Hash code for this polynomial. <b>Note: </b> only up to truncate. 1081 * @see java.lang.Object#hashCode() 1082 */ 1083 @Override 1084 public int hashCode() { 1085 int h = 0; 1086 //h = ( ring.hashCode() << 23 ); 1087 //h += val.hashCode(); 1088 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, truncate)) { 1089 C c = coefficient(i); 1090 if (!c.isZERO()) { 1091 h += i.hashCode(); 1092 h = (h << 23); 1093 } 1094 h += c.hashCode(); 1095 h = (h << 23); 1096 } 1097 return h; 1098 } 1099 1100 1101 /** 1102 * Is unit. 1103 * @return true, if this power series is invertible, else false. 1104 */ 1105 public boolean isUnit() { 1106 return leadingCoefficient().isUnit(); 1107 } 1108 1109 1110 /** 1111 * Sum a another power series. 1112 * @param ps other power series. 1113 * @return this + ps. 1114 */ 1115 public MultiVarPowerSeries<C> sum(final MultiVarPowerSeries<C> ps) { 1116 //final MultiVarPowerSeries<C> ps1 = this; // method name was ambiguous in generate 1117 int nt = Math.min(ring.truncate, Math.max(truncate(), ps.truncate())); 1118 1119 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1120 1121 1122 @Override 1123 public C generate(ExpVector e) { 1124 long tdeg = e.totalDeg(); 1125 if (lazyCoeffs.homCheck.get((int) tdeg)) { 1126 // generate respective homogeneous polynomial 1127 GenPolynomial<C> p = homogeneousPart(tdeg).sum(ps.homogeneousPart(tdeg)); 1128 coeffCache.put(tdeg, p); // overwrite 1129 homCheck.set((int) tdeg); 1130 C c = p.coefficient(e); 1131 //System.out.println("c = " + c + ", e = " + e + ", tdeg = " + tdeg); 1132 return c; 1133 } 1134 return coefficient(e).sum(ps.coefficient(e)); 1135 } 1136 }, nt); 1137 } 1138 1139 1140 /** 1141 * Subtract a another power series. 1142 * @param ps other power series. 1143 * @return this - ps. 1144 */ 1145 public MultiVarPowerSeries<C> subtract(final MultiVarPowerSeries<C> ps) { 1146 //final MultiVarPowerSeries<C> ps1 = this; // method name was ambiguous in generate 1147 int nt = Math.min(ring.truncate, Math.max(truncate(), ps.truncate())); 1148 1149 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1150 1151 1152 @Override 1153 public C generate(ExpVector e) { 1154 long tdeg = e.totalDeg(); 1155 if (lazyCoeffs.homCheck.get((int) tdeg)) { 1156 // generate respective homogeneous polynomial 1157 GenPolynomial<C> p = homogeneousPart(tdeg).subtract(ps.homogeneousPart(tdeg)); 1158 coeffCache.put(tdeg, p); // overwrite 1159 homCheck.set((int) tdeg); 1160 C c = p.coefficient(e); 1161 //System.out.println("p = " + p + ", e = " + e + ", tdeg = " + tdeg); 1162 return c; 1163 } 1164 return coefficient(e).subtract(ps.coefficient(e)); 1165 } 1166 }, nt); 1167 } 1168 1169 1170 /** 1171 * Multiply by another power series. 1172 * @param ps other power series. 1173 * @return this * ps. 1174 */ 1175 public MultiVarPowerSeries<C> multiply(final MultiVarPowerSeries<C> ps) { 1176 //final MultiVarPowerSeries<C> ps1 = this; // method name was ambiguous in generate 1177 int nt = Math.min(ring.truncate, truncate() + ps.truncate()); 1178 1179 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1180 1181 1182 @Override 1183 public C generate(ExpVector e) { 1184 long tdeg = e.totalDeg(); 1185 // generate respective homogeneous polynomial 1186 GenPolynomial<C> p = null; //fac.getZERO(); 1187 for (int k = 0; k <= tdeg; k++) { 1188 GenPolynomial<C> m = homogeneousPart(k).multiply(ps.homogeneousPart(tdeg - k)); 1189 if (p == null) { 1190 p = m; 1191 } else { 1192 p = p.sum(m); 1193 } 1194 } 1195 coeffCache.put(tdeg, p); // overwrite 1196 homCheck.set((int) tdeg); 1197 C c = p.coefficient(e); 1198 return c; 1199 } 1200 }, nt); 1201 } 1202 1203 1204 /** 1205 * Inverse power series. 1206 * @return ps with this * ps = 1. 1207 */ 1208 public MultiVarPowerSeries<C> inverse() { 1209 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1210 1211 1212 @Override 1213 public C generate(ExpVector e) { 1214 long tdeg = e.totalDeg(); 1215 C d = leadingCoefficient().inverse(); // may fail 1216 if (tdeg == 0) { 1217 return d; 1218 } 1219 GenPolynomial<C> p = null; //fac.getZERO(); 1220 for (int k = 0; k < tdeg; k++) { 1221 GenPolynomial<C> m = getHomPart(k).multiply(homogeneousPart(tdeg - k)); 1222 if (p == null) { 1223 p = m; 1224 } else { 1225 p = p.sum(m); 1226 } 1227 } 1228 p = p.multiply(d.negate()); 1229 //System.out.println("tdeg = " + tdeg + ", p = " + p); 1230 coeffCache.put(tdeg, p); // overwrite 1231 homCheck.set((int) tdeg); 1232 C c = p.coefficient(e); 1233 return c; 1234 } 1235 }); 1236 } 1237 1238 1239 /** 1240 * Divide by another power series. 1241 * @param ps nonzero power series with invertible coefficient. 1242 * @return this / ps. 1243 */ 1244 public MultiVarPowerSeries<C> divide(MultiVarPowerSeries<C> ps) { 1245 if (ps.isUnit()) { 1246 return multiply(ps.inverse()); 1247 } 1248 int m = order(); 1249 int n = ps.order(); 1250 if (m < n) { 1251 return ring.getZERO(); 1252 } 1253 ExpVector em = orderExpVector(); 1254 ExpVector en = ps.orderExpVector(); 1255 if (!ps.coefficient(en).isUnit()) { 1256 throw new ArithmeticException("division by non unit coefficient " + ps.coefficient(ps.evorder) 1257 + ", evorder = " + ps.evorder); 1258 } 1259 // now m >= n 1260 MultiVarPowerSeries<C> st, sps, q, sq; 1261 if (m == 0) { 1262 st = this; 1263 } else { 1264 st = this.shift(em.negate()); 1265 } 1266 if (n == 0) { 1267 sps = ps; 1268 } else { 1269 sps = ps.shift(en.negate()); 1270 } 1271 q = st.multiply(sps.inverse()); 1272 sq = q.shift(em.subtract(en)); 1273 return sq; 1274 } 1275 1276 1277 /** 1278 * Power series remainder. 1279 * @param ps nonzero power series with invertible leading coefficient. 1280 * @return remainder with this = quotient * ps + remainder. 1281 */ 1282 public MultiVarPowerSeries<C> remainder(MultiVarPowerSeries<C> ps) { 1283 int m = order(); 1284 int n = ps.order(); 1285 if (m >= n) { 1286 return ring.getZERO(); 1287 } 1288 return this; 1289 } 1290 1291 1292 /** 1293 * Differentiate with respect to variable r. 1294 * @param r variable for the direction. 1295 * @return differentiate(this). 1296 */ 1297 public MultiVarPowerSeries<C> differentiate(final int r) { 1298 if (r < 0 || ring.nvar < r) { 1299 throw new IllegalArgumentException("variable index out of bound"); 1300 } 1301 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1302 1303 1304 @Override 1305 public C generate(ExpVector i) { 1306 long d = i.getVal(r); 1307 ExpVector e = i.subst(r, d + 1); 1308 C v = coefficient(e); 1309 v = v.multiply(ring.coFac.fromInteger(d + 1)); 1310 return v; 1311 } 1312 }); 1313 } 1314 1315 1316 /** 1317 * Integrate with respect to variable r and with given constant. 1318 * @param c integration constant. 1319 * @param r variable for the direction. 1320 * @return integrate(this). 1321 */ 1322 public MultiVarPowerSeries<C> integrate(final C c, final int r) { 1323 if (r < 0 || ring.nvar < r) { 1324 throw new IllegalArgumentException("variable index out of bound"); 1325 } 1326 int nt = Math.min(ring.truncate, truncate + 1); 1327 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1328 1329 1330 @Override 1331 public C generate(ExpVector i) { 1332 if (i.isZERO()) { 1333 return c; 1334 } 1335 long d = i.getVal(r); 1336 if (d > 0) { 1337 ExpVector e = i.subst(r, d - 1); 1338 C v = coefficient(e); 1339 v = v.divide(ring.coFac.fromInteger(d)); 1340 return v; 1341 } 1342 return ring.coFac.getZERO(); 1343 } 1344 }, nt); 1345 } 1346 1347 1348 /** 1349 * Power series greatest common divisor. 1350 * @param ps power series. 1351 * @return gcd(this,ps). 1352 */ 1353 public MultiVarPowerSeries<C> gcd(MultiVarPowerSeries<C> ps) { 1354 if (ps.isZERO()) { 1355 return this; 1356 } 1357 if (this.isZERO()) { 1358 return ps; 1359 } 1360 ExpVector em = orderExpVector(); 1361 ExpVector en = ps.orderExpVector(); 1362 return ring.getONE().shift(em.gcd(en)); 1363 } 1364 1365 1366 /** 1367 * Power series extended greatest common divisor. <b>Note:</b> not 1368 * implemented. 1369 * @param S power series. 1370 * @return [ gcd(this,S), a, b ] with a*this + b*S = gcd(this,S). 1371 */ 1372 //SuppressWarnings("unchecked") 1373 public MultiVarPowerSeries<C>[] egcd(MultiVarPowerSeries<C> S) { 1374 throw new UnsupportedOperationException("egcd for power series not implemented"); 1375 } 1376 1377}