#1956: implement multivariate power series arithmetic
-------------------------------------------+--------------------------------
Reporter: was | Owner: malb
Type: enhancement | Status: needs_review
Priority: major | Milestone: sage-4.6
Component: commutative algebra | Keywords: multivariate power
series
Author: Niles Johnson | Upstream: N/A
Reviewer: Martin Albrecht, Simon King | Merged:
Work_issues: |
-------------------------------------------+--------------------------------
Comment(by niles):
Replying to [comment:47 niles]:
Mike, would you be able to take a look at this some time?
Here are some further comparisons with polynomial arithmetic, computing
``f**k`` for various values of ``k``. Obviously the power series
computation has a greater and greater advantage over the polynomial
computation as ``k`` grows, since the polynomial arithmetic computes many
more (irrelevant) terms. So the time tests below aren't exactly fair --
one could compute the necessary terms by computing small polynomial
powers, truncating, and repeating -- but part of the point of the power
series code is to automate this for the user. One side effect is that
arithmetic for small powers is slower, but for large powers the power
series arithmetic is fairly good: below, I have a comparison in a simple
special case: ``k = 2**n`` and non-zero constant coefficient.
In any case, the point of this ticket is to implement the basic
arithmetic. The code should be fast enough for basic use, and speed
improvement can be part of a separate ticket.
== Timings for ``f**k`` v.s. ``fp**k`` ==
{{{
sage: BaseRingList = [ZZ,QQ,GF(11),GF(65537)]
sage: for B in BaseRingList:
R = PowerSeriesRing(B,4,'t'); R
f = R.random_element(prec=10,bound=20)
fp = f.polynomial()
for k in [2,5,8,11]:
print 'f**'+str(k)+': (power series)'
timeit('f**k',number=5)
print 'fp**'+str(k)+': (polynomial)'
timeit('fp**k',number=5)
print '--'
....:
Multivariate Power Series Ring in t0, t1, t2, t3 over Integer Ring
f**2: (power series)
5 loops, best of 3: 879 µs per loop
fp**2: (polynomial)
5 loops, best of 3: 33.2 µs per loop
--
f**5: (power series)
5 loops, best of 3: 4.12 ms per loop
fp**5: (polynomial)
5 loops, best of 3: 2.52 ms per loop
--
f**8: (power series)
5 loops, best of 3: 6.47 ms per loop
fp**8: (polynomial)
5 loops, best of 3: 51.5 ms per loop
--
f**11: (power series)
5 loops, best of 3: 10.8 ms per loop
fp**11: (polynomial)
5 loops, best of 3: 652 ms per loop
--
Multivariate Power Series Ring in t0, t1, t2, t3 over Rational Field
f**2: (power series)
5 loops, best of 3: 1 ms per loop
fp**2: (polynomial)
5 loops, best of 3: 44.8 µs per loop
--
f**5: (power series)
5 loops, best of 3: 8.68 ms per loop
fp**5: (polynomial)
5 loops, best of 3: 4.12 ms per loop
--
f**8: (power series)
5 loops, best of 3: 26.4 ms per loop
fp**8: (polynomial)
5 loops, best of 3: 119 ms per loop
--
f**11: (power series)
5 loops, best of 3: 58.4 ms per loop
fp**11: (polynomial)
5 loops, best of 3: 2.28 s per loop
--
Multivariate Power Series Ring in t0, t1, t2, t3 over Finite Field of size
11
f**2: (power series)
5 loops, best of 3: 480 µs per loop
fp**2: (polynomial)
5 loops, best of 3: 10.4 µs per loop
--
f**5: (power series)
5 loops, best of 3: 2.49 ms per loop
fp**5: (polynomial)
5 loops, best of 3: 547 µs per loop
--
f**8: (power series)
5 loops, best of 3: 4.23 ms per loop
fp**8: (polynomial)
5 loops, best of 3: 11 ms per loop
--
f**11: (power series)
5 loops, best of 3: 8.6 ms per loop
fp**11: (polynomial)
5 loops, best of 3: 190 ms per loop
--
Multivariate Power Series Ring in t0, t1, t2, t3 over Finite Field of size
65537
f**2: (power series)
5 loops, best of 3: 614 µs per loop
fp**2: (polynomial)
5 loops, best of 3: 18.8 µs per loop
--
f**5: (power series)
5 loops, best of 3: 4.37 ms per loop
fp**5: (polynomial)
5 loops, best of 3: 3.13 ms per loop
--
f**8: (power series)
5 loops, best of 3: 16.5 ms per loop
fp**8: (polynomial)
5 loops, best of 3: 231 ms per loop
--
f**11: (power series)
5 loops, best of 3: 29.7 ms per loop
fp**11: (polynomial)
5 loops, best of 3: 4.85 s per loop
--
}}}
== Timings for ``f**(2**n)`` v.s. iterated square-and-truncate for
``fp**(2**n)`` ==
Here is a simple special case: for a power series `f` with non-zero
constant coefficient, the total-degree precision of `f**k` is equal to the
total-degree precision of `f` for all `k`. Below is a function which uses
polynomial arithmetic to compute coefficients of `fp**k` up to given total
degree for `k = 2^n`, where `fp = f.polynomial()`.
{{{
sage: def tot_deg_truncate(g,d):
return sum(g.monomial_coefficient(m)*m for m in [x for x in
g.monomials() if x.degree() < d])
....:
sage: def pow_recursive(g,k,d):
if k == 2:
return tot_deg_truncate(g**2,d)
else:
return pow_recursive(tot_deg_truncate(g**2,d),k/2,d)
....:
sage: R = PowerSeriesRing(ZZ,4,'t')
sage: f = 1 + R.random_element(prec=10,bound=20)
sage: f.constant_coefficient()
1
sage: fp = f.polynomial()
}}}
Check that ``pow_recursive`` gets the right answer:
{{{
sage: f**16 - pow_recursive(fp,16,10)
0 + O(t0, t1, t2, t3)^10
sage: f**64 - pow_recursive(fp,64,10)
0 + O(t0, t1, t2, t3)^10
}}}
Compare timings; although the polynomial arithmetic is faster, they are on
roughly the same order of magnitude.
{{{
sage: %timeit f**16
125 loops, best of 3: 4.8 ms per loop
sage: %timeit pow_recursive(fp,16,10)
625 loops, best of 3: 1.18 ms per loop
sage: %timeit f**64
125 loops, best of 3: 7.21 ms per loop
sage: %timeit pow_recursive(fp,64,10)
125 loops, best of 3: 1.76 ms per loop
}}}
Increasing the precision exaggerates the difference:
{{{
sage: f = 1 + R.random_element(prec=30,bound=20)
sage: f.constant_coefficient()
1
sage: fp = f.polynomial()
sage: f**16 - pow_recursive(fp,16,30)
0 + O(t0, t1, t2, t3)^30
sage: f**64 - pow_recursive(fp,64,30)
0 + O(t0, t1, t2, t3)^30
sage: %timeit f**16
25 loops, best of 3: 15.9 ms per loop
sage: %timeit pow_recursive(fp,16,30)
125 loops, best of 3: 1.81 ms per loop
sage: %timeit f**64
25 loops, best of 3: 23.9 ms per loop
sage: %timeit pow_recursive(fp,64,30)
125 loops, best of 3: 2.65 ms per loop
}}}
--
Ticket URL: <http://trac.sagemath.org/sage_trac/ticket/1956#comment:49>
Sage <http://www.sagemath.org>
Sage: Creating a Viable Open Source Alternative to Magma, Maple, Mathematica,
and MATLAB
--
You received this message because you are subscribed to the Google Groups
"sage-trac" group.
To post to this group, send email to [email protected].
To unsubscribe from this group, send email to
[email protected].
For more options, visit this group at
http://groups.google.com/group/sage-trac?hl=en.
