Coverage for /usr/lib/python3/dist-packages/sympy/series/order.py: 12%

1from sympy.core import S, sympify, Expr, Dummy, Add, Mul

2from sympy.core.cache import cacheit

3from sympy.core.containers import Tuple

4from sympy.core.function import Function, PoleError, expand_power_base, expand_log

5from sympy.core.sorting import default_sort_key

6from sympy.functions.elementary.exponential import exp, log

7from sympy.sets.sets import Complement

8from sympy.utilities.iterables import uniq, is_sequence

11class Order(Expr):

12 r""" Represents the limiting behavior of some function.

14 Explanation

15 ===========

17 The order of a function characterizes the function based on the limiting

18 behavior of the function as it goes to some limit. Only taking the limit

19 point to be a number is currently supported. This is expressed in

20 big O notation [1]_.

22 The formal definition for the order of a function `g(x)` about a point `a`

23 is such that `g(x) = O(f(x))` as `x \rightarrow a` if and only if there

24 exists a `\delta > 0` and an `M > 0` such that `|g(x)| \leq M|f(x)|` for

25 `|x-a| < \delta`. This is equivalent to `\limsup_{x \rightarrow a}

26 |g(x)/f(x)| < \infty`.

28 Let's illustrate it on the following example by taking the expansion of

29 `\sin(x)` about 0:

31 .. math ::

32 \sin(x) = x - x^3/3! + O(x^5)

34 where in this case `O(x^5) = x^5/5! - x^7/7! + \cdots`. By the definition

35 of `O`, there is a `\delta > 0` and an `M` such that:

37 .. math ::

38 |x^5/5! - x^7/7! + ....| <= M|x^5| \text{ for } |x| < \delta

40 or by the alternate definition:

42 .. math ::

43 \lim_{x \rightarrow 0} | (x^5/5! - x^7/7! + ....) / x^5| < \infty

45 which surely is true, because

47 .. math ::

48 \lim_{x \rightarrow 0} | (x^5/5! - x^7/7! + ....) / x^5| = 1/5!

51 As it is usually used, the order of a function can be intuitively thought

52 of representing all terms of powers greater than the one specified. For

53 example, `O(x^3)` corresponds to any terms proportional to `x^3,

54 x^4,\ldots` and any higher power. For a polynomial, this leaves terms

55 proportional to `x^2`, `x` and constants.

57 Examples

58 ========

60 >>> from sympy import O, oo, cos, pi

61 >>> from sympy.abc import x, y

63 >>> O(x + x**2)

64 O(x)

65 >>> O(x + x**2, (x, 0))

66 O(x)

67 >>> O(x + x**2, (x, oo))

68 O(x**2, (x, oo))

70 >>> O(1 + x*y)

71 O(1, x, y)

72 >>> O(1 + x*y, (x, 0), (y, 0))

73 O(1, x, y)

74 >>> O(1 + x*y, (x, oo), (y, oo))

75 O(x*y, (x, oo), (y, oo))

77 >>> O(1) in O(1, x)

78 True

79 >>> O(1, x) in O(1)

80 False

81 >>> O(x) in O(1, x)

82 True

83 >>> O(x**2) in O(x)

84 True

86 >>> O(x)*x

87 O(x**2)

88 >>> O(x) - O(x)

89 O(x)

90 >>> O(cos(x))

91 O(1)

92 >>> O(cos(x), (x, pi/2))

93 O(x - pi/2, (x, pi/2))

95 References

96 ==========

98 .. [1] `Big O notation <https://en.wikipedia.org/wiki/Big_O_notation>`_

100 Notes

101 =====

102

103 In ``O(f(x), x)`` the expression ``f(x)`` is assumed to have a leading

104 term. ``O(f(x), x)`` is automatically transformed to

105 ``O(f(x).as_leading_term(x),x)``.

106

107 ``O(expr*f(x), x)`` is ``O(f(x), x)``

108

109 ``O(expr, x)`` is ``O(1)``

110

111 ``O(0, x)`` is 0.

112

113 Multivariate O is also supported:

114

115 ``O(f(x, y), x, y)`` is transformed to

116 ``O(f(x, y).as_leading_term(x,y).as_leading_term(y), x, y)``

117

118 In the multivariate case, it is assumed the limits w.r.t. the various

119 symbols commute.

120

121 If no symbols are passed then all symbols in the expression are used

122 and the limit point is assumed to be zero.

124 """

126 is_Order = True

128 __slots__ = ()

130 @cacheit

131 def __new__(cls, expr, *args, **kwargs):

132 expr = sympify(expr)

133

134 if not args:

135 if expr.is_Order:

136 variables = expr.variables

137 point = expr.point

138 else:

139 variables = list(expr.free_symbols)

140 point = [S.Zero]*len(variables)

141 else:

142 args = list(args if is_sequence(args) else [args])

143 variables, point = [], []

144 if is_sequence(args[0]):

145 for a in args:

146 v, p = list(map(sympify, a))

147 variables.append(v)

148 point.append(p)

149 else:

150 variables = list(map(sympify, args))

151 point = [S.Zero]*len(variables)

152

153 if not all(v.is_symbol for v in variables):

154 raise TypeError('Variables are not symbols, got %s' % variables)

155

156 if len(list(uniq(variables))) != len(variables):

157 raise ValueError('Variables are supposed to be unique symbols, got %s' % variables)

158

159 if expr.is_Order:

160 expr_vp = dict(expr.args[1:])

161 new_vp = dict(expr_vp)

162 vp = dict(zip(variables, point))

163 for v, p in vp.items():

164 if v in new_vp.keys():

165 if p != new_vp[v]:

166 raise NotImplementedError(

167 "Mixing Order at different points is not supported.")

168 else:

169 new_vp[v] = p

170 if set(expr_vp.keys()) == set(new_vp.keys()):

171 return expr

172 else:

173 variables = list(new_vp.keys())

174 point = [new_vp[v] for v in variables]

175

176 if expr is S.NaN:

177 return S.NaN

178

179 if any(x in p.free_symbols for x in variables for p in point):

180 raise ValueError('Got %s as a point.' % point)

181

182 if variables:

183 if any(p != point[0] for p in point):

184 raise NotImplementedError(

185 "Multivariable orders at different points are not supported.")

186 if point[0] in (S.Infinity, S.Infinity*S.ImaginaryUnit):

187 s = {k: 1/Dummy() for k in variables}

188 rs = {1/v: 1/k for k, v in s.items()}

189 ps = [S.Zero for p in point]

190 elif point[0] in (S.NegativeInfinity, S.NegativeInfinity*S.ImaginaryUnit):

191 s = {k: -1/Dummy() for k in variables}

192 rs = {-1/v: -1/k for k, v in s.items()}

193 ps = [S.Zero for p in point]

194 elif point[0] is not S.Zero:

195 s = {k: Dummy() + point[0] for k in variables}

196 rs = {(v - point[0]).together(): k - point[0] for k, v in s.items()}

197 ps = [S.Zero for p in point]

198 else:

199 s = ()

200 rs = ()

201 ps = list(point)

202

203 expr = expr.subs(s)

204

205 if expr.is_Add:

206 expr = expr.factor()

207

208 if s:

209 args = tuple([r[0] for r in rs.items()])

210 else:

211 args = tuple(variables)

212

213 if len(variables) > 1:

214 # XXX: better way? We need this expand() to

215 # workaround e.g: expr = x*(x + y).

216 # (x*(x + y)).as_leading_term(x, y) currently returns

217 # x*y (wrong order term!). That's why we want to deal with

218 # expand()'ed expr (handled in "if expr.is_Add" branch below).

219 expr = expr.expand()

220

221 old_expr = None

222 while old_expr != expr:

223 old_expr = expr

224 if expr.is_Add:

225 lst = expr.extract_leading_order(args)

226 expr = Add(*[f.expr for (e, f) in lst])

227

228 elif expr:

229 try:

230 expr = expr.as_leading_term(*args)

231 except PoleError:

232 if isinstance(expr, Function) or\

233 all(isinstance(arg, Function) for arg in expr.args):

234 # It is not possible to simplify an expression

235 # containing only functions (which raise error on

236 # call to leading term) further

237 pass

238 else:

239 orders = []

240 pts = tuple(zip(args, ps))

241 for arg in expr.args:

242 try:

243 lt = arg.as_leading_term(*args)

244 except PoleError:

245 lt = arg

246 if lt not in args:

247 order = Order(lt)

248 else:

249 order = Order(lt, *pts)

250 orders.append(order)

251 if expr.is_Add:

252 new_expr = Order(Add(*orders), *pts)

253 if new_expr.is_Add:

254 new_expr = Order(Add(*[a.expr for a in new_expr.args]), *pts)

255 expr = new_expr.expr

256 elif expr.is_Mul:

257 expr = Mul(*[a.expr for a in orders])

258 elif expr.is_Pow:

259 e = expr.exp

260 b = expr.base

261 expr = exp(e * log(b))

262

263 # It would probably be better to handle this somewhere

264 # else. This is needed for a testcase in which there is a

265 # symbol with the assumptions zero=True.

266 if expr.is_zero:

267 expr = S.Zero

268 else:

269 expr = expr.as_independent(*args, as_Add=False)[1]

270

271 expr = expand_power_base(expr)

272 expr = expand_log(expr)

273

274 if len(args) == 1:

275 # The definition of O(f(x)) symbol explicitly stated that

276 # the argument of f(x) is irrelevant. That's why we can

277 # combine some power exponents (only "on top" of the

278 # expression tree for f(x)), e.g.:

279 # x**p * (-x)**q -> x**(p+q) for real p, q.

280 x = args[0]

281 margs = list(Mul.make_args(

282 expr.as_independent(x, as_Add=False)[1]))

283

284 for i, t in enumerate(margs):

285 if t.is_Pow:

286 b, q = t.args

287 if b in (x, -x) and q.is_real and not q.has(x):

288 margs[i] = x**q

289 elif b.is_Pow and not b.exp.has(x):

290 b, r = b.args

291 if b in (x, -x) and r.is_real:

292 margs[i] = x**(r*q)

293 elif b.is_Mul and b.args[0] is S.NegativeOne:

294 b = -b

295 if b.is_Pow and not b.exp.has(x):

296 b, r = b.args

297 if b in (x, -x) and r.is_real:

298 margs[i] = x**(r*q)

299

300 expr = Mul(*margs)

301

302 expr = expr.subs(rs)

303

304 if expr.is_Order:

305 expr = expr.expr

306

307 if not expr.has(*variables) and not expr.is_zero:

308 expr = S.One

309

310 # create Order instance:

311 vp = dict(zip(variables, point))

312 variables.sort(key=default_sort_key)

313 point = [vp[v] for v in variables]

314 args = (expr,) + Tuple(*zip(variables, point))

315 obj = Expr.__new__(cls, *args)

316 return obj

317

318 def _eval_nseries(self, x, n, logx, cdir=0):

319 return self

320

321 @property

322 def expr(self):

323 return self.args[0]

324

325 @property

326 def variables(self):

327 if self.args[1:]:

328 return tuple(x[0] for x in self.args[1:])

329 else:

330 return ()

331

332 @property

333 def point(self):

334 if self.args[1:]:

335 return tuple(x[1] for x in self.args[1:])

336 else:

337 return ()

338

339 @property

340 def free_symbols(self):

341 return self.expr.free_symbols | set(self.variables)

342

343 def _eval_power(b, e):

344 if e.is_Number and e.is_nonnegative:

345 return b.func(b.expr ** e, *b.args[1:])

346 if e == O(1):

347 return b

348 return

349

350 def as_expr_variables(self, order_symbols):

351 if order_symbols is None:

352 order_symbols = self.args[1:]

353 else:

354 if (not all(o[1] == order_symbols[0][1] for o in order_symbols) and

355 not all(p == self.point[0] for p in self.point)): # pragma: no cover

356 raise NotImplementedError('Order at points other than 0 '

357 'or oo not supported, got %s as a point.' % self.point)

358 if order_symbols and order_symbols[0][1] != self.point[0]:

359 raise NotImplementedError(

360 "Multiplying Order at different points is not supported.")

361 order_symbols = dict(order_symbols)

362 for s, p in dict(self.args[1:]).items():

363 if s not in order_symbols.keys():

364 order_symbols[s] = p

365 order_symbols = sorted(order_symbols.items(), key=lambda x: default_sort_key(x[0]))

366 return self.expr, tuple(order_symbols)

367

368 def removeO(self):

369 return S.Zero

370

371 def getO(self):

372 return self

373

374 @cacheit

375 def contains(self, expr):

376 r"""

377 Return True if expr belongs to Order(self.expr, \*self.variables).

378 Return False if self belongs to expr.

379 Return None if the inclusion relation cannot be determined

380 (e.g. when self and expr have different symbols).

381 """

382 expr = sympify(expr)

383 if expr.is_zero:

384 return True

385 if expr is S.NaN:

386 return False

387 point = self.point[0] if self.point else S.Zero

388 if expr.is_Order:

389 if (any(p != point for p in expr.point) or

390 any(p != point for p in self.point)):

391 return None

392 if expr.expr == self.expr:

393 # O(1) + O(1), O(1) + O(1, x), etc.

394 return all(x in self.args[1:] for x in expr.args[1:])

395 if expr.expr.is_Add:

396 return all(self.contains(x) for x in expr.expr.args)

397 if self.expr.is_Add and point.is_zero:

398 return any(self.func(x, *self.args[1:]).contains(expr)

399 for x in self.expr.args)

400 if self.variables and expr.variables:

401 common_symbols = tuple(

402 [s for s in self.variables if s in expr.variables])

403 elif self.variables:

404 common_symbols = self.variables

405 else:

406 common_symbols = expr.variables

407 if not common_symbols:

408 return None

409 if (self.expr.is_Pow and len(self.variables) == 1

410 and self.variables == expr.variables):

411 symbol = self.variables[0]

412 other = expr.expr.as_independent(symbol, as_Add=False)[1]

413 if (other.is_Pow and other.base == symbol and

414 self.expr.base == symbol):

415 if point.is_zero:

416 rv = (self.expr.exp - other.exp).is_nonpositive

417 if point.is_infinite:

418 rv = (self.expr.exp - other.exp).is_nonnegative

419 if rv is not None:

420 return rv

421

422 from sympy.simplify.powsimp import powsimp

423 r = None

424 ratio = self.expr/expr.expr

425 ratio = powsimp(ratio, deep=True, combine='exp')

426 for s in common_symbols:

427 from sympy.series.limits import Limit

428 l = Limit(ratio, s, point).doit(heuristics=False)

429 if not isinstance(l, Limit):

430 l = l != 0

431 else:

432 l = None

433 if r is None:

434 r = l

435 else:

436 if r != l:

437 return

438 return r

439

440 if self.expr.is_Pow and len(self.variables) == 1:

441 symbol = self.variables[0]

442 other = expr.as_independent(symbol, as_Add=False)[1]

443 if (other.is_Pow and other.base == symbol and

444 self.expr.base == symbol):

445 if point.is_zero:

446 rv = (self.expr.exp - other.exp).is_nonpositive

447 if point.is_infinite:

448 rv = (self.expr.exp - other.exp).is_nonnegative

449 if rv is not None:

450 return rv

451

452 obj = self.func(expr, *self.args[1:])

453 return self.contains(obj)

454

455 def __contains__(self, other):

456 result = self.contains(other)

457 if result is None:

458 raise TypeError('contains did not evaluate to a bool')

459 return result

460

461 def _eval_subs(self, old, new):

462 if old in self.variables:

463 newexpr = self.expr.subs(old, new)

464 i = self.variables.index(old)

465 newvars = list(self.variables)

466 newpt = list(self.point)

467 if new.is_symbol:

468 newvars[i] = new

469 else:

470 syms = new.free_symbols

471 if len(syms) == 1 or old in syms:

472 if old in syms:

473 var = self.variables[i]

474 else:

475 var = syms.pop()

476 # First, try to substitute self.point in the "new"

477 # expr to see if this is a fixed point.

478 # E.g. O(y).subs(y, sin(x))

479 point = new.subs(var, self.point[i])

480 if point != self.point[i]:

481 from sympy.solvers.solveset import solveset

482 d = Dummy()

483 sol = solveset(old - new.subs(var, d), d)

484 if isinstance(sol, Complement):

485 e1 = sol.args[0]

486 e2 = sol.args[1]

487 sol = set(e1) - set(e2)

488 res = [dict(zip((d, ), sol))]

489 point = d.subs(res[0]).limit(old, self.point[i])

490 newvars[i] = var

491 newpt[i] = point

492 elif old not in syms:

493 del newvars[i], newpt[i]

494 if not syms and new == self.point[i]:

495 newvars.extend(syms)

496 newpt.extend([S.Zero]*len(syms))

497 else:

498 return

499 return Order(newexpr, *zip(newvars, newpt))

500

501 def _eval_conjugate(self):

502 expr = self.expr._eval_conjugate()

503 if expr is not None:

504 return self.func(expr, *self.args[1:])

505

506 def _eval_derivative(self, x):

507 return self.func(self.expr.diff(x), *self.args[1:]) or self

508

509 def _eval_transpose(self):

510 expr = self.expr._eval_transpose()

511 if expr is not None:

512 return self.func(expr, *self.args[1:])

513

514 def __neg__(self):

515 return self

516

517O = Order