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