Coverage for /usr/lib/python3/dist-packages/sympy/series/limitseq.py: 16%

1"""Limits of sequences""" 

2 

3from sympy.calculus.accumulationbounds import AccumulationBounds 

4from sympy.core.add import Add 

5from sympy.core.function import PoleError 

6from sympy.core.power import Pow 

7from sympy.core.singleton import S 

8from sympy.core.symbol import Dummy 

9from sympy.core.sympify import sympify 

10from sympy.functions.combinatorial.numbers import fibonacci 

11from sympy.functions.combinatorial.factorials import factorial, subfactorial 

12from sympy.functions.special.gamma_functions import gamma 

13from sympy.functions.elementary.complexes import Abs 

14from sympy.functions.elementary.miscellaneous import Max, Min 

15from sympy.functions.elementary.trigonometric import cos, sin 

16from sympy.series.limits import Limit 

17 

18 

19def difference_delta(expr, n=None, step=1): 

20 """Difference Operator. 

21 

22 Explanation 

23 =========== 

24 

25 Discrete analog of differential operator. Given a sequence x[n], 

26 returns the sequence x[n + step] - x[n]. 

27 

28 Examples 

29 ======== 

30 

31 >>> from sympy import difference_delta as dd 

32 >>> from sympy.abc import n 

33 >>> dd(n*(n + 1), n) 

34 2*n + 2 

35 >>> dd(n*(n + 1), n, 2) 

36 4*n + 6 

37 

38 References 

39 ========== 

40 

41 .. [1] https://reference.wolfram.com/language/ref/DifferenceDelta.html 

42 """ 

43 expr = sympify(expr) 

44 

45 if n is None: 

46 f = expr.free_symbols 

47 if len(f) == 1: 

48 n = f.pop() 

49 elif len(f) == 0: 

50 return S.Zero 

51 else: 

52 raise ValueError("Since there is more than one variable in the" 

53 " expression, a variable must be supplied to" 

54 " take the difference of %s" % expr) 

55 step = sympify(step) 

56 if step.is_number is False or step.is_finite is False: 

57 raise ValueError("Step should be a finite number.") 

58 

59 if hasattr(expr, '_eval_difference_delta'): 

60 result = expr._eval_difference_delta(n, step) 

61 if result: 

62 return result 

63 

64 return expr.subs(n, n + step) - expr 

65 

66 

67def dominant(expr, n): 

68 """Finds the dominant term in a sum, that is a term that dominates 

69 every other term. 

70 

71 Explanation 

72 =========== 

73 

74 If limit(a/b, n, oo) is oo then a dominates b. 

75 If limit(a/b, n, oo) is 0 then b dominates a. 

76 Otherwise, a and b are comparable. 

77 

78 If there is no unique dominant term, then returns ``None``. 

79 

80 Examples 

81 ======== 

82 

83 >>> from sympy import Sum 

84 >>> from sympy.series.limitseq import dominant 

85 >>> from sympy.abc import n, k 

86 >>> dominant(5*n**3 + 4*n**2 + n + 1, n) 

87 5*n**3 

88 >>> dominant(2**n + Sum(k, (k, 0, n)), n) 

89 2**n 

90 

91 See Also 

92 ======== 

93 

94 sympy.series.limitseq.dominant 

95 """ 

96 terms = Add.make_args(expr.expand(func=True)) 

97 term0 = terms[-1] 

98 comp = [term0] # comparable terms 

99 for t in terms[:-1]: 

100 r = term0/t 

101 e = r.gammasimp() 

102 if e == r: 

103 e = r.factor() 

104 l = limit_seq(e, n) 

105 if l is None: 

106 return None 

107 elif l.is_zero: 

108 term0 = t 

109 comp = [term0] 

110 elif l not in [S.Infinity, S.NegativeInfinity]: 

111 comp.append(t) 

112 if len(comp) > 1: 

113 return None 

114 return term0 

115 

116 

117def _limit_inf(expr, n): 

118 try: 

119 return Limit(expr, n, S.Infinity).doit(deep=False) 

120 except (NotImplementedError, PoleError): 

121 return None 

122 

123 

124def _limit_seq(expr, n, trials): 

125 from sympy.concrete.summations import Sum 

126 

127 for i in range(trials): 

128 if not expr.has(Sum): 

129 result = _limit_inf(expr, n) 

130 if result is not None: 

131 return result 

132 

133 num, den = expr.as_numer_denom() 

134 if not den.has(n) or not num.has(n): 

135 result = _limit_inf(expr.doit(), n) 

136 if result is not None: 

137 return result 

138 return None 

139 

140 num, den = (difference_delta(t.expand(), n) for t in [num, den]) 

141 expr = (num / den).gammasimp() 

142 

143 if not expr.has(Sum): 

144 result = _limit_inf(expr, n) 

145 if result is not None: 

146 return result 

147 

148 num, den = expr.as_numer_denom() 

149 

150 num = dominant(num, n) 

151 if num is None: 

152 return None 

153 

154 den = dominant(den, n) 

155 if den is None: 

156 return None 

157 

158 expr = (num / den).gammasimp() 

159 

160 

161def limit_seq(expr, n=None, trials=5): 

162 """Finds the limit of a sequence as index ``n`` tends to infinity. 

163 

164 Parameters 

165 ========== 

166 

167 expr : Expr 

168 SymPy expression for the ``n-th`` term of the sequence 

169 n : Symbol, optional 

170 The index of the sequence, an integer that tends to positive 

171 infinity. If None, inferred from the expression unless it has 

172 multiple symbols. 

173 trials: int, optional 

174 The algorithm is highly recursive. ``trials`` is a safeguard from 

175 infinite recursion in case the limit is not easily computed by the 

176 algorithm. Try increasing ``trials`` if the algorithm returns ``None``. 

177 

178 Admissible Terms 

179 ================ 

180 

181 The algorithm is designed for sequences built from rational functions, 

182 indefinite sums, and indefinite products over an indeterminate n. Terms of 

183 alternating sign are also allowed, but more complex oscillatory behavior is 

184 not supported. 

185 

186 Examples 

187 ======== 

188 

189 >>> from sympy import limit_seq, Sum, binomial 

190 >>> from sympy.abc import n, k, m 

191 >>> limit_seq((5*n**3 + 3*n**2 + 4) / (3*n**3 + 4*n - 5), n) 

192 5/3 

193 >>> limit_seq(binomial(2*n, n) / Sum(binomial(2*k, k), (k, 1, n)), n) 

194 3/4 

195 >>> limit_seq(Sum(k**2 * Sum(2**m/m, (m, 1, k)), (k, 1, n)) / (2**n*n), n) 

196 4 

197 

198 See Also 

199 ======== 

200 

201 sympy.series.limitseq.dominant 

202 

203 References 

204 ========== 

205 

206 .. [1] Computing Limits of Sequences - Manuel Kauers 

207 """ 

208 

209 from sympy.concrete.summations import Sum 

210 if n is None: 

211 free = expr.free_symbols 

212 if len(free) == 1: 

213 n = free.pop() 

214 elif not free: 

215 return expr 

216 else: 

217 raise ValueError("Expression has more than one variable. " 

218 "Please specify a variable.") 

219 elif n not in expr.free_symbols: 

220 return expr 

221 

222 expr = expr.rewrite(fibonacci, S.GoldenRatio) 

223 expr = expr.rewrite(factorial, subfactorial, gamma) 

224 n_ = Dummy("n", integer=True, positive=True) 

225 n1 = Dummy("n", odd=True, positive=True) 

226 n2 = Dummy("n", even=True, positive=True) 

227 

228 # If there is a negative term raised to a power involving n, or a 

229 # trigonometric function, then consider even and odd n separately. 

230 powers = (p.as_base_exp() for p in expr.atoms(Pow)) 

231 if (any(b.is_negative and e.has(n) for b, e in powers) or 

232 expr.has(cos, sin)): 

233 L1 = _limit_seq(expr.xreplace({n: n1}), n1, trials) 

234 if L1 is not None: 

235 L2 = _limit_seq(expr.xreplace({n: n2}), n2, trials) 

236 if L1 != L2: 

237 if L1.is_comparable and L2.is_comparable: 

238 return AccumulationBounds(Min(L1, L2), Max(L1, L2)) 

239 else: 

240 return None 

241 else: 

242 L1 = _limit_seq(expr.xreplace({n: n_}), n_, trials) 

243 if L1 is not None: 

244 return L1 

245 else: 

246 if expr.is_Add: 

247 limits = [limit_seq(term, n, trials) for term in expr.args] 

248 if any(result is None for result in limits): 

249 return None 

250 else: 

251 return Add(*limits) 

252 # Maybe the absolute value is easier to deal with (though not if 

253 # it has a Sum). If it tends to 0, the limit is 0. 

254 elif not expr.has(Sum): 

255 lim = _limit_seq(Abs(expr.xreplace({n: n_})), n_, trials) 

256 if lim is not None and lim.is_zero: 

257 return S.Zero