Coverage for /usr/lib/python3/dist-packages/sympy/ntheory/continued_fraction.py: 10%

1from __future__ import annotations 

2from sympy.core.exprtools import factor_terms 

3from sympy.core.numbers import Integer, Rational 

4from sympy.core.singleton import S 

5from sympy.core.symbol import Dummy 

6from sympy.core.sympify import _sympify 

7from sympy.utilities.misc import as_int 

8 

9 

10def continued_fraction(a) -> list: 

11 """Return the continued fraction representation of a Rational or 

12 quadratic irrational. 

13 

14 Examples 

15 ======== 

16 

17 >>> from sympy.ntheory.continued_fraction import continued_fraction 

18 >>> from sympy import sqrt 

19 >>> continued_fraction((1 + 2*sqrt(3))/5) 

20 [0, 1, [8, 3, 34, 3]] 

21 

22 See Also 

23 ======== 

24 continued_fraction_periodic, continued_fraction_reduce, continued_fraction_convergents 

25 """ 

26 e = _sympify(a) 

27 if all(i.is_Rational for i in e.atoms()): 

28 if e.is_Integer: 

29 return continued_fraction_periodic(e, 1, 0) 

30 elif e.is_Rational: 

31 return continued_fraction_periodic(e.p, e.q, 0) 

32 elif e.is_Pow and e.exp is S.Half and e.base.is_Integer: 

33 return continued_fraction_periodic(0, 1, e.base) 

34 elif e.is_Mul and len(e.args) == 2 and ( 

35 e.args[0].is_Rational and 

36 e.args[1].is_Pow and 

37 e.args[1].base.is_Integer and 

38 e.args[1].exp is S.Half): 

39 a, b = e.args 

40 return continued_fraction_periodic(0, a.q, b.base, a.p) 

41 else: 

42 # this should not have to work very hard- no 

43 # simplification, cancel, etc... which should be 

44 # done by the user. e.g. This is a fancy 1 but 

45 # the user should simplify it first: 

46 # sqrt(2)*(1 + sqrt(2))/(sqrt(2) + 2) 

47 p, d = e.expand().as_numer_denom() 

48 if d.is_Integer: 

49 if p.is_Rational: 

50 return continued_fraction_periodic(p, d) 

51 # look for a + b*c 

52 # with c = sqrt(s) 

53 if p.is_Add and len(p.args) == 2: 

54 a, bc = p.args 

55 else: 

56 a = S.Zero 

57 bc = p 

58 if a.is_Integer: 

59 b = S.NaN 

60 if bc.is_Mul and len(bc.args) == 2: 

61 b, c = bc.args 

62 elif bc.is_Pow: 

63 b = Integer(1) 

64 c = bc 

65 if b.is_Integer and ( 

66 c.is_Pow and c.exp is S.Half and 

67 c.base.is_Integer): 

68 # (a + b*sqrt(c))/d 

69 c = c.base 

70 return continued_fraction_periodic(a, d, c, b) 

71 raise ValueError( 

72 'expecting a rational or quadratic irrational, not %s' % e) 

73 

74 

75def continued_fraction_periodic(p, q, d=0, s=1) -> list: 

76 r""" 

77 Find the periodic continued fraction expansion of a quadratic irrational. 

78 

79 Compute the continued fraction expansion of a rational or a 

80 quadratic irrational number, i.e. `\frac{p + s\sqrt{d}}{q}`, where 

81 `p`, `q \ne 0` and `d \ge 0` are integers. 

82 

83 Returns the continued fraction representation (canonical form) as 

84 a list of integers, optionally ending (for quadratic irrationals) 

85 with list of integers representing the repeating digits. 

86 

87 Parameters 

88 ========== 

89 

90 p : int 

91 the rational part of the number's numerator 

92 q : int 

93 the denominator of the number 

94 d : int, optional 

95 the irrational part (discriminator) of the number's numerator 

96 s : int, optional 

97 the coefficient of the irrational part 

98 

99 Examples 

100 ======== 

101 

102 >>> from sympy.ntheory.continued_fraction import continued_fraction_periodic 

103 >>> continued_fraction_periodic(3, 2, 7) 

104 [2, [1, 4, 1, 1]] 

105 

106 Golden ratio has the simplest continued fraction expansion: 

107 

108 >>> continued_fraction_periodic(1, 2, 5) 

109 [[1]] 

110 

111 If the discriminator is zero or a perfect square then the number will be a 

112 rational number: 

113 

114 >>> continued_fraction_periodic(4, 3, 0) 

115 [1, 3] 

116 >>> continued_fraction_periodic(4, 3, 49) 

117 [3, 1, 2] 

118 

119 See Also 

120 ======== 

121 

122 continued_fraction_iterator, continued_fraction_reduce 

123 

124 References 

125 ========== 

126 

127 .. [1] https://en.wikipedia.org/wiki/Periodic_continued_fraction 

128 .. [2] K. Rosen. Elementary Number theory and its applications. 

129 Addison-Wesley, 3 Sub edition, pages 379-381, January 1992. 

130 

131 """ 

132 from sympy.functions import sqrt, floor 

133 

134 p, q, d, s = list(map(as_int, [p, q, d, s])) 

135 

136 if d < 0: 

137 raise ValueError("expected non-negative for `d` but got %s" % d) 

138 

139 if q == 0: 

140 raise ValueError("The denominator cannot be 0.") 

141 

142 if not s: 

143 d = 0 

144 

145 # check for rational case 

146 sd = sqrt(d) 

147 if sd.is_Integer: 

148 return list(continued_fraction_iterator(Rational(p + s*sd, q))) 

149 

150 # irrational case with sd != Integer 

151 if q < 0: 

152 p, q, s = -p, -q, -s 

153 

154 n = (p + s*sd)/q 

155 if n < 0: 

156 w = floor(-n) 

157 f = -n - w 

158 one_f = continued_fraction(1 - f) # 1-f < 1 so cf is [0 ... [...]] 

159 one_f[0] -= w + 1 

160 return one_f 

161 

162 d *= s**2 

163 sd *= s 

164 

165 if (d - p**2)%q: 

166 d *= q**2 

167 sd *= q 

168 p *= q 

169 q *= q 

170 

171 terms: list[int] = [] 

172 pq = {} 

173 

174 while (p, q) not in pq: 

175 pq[(p, q)] = len(terms) 

176 terms.append((p + sd)//q) 

177 p = terms[-1]*q - p 

178 q = (d - p**2)//q 

179 

180 i = pq[(p, q)] 

181 return terms[:i] + [terms[i:]] # type: ignore 

182 

183 

184def continued_fraction_reduce(cf): 

185 """ 

186 Reduce a continued fraction to a rational or quadratic irrational. 

187 

188 Compute the rational or quadratic irrational number from its 

189 terminating or periodic continued fraction expansion. The 

190 continued fraction expansion (cf) should be supplied as a 

191 terminating iterator supplying the terms of the expansion. For 

192 terminating continued fractions, this is equivalent to 

193 ``list(continued_fraction_convergents(cf))[-1]``, only a little more 

194 efficient. If the expansion has a repeating part, a list of the 

195 repeating terms should be returned as the last element from the 

196 iterator. This is the format returned by 

197 continued_fraction_periodic. 

198 

199 For quadratic irrationals, returns the largest solution found, 

200 which is generally the one sought, if the fraction is in canonical 

201 form (all terms positive except possibly the first). 

202 

203 Examples 

204 ======== 

205 

206 >>> from sympy.ntheory.continued_fraction import continued_fraction_reduce 

207 >>> continued_fraction_reduce([1, 2, 3, 4, 5]) 

208 225/157 

209 >>> continued_fraction_reduce([-2, 1, 9, 7, 1, 2]) 

210 -256/233 

211 >>> continued_fraction_reduce([2, 1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8]).n(10) 

212 2.718281835 

213 >>> continued_fraction_reduce([1, 4, 2, [3, 1]]) 

214 (sqrt(21) + 287)/238 

215 >>> continued_fraction_reduce([[1]]) 

216 (1 + sqrt(5))/2 

217 >>> from sympy.ntheory.continued_fraction import continued_fraction_periodic 

218 >>> continued_fraction_reduce(continued_fraction_periodic(8, 5, 13)) 

219 (sqrt(13) + 8)/5 

220 

221 See Also 

222 ======== 

223 

224 continued_fraction_periodic 

225 

226 """ 

227 from sympy.solvers import solve 

228 

229 period = [] 

230 x = Dummy('x') 

231 

232 def untillist(cf): 

233 for nxt in cf: 

234 if isinstance(nxt, list): 

235 period.extend(nxt) 

236 yield x 

237 break 

238 yield nxt 

239 

240 a = S.Zero 

241 for a in continued_fraction_convergents(untillist(cf)): 

242 pass 

243 

244 if period: 

245 y = Dummy('y') 

246 solns = solve(continued_fraction_reduce(period + [y]) - y, y) 

247 solns.sort() 

248 pure = solns[-1] 

249 rv = a.subs(x, pure).radsimp() 

250 else: 

251 rv = a 

252 if rv.is_Add: 

253 rv = factor_terms(rv) 

254 if rv.is_Mul and rv.args[0] == -1: 

255 rv = rv.func(*rv.args) 

256 return rv 

257 

258 

259def continued_fraction_iterator(x): 

260 """ 

261 Return continued fraction expansion of x as iterator. 

262 

263 Examples 

264 ======== 

265 

266 >>> from sympy import Rational, pi 

267 >>> from sympy.ntheory.continued_fraction import continued_fraction_iterator 

268 

269 >>> list(continued_fraction_iterator(Rational(3, 8))) 

270 [0, 2, 1, 2] 

271 >>> list(continued_fraction_iterator(Rational(-3, 8))) 

272 [-1, 1, 1, 1, 2] 

273 

274 >>> for i, v in enumerate(continued_fraction_iterator(pi)): 

275 ... if i > 7: 

276 ... break 

277 ... print(v) 

278 3 

279 7 

280 15 

281 1 

282 292 

283 1 

284 1 

285 1 

286 

287 References 

288 ========== 

289 

290 .. [1] https://en.wikipedia.org/wiki/Continued_fraction 

291 

292 """ 

293 from sympy.functions import floor 

294 while True: 

295 i = floor(x) 

296 yield i 

297 x -= i 

298 if not x: 

299 break 

300 x = 1/x 

301 

302 

303def continued_fraction_convergents(cf): 

304 """ 

305 Return an iterator over the convergents of a continued fraction (cf). 

306 

307 The parameter should be an iterable returning successive 

308 partial quotients of the continued fraction, such as might be 

309 returned by continued_fraction_iterator. In computing the 

310 convergents, the continued fraction need not be strictly in 

311 canonical form (all integers, all but the first positive). 

312 Rational and negative elements may be present in the expansion. 

313 

314 Examples 

315 ======== 

316 

317 >>> from sympy.core import pi 

318 >>> from sympy import S 

319 >>> from sympy.ntheory.continued_fraction import \ 

320 continued_fraction_convergents, continued_fraction_iterator 

321 

322 >>> list(continued_fraction_convergents([0, 2, 1, 2])) 

323 [0, 1/2, 1/3, 3/8] 

324 

325 >>> list(continued_fraction_convergents([1, S('1/2'), -7, S('1/4')])) 

326 [1, 3, 19/5, 7] 

327 

328 >>> it = continued_fraction_convergents(continued_fraction_iterator(pi)) 

329 >>> for n in range(7): 

330 ... print(next(it)) 

331 3 

332 22/7 

333 333/106 

334 355/113 

335 103993/33102 

336 104348/33215 

337 208341/66317 

338 

339 See Also 

340 ======== 

341 

342 continued_fraction_iterator 

343 

344 """ 

345 p_2, q_2 = S.Zero, S.One 

346 p_1, q_1 = S.One, S.Zero 

347 for a in cf: 

348 p, q = a*p_1 + p_2, a*q_1 + q_2 

349 p_2, q_2 = p_1, q_1 

350 p_1, q_1 = p, q 

351 yield p/q