Opened 9 years ago

Closed 9 years ago

#12241 closed enhancement (fixed)

exp, log, derivative of multivariate power series

Reported by: vbraun Owned by: malb
Priority: major Milestone: sage-5.0
Component: commutative algebra Keywords: multivariate power series
Cc: niles Merged in: sage-5.0.beta6
Authors: Volker Braun Reviewers: Niles Johnson
Report Upstream: N/A Work issues:
Branch: Commit:
Dependencies: Stopgaps:

Description

This ticket implements exp, log, and derivative of multivariate power series:

sage: T.<a,b> = PowerSeriesRing(ZZ,2)
sage: f = a + b + a^2*b
sage: f = a + b + a^2*b + T.O(4)
sage: f.derivative(a)
1 + 2*a*b + O(a, b)^3
sage: exp(f)
1 + a + b + 1/2*a^2 + a*b + 1/2*b^2 + 1/6*a^3 + 3/2*a^2*b + 1/2*a*b^2 + 1/6*b^3 + O(a, b)^4
sage: log(1+f)
a + b - 1/2*a^2 - a*b - 1/2*b^2 + 1/3*a^3 + 2*a^2*b + a*b^2 + 1/3*b^3 + O(a, b)^4

Attachments (1)

trac_12241_more_multivariate_powerseries.patch (10.5 KB) - added by vbraun 9 years ago.
Updated patch

Download all attachments as: .zip

Change History (7)

comment:1 Changed 9 years ago by vbraun

  • Cc niles added
  • Keywords multivariate power series added
  • Status changed from new to needs_review

comment:2 Changed 9 years ago by vbraun

Fixed bug where valuation >= 2 was exponentiated with the wrong precision.

comment:3 Changed 9 years ago by niles

  • Reviewers set to Niles Johnson
  • Status changed from needs_review to needs_work
  • Work issues set to efficiency

This is good -- thanks for working on it! Here are some comments:

derivative

Looks great!

pow

  • I think the default behavior should be same as repeated multiplication by self -- is there a good reason not to do this?
  • Also, when precision is much larger than valuation, pow is slower than repeated multiplication anyway:
    sage: R.<a,b,c> = PowerSeriesRing(ZZ)
    sage: f = a+R.random_element(100)
    sage: f
    a - a^23*c^11 + a^29*b^8*c^55 + 2*a^36*b^54*c^6 + O(a, b, c)^100
    sage: %timeit f^10
    5 loops, best of 3: 95.8 ms per loop
    sage: %timeit f.pow(10)
    5 loops, best of 3: 215 ms per loop
    

I think this should be addressed by rewriting pow to use _bg_value, as has been done with _mul_.

  • Unexpected behavior when requesting larger precision from pow:
    sage: f = a*R.random_element(); f
    2*a^3*b - 7*a^3*b^4*c^3 - 3*a^6*b^6 - 109*a^3*b^3*c^6 - 2*a*b^9*c^2 + O(a, b, c)^13
    sage: f^2
    4*a^6*b^2 - 28*a^6*b^5*c^3 - 12*a^9*b^7 - 436*a^6*b^4*c^6 - 8*a^4*b^10*c^2 + O(a, b, c)^17
    sage: f.pow(2,prec=17)
    4*a^6*b^2 + O(a, b, c)^13
    

exp / log

  • The warning about constant coefficients is confusing and vague. I would prefer something more straightforward like: "If f has constant coefficient c, and exp(c) is transcendental, then exp(f) would have to be a power series over the Symbolic Ring. These are not yet implemented and therefore such cases raise an error:"
  • Another workaround for this limitation is to change base ring to one which is closed under exponentiation, such as RR or CC:
    sage: R.<a,b,c> = PowerSeriesRing(ZZ)
    sage: f = 2+R.random_element()
    sage: exp(f)
     . . .
    TypeError: unsupported operand parent(s) for '*': 'Symbolic Ring' and 'Multivariate Power Series Ring in a, b, c over Rational Field'
    sage: S = R.change_ring(CC)
    sage: f = S(f)
    sage: exp(f)
    7.38905609893065 + (-7.38905609893065)*a*b^3*c + (-22.1671682967919)*a^6*b^2*c + (-7.38905609893065)*a^5*b^2*c^2 + (-7.38905609893065)*b^2*c^7 + 3.69452804946533*a^2*b^6*c^2 + (-7.38905609893065)*a*b^10 + O(a, b, c)^12
    
  • For elements with infinite precision, consider using default_prec if no precision is specified:
    sage: R.<a,b,c> = PowerSeriesRing(ZZ)
    sage: R.default_prec()
    12
    
  • The implementation seems unnecessarily slow, because of repeated calls to pow: For precision n, calling pow for each term will require a total of something like n^2 multiplications of the input with itself. But these functions should require only n. Here is a suggested alternate implementation for log:
    def fast_log(f):
        prec = f.prec()
        x = 1-f
        r = 1
        logx = 0
        for i in range(1,prec):
            r = r*x + R.O(prec)
            log += r/i
        return -logx
    
    sage: f = 1+R.random_element(30); f
    1 - 6*a^5*b^2*c^7 + a^14*b^13 + O(a, b, c)^30
    sage: %timeit log(f)
    5 loops, best of 3: 224 ms per loop
    sage: %timeit fast_log(f)
    25 loops, best of 3: 26.4 ms per loop
    sage: fast_log(f) == log(f)
    True
    

And this might be faster if you use _bg_value to do the multiplication here too.

comment:4 Changed 9 years ago by vbraun

  • Status changed from needs_work to needs_review

I've removed the pow() method and sped up log/exp:

sage: R.<a,b,c> = PowerSeriesRing(ZZ)
sage: f = 1+R.random_element(30); f
1 - a^8*b^2*c^3 + 2*a*b^6*c^8 - 3*a^2*b^16*c - 17*a^13*b^3*c^8 + 2*a^2*b^8*c^15 + O(a, b, c)^30
sage: %timeit log(f)
125 loops, best of 3: 3.6 ms per loop

Changed 9 years ago by vbraun

Updated patch

comment:5 Changed 9 years ago by niles

  • Status changed from needs_review to positive_review
  • Work issues efficiency deleted

Looks good, passes all tests, positive review!

comment:6 Changed 9 years ago by jdemeyer

  • Merged in set to sage-5.0.beta6
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.