Opened 7 years ago

Closed 5 years ago

#16198 closed defect (fixed)

allow constant != 1 in log(power series)

Reported by: rws Owned by:
Priority: minor Milestone: sage-7.0
Component: commutative algebra Keywords: log, function, series expansion
Cc: Merged in:
Authors: Ralf Stephan Reviewers: Volker Braun
Report Upstream: N/A Work issues:
Branch: ee40e93 (Commits, GitHub, GitLab) Commit: ee40e93a8b09aa20721240e591c6cf4332abb06d
Dependencies: Stopgaps:

Status badges

Description (last modified by jdemeyer)

log(power series) fails when the constant coefficient is not 1:

sage: R.<x> = PowerSeriesRing(RR)
sage: log(2+x)
---------------------------------------------------------------------------
ArithmeticError                           Traceback (most recent call last)
<ipython-input-4-49ca3295e0d6> in <module>()
----> 1 log(Integer(2)+x)

/usr/local/src/sage-git/local/lib/python2.7/site-packages/sage/functions/log.pyc in __call__(self, *args, **kwds)
    311         if base is None:
    312             if len(args) == 1:
--> 313                 return GinacFunction.__call__(self, *args, **kwds)
    314             # second argument is base
    315             base = args[1]

/usr/local/src/sage-git/local/lib/python2.7/site-packages/sage/symbolic/function.so in sage.symbolic.function.BuiltinFunction.__call__ (build/cythonized/sage/symbolic/function.cpp:9120)()

/usr/local/src/sage-git/local/lib/python2.7/site-packages/sage/rings/power_series_ring_element.so in sage.rings.power_series_ring_element.PowerSeries.log (build/cythonized/sage/rings/power_series_ring_element.c:15894)()

ArithmeticError: constant term of power series is not 1

Pari has no problems:

? log(2-x)
%1 = 0.69314718055994530941723212145817656807 - 1/2*x - 1/8*x^2 - 1/24*x^3 - 1/64*x^4 - 1/160*x^5 - 1/384*x^6 - 1/896*x^7 - 1/2048*x^8 - 1/4608*x^9 - 1/10240*x^10 - 1/22528*x^11 - 1/49152*x^12 - 1/106496*x^13 - 1/229376*x^14 - 1/491520*x^15 + O(x^16)

Effectively the log of the constant is taken but this simply changes the base ring from whatever to SR if the constant is not equal to one.

Change History (12)

comment:1 Changed 7 years ago by kcrisman

I don't think Ginac is actually involved. Take a look at the last piece of the traceback, which points to here. In particular, in sage/rings/power_series_ring_element.pyx,

        if prec is None:
            prec = self._parent.default_prec()

        if not self[0].is_one():
            raise ArithmeticError("constant term of power series is not 1")

        zero = self.parent().zero()
        t = zero.solve_linear_de(prec,b=self.derivative()/self,f0=0)
        return t

So this is actually hard-coded in for some reason, referring eventually to _solve_linear_de which was factored out of this stuff a very long time ago.

It could be worth asking around why/whether this is really necessary. I don't know that I want the above instead of

sage: log(2-x)
log(-x + 2)
sage: _.taylor(x,0,16)
-1/1048576*x^16 - 1/491520*x^15 - 1/229376*x^14 - 1/106496*x^13 - 1/49152*x^12 - 1/22528*x^11 - 1/10240*x^10 - 1/4608*x^9 - 1/2048*x^8 - 1/896*x^7 - 1/384*x^6 - 1/160*x^5 - 1/64*x^4 - 1/24*x^3 - 1/8*x^2 - 1/2*x + log(2)

note the log(2) instead of an approximation.

comment:2 follow-up: Changed 7 years ago by rws

  • Description modified (diff)
  • Priority changed from major to minor
  • Summary changed from replace default algorithm for log(power series) to allow constant <> 1 in log(power series)

So it looks like any constant from exp(n), n in NN0 will leave this in PowerSeries(ZZ) but anything else will make it necessary to rebase in PowerSeries(SR).

comment:3 Changed 7 years ago by vbraun_spam

  • Milestone changed from sage-6.2 to sage-6.3

comment:4 Changed 7 years ago by vbraun_spam

  • Milestone changed from sage-6.3 to sage-6.4

comment:5 Changed 6 years ago by rws

  • Milestone changed from sage-6.4 to sage-duplicate/invalid/wontfix
  • Status changed from new to needs_review

This was ill-conceived from the beginning, because only SR is a superset of the coefficients in the Pari result (floating point and rational). That's why Pari is questionable, even if it's convenient.

As I find it not useful to duplicate symbolics functionality in the power series ring, the best way to get such a power series is via SR and the newly available #16203:

            sage: R.<x> = PowerSeriesRing(SR)
            sage: ex=(log(2-y)).series(y,4); R(ex)
            log(2) - 1/2*x - 1/8*x^2 - 1/24*x^3 + O(x^4)
            sage: ex=(gamma(1-y)).series(y,3); R(ex)
            1 + euler_gamma*x + (1/2*euler_gamma^2 + 1/12*pi^2)*x^2 + O(x^3)

comment:6 Changed 6 years ago by rws

  • Status changed from needs_review to positive_review

comment:7 in reply to: ↑ 2 Changed 6 years ago by jdemeyer

  • Milestone changed from sage-duplicate/invalid/wontfix to sage-6.5
  • Status changed from positive_review to needs_work
  • Summary changed from allow constant <> 1 in log(power series) to allow constant != 1 in log(power series)

Replying to rws:

So it looks like any constant from exp(n), n in NN0 will leave this in PowerSeries(ZZ) but anything else will make it necessary to rebase in PowerSeries(SR).

If the base ring is RR, then log() is perfectly defined for non-one constant terms, so that should certainly work.

comment:8 Changed 6 years ago by jdemeyer

  • Description modified (diff)

comment:9 Changed 5 years ago by rws

  • Branch set to u/rws/allow_constant____1_in_log_power_series_

comment:10 Changed 5 years ago by rws

  • Authors set to Ralf Stephan
  • Commit set to ee40e93a8b09aa20721240e591c6cf4332abb06d
  • Component changed from calculus to commutative algebra
  • Milestone changed from sage-6.5 to sage-7.0
  • Status changed from needs_work to needs_review

New commits:

ee40e9316198: allow constant != 1 in log(power series)

comment:11 Changed 5 years ago by vbraun

  • Reviewers set to Volker Braun
  • Status changed from needs_review to positive_review

comment:12 Changed 5 years ago by vbraun

  • Branch changed from u/rws/allow_constant____1_in_log_power_series_ to ee40e93a8b09aa20721240e591c6cf4332abb06d
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.