Skip to content

Commit

Permalink
Merge pull request #20002 from sachin-4099/gsoc#17
Browse files Browse the repository at this point in the history
[GSoC] Functions: Fixes _eval_nseries() of log
  • Loading branch information
jksuom committed Aug 26, 2020
2 parents 0dcb88d + d48ecfe commit d20fcb2
Show file tree
Hide file tree
Showing 3 changed files with 58 additions and 10 deletions.
58 changes: 50 additions & 8 deletions sympy/functions/elementary/exponential.py
Expand Up @@ -919,6 +919,7 @@ def _eval_nseries(self, x, n, logx, cdir=0):
# NOTE Please see the comment at the beginning of this file, labelled
# IMPORTANT.
from sympy import im, cancel, I, Order, logcombine
from itertools import product
if not logx:
logx = log(x)
if self.args[0] == x:
Expand All @@ -932,6 +933,20 @@ def _eval_nseries(self, x, n, logx, cdir=0):
r = log(k) + l*logx # XXX true regardless of assumptions?
return r

def coeff_exp(term, x):
coeff, exp = S.One, S.Zero
for factor in Mul.make_args(term):
if factor.has(x):
base, exp = factor.as_base_exp()
if base != x:
try:
return term.leadterm(x)
except ValueError:
return term, S.Zero
else:
coeff *= factor
return coeff, exp

# TODO new and probably slow
try:
a, b = arg.leadterm(x)
Expand All @@ -942,22 +957,49 @@ def _eval_nseries(self, x, n, logx, cdir=0):
n += 1
s = arg.nseries(x, n=n, logx=logx)
a, b = s.removeO().leadterm(x)
p = cancel(s/(a*x**b) - 1)
p = cancel(s/(a*x**b) - 1).expand().powsimp()
if p.has(exp):
p = logcombine(p)
g = None
l = []
for i in range(n + 2):
g = log.taylor_term(i, p, g)
g = g.nseries(x, n=n, logx=logx)
l.append(g)
if isinstance(p, Order):
n = p.getn()
_, d = coeff_exp(p, x)
if not d.is_positive:
return log(a) + b*logx + Order(x**n, x)

def mul(d1, d2):
res = {}
for e1, e2 in product(d1, d2):
ex = e1 + e2
if ex < n:
res[ex] = res.get(ex, S.Zero) + d1[e1]*d2[e2]
return res

pterms = {}

for term in Add.make_args(p):
co1, e1 = coeff_exp(term, x)
pterms[e1] = pterms.get(e1, S.Zero) + co1.removeO()

k = S.One
terms = {}
pk = pterms

while k*d < n:
coeff = -(-1)**k/k
for ex in pk:
terms[ex] = terms.get(ex, S.Zero) + coeff*pk[ex]
pk = mul(pk, pterms)
k += S.One

res = log(a) + b*logx
for ex in terms:
res += terms[ex]*x**(ex)

if cdir != 0:
cdir = self.args[0].dir(x, cdir)
if a.is_real and a.is_negative and im(cdir) < 0:
res -= 2*I*S.Pi
return res + Add(*l) + Order(p**n, x)
return res + Order(x**n, x)

def _eval_as_leading_term(self, x, cdir=0):
from sympy import I, im
Expand Down
4 changes: 2 additions & 2 deletions sympy/series/tests/test_nseries.py
Expand Up @@ -277,8 +277,8 @@ def test_issue_3463():
r = log(5)/log(3)
p = w**(-1 + r)
e = 1/x*(-log(w**(1 + r)) + log(w + w**r))
e_ser = -r*log(w)/x + p/x - p**2/(2*x) + O(p**3)
assert e.nseries(w, n=3) == e_ser
e_ser = -r*log(w)/x + p/x - p**2/(2*x) + O(w)
assert e.nseries(w, n=1) == e_ser


def test_sin():
Expand Down
6 changes: 6 additions & 0 deletions sympy/series/tests/test_series.py
Expand Up @@ -257,6 +257,12 @@ def test_issue_18008():
O(x**(-4), (x, oo))


def test_issue_18842():
f = log(x/(1 - x))
assert f.series(x, 0.491, n=1).removeO().nsimplify() == \
-S(180019443780011)/5000000000000000


def test_issue_19534():
dt = symbols('dt', real=True)
expr = 16*dt*(0.125*dt*(2.0*dt + 1.0) + 0.875*dt + 1.0)/45 + \
Expand Down

0 comments on commit d20fcb2

Please sign in to comment.