Coverage for /usr/lib/python3/dist-packages/mpmath/ctx_mp_python.py: 30%

725 statements  

« prev     ^ index     » next       coverage.py v7.9.1, created at 2025-06-14 15:55 +0200

1#from ctx_base import StandardBaseContext 

2 

3from .libmp.backend import basestring, exec_ 

4 

5from .libmp import (MPZ, MPZ_ZERO, MPZ_ONE, int_types, repr_dps, 

6 round_floor, round_ceiling, dps_to_prec, round_nearest, prec_to_dps, 

7 ComplexResult, to_pickable, from_pickable, normalize, 

8 from_int, from_float, from_npfloat, from_Decimal, from_str, to_int, to_float, to_str, 

9 from_rational, from_man_exp, 

10 fone, fzero, finf, fninf, fnan, 

11 mpf_abs, mpf_pos, mpf_neg, mpf_add, mpf_sub, mpf_mul, mpf_mul_int, 

12 mpf_div, mpf_rdiv_int, mpf_pow_int, mpf_mod, 

13 mpf_eq, mpf_cmp, mpf_lt, mpf_gt, mpf_le, mpf_ge, 

14 mpf_hash, mpf_rand, 

15 mpf_sum, 

16 bitcount, to_fixed, 

17 mpc_to_str, 

18 mpc_to_complex, mpc_hash, mpc_pos, mpc_is_nonzero, mpc_neg, mpc_conjugate, 

19 mpc_abs, mpc_add, mpc_add_mpf, mpc_sub, mpc_sub_mpf, mpc_mul, mpc_mul_mpf, 

20 mpc_mul_int, mpc_div, mpc_div_mpf, mpc_pow, mpc_pow_mpf, mpc_pow_int, 

21 mpc_mpf_div, 

22 mpf_pow, 

23 mpf_pi, mpf_degree, mpf_e, mpf_phi, mpf_ln2, mpf_ln10, 

24 mpf_euler, mpf_catalan, mpf_apery, mpf_khinchin, 

25 mpf_glaisher, mpf_twinprime, mpf_mertens, 

26 int_types) 

27 

28from . import rational 

29from . import function_docs 

30 

31new = object.__new__ 

32 

33class mpnumeric(object): 

34 """Base class for mpf and mpc.""" 

35 __slots__ = [] 

36 def __new__(cls, val): 

37 raise NotImplementedError 

38 

39class _mpf(mpnumeric): 

40 """ 

41 An mpf instance holds a real-valued floating-point number. mpf:s 

42 work analogously to Python floats, but support arbitrary-precision 

43 arithmetic. 

44 """ 

45 __slots__ = ['_mpf_'] 

46 

47 def __new__(cls, val=fzero, **kwargs): 

48 """A new mpf can be created from a Python float, an int, a 

49 or a decimal string representing a number in floating-point 

50 format.""" 

51 prec, rounding = cls.context._prec_rounding 

52 if kwargs: 

53 prec = kwargs.get('prec', prec) 

54 if 'dps' in kwargs: 

55 prec = dps_to_prec(kwargs['dps']) 

56 rounding = kwargs.get('rounding', rounding) 

57 if type(val) is cls: 

58 sign, man, exp, bc = val._mpf_ 

59 if (not man) and exp: 

60 return val 

61 v = new(cls) 

62 v._mpf_ = normalize(sign, man, exp, bc, prec, rounding) 

63 return v 

64 elif type(val) is tuple: 

65 if len(val) == 2: 

66 v = new(cls) 

67 v._mpf_ = from_man_exp(val[0], val[1], prec, rounding) 

68 return v 

69 if len(val) == 4: 

70 sign, man, exp, bc = val 

71 v = new(cls) 

72 v._mpf_ = normalize(sign, MPZ(man), exp, bc, prec, rounding) 

73 return v 

74 raise ValueError 

75 else: 

76 v = new(cls) 

77 v._mpf_ = mpf_pos(cls.mpf_convert_arg(val, prec, rounding), prec, rounding) 

78 return v 

79 

80 @classmethod 

81 def mpf_convert_arg(cls, x, prec, rounding): 

82 if isinstance(x, int_types): return from_int(x) 

83 if isinstance(x, float): return from_float(x) 

84 if isinstance(x, basestring): return from_str(x, prec, rounding) 

85 if isinstance(x, cls.context.constant): return x.func(prec, rounding) 

86 if hasattr(x, '_mpf_'): return x._mpf_ 

87 if hasattr(x, '_mpmath_'): 

88 t = cls.context.convert(x._mpmath_(prec, rounding)) 

89 if hasattr(t, '_mpf_'): 

90 return t._mpf_ 

91 if hasattr(x, '_mpi_'): 

92 a, b = x._mpi_ 

93 if a == b: 

94 return a 

95 raise ValueError("can only create mpf from zero-width interval") 

96 raise TypeError("cannot create mpf from " + repr(x)) 

97 

98 @classmethod 

99 def mpf_convert_rhs(cls, x): 

100 if isinstance(x, int_types): return from_int(x) 

101 if isinstance(x, float): return from_float(x) 

102 if isinstance(x, complex_types): return cls.context.mpc(x) 

103 if isinstance(x, rational.mpq): 

104 p, q = x._mpq_ 

105 return from_rational(p, q, cls.context.prec) 

106 if hasattr(x, '_mpf_'): return x._mpf_ 

107 if hasattr(x, '_mpmath_'): 

108 t = cls.context.convert(x._mpmath_(*cls.context._prec_rounding)) 

109 if hasattr(t, '_mpf_'): 

110 return t._mpf_ 

111 return t 

112 return NotImplemented 

113 

114 @classmethod 

115 def mpf_convert_lhs(cls, x): 

116 x = cls.mpf_convert_rhs(x) 

117 if type(x) is tuple: 

118 return cls.context.make_mpf(x) 

119 return x 

120 

121 man_exp = property(lambda self: self._mpf_[1:3]) 

122 man = property(lambda self: self._mpf_[1]) 

123 exp = property(lambda self: self._mpf_[2]) 

124 bc = property(lambda self: self._mpf_[3]) 

125 

126 real = property(lambda self: self) 

127 imag = property(lambda self: self.context.zero) 

128 

129 conjugate = lambda self: self 

130 

131 def __getstate__(self): return to_pickable(self._mpf_) 

132 def __setstate__(self, val): self._mpf_ = from_pickable(val) 

133 

134 def __repr__(s): 

135 if s.context.pretty: 

136 return str(s) 

137 return "mpf('%s')" % to_str(s._mpf_, s.context._repr_digits) 

138 

139 def __str__(s): return to_str(s._mpf_, s.context._str_digits) 

140 def __hash__(s): return mpf_hash(s._mpf_) 

141 def __int__(s): return int(to_int(s._mpf_)) 

142 def __long__(s): return long(to_int(s._mpf_)) 

143 def __float__(s): return to_float(s._mpf_, rnd=s.context._prec_rounding[1]) 

144 def __complex__(s): return complex(float(s)) 

145 def __nonzero__(s): return s._mpf_ != fzero 

146 

147 __bool__ = __nonzero__ 

148 

149 def __abs__(s): 

150 cls, new, (prec, rounding) = s._ctxdata 

151 v = new(cls) 

152 v._mpf_ = mpf_abs(s._mpf_, prec, rounding) 

153 return v 

154 

155 def __pos__(s): 

156 cls, new, (prec, rounding) = s._ctxdata 

157 v = new(cls) 

158 v._mpf_ = mpf_pos(s._mpf_, prec, rounding) 

159 return v 

160 

161 def __neg__(s): 

162 cls, new, (prec, rounding) = s._ctxdata 

163 v = new(cls) 

164 v._mpf_ = mpf_neg(s._mpf_, prec, rounding) 

165 return v 

166 

167 def _cmp(s, t, func): 

168 if hasattr(t, '_mpf_'): 

169 t = t._mpf_ 

170 else: 

171 t = s.mpf_convert_rhs(t) 

172 if t is NotImplemented: 

173 return t 

174 return func(s._mpf_, t) 

175 

176 def __cmp__(s, t): return s._cmp(t, mpf_cmp) 

177 def __lt__(s, t): return s._cmp(t, mpf_lt) 

178 def __gt__(s, t): return s._cmp(t, mpf_gt) 

179 def __le__(s, t): return s._cmp(t, mpf_le) 

180 def __ge__(s, t): return s._cmp(t, mpf_ge) 

181 

182 def __ne__(s, t): 

183 v = s.__eq__(t) 

184 if v is NotImplemented: 

185 return v 

186 return not v 

187 

188 def __rsub__(s, t): 

189 cls, new, (prec, rounding) = s._ctxdata 

190 if type(t) in int_types: 

191 v = new(cls) 

192 v._mpf_ = mpf_sub(from_int(t), s._mpf_, prec, rounding) 

193 return v 

194 t = s.mpf_convert_lhs(t) 

195 if t is NotImplemented: 

196 return t 

197 return t - s 

198 

199 def __rdiv__(s, t): 

200 cls, new, (prec, rounding) = s._ctxdata 

201 if isinstance(t, int_types): 

202 v = new(cls) 

203 v._mpf_ = mpf_rdiv_int(t, s._mpf_, prec, rounding) 

204 return v 

205 t = s.mpf_convert_lhs(t) 

206 if t is NotImplemented: 

207 return t 

208 return t / s 

209 

210 def __rpow__(s, t): 

211 t = s.mpf_convert_lhs(t) 

212 if t is NotImplemented: 

213 return t 

214 return t ** s 

215 

216 def __rmod__(s, t): 

217 t = s.mpf_convert_lhs(t) 

218 if t is NotImplemented: 

219 return t 

220 return t % s 

221 

222 def sqrt(s): 

223 return s.context.sqrt(s) 

224 

225 def ae(s, t, rel_eps=None, abs_eps=None): 

226 return s.context.almosteq(s, t, rel_eps, abs_eps) 

227 

228 def to_fixed(self, prec): 

229 return to_fixed(self._mpf_, prec) 

230 

231 def __round__(self, *args): 

232 return round(float(self), *args) 

233 

234mpf_binary_op = """ 

235def %NAME%(self, other): 

236 mpf, new, (prec, rounding) = self._ctxdata 

237 sval = self._mpf_ 

238 if hasattr(other, '_mpf_'): 

239 tval = other._mpf_ 

240 %WITH_MPF% 

241 ttype = type(other) 

242 if ttype in int_types: 

243 %WITH_INT% 

244 elif ttype is float: 

245 tval = from_float(other) 

246 %WITH_MPF% 

247 elif hasattr(other, '_mpc_'): 

248 tval = other._mpc_ 

249 mpc = type(other) 

250 %WITH_MPC% 

251 elif ttype is complex: 

252 tval = from_float(other.real), from_float(other.imag) 

253 mpc = self.context.mpc 

254 %WITH_MPC% 

255 if isinstance(other, mpnumeric): 

256 return NotImplemented 

257 try: 

258 other = mpf.context.convert(other, strings=False) 

259 except TypeError: 

260 return NotImplemented 

261 return self.%NAME%(other) 

262""" 

263 

264return_mpf = "; obj = new(mpf); obj._mpf_ = val; return obj" 

265return_mpc = "; obj = new(mpc); obj._mpc_ = val; return obj" 

266 

267mpf_pow_same = """ 

268 try: 

269 val = mpf_pow(sval, tval, prec, rounding) %s 

270 except ComplexResult: 

271 if mpf.context.trap_complex: 

272 raise 

273 mpc = mpf.context.mpc 

274 val = mpc_pow((sval, fzero), (tval, fzero), prec, rounding) %s 

275""" % (return_mpf, return_mpc) 

276 

277def binary_op(name, with_mpf='', with_int='', with_mpc=''): 

278 code = mpf_binary_op 

279 code = code.replace("%WITH_INT%", with_int) 

280 code = code.replace("%WITH_MPC%", with_mpc) 

281 code = code.replace("%WITH_MPF%", with_mpf) 

282 code = code.replace("%NAME%", name) 

283 np = {} 

284 exec_(code, globals(), np) 

285 return np[name] 

286 

287_mpf.__eq__ = binary_op('__eq__', 

288 'return mpf_eq(sval, tval)', 

289 'return mpf_eq(sval, from_int(other))', 

290 'return (tval[1] == fzero) and mpf_eq(tval[0], sval)') 

291 

292_mpf.__add__ = binary_op('__add__', 

293 'val = mpf_add(sval, tval, prec, rounding)' + return_mpf, 

294 'val = mpf_add(sval, from_int(other), prec, rounding)' + return_mpf, 

295 'val = mpc_add_mpf(tval, sval, prec, rounding)' + return_mpc) 

296 

297_mpf.__sub__ = binary_op('__sub__', 

298 'val = mpf_sub(sval, tval, prec, rounding)' + return_mpf, 

299 'val = mpf_sub(sval, from_int(other), prec, rounding)' + return_mpf, 

300 'val = mpc_sub((sval, fzero), tval, prec, rounding)' + return_mpc) 

301 

302_mpf.__mul__ = binary_op('__mul__', 

303 'val = mpf_mul(sval, tval, prec, rounding)' + return_mpf, 

304 'val = mpf_mul_int(sval, other, prec, rounding)' + return_mpf, 

305 'val = mpc_mul_mpf(tval, sval, prec, rounding)' + return_mpc) 

306 

307_mpf.__div__ = binary_op('__div__', 

308 'val = mpf_div(sval, tval, prec, rounding)' + return_mpf, 

309 'val = mpf_div(sval, from_int(other), prec, rounding)' + return_mpf, 

310 'val = mpc_mpf_div(sval, tval, prec, rounding)' + return_mpc) 

311 

312_mpf.__mod__ = binary_op('__mod__', 

313 'val = mpf_mod(sval, tval, prec, rounding)' + return_mpf, 

314 'val = mpf_mod(sval, from_int(other), prec, rounding)' + return_mpf, 

315 'raise NotImplementedError("complex modulo")') 

316 

317_mpf.__pow__ = binary_op('__pow__', 

318 mpf_pow_same, 

319 'val = mpf_pow_int(sval, other, prec, rounding)' + return_mpf, 

320 'val = mpc_pow((sval, fzero), tval, prec, rounding)' + return_mpc) 

321 

322_mpf.__radd__ = _mpf.__add__ 

323_mpf.__rmul__ = _mpf.__mul__ 

324_mpf.__truediv__ = _mpf.__div__ 

325_mpf.__rtruediv__ = _mpf.__rdiv__ 

326 

327 

328class _constant(_mpf): 

329 """Represents a mathematical constant with dynamic precision. 

330 When printed or used in an arithmetic operation, a constant 

331 is converted to a regular mpf at the working precision. A 

332 regular mpf can also be obtained using the operation +x.""" 

333 

334 def __new__(cls, func, name, docname=''): 

335 a = object.__new__(cls) 

336 a.name = name 

337 a.func = func 

338 a.__doc__ = getattr(function_docs, docname, '') 

339 return a 

340 

341 def __call__(self, prec=None, dps=None, rounding=None): 

342 prec2, rounding2 = self.context._prec_rounding 

343 if not prec: prec = prec2 

344 if not rounding: rounding = rounding2 

345 if dps: prec = dps_to_prec(dps) 

346 return self.context.make_mpf(self.func(prec, rounding)) 

347 

348 @property 

349 def _mpf_(self): 

350 prec, rounding = self.context._prec_rounding 

351 return self.func(prec, rounding) 

352 

353 def __repr__(self): 

354 return "<%s: %s~>" % (self.name, self.context.nstr(self(dps=15))) 

355 

356 

357class _mpc(mpnumeric): 

358 """ 

359 An mpc represents a complex number using a pair of mpf:s (one 

360 for the real part and another for the imaginary part.) The mpc 

361 class behaves fairly similarly to Python's complex type. 

362 """ 

363 

364 __slots__ = ['_mpc_'] 

365 

366 def __new__(cls, real=0, imag=0): 

367 s = object.__new__(cls) 

368 if isinstance(real, complex_types): 

369 real, imag = real.real, real.imag 

370 elif hasattr(real, '_mpc_'): 

371 s._mpc_ = real._mpc_ 

372 return s 

373 real = cls.context.mpf(real) 

374 imag = cls.context.mpf(imag) 

375 s._mpc_ = (real._mpf_, imag._mpf_) 

376 return s 

377 

378 real = property(lambda self: self.context.make_mpf(self._mpc_[0])) 

379 imag = property(lambda self: self.context.make_mpf(self._mpc_[1])) 

380 

381 def __getstate__(self): 

382 return to_pickable(self._mpc_[0]), to_pickable(self._mpc_[1]) 

383 

384 def __setstate__(self, val): 

385 self._mpc_ = from_pickable(val[0]), from_pickable(val[1]) 

386 

387 def __repr__(s): 

388 if s.context.pretty: 

389 return str(s) 

390 r = repr(s.real)[4:-1] 

391 i = repr(s.imag)[4:-1] 

392 return "%s(real=%s, imag=%s)" % (type(s).__name__, r, i) 

393 

394 def __str__(s): 

395 return "(%s)" % mpc_to_str(s._mpc_, s.context._str_digits) 

396 

397 def __complex__(s): 

398 return mpc_to_complex(s._mpc_, rnd=s.context._prec_rounding[1]) 

399 

400 def __pos__(s): 

401 cls, new, (prec, rounding) = s._ctxdata 

402 v = new(cls) 

403 v._mpc_ = mpc_pos(s._mpc_, prec, rounding) 

404 return v 

405 

406 def __abs__(s): 

407 prec, rounding = s.context._prec_rounding 

408 v = new(s.context.mpf) 

409 v._mpf_ = mpc_abs(s._mpc_, prec, rounding) 

410 return v 

411 

412 def __neg__(s): 

413 cls, new, (prec, rounding) = s._ctxdata 

414 v = new(cls) 

415 v._mpc_ = mpc_neg(s._mpc_, prec, rounding) 

416 return v 

417 

418 def conjugate(s): 

419 cls, new, (prec, rounding) = s._ctxdata 

420 v = new(cls) 

421 v._mpc_ = mpc_conjugate(s._mpc_, prec, rounding) 

422 return v 

423 

424 def __nonzero__(s): 

425 return mpc_is_nonzero(s._mpc_) 

426 

427 __bool__ = __nonzero__ 

428 

429 def __hash__(s): 

430 return mpc_hash(s._mpc_) 

431 

432 @classmethod 

433 def mpc_convert_lhs(cls, x): 

434 try: 

435 y = cls.context.convert(x) 

436 return y 

437 except TypeError: 

438 return NotImplemented 

439 

440 def __eq__(s, t): 

441 if not hasattr(t, '_mpc_'): 

442 if isinstance(t, str): 

443 return False 

444 t = s.mpc_convert_lhs(t) 

445 if t is NotImplemented: 

446 return t 

447 return s.real == t.real and s.imag == t.imag 

448 

449 def __ne__(s, t): 

450 b = s.__eq__(t) 

451 if b is NotImplemented: 

452 return b 

453 return not b 

454 

455 def _compare(*args): 

456 raise TypeError("no ordering relation is defined for complex numbers") 

457 

458 __gt__ = _compare 

459 __le__ = _compare 

460 __gt__ = _compare 

461 __ge__ = _compare 

462 

463 def __add__(s, t): 

464 cls, new, (prec, rounding) = s._ctxdata 

465 if not hasattr(t, '_mpc_'): 

466 t = s.mpc_convert_lhs(t) 

467 if t is NotImplemented: 

468 return t 

469 if hasattr(t, '_mpf_'): 

470 v = new(cls) 

471 v._mpc_ = mpc_add_mpf(s._mpc_, t._mpf_, prec, rounding) 

472 return v 

473 v = new(cls) 

474 v._mpc_ = mpc_add(s._mpc_, t._mpc_, prec, rounding) 

475 return v 

476 

477 def __sub__(s, t): 

478 cls, new, (prec, rounding) = s._ctxdata 

479 if not hasattr(t, '_mpc_'): 

480 t = s.mpc_convert_lhs(t) 

481 if t is NotImplemented: 

482 return t 

483 if hasattr(t, '_mpf_'): 

484 v = new(cls) 

485 v._mpc_ = mpc_sub_mpf(s._mpc_, t._mpf_, prec, rounding) 

486 return v 

487 v = new(cls) 

488 v._mpc_ = mpc_sub(s._mpc_, t._mpc_, prec, rounding) 

489 return v 

490 

491 def __mul__(s, t): 

492 cls, new, (prec, rounding) = s._ctxdata 

493 if not hasattr(t, '_mpc_'): 

494 if isinstance(t, int_types): 

495 v = new(cls) 

496 v._mpc_ = mpc_mul_int(s._mpc_, t, prec, rounding) 

497 return v 

498 t = s.mpc_convert_lhs(t) 

499 if t is NotImplemented: 

500 return t 

501 if hasattr(t, '_mpf_'): 

502 v = new(cls) 

503 v._mpc_ = mpc_mul_mpf(s._mpc_, t._mpf_, prec, rounding) 

504 return v 

505 t = s.mpc_convert_lhs(t) 

506 v = new(cls) 

507 v._mpc_ = mpc_mul(s._mpc_, t._mpc_, prec, rounding) 

508 return v 

509 

510 def __div__(s, t): 

511 cls, new, (prec, rounding) = s._ctxdata 

512 if not hasattr(t, '_mpc_'): 

513 t = s.mpc_convert_lhs(t) 

514 if t is NotImplemented: 

515 return t 

516 if hasattr(t, '_mpf_'): 

517 v = new(cls) 

518 v._mpc_ = mpc_div_mpf(s._mpc_, t._mpf_, prec, rounding) 

519 return v 

520 v = new(cls) 

521 v._mpc_ = mpc_div(s._mpc_, t._mpc_, prec, rounding) 

522 return v 

523 

524 def __pow__(s, t): 

525 cls, new, (prec, rounding) = s._ctxdata 

526 if isinstance(t, int_types): 

527 v = new(cls) 

528 v._mpc_ = mpc_pow_int(s._mpc_, t, prec, rounding) 

529 return v 

530 t = s.mpc_convert_lhs(t) 

531 if t is NotImplemented: 

532 return t 

533 v = new(cls) 

534 if hasattr(t, '_mpf_'): 

535 v._mpc_ = mpc_pow_mpf(s._mpc_, t._mpf_, prec, rounding) 

536 else: 

537 v._mpc_ = mpc_pow(s._mpc_, t._mpc_, prec, rounding) 

538 return v 

539 

540 __radd__ = __add__ 

541 

542 def __rsub__(s, t): 

543 t = s.mpc_convert_lhs(t) 

544 if t is NotImplemented: 

545 return t 

546 return t - s 

547 

548 def __rmul__(s, t): 

549 cls, new, (prec, rounding) = s._ctxdata 

550 if isinstance(t, int_types): 

551 v = new(cls) 

552 v._mpc_ = mpc_mul_int(s._mpc_, t, prec, rounding) 

553 return v 

554 t = s.mpc_convert_lhs(t) 

555 if t is NotImplemented: 

556 return t 

557 return t * s 

558 

559 def __rdiv__(s, t): 

560 t = s.mpc_convert_lhs(t) 

561 if t is NotImplemented: 

562 return t 

563 return t / s 

564 

565 def __rpow__(s, t): 

566 t = s.mpc_convert_lhs(t) 

567 if t is NotImplemented: 

568 return t 

569 return t ** s 

570 

571 __truediv__ = __div__ 

572 __rtruediv__ = __rdiv__ 

573 

574 def ae(s, t, rel_eps=None, abs_eps=None): 

575 return s.context.almosteq(s, t, rel_eps, abs_eps) 

576 

577 

578complex_types = (complex, _mpc) 

579 

580 

581class PythonMPContext(object): 

582 

583 def __init__(ctx): 

584 ctx._prec_rounding = [53, round_nearest] 

585 ctx.mpf = type('mpf', (_mpf,), {}) 

586 ctx.mpc = type('mpc', (_mpc,), {}) 

587 ctx.mpf._ctxdata = [ctx.mpf, new, ctx._prec_rounding] 

588 ctx.mpc._ctxdata = [ctx.mpc, new, ctx._prec_rounding] 

589 ctx.mpf.context = ctx 

590 ctx.mpc.context = ctx 

591 ctx.constant = type('constant', (_constant,), {}) 

592 ctx.constant._ctxdata = [ctx.mpf, new, ctx._prec_rounding] 

593 ctx.constant.context = ctx 

594 

595 def make_mpf(ctx, v): 

596 a = new(ctx.mpf) 

597 a._mpf_ = v 

598 return a 

599 

600 def make_mpc(ctx, v): 

601 a = new(ctx.mpc) 

602 a._mpc_ = v 

603 return a 

604 

605 def default(ctx): 

606 ctx._prec = ctx._prec_rounding[0] = 53 

607 ctx._dps = 15 

608 ctx.trap_complex = False 

609 

610 def _set_prec(ctx, n): 

611 ctx._prec = ctx._prec_rounding[0] = max(1, int(n)) 

612 ctx._dps = prec_to_dps(n) 

613 

614 def _set_dps(ctx, n): 

615 ctx._prec = ctx._prec_rounding[0] = dps_to_prec(n) 

616 ctx._dps = max(1, int(n)) 

617 

618 prec = property(lambda ctx: ctx._prec, _set_prec) 

619 dps = property(lambda ctx: ctx._dps, _set_dps) 

620 

621 def convert(ctx, x, strings=True): 

622 """ 

623 Converts *x* to an ``mpf`` or ``mpc``. If *x* is of type ``mpf``, 

624 ``mpc``, ``int``, ``float``, ``complex``, the conversion 

625 will be performed losslessly. 

626 

627 If *x* is a string, the result will be rounded to the present 

628 working precision. Strings representing fractions or complex 

629 numbers are permitted. 

630 

631 >>> from mpmath import * 

632 >>> mp.dps = 15; mp.pretty = False 

633 >>> mpmathify(3.5) 

634 mpf('3.5') 

635 >>> mpmathify('2.1') 

636 mpf('2.1000000000000001') 

637 >>> mpmathify('3/4') 

638 mpf('0.75') 

639 >>> mpmathify('2+3j') 

640 mpc(real='2.0', imag='3.0') 

641 

642 """ 

643 if type(x) in ctx.types: return x 

644 if isinstance(x, int_types): return ctx.make_mpf(from_int(x)) 

645 if isinstance(x, float): return ctx.make_mpf(from_float(x)) 

646 if isinstance(x, complex): 

647 return ctx.make_mpc((from_float(x.real), from_float(x.imag))) 

648 if type(x).__module__ == 'numpy': return ctx.npconvert(x) 

649 if isinstance(x, numbers.Rational): # e.g. Fraction 

650 try: x = rational.mpq(int(x.numerator), int(x.denominator)) 

651 except: pass 

652 prec, rounding = ctx._prec_rounding 

653 if isinstance(x, rational.mpq): 

654 p, q = x._mpq_ 

655 return ctx.make_mpf(from_rational(p, q, prec)) 

656 if strings and isinstance(x, basestring): 

657 try: 

658 _mpf_ = from_str(x, prec, rounding) 

659 return ctx.make_mpf(_mpf_) 

660 except ValueError: 

661 pass 

662 if hasattr(x, '_mpf_'): return ctx.make_mpf(x._mpf_) 

663 if hasattr(x, '_mpc_'): return ctx.make_mpc(x._mpc_) 

664 if hasattr(x, '_mpmath_'): 

665 return ctx.convert(x._mpmath_(prec, rounding)) 

666 if type(x).__module__ == 'decimal': 

667 try: return ctx.make_mpf(from_Decimal(x, prec, rounding)) 

668 except: pass 

669 return ctx._convert_fallback(x, strings) 

670 

671 def npconvert(ctx, x): 

672 """ 

673 Converts *x* to an ``mpf`` or ``mpc``. *x* should be a numpy 

674 scalar. 

675 """ 

676 import numpy as np 

677 if isinstance(x, np.integer): return ctx.make_mpf(from_int(int(x))) 

678 if isinstance(x, np.floating): return ctx.make_mpf(from_npfloat(x)) 

679 if isinstance(x, np.complexfloating): 

680 return ctx.make_mpc((from_npfloat(x.real), from_npfloat(x.imag))) 

681 raise TypeError("cannot create mpf from " + repr(x)) 

682 

683 def isnan(ctx, x): 

684 """ 

685 Return *True* if *x* is a NaN (not-a-number), or for a complex 

686 number, whether either the real or complex part is NaN; 

687 otherwise return *False*:: 

688 

689 >>> from mpmath import * 

690 >>> isnan(3.14) 

691 False 

692 >>> isnan(nan) 

693 True 

694 >>> isnan(mpc(3.14,2.72)) 

695 False 

696 >>> isnan(mpc(3.14,nan)) 

697 True 

698 

699 """ 

700 if hasattr(x, "_mpf_"): 

701 return x._mpf_ == fnan 

702 if hasattr(x, "_mpc_"): 

703 return fnan in x._mpc_ 

704 if isinstance(x, int_types) or isinstance(x, rational.mpq): 

705 return False 

706 x = ctx.convert(x) 

707 if hasattr(x, '_mpf_') or hasattr(x, '_mpc_'): 

708 return ctx.isnan(x) 

709 raise TypeError("isnan() needs a number as input") 

710 

711 def isinf(ctx, x): 

712 """ 

713 Return *True* if the absolute value of *x* is infinite; 

714 otherwise return *False*:: 

715 

716 >>> from mpmath import * 

717 >>> isinf(inf) 

718 True 

719 >>> isinf(-inf) 

720 True 

721 >>> isinf(3) 

722 False 

723 >>> isinf(3+4j) 

724 False 

725 >>> isinf(mpc(3,inf)) 

726 True 

727 >>> isinf(mpc(inf,3)) 

728 True 

729 

730 """ 

731 if hasattr(x, "_mpf_"): 

732 return x._mpf_ in (finf, fninf) 

733 if hasattr(x, "_mpc_"): 

734 re, im = x._mpc_ 

735 return re in (finf, fninf) or im in (finf, fninf) 

736 if isinstance(x, int_types) or isinstance(x, rational.mpq): 

737 return False 

738 x = ctx.convert(x) 

739 if hasattr(x, '_mpf_') or hasattr(x, '_mpc_'): 

740 return ctx.isinf(x) 

741 raise TypeError("isinf() needs a number as input") 

742 

743 def isnormal(ctx, x): 

744 """ 

745 Determine whether *x* is "normal" in the sense of floating-point 

746 representation; that is, return *False* if *x* is zero, an 

747 infinity or NaN; otherwise return *True*. By extension, a 

748 complex number *x* is considered "normal" if its magnitude is 

749 normal:: 

750 

751 >>> from mpmath import * 

752 >>> isnormal(3) 

753 True 

754 >>> isnormal(0) 

755 False 

756 >>> isnormal(inf); isnormal(-inf); isnormal(nan) 

757 False 

758 False 

759 False 

760 >>> isnormal(0+0j) 

761 False 

762 >>> isnormal(0+3j) 

763 True 

764 >>> isnormal(mpc(2,nan)) 

765 False 

766 """ 

767 if hasattr(x, "_mpf_"): 

768 return bool(x._mpf_[1]) 

769 if hasattr(x, "_mpc_"): 

770 re, im = x._mpc_ 

771 re_normal = bool(re[1]) 

772 im_normal = bool(im[1]) 

773 if re == fzero: return im_normal 

774 if im == fzero: return re_normal 

775 return re_normal and im_normal 

776 if isinstance(x, int_types) or isinstance(x, rational.mpq): 

777 return bool(x) 

778 x = ctx.convert(x) 

779 if hasattr(x, '_mpf_') or hasattr(x, '_mpc_'): 

780 return ctx.isnormal(x) 

781 raise TypeError("isnormal() needs a number as input") 

782 

783 def isint(ctx, x, gaussian=False): 

784 """ 

785 Return *True* if *x* is integer-valued; otherwise return 

786 *False*:: 

787 

788 >>> from mpmath import * 

789 >>> isint(3) 

790 True 

791 >>> isint(mpf(3)) 

792 True 

793 >>> isint(3.2) 

794 False 

795 >>> isint(inf) 

796 False 

797 

798 Optionally, Gaussian integers can be checked for:: 

799 

800 >>> isint(3+0j) 

801 True 

802 >>> isint(3+2j) 

803 False 

804 >>> isint(3+2j, gaussian=True) 

805 True 

806 

807 """ 

808 if isinstance(x, int_types): 

809 return True 

810 if hasattr(x, "_mpf_"): 

811 sign, man, exp, bc = xval = x._mpf_ 

812 return bool((man and exp >= 0) or xval == fzero) 

813 if hasattr(x, "_mpc_"): 

814 re, im = x._mpc_ 

815 rsign, rman, rexp, rbc = re 

816 isign, iman, iexp, ibc = im 

817 re_isint = (rman and rexp >= 0) or re == fzero 

818 if gaussian: 

819 im_isint = (iman and iexp >= 0) or im == fzero 

820 return re_isint and im_isint 

821 return re_isint and im == fzero 

822 if isinstance(x, rational.mpq): 

823 p, q = x._mpq_ 

824 return p % q == 0 

825 x = ctx.convert(x) 

826 if hasattr(x, '_mpf_') or hasattr(x, '_mpc_'): 

827 return ctx.isint(x, gaussian) 

828 raise TypeError("isint() needs a number as input") 

829 

830 def fsum(ctx, terms, absolute=False, squared=False): 

831 """ 

832 Calculates a sum containing a finite number of terms (for infinite 

833 series, see :func:`~mpmath.nsum`). The terms will be converted to 

834 mpmath numbers. For len(terms) > 2, this function is generally 

835 faster and produces more accurate results than the builtin 

836 Python function :func:`sum`. 

837 

838 >>> from mpmath import * 

839 >>> mp.dps = 15; mp.pretty = False 

840 >>> fsum([1, 2, 0.5, 7]) 

841 mpf('10.5') 

842 

843 With squared=True each term is squared, and with absolute=True 

844 the absolute value of each term is used. 

845 """ 

846 prec, rnd = ctx._prec_rounding 

847 real = [] 

848 imag = [] 

849 for term in terms: 

850 reval = imval = 0 

851 if hasattr(term, "_mpf_"): 

852 reval = term._mpf_ 

853 elif hasattr(term, "_mpc_"): 

854 reval, imval = term._mpc_ 

855 else: 

856 term = ctx.convert(term) 

857 if hasattr(term, "_mpf_"): 

858 reval = term._mpf_ 

859 elif hasattr(term, "_mpc_"): 

860 reval, imval = term._mpc_ 

861 else: 

862 raise NotImplementedError 

863 if imval: 

864 if squared: 

865 if absolute: 

866 real.append(mpf_mul(reval,reval)) 

867 real.append(mpf_mul(imval,imval)) 

868 else: 

869 reval, imval = mpc_pow_int((reval,imval),2,prec+10) 

870 real.append(reval) 

871 imag.append(imval) 

872 elif absolute: 

873 real.append(mpc_abs((reval,imval), prec)) 

874 else: 

875 real.append(reval) 

876 imag.append(imval) 

877 else: 

878 if squared: 

879 reval = mpf_mul(reval, reval) 

880 elif absolute: 

881 reval = mpf_abs(reval) 

882 real.append(reval) 

883 s = mpf_sum(real, prec, rnd, absolute) 

884 if imag: 

885 s = ctx.make_mpc((s, mpf_sum(imag, prec, rnd))) 

886 else: 

887 s = ctx.make_mpf(s) 

888 return s 

889 

890 def fdot(ctx, A, B=None, conjugate=False): 

891 r""" 

892 Computes the dot product of the iterables `A` and `B`, 

893 

894 .. math :: 

895 

896 \sum_{k=0} A_k B_k. 

897 

898 Alternatively, :func:`~mpmath.fdot` accepts a single iterable of pairs. 

899 In other words, ``fdot(A,B)`` and ``fdot(zip(A,B))`` are equivalent. 

900 The elements are automatically converted to mpmath numbers. 

901 

902 With ``conjugate=True``, the elements in the second vector 

903 will be conjugated: 

904 

905 .. math :: 

906 

907 \sum_{k=0} A_k \overline{B_k} 

908 

909 **Examples** 

910 

911 >>> from mpmath import * 

912 >>> mp.dps = 15; mp.pretty = False 

913 >>> A = [2, 1.5, 3] 

914 >>> B = [1, -1, 2] 

915 >>> fdot(A, B) 

916 mpf('6.5') 

917 >>> list(zip(A, B)) 

918 [(2, 1), (1.5, -1), (3, 2)] 

919 >>> fdot(_) 

920 mpf('6.5') 

921 >>> A = [2, 1.5, 3j] 

922 >>> B = [1+j, 3, -1-j] 

923 >>> fdot(A, B) 

924 mpc(real='9.5', imag='-1.0') 

925 >>> fdot(A, B, conjugate=True) 

926 mpc(real='3.5', imag='-5.0') 

927 

928 """ 

929 if B is not None: 

930 A = zip(A, B) 

931 prec, rnd = ctx._prec_rounding 

932 real = [] 

933 imag = [] 

934 hasattr_ = hasattr 

935 types = (ctx.mpf, ctx.mpc) 

936 for a, b in A: 

937 if type(a) not in types: a = ctx.convert(a) 

938 if type(b) not in types: b = ctx.convert(b) 

939 a_real = hasattr_(a, "_mpf_") 

940 b_real = hasattr_(b, "_mpf_") 

941 if a_real and b_real: 

942 real.append(mpf_mul(a._mpf_, b._mpf_)) 

943 continue 

944 a_complex = hasattr_(a, "_mpc_") 

945 b_complex = hasattr_(b, "_mpc_") 

946 if a_real and b_complex: 

947 aval = a._mpf_ 

948 bre, bim = b._mpc_ 

949 if conjugate: 

950 bim = mpf_neg(bim) 

951 real.append(mpf_mul(aval, bre)) 

952 imag.append(mpf_mul(aval, bim)) 

953 elif b_real and a_complex: 

954 are, aim = a._mpc_ 

955 bval = b._mpf_ 

956 real.append(mpf_mul(are, bval)) 

957 imag.append(mpf_mul(aim, bval)) 

958 elif a_complex and b_complex: 

959 #re, im = mpc_mul(a._mpc_, b._mpc_, prec+20) 

960 are, aim = a._mpc_ 

961 bre, bim = b._mpc_ 

962 if conjugate: 

963 bim = mpf_neg(bim) 

964 real.append(mpf_mul(are, bre)) 

965 real.append(mpf_neg(mpf_mul(aim, bim))) 

966 imag.append(mpf_mul(are, bim)) 

967 imag.append(mpf_mul(aim, bre)) 

968 else: 

969 raise NotImplementedError 

970 s = mpf_sum(real, prec, rnd) 

971 if imag: 

972 s = ctx.make_mpc((s, mpf_sum(imag, prec, rnd))) 

973 else: 

974 s = ctx.make_mpf(s) 

975 return s 

976 

977 def _wrap_libmp_function(ctx, mpf_f, mpc_f=None, mpi_f=None, doc="<no doc>"): 

978 """ 

979 Given a low-level mpf_ function, and optionally similar functions 

980 for mpc_ and mpi_, defines the function as a context method. 

981 

982 It is assumed that the return type is the same as that of 

983 the input; the exception is that propagation from mpf to mpc is possible 

984 by raising ComplexResult. 

985 

986 """ 

987 def f(x, **kwargs): 

988 if type(x) not in ctx.types: 

989 x = ctx.convert(x) 

990 prec, rounding = ctx._prec_rounding 

991 if kwargs: 

992 prec = kwargs.get('prec', prec) 

993 if 'dps' in kwargs: 

994 prec = dps_to_prec(kwargs['dps']) 

995 rounding = kwargs.get('rounding', rounding) 

996 if hasattr(x, '_mpf_'): 

997 try: 

998 return ctx.make_mpf(mpf_f(x._mpf_, prec, rounding)) 

999 except ComplexResult: 

1000 # Handle propagation to complex 

1001 if ctx.trap_complex: 

1002 raise 

1003 return ctx.make_mpc(mpc_f((x._mpf_, fzero), prec, rounding)) 

1004 elif hasattr(x, '_mpc_'): 

1005 return ctx.make_mpc(mpc_f(x._mpc_, prec, rounding)) 

1006 raise NotImplementedError("%s of a %s" % (name, type(x))) 

1007 name = mpf_f.__name__[4:] 

1008 f.__doc__ = function_docs.__dict__.get(name, "Computes the %s of x" % doc) 

1009 return f 

1010 

1011 # Called by SpecialFunctions.__init__() 

1012 @classmethod 

1013 def _wrap_specfun(cls, name, f, wrap): 

1014 if wrap: 

1015 def f_wrapped(ctx, *args, **kwargs): 

1016 convert = ctx.convert 

1017 args = [convert(a) for a in args] 

1018 prec = ctx.prec 

1019 try: 

1020 ctx.prec += 10 

1021 retval = f(ctx, *args, **kwargs) 

1022 finally: 

1023 ctx.prec = prec 

1024 return +retval 

1025 else: 

1026 f_wrapped = f 

1027 f_wrapped.__doc__ = function_docs.__dict__.get(name, f.__doc__) 

1028 setattr(cls, name, f_wrapped) 

1029 

1030 def _convert_param(ctx, x): 

1031 if hasattr(x, "_mpc_"): 

1032 v, im = x._mpc_ 

1033 if im != fzero: 

1034 return x, 'C' 

1035 elif hasattr(x, "_mpf_"): 

1036 v = x._mpf_ 

1037 else: 

1038 if type(x) in int_types: 

1039 return int(x), 'Z' 

1040 p = None 

1041 if isinstance(x, tuple): 

1042 p, q = x 

1043 elif hasattr(x, '_mpq_'): 

1044 p, q = x._mpq_ 

1045 elif isinstance(x, basestring) and '/' in x: 

1046 p, q = x.split('/') 

1047 p = int(p) 

1048 q = int(q) 

1049 if p is not None: 

1050 if not p % q: 

1051 return p // q, 'Z' 

1052 return ctx.mpq(p,q), 'Q' 

1053 x = ctx.convert(x) 

1054 if hasattr(x, "_mpc_"): 

1055 v, im = x._mpc_ 

1056 if im != fzero: 

1057 return x, 'C' 

1058 elif hasattr(x, "_mpf_"): 

1059 v = x._mpf_ 

1060 else: 

1061 return x, 'U' 

1062 sign, man, exp, bc = v 

1063 if man: 

1064 if exp >= -4: 

1065 if sign: 

1066 man = -man 

1067 if exp >= 0: 

1068 return int(man) << exp, 'Z' 

1069 if exp >= -4: 

1070 p, q = int(man), (1<<(-exp)) 

1071 return ctx.mpq(p,q), 'Q' 

1072 x = ctx.make_mpf(v) 

1073 return x, 'R' 

1074 elif not exp: 

1075 return 0, 'Z' 

1076 else: 

1077 return x, 'U' 

1078 

1079 def _mpf_mag(ctx, x): 

1080 sign, man, exp, bc = x 

1081 if man: 

1082 return exp+bc 

1083 if x == fzero: 

1084 return ctx.ninf 

1085 if x == finf or x == fninf: 

1086 return ctx.inf 

1087 return ctx.nan 

1088 

1089 def mag(ctx, x): 

1090 """ 

1091 Quick logarithmic magnitude estimate of a number. Returns an 

1092 integer or infinity `m` such that `|x| <= 2^m`. It is not 

1093 guaranteed that `m` is an optimal bound, but it will never 

1094 be too large by more than 2 (and probably not more than 1). 

1095 

1096 **Examples** 

1097 

1098 >>> from mpmath import * 

1099 >>> mp.pretty = True 

1100 >>> mag(10), mag(10.0), mag(mpf(10)), int(ceil(log(10,2))) 

1101 (4, 4, 4, 4) 

1102 >>> mag(10j), mag(10+10j) 

1103 (4, 5) 

1104 >>> mag(0.01), int(ceil(log(0.01,2))) 

1105 (-6, -6) 

1106 >>> mag(0), mag(inf), mag(-inf), mag(nan) 

1107 (-inf, +inf, +inf, nan) 

1108 

1109 """ 

1110 if hasattr(x, "_mpf_"): 

1111 return ctx._mpf_mag(x._mpf_) 

1112 elif hasattr(x, "_mpc_"): 

1113 r, i = x._mpc_ 

1114 if r == fzero: 

1115 return ctx._mpf_mag(i) 

1116 if i == fzero: 

1117 return ctx._mpf_mag(r) 

1118 return 1+max(ctx._mpf_mag(r), ctx._mpf_mag(i)) 

1119 elif isinstance(x, int_types): 

1120 if x: 

1121 return bitcount(abs(x)) 

1122 return ctx.ninf 

1123 elif isinstance(x, rational.mpq): 

1124 p, q = x._mpq_ 

1125 if p: 

1126 return 1 + bitcount(abs(p)) - bitcount(q) 

1127 return ctx.ninf 

1128 else: 

1129 x = ctx.convert(x) 

1130 if hasattr(x, "_mpf_") or hasattr(x, "_mpc_"): 

1131 return ctx.mag(x) 

1132 else: 

1133 raise TypeError("requires an mpf/mpc") 

1134 

1135 

1136# Register with "numbers" ABC 

1137# We do not subclass, hence we do not use the @abstractmethod checks. While 

1138# this is less invasive it may turn out that we do not actually support 

1139# parts of the expected interfaces. See 

1140# http://docs.python.org/2/library/numbers.html for list of abstract 

1141# methods. 

1142try: 

1143 import numbers 

1144 numbers.Complex.register(_mpc) 

1145 numbers.Real.register(_mpf) 

1146except ImportError: 

1147 pass