Coverage for /usr/lib/python3/dist-packages/sympy/series/sequences.py: 29%

1from sympy.core.basic import Basic 

2from sympy.core.cache import cacheit 

3from sympy.core.containers import Tuple 

4from sympy.core.decorators import call_highest_priority 

5from sympy.core.parameters import global_parameters 

6from sympy.core.function import AppliedUndef, expand 

7from sympy.core.mul import Mul 

8from sympy.core.numbers import Integer 

9from sympy.core.relational import Eq 

10from sympy.core.singleton import S, Singleton 

11from sympy.core.sorting import ordered 

12from sympy.core.symbol import Dummy, Symbol, Wild 

13from sympy.core.sympify import sympify 

14from sympy.matrices import Matrix 

15from sympy.polys import lcm, factor 

16from sympy.sets.sets import Interval, Intersection 

17from sympy.tensor.indexed import Idx 

18from sympy.utilities.iterables import flatten, is_sequence, iterable 

19 

20 

21############################################################################### 

22# SEQUENCES # 

23############################################################################### 

24 

25 

26class SeqBase(Basic): 

27 """Base class for sequences""" 

28 

29 is_commutative = True 

30 _op_priority = 15 

31 

32 @staticmethod 

33 def _start_key(expr): 

34 """Return start (if possible) else S.Infinity. 

35 

36 adapted from Set._infimum_key 

37 """ 

38 try: 

39 start = expr.start 

40 except (NotImplementedError, 

41 AttributeError, ValueError): 

42 start = S.Infinity 

43 return start 

44 

45 def _intersect_interval(self, other): 

46 """Returns start and stop. 

47 

48 Takes intersection over the two intervals. 

49 """ 

50 interval = Intersection(self.interval, other.interval) 

51 return interval.inf, interval.sup 

52 

53 @property 

54 def gen(self): 

55 """Returns the generator for the sequence""" 

56 raise NotImplementedError("(%s).gen" % self) 

57 

58 @property 

59 def interval(self): 

60 """The interval on which the sequence is defined""" 

61 raise NotImplementedError("(%s).interval" % self) 

62 

63 @property 

64 def start(self): 

65 """The starting point of the sequence. This point is included""" 

66 raise NotImplementedError("(%s).start" % self) 

67 

68 @property 

69 def stop(self): 

70 """The ending point of the sequence. This point is included""" 

71 raise NotImplementedError("(%s).stop" % self) 

72 

73 @property 

74 def length(self): 

75 """Length of the sequence""" 

76 raise NotImplementedError("(%s).length" % self) 

77 

78 @property 

79 def variables(self): 

80 """Returns a tuple of variables that are bounded""" 

81 return () 

82 

83 @property 

84 def free_symbols(self): 

85 """ 

86 This method returns the symbols in the object, excluding those 

87 that take on a specific value (i.e. the dummy symbols). 

88 

89 Examples 

90 ======== 

91 

92 >>> from sympy import SeqFormula 

93 >>> from sympy.abc import n, m 

94 >>> SeqFormula(m*n**2, (n, 0, 5)).free_symbols 

95 {m} 

96 """ 

97 return ({j for i in self.args for j in i.free_symbols 

98 .difference(self.variables)}) 

99 

100 @cacheit 

101 def coeff(self, pt): 

102 """Returns the coefficient at point pt""" 

103 if pt < self.start or pt > self.stop: 

104 raise IndexError("Index %s out of bounds %s" % (pt, self.interval)) 

105 return self._eval_coeff(pt) 

106 

107 def _eval_coeff(self, pt): 

108 raise NotImplementedError("The _eval_coeff method should be added to" 

109 "%s to return coefficient so it is available" 

110 "when coeff calls it." 

111 % self.func) 

112 

113 def _ith_point(self, i): 

114 """Returns the i'th point of a sequence. 

115 

116 Explanation 

117 =========== 

118 

119 If start point is negative infinity, point is returned from the end. 

120 Assumes the first point to be indexed zero. 

121 

122 Examples 

123 ========= 

124 

125 >>> from sympy import oo 

126 >>> from sympy.series.sequences import SeqPer 

127 

128 bounded 

129 

130 >>> SeqPer((1, 2, 3), (-10, 10))._ith_point(0) 

131 -10 

132 >>> SeqPer((1, 2, 3), (-10, 10))._ith_point(5) 

133 -5 

134 

135 End is at infinity 

136 

137 >>> SeqPer((1, 2, 3), (0, oo))._ith_point(5) 

138 5 

139 

140 Starts at negative infinity 

141 

142 >>> SeqPer((1, 2, 3), (-oo, 0))._ith_point(5) 

143 -5 

144 """ 

145 if self.start is S.NegativeInfinity: 

146 initial = self.stop 

147 else: 

148 initial = self.start 

149 

150 if self.start is S.NegativeInfinity: 

151 step = -1 

152 else: 

153 step = 1 

154 

155 return initial + i*step 

156 

157 def _add(self, other): 

158 """ 

159 Should only be used internally. 

160 

161 Explanation 

162 =========== 

163 

164 self._add(other) returns a new, term-wise added sequence if self 

165 knows how to add with other, otherwise it returns ``None``. 

166 

167 ``other`` should only be a sequence object. 

168 

169 Used within :class:`SeqAdd` class. 

170 """ 

171 return None 

172 

173 def _mul(self, other): 

174 """ 

175 Should only be used internally. 

176 

177 Explanation 

178 =========== 

179 

180 self._mul(other) returns a new, term-wise multiplied sequence if self 

181 knows how to multiply with other, otherwise it returns ``None``. 

182 

183 ``other`` should only be a sequence object. 

184 

185 Used within :class:`SeqMul` class. 

186 """ 

187 return None 

188 

189 def coeff_mul(self, other): 

190 """ 

191 Should be used when ``other`` is not a sequence. Should be 

192 defined to define custom behaviour. 

193 

194 Examples 

195 ======== 

196 

197 >>> from sympy import SeqFormula 

198 >>> from sympy.abc import n 

199 >>> SeqFormula(n**2).coeff_mul(2) 

200 SeqFormula(2*n**2, (n, 0, oo)) 

201 

202 Notes 

203 ===== 

204 

205 '*' defines multiplication of sequences with sequences only. 

206 """ 

207 return Mul(self, other) 

208 

209 def __add__(self, other): 

210 """Returns the term-wise addition of 'self' and 'other'. 

211 

212 ``other`` should be a sequence. 

213 

214 Examples 

215 ======== 

216 

217 >>> from sympy import SeqFormula 

218 >>> from sympy.abc import n 

219 >>> SeqFormula(n**2) + SeqFormula(n**3) 

220 SeqFormula(n**3 + n**2, (n, 0, oo)) 

221 """ 

222 if not isinstance(other, SeqBase): 

223 raise TypeError('cannot add sequence and %s' % type(other)) 

224 return SeqAdd(self, other) 

225 

226 @call_highest_priority('__add__') 

227 def __radd__(self, other): 

228 return self + other 

229 

230 def __sub__(self, other): 

231 """Returns the term-wise subtraction of ``self`` and ``other``. 

232 

233 ``other`` should be a sequence. 

234 

235 Examples 

236 ======== 

237 

238 >>> from sympy import SeqFormula 

239 >>> from sympy.abc import n 

240 >>> SeqFormula(n**2) - (SeqFormula(n)) 

241 SeqFormula(n**2 - n, (n, 0, oo)) 

242 """ 

243 if not isinstance(other, SeqBase): 

244 raise TypeError('cannot subtract sequence and %s' % type(other)) 

245 return SeqAdd(self, -other) 

246 

247 @call_highest_priority('__sub__') 

248 def __rsub__(self, other): 

249 return (-self) + other 

250 

251 def __neg__(self): 

252 """Negates the sequence. 

253 

254 Examples 

255 ======== 

256 

257 >>> from sympy import SeqFormula 

258 >>> from sympy.abc import n 

259 >>> -SeqFormula(n**2) 

260 SeqFormula(-n**2, (n, 0, oo)) 

261 """ 

262 return self.coeff_mul(-1) 

263 

264 def __mul__(self, other): 

265 """Returns the term-wise multiplication of 'self' and 'other'. 

266 

267 ``other`` should be a sequence. For ``other`` not being a 

268 sequence see :func:`coeff_mul` method. 

269 

270 Examples 

271 ======== 

272 

273 >>> from sympy import SeqFormula 

274 >>> from sympy.abc import n 

275 >>> SeqFormula(n**2) * (SeqFormula(n)) 

276 SeqFormula(n**3, (n, 0, oo)) 

277 """ 

278 if not isinstance(other, SeqBase): 

279 raise TypeError('cannot multiply sequence and %s' % type(other)) 

280 return SeqMul(self, other) 

281 

282 @call_highest_priority('__mul__') 

283 def __rmul__(self, other): 

284 return self * other 

285 

286 def __iter__(self): 

287 for i in range(self.length): 

288 pt = self._ith_point(i) 

289 yield self.coeff(pt) 

290 

291 def __getitem__(self, index): 

292 if isinstance(index, int): 

293 index = self._ith_point(index) 

294 return self.coeff(index) 

295 elif isinstance(index, slice): 

296 start, stop = index.start, index.stop 

297 if start is None: 

298 start = 0 

299 if stop is None: 

300 stop = self.length 

301 return [self.coeff(self._ith_point(i)) for i in 

302 range(start, stop, index.step or 1)] 

303 

304 def find_linear_recurrence(self,n,d=None,gfvar=None): 

305 r""" 

306 Finds the shortest linear recurrence that satisfies the first n 

307 terms of sequence of order `\leq` ``n/2`` if possible. 

308 If ``d`` is specified, find shortest linear recurrence of order 

309 `\leq` min(d, n/2) if possible. 

310 Returns list of coefficients ``[b(1), b(2), ...]`` corresponding to the 

311 recurrence relation ``x(n) = b(1)*x(n-1) + b(2)*x(n-2) + ...`` 

312 Returns ``[]`` if no recurrence is found. 

313 If gfvar is specified, also returns ordinary generating function as a 

314 function of gfvar. 

315 

316 Examples 

317 ======== 

318 

319 >>> from sympy import sequence, sqrt, oo, lucas 

320 >>> from sympy.abc import n, x, y 

321 >>> sequence(n**2).find_linear_recurrence(10, 2) 

322 [] 

323 >>> sequence(n**2).find_linear_recurrence(10) 

324 [3, -3, 1] 

325 >>> sequence(2**n).find_linear_recurrence(10) 

326 [2] 

327 >>> sequence(23*n**4+91*n**2).find_linear_recurrence(10) 

328 [5, -10, 10, -5, 1] 

329 >>> sequence(sqrt(5)*(((1 + sqrt(5))/2)**n - (-(1 + sqrt(5))/2)**(-n))/5).find_linear_recurrence(10) 

330 [1, 1] 

331 >>> sequence(x+y*(-2)**(-n), (n, 0, oo)).find_linear_recurrence(30) 

332 [1/2, 1/2] 

333 >>> sequence(3*5**n + 12).find_linear_recurrence(20,gfvar=x) 

334 ([6, -5], 3*(5 - 21*x)/((x - 1)*(5*x - 1))) 

335 >>> sequence(lucas(n)).find_linear_recurrence(15,gfvar=x) 

336 ([1, 1], (x - 2)/(x**2 + x - 1)) 

337 """ 

338 from sympy.simplify import simplify 

339 x = [simplify(expand(t)) for t in self[:n]] 

340 lx = len(x) 

341 if d is None: 

342 r = lx//2 

343 else: 

344 r = min(d,lx//2) 

345 coeffs = [] 

346 for l in range(1, r+1): 

347 l2 = 2*l 

348 mlist = [] 

349 for k in range(l): 

350 mlist.append(x[k:k+l]) 

351 m = Matrix(mlist) 

352 if m.det() != 0: 

353 y = simplify(m.LUsolve(Matrix(x[l:l2]))) 

354 if lx == l2: 

355 coeffs = flatten(y[::-1]) 

356 break 

357 mlist = [] 

358 for k in range(l,lx-l): 

359 mlist.append(x[k:k+l]) 

360 m = Matrix(mlist) 

361 if m*y == Matrix(x[l2:]): 

362 coeffs = flatten(y[::-1]) 

363 break 

364 if gfvar is None: 

365 return coeffs 

366 else: 

367 l = len(coeffs) 

368 if l == 0: 

369 return [], None 

370 else: 

371 n, d = x[l-1]*gfvar**(l-1), 1 - coeffs[l-1]*gfvar**l 

372 for i in range(l-1): 

373 n += x[i]*gfvar**i 

374 for j in range(l-i-1): 

375 n -= coeffs[i]*x[j]*gfvar**(i+j+1) 

376 d -= coeffs[i]*gfvar**(i+1) 

377 return coeffs, simplify(factor(n)/factor(d)) 

378 

379class EmptySequence(SeqBase, metaclass=Singleton): 

380 """Represents an empty sequence. 

381 

382 The empty sequence is also available as a singleton as 

383 ``S.EmptySequence``. 

384 

385 Examples 

386 ======== 

387 

388 >>> from sympy import EmptySequence, SeqPer 

389 >>> from sympy.abc import x 

390 >>> EmptySequence 

391 EmptySequence 

392 >>> SeqPer((1, 2), (x, 0, 10)) + EmptySequence 

393 SeqPer((1, 2), (x, 0, 10)) 

394 >>> SeqPer((1, 2)) * EmptySequence 

395 EmptySequence 

396 >>> EmptySequence.coeff_mul(-1) 

397 EmptySequence 

398 """ 

399 

400 @property 

401 def interval(self): 

402 return S.EmptySet 

403 

404 @property 

405 def length(self): 

406 return S.Zero 

407 

408 def coeff_mul(self, coeff): 

409 """See docstring of SeqBase.coeff_mul""" 

410 return self 

411 

412 def __iter__(self): 

413 return iter([]) 

414 

415 

416class SeqExpr(SeqBase): 

417 """Sequence expression class. 

418 

419 Various sequences should inherit from this class. 

420 

421 Examples 

422 ======== 

423 

424 >>> from sympy.series.sequences import SeqExpr 

425 >>> from sympy.abc import x 

426 >>> from sympy import Tuple 

427 >>> s = SeqExpr(Tuple(1, 2, 3), Tuple(x, 0, 10)) 

428 >>> s.gen 

429 (1, 2, 3) 

430 >>> s.interval 

431 Interval(0, 10) 

432 >>> s.length 

433 11 

434 

435 See Also 

436 ======== 

437 

438 sympy.series.sequences.SeqPer 

439 sympy.series.sequences.SeqFormula 

440 """ 

441 

442 @property 

443 def gen(self): 

444 return self.args[0] 

445 

446 @property 

447 def interval(self): 

448 return Interval(self.args[1][1], self.args[1][2]) 

449 

450 @property 

451 def start(self): 

452 return self.interval.inf 

453 

454 @property 

455 def stop(self): 

456 return self.interval.sup 

457 

458 @property 

459 def length(self): 

460 return self.stop - self.start + 1 

461 

462 @property 

463 def variables(self): 

464 return (self.args[1][0],) 

465 

466 

467class SeqPer(SeqExpr): 

468 """ 

469 Represents a periodic sequence. 

470 

471 The elements are repeated after a given period. 

472 

473 Examples 

474 ======== 

475 

476 >>> from sympy import SeqPer, oo 

477 >>> from sympy.abc import k 

478 

479 >>> s = SeqPer((1, 2, 3), (0, 5)) 

480 >>> s.periodical 

481 (1, 2, 3) 

482 >>> s.period 

483 3 

484 

485 For value at a particular point 

486 

487 >>> s.coeff(3) 

488 1 

489 

490 supports slicing 

491 

492 >>> s[:] 

493 [1, 2, 3, 1, 2, 3] 

494 

495 iterable 

496 

497 >>> list(s) 

498 [1, 2, 3, 1, 2, 3] 

499 

500 sequence starts from negative infinity 

501 

502 >>> SeqPer((1, 2, 3), (-oo, 0))[0:6] 

503 [1, 2, 3, 1, 2, 3] 

504 

505 Periodic formulas 

506 

507 >>> SeqPer((k, k**2, k**3), (k, 0, oo))[0:6] 

508 [0, 1, 8, 3, 16, 125] 

509 

510 See Also 

511 ======== 

512 

513 sympy.series.sequences.SeqFormula 

514 """ 

515 

516 def __new__(cls, periodical, limits=None): 

517 periodical = sympify(periodical) 

518 

519 def _find_x(periodical): 

520 free = periodical.free_symbols 

521 if len(periodical.free_symbols) == 1: 

522 return free.pop() 

523 else: 

524 return Dummy('k') 

525 

526 x, start, stop = None, None, None 

527 if limits is None: 

528 x, start, stop = _find_x(periodical), 0, S.Infinity 

529 if is_sequence(limits, Tuple): 

530 if len(limits) == 3: 

531 x, start, stop = limits 

532 elif len(limits) == 2: 

533 x = _find_x(periodical) 

534 start, stop = limits 

535 

536 if not isinstance(x, (Symbol, Idx)) or start is None or stop is None: 

537 raise ValueError('Invalid limits given: %s' % str(limits)) 

538 

539 if start is S.NegativeInfinity and stop is S.Infinity: 

540 raise ValueError("Both the start and end value" 

541 "cannot be unbounded") 

542 

543 limits = sympify((x, start, stop)) 

544 

545 if is_sequence(periodical, Tuple): 

546 periodical = sympify(tuple(flatten(periodical))) 

547 else: 

548 raise ValueError("invalid period %s should be something " 

549 "like e.g (1, 2) " % periodical) 

550 

551 if Interval(limits[1], limits[2]) is S.EmptySet: 

552 return S.EmptySequence 

553 

554 return Basic.__new__(cls, periodical, limits) 

555 

556 @property 

557 def period(self): 

558 return len(self.gen) 

559 

560 @property 

561 def periodical(self): 

562 return self.gen 

563 

564 def _eval_coeff(self, pt): 

565 if self.start is S.NegativeInfinity: 

566 idx = (self.stop - pt) % self.period 

567 else: 

568 idx = (pt - self.start) % self.period 

569 return self.periodical[idx].subs(self.variables[0], pt) 

570 

571 def _add(self, other): 

572 """See docstring of SeqBase._add""" 

573 if isinstance(other, SeqPer): 

574 per1, lper1 = self.periodical, self.period 

575 per2, lper2 = other.periodical, other.period 

576 

577 per_length = lcm(lper1, lper2) 

578 

579 new_per = [] 

580 for x in range(per_length): 

581 ele1 = per1[x % lper1] 

582 ele2 = per2[x % lper2] 

583 new_per.append(ele1 + ele2) 

584 

585 start, stop = self._intersect_interval(other) 

586 return SeqPer(new_per, (self.variables[0], start, stop)) 

587 

588 def _mul(self, other): 

589 """See docstring of SeqBase._mul""" 

590 if isinstance(other, SeqPer): 

591 per1, lper1 = self.periodical, self.period 

592 per2, lper2 = other.periodical, other.period 

593 

594 per_length = lcm(lper1, lper2) 

595 

596 new_per = [] 

597 for x in range(per_length): 

598 ele1 = per1[x % lper1] 

599 ele2 = per2[x % lper2] 

600 new_per.append(ele1 * ele2) 

601 

602 start, stop = self._intersect_interval(other) 

603 return SeqPer(new_per, (self.variables[0], start, stop)) 

604 

605 def coeff_mul(self, coeff): 

606 """See docstring of SeqBase.coeff_mul""" 

607 coeff = sympify(coeff) 

608 per = [x * coeff for x in self.periodical] 

609 return SeqPer(per, self.args[1]) 

610 

611 

612class SeqFormula(SeqExpr): 

613 """ 

614 Represents sequence based on a formula. 

615 

616 Elements are generated using a formula. 

617 

618 Examples 

619 ======== 

620 

621 >>> from sympy import SeqFormula, oo, Symbol 

622 >>> n = Symbol('n') 

623 >>> s = SeqFormula(n**2, (n, 0, 5)) 

624 >>> s.formula 

625 n**2 

626 

627 For value at a particular point 

628 

629 >>> s.coeff(3) 

630 9 

631 

632 supports slicing 

633 

634 >>> s[:] 

635 [0, 1, 4, 9, 16, 25] 

636 

637 iterable 

638 

639 >>> list(s) 

640 [0, 1, 4, 9, 16, 25] 

641 

642 sequence starts from negative infinity 

643 

644 >>> SeqFormula(n**2, (-oo, 0))[0:6] 

645 [0, 1, 4, 9, 16, 25] 

646 

647 See Also 

648 ======== 

649 

650 sympy.series.sequences.SeqPer 

651 """ 

652 

653 def __new__(cls, formula, limits=None): 

654 formula = sympify(formula) 

655 

656 def _find_x(formula): 

657 free = formula.free_symbols 

658 if len(free) == 1: 

659 return free.pop() 

660 elif not free: 

661 return Dummy('k') 

662 else: 

663 raise ValueError( 

664 " specify dummy variables for %s. If the formula contains" 

665 " more than one free symbol, a dummy variable should be" 

666 " supplied explicitly e.g., SeqFormula(m*n**2, (n, 0, 5))" 

667 % formula) 

668 

669 x, start, stop = None, None, None 

670 if limits is None: 

671 x, start, stop = _find_x(formula), 0, S.Infinity 

672 if is_sequence(limits, Tuple): 

673 if len(limits) == 3: 

674 x, start, stop = limits 

675 elif len(limits) == 2: 

676 x = _find_x(formula) 

677 start, stop = limits 

678 

679 if not isinstance(x, (Symbol, Idx)) or start is None or stop is None: 

680 raise ValueError('Invalid limits given: %s' % str(limits)) 

681 

682 if start is S.NegativeInfinity and stop is S.Infinity: 

683 raise ValueError("Both the start and end value " 

684 "cannot be unbounded") 

685 limits = sympify((x, start, stop)) 

686 

687 if Interval(limits[1], limits[2]) is S.EmptySet: 

688 return S.EmptySequence 

689 

690 return Basic.__new__(cls, formula, limits) 

691 

692 @property 

693 def formula(self): 

694 return self.gen 

695 

696 def _eval_coeff(self, pt): 

697 d = self.variables[0] 

698 return self.formula.subs(d, pt) 

699 

700 def _add(self, other): 

701 """See docstring of SeqBase._add""" 

702 if isinstance(other, SeqFormula): 

703 form1, v1 = self.formula, self.variables[0] 

704 form2, v2 = other.formula, other.variables[0] 

705 formula = form1 + form2.subs(v2, v1) 

706 start, stop = self._intersect_interval(other) 

707 return SeqFormula(formula, (v1, start, stop)) 

708 

709 def _mul(self, other): 

710 """See docstring of SeqBase._mul""" 

711 if isinstance(other, SeqFormula): 

712 form1, v1 = self.formula, self.variables[0] 

713 form2, v2 = other.formula, other.variables[0] 

714 formula = form1 * form2.subs(v2, v1) 

715 start, stop = self._intersect_interval(other) 

716 return SeqFormula(formula, (v1, start, stop)) 

717 

718 def coeff_mul(self, coeff): 

719 """See docstring of SeqBase.coeff_mul""" 

720 coeff = sympify(coeff) 

721 formula = self.formula * coeff 

722 return SeqFormula(formula, self.args[1]) 

723 

724 def expand(self, *args, **kwargs): 

725 return SeqFormula(expand(self.formula, *args, **kwargs), self.args[1]) 

726 

727class RecursiveSeq(SeqBase): 

728 """ 

729 A finite degree recursive sequence. 

730 

731 Explanation 

732 =========== 

733 

734 That is, a sequence a(n) that depends on a fixed, finite number of its 

735 previous values. The general form is 

736 

737 a(n) = f(a(n - 1), a(n - 2), ..., a(n - d)) 

738 

739 for some fixed, positive integer d, where f is some function defined by a 

740 SymPy expression. 

741 

742 Parameters 

743 ========== 

744 

745 recurrence : SymPy expression defining recurrence 

746 This is *not* an equality, only the expression that the nth term is 

747 equal to. For example, if :code:`a(n) = f(a(n - 1), ..., a(n - d))`, 

748 then the expression should be :code:`f(a(n - 1), ..., a(n - d))`. 

749 

750 yn : applied undefined function 

751 Represents the nth term of the sequence as e.g. :code:`y(n)` where 

752 :code:`y` is an undefined function and `n` is the sequence index. 

753 

754 n : symbolic argument 

755 The name of the variable that the recurrence is in, e.g., :code:`n` if 

756 the recurrence function is :code:`y(n)`. 

757 

758 initial : iterable with length equal to the degree of the recurrence 

759 The initial values of the recurrence. 

760 

761 start : start value of sequence (inclusive) 

762 

763 Examples 

764 ======== 

765 

766 >>> from sympy import Function, symbols 

767 >>> from sympy.series.sequences import RecursiveSeq 

768 >>> y = Function("y") 

769 >>> n = symbols("n") 

770 >>> fib = RecursiveSeq(y(n - 1) + y(n - 2), y(n), n, [0, 1]) 

771 

772 >>> fib.coeff(3) # Value at a particular point 

773 2 

774 

775 >>> fib[:6] # supports slicing 

776 [0, 1, 1, 2, 3, 5] 

777 

778 >>> fib.recurrence # inspect recurrence 

779 Eq(y(n), y(n - 2) + y(n - 1)) 

780 

781 >>> fib.degree # automatically determine degree 

782 2 

783 

784 >>> for x in zip(range(10), fib): # supports iteration 

785 ... print(x) 

786 (0, 0) 

787 (1, 1) 

788 (2, 1) 

789 (3, 2) 

790 (4, 3) 

791 (5, 5) 

792 (6, 8) 

793 (7, 13) 

794 (8, 21) 

795 (9, 34) 

796 

797 See Also 

798 ======== 

799 

800 sympy.series.sequences.SeqFormula 

801 

802 """ 

803 

804 def __new__(cls, recurrence, yn, n, initial=None, start=0): 

805 if not isinstance(yn, AppliedUndef): 

806 raise TypeError("recurrence sequence must be an applied undefined function" 

807 ", found `{}`".format(yn)) 

808 

809 if not isinstance(n, Basic) or not n.is_symbol: 

810 raise TypeError("recurrence variable must be a symbol" 

811 ", found `{}`".format(n)) 

812 

813 if yn.args != (n,): 

814 raise TypeError("recurrence sequence does not match symbol") 

815 

816 y = yn.func 

817 

818 k = Wild("k", exclude=(n,)) 

819 degree = 0 

820 

821 # Find all applications of y in the recurrence and check that: 

822 # 1. The function y is only being used with a single argument; and 

823 # 2. All arguments are n + k for constant negative integers k. 

824 

825 prev_ys = recurrence.find(y) 

826 for prev_y in prev_ys: 

827 if len(prev_y.args) != 1: 

828 raise TypeError("Recurrence should be in a single variable") 

829 

830 shift = prev_y.args[0].match(n + k)[k] 

831 if not (shift.is_constant() and shift.is_integer and shift < 0): 

832 raise TypeError("Recurrence should have constant," 

833 " negative, integer shifts" 

834 " (found {})".format(prev_y)) 

835 

836 if -shift > degree: 

837 degree = -shift 

838 

839 if not initial: 

840 initial = [Dummy("c_{}".format(k)) for k in range(degree)] 

841 

842 if len(initial) != degree: 

843 raise ValueError("Number of initial terms must equal degree") 

844 

845 degree = Integer(degree) 

846 start = sympify(start) 

847 

848 initial = Tuple(*(sympify(x) for x in initial)) 

849 

850 seq = Basic.__new__(cls, recurrence, yn, n, initial, start) 

851 

852 seq.cache = {y(start + k): init for k, init in enumerate(initial)} 

853 seq.degree = degree 

854 

855 return seq 

856 

857 @property 

858 def _recurrence(self): 

859 """Equation defining recurrence.""" 

860 return self.args[0] 

861 

862 @property 

863 def recurrence(self): 

864 """Equation defining recurrence.""" 

865 return Eq(self.yn, self.args[0]) 

866 

867 @property 

868 def yn(self): 

869 """Applied function representing the nth term""" 

870 return self.args[1] 

871 

872 @property 

873 def y(self): 

874 """Undefined function for the nth term of the sequence""" 

875 return self.yn.func 

876 

877 @property 

878 def n(self): 

879 """Sequence index symbol""" 

880 return self.args[2] 

881 

882 @property 

883 def initial(self): 

884 """The initial values of the sequence""" 

885 return self.args[3] 

886 

887 @property 

888 def start(self): 

889 """The starting point of the sequence. This point is included""" 

890 return self.args[4] 

891 

892 @property 

893 def stop(self): 

894 """The ending point of the sequence. (oo)""" 

895 return S.Infinity 

896 

897 @property 

898 def interval(self): 

899 """Interval on which sequence is defined.""" 

900 return (self.start, S.Infinity) 

901 

902 def _eval_coeff(self, index): 

903 if index - self.start < len(self.cache): 

904 return self.cache[self.y(index)] 

905 

906 for current in range(len(self.cache), index + 1): 

907 # Use xreplace over subs for performance. 

908 # See issue #10697. 

909 seq_index = self.start + current 

910 current_recurrence = self._recurrence.xreplace({self.n: seq_index}) 

911 new_term = current_recurrence.xreplace(self.cache) 

912 

913 self.cache[self.y(seq_index)] = new_term 

914 

915 return self.cache[self.y(self.start + current)] 

916 

917 def __iter__(self): 

918 index = self.start 

919 while True: 

920 yield self._eval_coeff(index) 

921 index += 1 

922 

923 

924def sequence(seq, limits=None): 

925 """ 

926 Returns appropriate sequence object. 

927 

928 Explanation