001/* 002 * $Id: MultiVarPowerSeries.java 5310 2015-10-25 18:35:31Z 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 @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 @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 if (m == null) { 402 return ring.getZERO(); 403 } 404 ExpVector e = m.getKey(); 405 long d = e.totalDeg(); 406 MultiVarCoefficients<C> mc = lazyCoeffs; 407 HashMap<Long, GenPolynomial<C>> cc = new HashMap<Long, GenPolynomial<C>>(mc.coeffCache); 408 GenPolynomial<C> p = cc.get(d); 409 if (p != null && !p.isZERO()) { 410 p = p.subtract(m.getValue(), e); // p contains this term after orderMonomial() 411 cc.put(d, p); 412 } 413 HashSet<ExpVector> z = new HashSet<ExpVector>(mc.zeroCache); 414 if (!mc.homCheck.get((int) d)) { 415 z.add(e); 416 //System.out.println("e = " + e); 417 } 418 419 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(mc.pfac, cc, z, mc.homCheck) { 420 421 422 @Override 423 public C generate(ExpVector i) { 424 return coefficient(i); 425 } 426 }); 427 } 428 429 430 /** 431 * Shift coefficients. Multiply by exponent vector. 432 * @param k shift ExpVector. 433 * @return new power series with coefficient(i) = old.coefficient(i-k). 434 */ 435 public MultiVarPowerSeries<C> shift(final ExpVector k) { 436 if (k == null) { 437 throw new IllegalArgumentException("null ExpVector not allowed"); 438 } 439 if (k.signum() == 0) { 440 return this; 441 } 442 int nt = Math.min(truncate + (int) k.totalDeg(), ring.truncate); 443 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 444 445 446 @Override 447 public C generate(ExpVector i) { 448 ExpVector d = i.subtract(k); 449 if (d.signum() < 0) { 450 return ring.coFac.getZERO(); 451 } 452 return coefficient(d); 453 } 454 }, nt); 455 } 456 457 458 /** 459 * Multiply by exponent vector and coefficient. 460 * @param k shift ExpVector. 461 * @param c coefficient multiplier. 462 * @return new power series with coefficient(i) = old.coefficient(i-k)*c. 463 */ 464 public MultiVarPowerSeries<C> multiply(final C c, final ExpVector k) { 465 if (k == null) { 466 throw new IllegalArgumentException("null ExpVector not allowed"); 467 } 468 if (k.signum() == 0) { 469 return this.multiply(c); 470 } 471 if (c.signum() == 0) { 472 return ring.getZERO(); 473 } 474 int nt = Math.min(ring.truncate, truncate + (int) k.totalDeg()); 475 476 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 477 478 479 @Override 480 public C generate(ExpVector i) { 481 ExpVector d = i.subtract(k); 482 if (d.signum() < 0) { 483 return ring.coFac.getZERO(); 484 } 485 long tdegd = d.totalDeg(); 486 if (lazyCoeffs.homCheck.get((int) tdegd)) { 487 GenPolynomial<C> p = homogeneousPart(tdegd).multiply(c, k); 488 long tdegi = i.totalDeg(); 489 coeffCache.put(tdegi, p); // overwrite 490 homCheck.set((int) tdegi); 491 C b = p.coefficient(i); 492 //System.out.println("b = " + b + ", i = " + i + ", tdegi = " + tdegi+ ", tdegd = " + tdegd); 493 //System.out.println("p = " + p + ", i = " + i); 494 return b; 495 } 496 return coefficient(d).multiply(c); 497 } 498 }, nt); 499 } 500 501 502 /** 503 * Sum monomial. 504 * @param m ExpVector , coeffcient pair 505 * @return this + ONE.multiply(m.coefficient,m.exponent). 506 */ 507 public MultiVarPowerSeries<C> sum(Map.Entry<ExpVector, C> m) { 508 if (m == null) { 509 throw new IllegalArgumentException("null Map.Entry not allowed"); 510 } 511 return sum(m.getValue(), m.getKey()); 512 } 513 514 515 /** 516 * Sum exponent vector and coefficient. 517 * @param k ExpVector. 518 * @param c coefficient. 519 * @return this + ONE.multiply(c,k). 520 */ 521 public MultiVarPowerSeries<C> sum(final C c, final ExpVector k) { 522 if (k == null) { 523 throw new IllegalArgumentException("null ExpVector not allowed"); 524 } 525 if (c.signum() == 0) { 526 return this; 527 } 528 long d = k.totalDeg(); 529 MultiVarCoefficients<C> mc = lazyCoeffs; 530 HashMap<Long, GenPolynomial<C>> cc = new HashMap<Long, GenPolynomial<C>>(mc.coeffCache); 531 GenPolynomial<C> p = cc.get(d); 532 if (p == null) { 533 p = mc.pfac.getZERO(); 534 } 535 p = p.sum(c, k); 536 //System.out.println("p = " + p); 537 cc.put(d, p); 538 HashSet<ExpVector> z = new HashSet<ExpVector>(mc.zeroCache); 539 //System.out.println("z = " + z); 540 if (p.coefficient(k).isZERO() && !mc.homCheck.get((int) d)) { 541 z.add(k); 542 } 543 544 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(mc.pfac, cc, z, mc.homCheck) { 545 546 547 @Override 548 public C generate(ExpVector i) { 549 return coefficient(i); 550 } 551 }); 552 } 553 554 555 /** 556 * Subtract exponent vector and coefficient. 557 * @param k ExpVector. 558 * @param c coefficient. 559 * @return this - ONE.multiply(c,k). 560 */ 561 public MultiVarPowerSeries<C> subtract(final C c, final ExpVector k) { 562 if (k == null) { 563 throw new IllegalArgumentException("null ExpVector not allowed"); 564 } 565 if (c.signum() == 0) { 566 return this; 567 } 568 long d = k.totalDeg(); 569 MultiVarCoefficients<C> mc = lazyCoeffs; 570 HashMap<Long, GenPolynomial<C>> cc = new HashMap<Long, GenPolynomial<C>>(mc.coeffCache); 571 GenPolynomial<C> p = cc.get(d); 572 if (p == null) { 573 p = mc.pfac.getZERO(); 574 } 575 p = p.subtract(c, k); 576 cc.put(d, p); 577 HashSet<ExpVector> z = new HashSet<ExpVector>(mc.zeroCache); 578 //System.out.println("z = " + z); 579 if (p.coefficient(k).isZERO() && !mc.homCheck.get((int) d)) { 580 z.add(k); 581 } 582 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(mc.pfac, cc, z, mc.homCheck) { 583 584 585 @Override 586 public C generate(ExpVector i) { 587 return coefficient(i); 588 } 589 }); 590 } 591 592 593 /** 594 * Sum exponent vector and coefficient. 595 * @param mvc cached coefficients. 596 * @return this + mvc. 597 */ 598 public MultiVarPowerSeries<C> sum(MultiVarCoefficients<C> mvc) { 599 MultiVarCoefficients<C> mc = lazyCoeffs; 600 TreeMap<Long, GenPolynomial<C>> cc = new TreeMap<Long, GenPolynomial<C>>(mc.coeffCache); 601 TreeMap<Long, GenPolynomial<C>> ccv = new TreeMap<Long, GenPolynomial<C>>(mvc.coeffCache); 602 long d1 = (cc.size() > 0 ? cc.lastKey() : 0); 603 long d2 = (ccv.size() > 0 ? ccv.lastKey() : 0); 604 HashSet<ExpVector> z = new HashSet<ExpVector>(mc.zeroCache); 605 z.addAll(mvc.zeroCache); 606 long d = Math.max(d1, d2); 607 BitSet hc = new BitSet((int) d); 608 for (long i = 0; i <= d; i++) { 609 GenPolynomial<C> p1 = cc.get(i); 610 GenPolynomial<C> p2 = mvc.coeffCache.get(i); 611 if (p1 == null) { 612 p1 = mc.pfac.getZERO(); 613 } 614 if (p2 == null) { 615 p2 = mc.pfac.getZERO(); 616 } 617 GenPolynomial<C> p = p1.sum(p2); 618 //System.out.println("p = " + p); 619 cc.put(i, p); 620 if (mc.homCheck.get((int) i) && mvc.homCheck.get((int) i)) { 621 hc.set((int) i); 622 } else { 623 Set<ExpVector> ev = new HashSet<ExpVector>(p1.getMap().keySet()); 624 ev.addAll(p2.getMap().keySet()); 625 ev.removeAll(p.getMap().keySet()); 626 z.addAll(ev); 627 } 628 } 629 //System.out.println("cc = " + cc); 630 631 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(mc.pfac, 632 new HashMap<Long, GenPolynomial<C>>(cc), z, hc) { 633 634 635 @Override 636 public C generate(ExpVector i) { 637 return coefficient(i); 638 } 639 }); 640 } 641 642 643 /** 644 * Select coefficients. 645 * @param sel selector functor. 646 * @return new power series with selected coefficients. 647 */ 648 public MultiVarPowerSeries<C> select(final Selector<? super C> sel) { 649 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 650 651 652 @Override 653 public C generate(ExpVector i) { 654 C c = coefficient(i); 655 if (sel.select(c)) { 656 return c; 657 } 658 return ring.coFac.getZERO(); 659 } 660 }); 661 } 662 663 664 /** 665 * Shift select coefficients. Not selected coefficients are removed from the 666 * result series. 667 * @param sel selector functor. 668 * @return new power series with shifted selected coefficients. 669 */ 670 public MultiVarPowerSeries<C> shiftSelect(final Selector<? super C> sel) { 671 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 672 673 674 ExpVectorIterable ib = new ExpVectorIterable(ring.nvar, true, truncate); 675 676 677 Iterator<ExpVector> pos = ib.iterator(); 678 679 680 @Override 681 public C generate(ExpVector i) { 682 C c; 683 if (i.signum() > 0) { 684 int[] deps = i.dependencyOnVariables(); 685 ExpVector x = i.subst(deps[0], i.getVal(deps[0]) - 1L); 686 c = get(x); // ensure all coefficients are generated 687 } 688 do { 689 c = null; 690 if (pos.hasNext()) { 691 ExpVector e = pos.next(); 692 c = coefficient(e); 693 } else { 694 break; 695 } 696 } while (!sel.select(c)); 697 if (c == null) { // not correct because not known 698 c = ring.coFac.getZERO(); 699 } 700 return c; 701 } 702 }); 703 } 704 705 706 /** 707 * Map a unary function to this power series. 708 * @param f evaluation functor. 709 * @return new power series with coefficients f(this(i)). 710 */ 711 public MultiVarPowerSeries<C> map(final UnaryFunctor<? super C, C> f) { 712 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 713 714 715 @Override 716 public C generate(ExpVector i) { 717 return f.eval(coefficient(i)); 718 } 719 }); 720 } 721 722 723 /** 724 * Map a binary function to this and another power series. 725 * @param f evaluation functor with coefficients f(this(i),other(i)). 726 * @param ps other power series. 727 * @return new power series. 728 */ 729 public MultiVarPowerSeries<C> zip(final BinaryFunctor<? super C, ? super C, C> f, 730 final MultiVarPowerSeries<C> ps) { 731 int m = Math.min(ring.truncate, Math.max(truncate, ps.truncate())); 732 733 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 734 735 736 @Override 737 public C generate(ExpVector i) { 738 return f.eval(coefficient(i), ps.coefficient(i)); 739 } 740 }, m); 741 } 742 743 744 /** 745 * Sum of two power series, using zip(). 746 * @param ps other power series. 747 * @return this + ps. 748 */ 749 public MultiVarPowerSeries<C> sumZip(MultiVarPowerSeries<C> ps) { 750 return zip(new BinaryFunctor<C, C, C>() { 751 752 753 @Override 754 public C eval(C c1, C c2) { 755 return c1.sum(c2); 756 } 757 }, ps); 758 } 759 760 761 /** 762 * Subtraction of two power series, using zip(). 763 * @param ps other power series. 764 * @return this - ps. 765 */ 766 public MultiVarPowerSeries<C> subtractZip(MultiVarPowerSeries<C> ps) { 767 return zip(new BinaryFunctor<C, C, C>() { 768 769 770 @Override 771 public C eval(C c1, C c2) { 772 return c1.subtract(c2); 773 } 774 }, ps); 775 } 776 777 778 /** 779 * Multiply by coefficient. 780 * @param a coefficient. 781 * @return this * a. 782 */ 783 public MultiVarPowerSeries<C> multiply(final C a) { 784 if (a.isZERO()) { 785 return ring.getZERO(); 786 } 787 if (a.isONE()) { 788 return this; 789 } 790 return map(new UnaryFunctor<C, C>() { 791 792 793 @Override 794 public C eval(C c) { 795 return c.multiply(a); 796 } 797 }); 798 } 799 800 801 /** 802 * Monic. 803 * @return 1/orderCoeff() * this. 804 */ 805 public MultiVarPowerSeries<C> monic() { 806 ExpVector e = orderExpVector(); 807 if (e == null) { 808 return this; 809 } 810 C a = coefficient(e); 811 if (a.isONE()) { 812 return this; 813 } 814 if (a.isZERO()) { // sic 815 return this; 816 } 817 final C b = a.inverse(); 818 return map(new UnaryFunctor<C, C>() { 819 820 821 @Override 822 public C eval(C c) { 823 return b.multiply(c); 824 } 825 }); 826 } 827 828 829 /** 830 * Negate. 831 * @return - this. 832 */ 833 public MultiVarPowerSeries<C> negate() { 834 return map(new UnaryFunctor<C, C>() { 835 836 837 @Override 838 public C eval(C c) { 839 return c.negate(); 840 } 841 }); 842 } 843 844 845 /** 846 * Absolute value. 847 * @return abs(this). 848 */ 849 public MultiVarPowerSeries<C> abs() { 850 if (signum() < 0) { 851 return negate(); 852 } 853 return this; 854 } 855 856 857 /** 858 * Evaluate at given point. 859 * @return ps(a). 860 */ 861 public C evaluate(List<C> a) { 862 C v = ring.coFac.getZERO(); 863 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, truncate)) { 864 C c = coefficient(i); 865 if (c.isZERO()) { 866 continue; 867 } 868 c = c.multiply(i.evaluate(ring.coFac, a)); 869 v = v.sum(c); 870 } 871 return v; 872 } 873 874 875 /** 876 * Order. 877 * @return index of first non zero coefficient. 878 */ 879 public int order() { 880 if (order >= 0) { 881 return order; 882 } 883 // must compute it 884 GenPolynomial<C> p = null; 885 int t = 0; 886 while (lazyCoeffs.homCheck.get(t)) { 887 p = lazyCoeffs.coeffCache.get((long)t); 888 if (p == null || p.isZERO()) { // ?? 889 t++; 890 continue; 891 } 892 order = t; 893 evorder = p.trailingExpVector(); 894 //System.out.println("order = " + t); 895 return order; 896 } 897 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, truncate)) { 898 if (!coefficient(i).isZERO()) { 899 order = (int) i.totalDeg(); //ord; 900 evorder = i; 901 //System.out.println("order = " + order + ", evorder = " + evorder); 902 return order; 903 } 904 } 905 order = truncate + 1; 906 // evorder is null 907 return order; 908 } 909 910 911 /** 912 * Order ExpVector. 913 * @return ExpVector of first non zero coefficient. 914 */ 915 public ExpVector orderExpVector() { 916 //int x = 917 order(); // ensure evorder is set 918 return evorder; 919 } 920 921 922 /** 923 * Order monomial. 924 * @return first map entry or null. 925 */ 926 public Map.Entry<ExpVector, C> orderMonomial() { 927 ExpVector e = orderExpVector(); 928 if (e == null) { 929 return null; 930 } 931 //JAVA6only: 932 //return new AbstractMap.SimpleImmutableEntry<ExpVector, C>(e, coefficient(e)); 933 return new MapEntry<ExpVector, C>(e, coefficient(e)); 934 } 935 936 937 /** 938 * Truncate. 939 * @return truncate index of power series. 940 */ 941 public int truncate() { 942 return truncate; 943 } 944 945 946 /** 947 * Set truncate. 948 * @param t new truncate index. 949 * @return old truncate index of power series. 950 */ 951 public int setTruncate(int t) { 952 if (t < 0) { 953 throw new IllegalArgumentException("negative truncate not allowed"); 954 } 955 int ot = truncate; 956 if (order >= 0) { 957 if (order > truncate) { 958 order = -1; // reset 959 evorder = null; 960 } 961 } 962 truncate = t; 963 return ot; 964 } 965 966 967 /** 968 * Ecart. 969 * @return ecart. 970 */ 971 public long ecart() { 972 ExpVector e = orderExpVector(); 973 if (e == null) { 974 return 0L; 975 } 976 long d = e.totalDeg(); 977 long hd = d; 978 for (long i = d + 1L; i <= truncate; i++) { 979 if (!homogeneousPart(i).isZERO()) { 980 hd = i; 981 } 982 } 983 //System.out.println("d = " + d + ", hd = " + hd + ", coeffCache = " + lazyCoeffs.coeffCache + ", this = " + this); 984 return hd - d; 985 } 986 987 988 /** 989 * Signum. 990 * @return sign of first non zero coefficient. 991 */ 992 public int signum() { 993 ExpVector ev = orderExpVector(); 994 if (ev != null) { 995 return coefficient(ev).signum(); 996 } 997 return 0; 998 } 999 1000 1001 /** 1002 * Compare to. <b>Note: </b> compare only up to max(truncates). 1003 * @return sign of first non zero coefficient of this-ps. 1004 */ 1005 @Override 1006 public int compareTo(MultiVarPowerSeries<C> ps) { 1007 final int m = truncate(); 1008 final int n = ps.truncate(); 1009 final int pos = Math.min(ring.truncate, Math.min(m, n)); 1010 int s = 0; 1011 //System.out.println("coeffCache_c1 = " + lazyCoeffs.coeffCache); 1012 //System.out.println("coeffCache_c2 = " + ps.lazyCoeffs.coeffCache); 1013 // test homogeneous parts first is slower 1014 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, pos)) { 1015 s = coefficient(i).compareTo(ps.coefficient(i)); 1016 if (s != 0) { 1017 //System.out.println("i = " + i + ", coeff = " + coefficient(i) + ", ps.coeff = " + ps.coefficient(i)); 1018 return s; 1019 } 1020 } 1021 for (int j = pos + 1; j <= Math.min(ring.truncate, Math.max(m, n)); j++) { 1022 for (ExpVector i : new ExpVectorIterable(ring.nvar, j)) { 1023 s = coefficient(i).compareTo(ps.coefficient(i)); 1024 //System.out.println("i = " + i + ", coeff = " + coefficient(i) + ", ps.coeff = " + ps.coefficient(i)); 1025 if (s != 0) { 1026 //System.out.println("i = " + i + ", coeff = " + coefficient(i) + ", ps.coeff = " + ps.coefficient(i)); 1027 return s; 1028 } 1029 } 1030 } 1031 return s; 1032 } 1033 1034 1035 /** 1036 * Is power series zero. <b>Note: </b> compare only up to truncate. 1037 * @return If this is 0 then true is returned, else false. 1038 * @see edu.jas.structure.RingElem#isZERO() 1039 */ 1040 public boolean isZERO() { 1041 return (signum() == 0); 1042 } 1043 1044 1045 /** 1046 * Is power series one. <b>Note: </b> compare only up to truncate. 1047 * @return If this is 1 then true is returned, else false. 1048 * @see edu.jas.structure.RingElem#isONE() 1049 */ 1050 public boolean isONE() { 1051 if (!leadingCoefficient().isONE()) { 1052 return false; 1053 } 1054 return (compareTo(ring.ONE) == 0); 1055 //return reductum().isZERO(); 1056 } 1057 1058 1059 /** 1060 * Comparison with any other object. <b>Note: </b> compare only up to 1061 * truncate. 1062 * @see java.lang.Object#equals(java.lang.Object) 1063 */ 1064 @Override 1065 @SuppressWarnings("unchecked") 1066 public boolean equals(Object B) { 1067 if (B == null) { 1068 return false; 1069 } 1070 if (!(B instanceof MultiVarPowerSeries)) { 1071 return false; 1072 } 1073 MultiVarPowerSeries<C> a = (MultiVarPowerSeries<C>) B; 1074 return compareTo(a) == 0; 1075 } 1076 1077 1078 /** 1079 * Hash code for this polynomial. <b>Note: </b> only up to truncate. 1080 * @see java.lang.Object#hashCode() 1081 */ 1082 @Override 1083 public int hashCode() { 1084 int h = 0; 1085 //h = ( ring.hashCode() << 23 ); 1086 //h += val.hashCode(); 1087 for (ExpVector i : new ExpVectorIterable(ring.nvar, true, truncate)) { 1088 C c = coefficient(i); 1089 if (!c.isZERO()) { 1090 h += i.hashCode(); 1091 h = (h << 23); 1092 } 1093 h += c.hashCode(); 1094 h = (h << 23); 1095 } 1096 return h; 1097 } 1098 1099 1100 /** 1101 * Is unit. 1102 * @return true, if this power series is invertible, else false. 1103 */ 1104 public boolean isUnit() { 1105 return leadingCoefficient().isUnit(); 1106 } 1107 1108 1109 /** 1110 * Sum a another power series. 1111 * @param ps other power series. 1112 * @return this + ps. 1113 */ 1114 public MultiVarPowerSeries<C> sum(final MultiVarPowerSeries<C> ps) { 1115 //final MultiVarPowerSeries<C> ps1 = this; // method name was ambiguous in generate 1116 int nt = Math.min(ring.truncate, Math.max(truncate(), ps.truncate())); 1117 1118 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1119 1120 1121 @Override 1122 public C generate(ExpVector e) { 1123 long tdeg = e.totalDeg(); 1124 if (lazyCoeffs.homCheck.get((int) tdeg)) { 1125 // generate respective homogeneous polynomial 1126 GenPolynomial<C> p = homogeneousPart(tdeg).sum(ps.homogeneousPart(tdeg)); 1127 coeffCache.put(tdeg, p); // overwrite 1128 homCheck.set((int) tdeg); 1129 C c = p.coefficient(e); 1130 //System.out.println("c = " + c + ", e = " + e + ", tdeg = " + tdeg); 1131 return c; 1132 } 1133 return coefficient(e).sum(ps.coefficient(e)); 1134 } 1135 }, nt); 1136 } 1137 1138 1139 /** 1140 * Subtract a another power series. 1141 * @param ps other power series. 1142 * @return this - ps. 1143 */ 1144 public MultiVarPowerSeries<C> subtract(final MultiVarPowerSeries<C> ps) { 1145 //final MultiVarPowerSeries<C> ps1 = this; // method name was ambiguous in generate 1146 int nt = Math.min(ring.truncate, Math.max(truncate(), ps.truncate())); 1147 1148 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1149 1150 1151 @Override 1152 public C generate(ExpVector e) { 1153 long tdeg = e.totalDeg(); 1154 if (lazyCoeffs.homCheck.get((int) tdeg)) { 1155 // generate respective homogeneous polynomial 1156 GenPolynomial<C> p = homogeneousPart(tdeg).subtract(ps.homogeneousPart(tdeg)); 1157 coeffCache.put(tdeg, p); // overwrite 1158 homCheck.set((int) tdeg); 1159 C c = p.coefficient(e); 1160 //System.out.println("p = " + p + ", e = " + e + ", tdeg = " + tdeg); 1161 return c; 1162 } 1163 return coefficient(e).subtract(ps.coefficient(e)); 1164 } 1165 }, nt); 1166 } 1167 1168 1169 /** 1170 * Multiply by another power series. 1171 * @param ps other power series. 1172 * @return this * ps. 1173 */ 1174 public MultiVarPowerSeries<C> multiply(final MultiVarPowerSeries<C> ps) { 1175 //final MultiVarPowerSeries<C> ps1 = this; // method name was ambiguous in generate 1176 int nt = Math.min(ring.truncate, truncate() + ps.truncate()); 1177 1178 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1179 1180 1181 @Override 1182 public C generate(ExpVector e) { 1183 long tdeg = e.totalDeg(); 1184 // generate respective homogeneous polynomial 1185 GenPolynomial<C> p = null; //fac.getZERO(); 1186 for (int k = 0; k <= tdeg; k++) { 1187 GenPolynomial<C> m = homogeneousPart(k).multiply(ps.homogeneousPart(tdeg - k)); 1188 if (p == null) { 1189 p = m; 1190 } else { 1191 p = p.sum(m); 1192 } 1193 } 1194 coeffCache.put(tdeg, p); // overwrite 1195 homCheck.set((int) tdeg); 1196 C c = p.coefficient(e); 1197 return c; 1198 } 1199 }, nt); 1200 } 1201 1202 1203 /** 1204 * Inverse power series. 1205 * @return ps with this * ps = 1. 1206 */ 1207 public MultiVarPowerSeries<C> inverse() { 1208 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1209 1210 1211 @Override 1212 public C generate(ExpVector e) { 1213 long tdeg = e.totalDeg(); 1214 C d = leadingCoefficient().inverse(); // may fail 1215 if (tdeg == 0) { 1216 return d; 1217 } 1218 GenPolynomial<C> p = null; //fac.getZERO(); 1219 for (int k = 0; k < tdeg; k++) { 1220 GenPolynomial<C> m = getHomPart(k).multiply(homogeneousPart(tdeg - k)); 1221 if (p == null) { 1222 p = m; 1223 } else { 1224 p = p.sum(m); 1225 } 1226 } 1227 p = p.multiply(d.negate()); 1228 //System.out.println("tdeg = " + tdeg + ", p = " + p); 1229 coeffCache.put(tdeg, p); // overwrite 1230 homCheck.set((int) tdeg); 1231 C c = p.coefficient(e); 1232 return c; 1233 } 1234 }); 1235 } 1236 1237 1238 /** 1239 * Divide by another power series. 1240 * @param ps nonzero power series with invertible coefficient. 1241 * @return this / ps. 1242 */ 1243 public MultiVarPowerSeries<C> divide(MultiVarPowerSeries<C> ps) { 1244 if (ps.isUnit()) { 1245 return multiply(ps.inverse()); 1246 } 1247 int m = order(); 1248 int n = ps.order(); 1249 if (m < n) { 1250 return ring.getZERO(); 1251 } 1252 ExpVector em = orderExpVector(); 1253 ExpVector en = ps.orderExpVector(); 1254 if (!ps.coefficient(en).isUnit()) { 1255 throw new ArithmeticException("division by non unit coefficient " + ps.coefficient(ps.evorder) 1256 + ", evorder = " + ps.evorder); 1257 } 1258 // now m >= n 1259 MultiVarPowerSeries<C> st, sps, q, sq; 1260 if (m == 0) { 1261 st = this; 1262 } else { 1263 st = this.shift(em.negate()); 1264 } 1265 if (n == 0) { 1266 sps = ps; 1267 } else { 1268 sps = ps.shift(en.negate()); 1269 } 1270 q = st.multiply(sps.inverse()); 1271 sq = q.shift(em.subtract(en)); 1272 return sq; 1273 } 1274 1275 1276 /** 1277 * Power series remainder. 1278 * @param ps nonzero power series with invertible leading coefficient. 1279 * @return remainder with this = quotient * ps + remainder. 1280 */ 1281 public MultiVarPowerSeries<C> remainder(MultiVarPowerSeries<C> ps) { 1282 int m = order(); 1283 int n = ps.order(); 1284 if (m >= n) { 1285 return ring.getZERO(); 1286 } 1287 return this; 1288 } 1289 1290 1291 /** 1292 * Quotient and remainder by division of this by S. 1293 * @param S a MultiVarPowerSeries 1294 * @return [this/S, this - (this/S)*S]. 1295 */ 1296 public MultiVarPowerSeries<C>[] quotientRemainder(MultiVarPowerSeries<C> S) { 1297 return new MultiVarPowerSeries[] { divide(S), remainder(S) }; 1298 } 1299 1300 1301 /** 1302 * Differentiate with respect to variable r. 1303 * @param r variable for the direction. 1304 * @return differentiate(this). 1305 */ 1306 public MultiVarPowerSeries<C> differentiate(final int r) { 1307 if (r < 0 || ring.nvar < r) { 1308 throw new IllegalArgumentException("variable index out of bound"); 1309 } 1310 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1311 1312 1313 @Override 1314 public C generate(ExpVector i) { 1315 long d = i.getVal(r); 1316 ExpVector e = i.subst(r, d + 1); 1317 C v = coefficient(e); 1318 v = v.multiply(ring.coFac.fromInteger(d + 1)); 1319 return v; 1320 } 1321 }); 1322 } 1323 1324 1325 /** 1326 * Integrate with respect to variable r and with given constant. 1327 * @param c integration constant. 1328 * @param r variable for the direction. 1329 * @return integrate(this). 1330 */ 1331 public MultiVarPowerSeries<C> integrate(final C c, final int r) { 1332 if (r < 0 || ring.nvar < r) { 1333 throw new IllegalArgumentException("variable index out of bound"); 1334 } 1335 int nt = Math.min(ring.truncate, truncate + 1); 1336 return new MultiVarPowerSeries<C>(ring, new MultiVarCoefficients<C>(ring) { 1337 1338 1339 @Override 1340 public C generate(ExpVector i) { 1341 if (i.isZERO()) { 1342 return c; 1343 } 1344 long d = i.getVal(r); 1345 if (d > 0) { 1346 ExpVector e = i.subst(r, d - 1); 1347 C v = coefficient(e); 1348 v = v.divide(ring.coFac.fromInteger(d)); 1349 return v; 1350 } 1351 return ring.coFac.getZERO(); 1352 } 1353 }, nt); 1354 } 1355 1356 1357 /** 1358 * Power series greatest common divisor. 1359 * @param ps power series. 1360 * @return gcd(this,ps). 1361 */ 1362 public MultiVarPowerSeries<C> gcd(MultiVarPowerSeries<C> ps) { 1363 if (ps.isZERO()) { 1364 return this; 1365 } 1366 if (this.isZERO()) { 1367 return ps; 1368 } 1369 ExpVector em = orderExpVector(); 1370 ExpVector en = ps.orderExpVector(); 1371 return ring.getONE().shift(em.gcd(en)); 1372 } 1373 1374 1375 /** 1376 * Power series extended greatest common divisor. <b>Note:</b> not 1377 * implemented. 1378 * @param S power series. 1379 * @return [ gcd(this,S), a, b ] with a*this + b*S = gcd(this,S). 1380 */ 1381 //SuppressWarnings("unchecked") 1382 public MultiVarPowerSeries<C>[] egcd(MultiVarPowerSeries<C> S) { 1383 throw new UnsupportedOperationException("egcd for power series not implemented"); 1384 } 1385 1386}