_IntegerGMP.py 27 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782
  1. # ===================================================================
  2. #
  3. # Copyright (c) 2014, Legrandin <helderijs@gmail.com>
  4. # All rights reserved.
  5. #
  6. # Redistribution and use in source and binary forms, with or without
  7. # modification, are permitted provided that the following conditions
  8. # are met:
  9. #
  10. # 1. Redistributions of source code must retain the above copyright
  11. # notice, this list of conditions and the following disclaimer.
  12. # 2. Redistributions in binary form must reproduce the above copyright
  13. # notice, this list of conditions and the following disclaimer in
  14. # the documentation and/or other materials provided with the
  15. # distribution.
  16. #
  17. # THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
  18. # "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
  19. # LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
  20. # FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
  21. # COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
  22. # INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
  23. # BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
  24. # LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
  25. # CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
  26. # LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
  27. # ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  28. # POSSIBILITY OF SUCH DAMAGE.
  29. # ===================================================================
  30. import sys
  31. from Crypto.Util.py3compat import tobytes, is_native_int
  32. from Crypto.Util._raw_api import (backend, load_lib,
  33. get_raw_buffer, get_c_string,
  34. null_pointer, create_string_buffer,
  35. c_ulong, c_size_t, c_uint8_ptr)
  36. from ._IntegerBase import IntegerBase
  37. gmp_defs = """typedef unsigned long UNIX_ULONG;
  38. typedef struct { int a; int b; void *c; } MPZ;
  39. typedef MPZ mpz_t[1];
  40. typedef UNIX_ULONG mp_bitcnt_t;
  41. void __gmpz_init (mpz_t x);
  42. void __gmpz_init_set (mpz_t rop, const mpz_t op);
  43. void __gmpz_init_set_ui (mpz_t rop, UNIX_ULONG op);
  44. UNIX_ULONG __gmpz_get_ui (const mpz_t op);
  45. void __gmpz_set (mpz_t rop, const mpz_t op);
  46. void __gmpz_set_ui (mpz_t rop, UNIX_ULONG op);
  47. void __gmpz_add (mpz_t rop, const mpz_t op1, const mpz_t op2);
  48. void __gmpz_add_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  49. void __gmpz_sub_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  50. void __gmpz_addmul (mpz_t rop, const mpz_t op1, const mpz_t op2);
  51. void __gmpz_addmul_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  52. void __gmpz_submul_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  53. void __gmpz_import (mpz_t rop, size_t count, int order, size_t size,
  54. int endian, size_t nails, const void *op);
  55. void * __gmpz_export (void *rop, size_t *countp, int order,
  56. size_t size,
  57. int endian, size_t nails, const mpz_t op);
  58. size_t __gmpz_sizeinbase (const mpz_t op, int base);
  59. void __gmpz_sub (mpz_t rop, const mpz_t op1, const mpz_t op2);
  60. void __gmpz_mul (mpz_t rop, const mpz_t op1, const mpz_t op2);
  61. void __gmpz_mul_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  62. int __gmpz_cmp (const mpz_t op1, const mpz_t op2);
  63. void __gmpz_powm (mpz_t rop, const mpz_t base, const mpz_t exp, const
  64. mpz_t mod);
  65. void __gmpz_powm_ui (mpz_t rop, const mpz_t base, UNIX_ULONG exp,
  66. const mpz_t mod);
  67. void __gmpz_pow_ui (mpz_t rop, const mpz_t base, UNIX_ULONG exp);
  68. void __gmpz_sqrt(mpz_t rop, const mpz_t op);
  69. void __gmpz_mod (mpz_t r, const mpz_t n, const mpz_t d);
  70. void __gmpz_neg (mpz_t rop, const mpz_t op);
  71. void __gmpz_abs (mpz_t rop, const mpz_t op);
  72. void __gmpz_and (mpz_t rop, const mpz_t op1, const mpz_t op2);
  73. void __gmpz_ior (mpz_t rop, const mpz_t op1, const mpz_t op2);
  74. void __gmpz_clear (mpz_t x);
  75. void __gmpz_tdiv_q_2exp (mpz_t q, const mpz_t n, mp_bitcnt_t b);
  76. void __gmpz_fdiv_q (mpz_t q, const mpz_t n, const mpz_t d);
  77. void __gmpz_mul_2exp (mpz_t rop, const mpz_t op1, mp_bitcnt_t op2);
  78. int __gmpz_tstbit (const mpz_t op, mp_bitcnt_t bit_index);
  79. int __gmpz_perfect_square_p (const mpz_t op);
  80. int __gmpz_jacobi (const mpz_t a, const mpz_t b);
  81. void __gmpz_gcd (mpz_t rop, const mpz_t op1, const mpz_t op2);
  82. UNIX_ULONG __gmpz_gcd_ui (mpz_t rop, const mpz_t op1,
  83. UNIX_ULONG op2);
  84. void __gmpz_lcm (mpz_t rop, const mpz_t op1, const mpz_t op2);
  85. int __gmpz_invert (mpz_t rop, const mpz_t op1, const mpz_t op2);
  86. int __gmpz_divisible_p (const mpz_t n, const mpz_t d);
  87. int __gmpz_divisible_ui_p (const mpz_t n, UNIX_ULONG d);
  88. """
  89. if sys.platform == "win32":
  90. raise ImportError("Not using GMP on Windows")
  91. lib = load_lib("gmp", gmp_defs)
  92. implementation = {"library": "gmp", "api": backend}
  93. if hasattr(lib, "__mpir_version"):
  94. raise ImportError("MPIR library detected")
  95. # In order to create a function that returns a pointer to
  96. # a new MPZ structure, we need to break the abstraction
  97. # and know exactly what ffi backend we have
  98. if implementation["api"] == "ctypes":
  99. from ctypes import Structure, c_int, c_void_p, byref
  100. class _MPZ(Structure):
  101. _fields_ = [('_mp_alloc', c_int),
  102. ('_mp_size', c_int),
  103. ('_mp_d', c_void_p)]
  104. def new_mpz():
  105. return byref(_MPZ())
  106. else:
  107. # We are using CFFI
  108. from Crypto.Util._raw_api import ffi
  109. def new_mpz():
  110. return ffi.new("MPZ*")
  111. # Lazy creation of GMP methods
  112. class _GMP(object):
  113. def __getattr__(self, name):
  114. if name.startswith("mpz_"):
  115. func_name = "__gmpz_" + name[4:]
  116. elif name.startswith("gmp_"):
  117. func_name = "__gmp_" + name[4:]
  118. else:
  119. raise AttributeError("Attribute %s is invalid" % name)
  120. func = getattr(lib, func_name)
  121. setattr(self, name, func)
  122. return func
  123. _gmp = _GMP()
  124. class IntegerGMP(IntegerBase):
  125. """A fast, arbitrary precision integer"""
  126. _zero_mpz_p = new_mpz()
  127. _gmp.mpz_init_set_ui(_zero_mpz_p, c_ulong(0))
  128. def __init__(self, value):
  129. """Initialize the integer to the given value."""
  130. self._mpz_p = new_mpz()
  131. self._initialized = False
  132. if isinstance(value, float):
  133. raise ValueError("A floating point type is not a natural number")
  134. if is_native_int(value):
  135. _gmp.mpz_init(self._mpz_p)
  136. self._initialized = True
  137. if value == 0:
  138. return
  139. tmp = new_mpz()
  140. _gmp.mpz_init(tmp)
  141. try:
  142. positive = value >= 0
  143. reduce = abs(value)
  144. slots = (reduce.bit_length() - 1) // 32 + 1
  145. while slots > 0:
  146. slots = slots - 1
  147. _gmp.mpz_set_ui(tmp,
  148. c_ulong(0xFFFFFFFF & (reduce >> (slots * 32))))
  149. _gmp.mpz_mul_2exp(tmp, tmp, c_ulong(slots * 32))
  150. _gmp.mpz_add(self._mpz_p, self._mpz_p, tmp)
  151. finally:
  152. _gmp.mpz_clear(tmp)
  153. if not positive:
  154. _gmp.mpz_neg(self._mpz_p, self._mpz_p)
  155. elif isinstance(value, IntegerGMP):
  156. _gmp.mpz_init_set(self._mpz_p, value._mpz_p)
  157. self._initialized = True
  158. else:
  159. raise NotImplementedError
  160. # Conversions
  161. def __int__(self):
  162. tmp = new_mpz()
  163. _gmp.mpz_init_set(tmp, self._mpz_p)
  164. try:
  165. value = 0
  166. slot = 0
  167. while _gmp.mpz_cmp(tmp, self._zero_mpz_p) != 0:
  168. lsb = _gmp.mpz_get_ui(tmp) & 0xFFFFFFFF
  169. value |= lsb << (slot * 32)
  170. _gmp.mpz_tdiv_q_2exp(tmp, tmp, c_ulong(32))
  171. slot = slot + 1
  172. finally:
  173. _gmp.mpz_clear(tmp)
  174. if self < 0:
  175. value = -value
  176. return int(value)
  177. def __str__(self):
  178. return str(int(self))
  179. def __repr__(self):
  180. return "Integer(%s)" % str(self)
  181. # Only Python 2.x
  182. def __hex__(self):
  183. return hex(int(self))
  184. # Only Python 3.x
  185. def __index__(self):
  186. return int(self)
  187. def to_bytes(self, block_size=0, byteorder='big'):
  188. """Convert the number into a byte string.
  189. This method encodes the number in network order and prepends
  190. as many zero bytes as required. It only works for non-negative
  191. values.
  192. :Parameters:
  193. block_size : integer
  194. The exact size the output byte string must have.
  195. If zero, the string has the minimal length.
  196. byteorder : string
  197. 'big' for big-endian integers (default), 'little' for litte-endian.
  198. :Returns:
  199. A byte string.
  200. :Raise ValueError:
  201. If the value is negative or if ``block_size`` is
  202. provided and the length of the byte string would exceed it.
  203. """
  204. if self < 0:
  205. raise ValueError("Conversion only valid for non-negative numbers")
  206. buf_len = (_gmp.mpz_sizeinbase(self._mpz_p, 2) + 7) // 8
  207. if buf_len > block_size > 0:
  208. raise ValueError("Number is too big to convert to byte string"
  209. " of prescribed length")
  210. buf = create_string_buffer(buf_len)
  211. _gmp.mpz_export(
  212. buf,
  213. null_pointer, # Ignore countp
  214. 1, # Big endian
  215. c_size_t(1), # Each word is 1 byte long
  216. 0, # Endianess within a word - not relevant
  217. c_size_t(0), # No nails
  218. self._mpz_p)
  219. result = b'\x00' * max(0, block_size - buf_len) + get_raw_buffer(buf)
  220. if byteorder == 'big':
  221. pass
  222. elif byteorder == 'little':
  223. result = bytearray(result)
  224. result.reverse()
  225. result = bytes(result)
  226. else:
  227. raise ValueError("Incorrect byteorder")
  228. return result
  229. @staticmethod
  230. def from_bytes(byte_string, byteorder='big'):
  231. """Convert a byte string into a number.
  232. :Parameters:
  233. byte_string : byte string
  234. The input number, encoded in network order.
  235. It can only be non-negative.
  236. byteorder : string
  237. 'big' for big-endian integers (default), 'little' for litte-endian.
  238. :Return:
  239. The ``Integer`` object carrying the same value as the input.
  240. """
  241. result = IntegerGMP(0)
  242. if byteorder == 'big':
  243. pass
  244. elif byteorder == 'little':
  245. byte_string = bytearray(byte_string)
  246. byte_string.reverse()
  247. else:
  248. raise ValueError("Incorrect byteorder")
  249. _gmp.mpz_import(
  250. result._mpz_p,
  251. c_size_t(len(byte_string)), # Amount of words to read
  252. 1, # Big endian
  253. c_size_t(1), # Each word is 1 byte long
  254. 0, # Endianess within a word - not relevant
  255. c_size_t(0), # No nails
  256. c_uint8_ptr(byte_string))
  257. return result
  258. # Relations
  259. def _apply_and_return(self, func, term):
  260. if not isinstance(term, IntegerGMP):
  261. term = IntegerGMP(term)
  262. return func(self._mpz_p, term._mpz_p)
  263. def __eq__(self, term):
  264. if not (isinstance(term, IntegerGMP) or is_native_int(term)):
  265. return False
  266. return self._apply_and_return(_gmp.mpz_cmp, term) == 0
  267. def __ne__(self, term):
  268. if not (isinstance(term, IntegerGMP) or is_native_int(term)):
  269. return True
  270. return self._apply_and_return(_gmp.mpz_cmp, term) != 0
  271. def __lt__(self, term):
  272. return self._apply_and_return(_gmp.mpz_cmp, term) < 0
  273. def __le__(self, term):
  274. return self._apply_and_return(_gmp.mpz_cmp, term) <= 0
  275. def __gt__(self, term):
  276. return self._apply_and_return(_gmp.mpz_cmp, term) > 0
  277. def __ge__(self, term):
  278. return self._apply_and_return(_gmp.mpz_cmp, term) >= 0
  279. def __nonzero__(self):
  280. return _gmp.mpz_cmp(self._mpz_p, self._zero_mpz_p) != 0
  281. __bool__ = __nonzero__
  282. def is_negative(self):
  283. return _gmp.mpz_cmp(self._mpz_p, self._zero_mpz_p) < 0
  284. # Arithmetic operations
  285. def __add__(self, term):
  286. result = IntegerGMP(0)
  287. if not isinstance(term, IntegerGMP):
  288. try:
  289. term = IntegerGMP(term)
  290. except NotImplementedError:
  291. return NotImplemented
  292. _gmp.mpz_add(result._mpz_p,
  293. self._mpz_p,
  294. term._mpz_p)
  295. return result
  296. def __sub__(self, term):
  297. result = IntegerGMP(0)
  298. if not isinstance(term, IntegerGMP):
  299. try:
  300. term = IntegerGMP(term)
  301. except NotImplementedError:
  302. return NotImplemented
  303. _gmp.mpz_sub(result._mpz_p,
  304. self._mpz_p,
  305. term._mpz_p)
  306. return result
  307. def __mul__(self, term):
  308. result = IntegerGMP(0)
  309. if not isinstance(term, IntegerGMP):
  310. try:
  311. term = IntegerGMP(term)
  312. except NotImplementedError:
  313. return NotImplemented
  314. _gmp.mpz_mul(result._mpz_p,
  315. self._mpz_p,
  316. term._mpz_p)
  317. return result
  318. def __floordiv__(self, divisor):
  319. if not isinstance(divisor, IntegerGMP):
  320. divisor = IntegerGMP(divisor)
  321. if _gmp.mpz_cmp(divisor._mpz_p,
  322. self._zero_mpz_p) == 0:
  323. raise ZeroDivisionError("Division by zero")
  324. result = IntegerGMP(0)
  325. _gmp.mpz_fdiv_q(result._mpz_p,
  326. self._mpz_p,
  327. divisor._mpz_p)
  328. return result
  329. def __mod__(self, divisor):
  330. if not isinstance(divisor, IntegerGMP):
  331. divisor = IntegerGMP(divisor)
  332. comp = _gmp.mpz_cmp(divisor._mpz_p,
  333. self._zero_mpz_p)
  334. if comp == 0:
  335. raise ZeroDivisionError("Division by zero")
  336. if comp < 0:
  337. raise ValueError("Modulus must be positive")
  338. result = IntegerGMP(0)
  339. _gmp.mpz_mod(result._mpz_p,
  340. self._mpz_p,
  341. divisor._mpz_p)
  342. return result
  343. def inplace_pow(self, exponent, modulus=None):
  344. if modulus is None:
  345. if exponent < 0:
  346. raise ValueError("Exponent must not be negative")
  347. # Normal exponentiation
  348. if exponent > 256:
  349. raise ValueError("Exponent is too big")
  350. _gmp.mpz_pow_ui(self._mpz_p,
  351. self._mpz_p, # Base
  352. c_ulong(int(exponent))
  353. )
  354. else:
  355. # Modular exponentiation
  356. if not isinstance(modulus, IntegerGMP):
  357. modulus = IntegerGMP(modulus)
  358. if not modulus:
  359. raise ZeroDivisionError("Division by zero")
  360. if modulus.is_negative():
  361. raise ValueError("Modulus must be positive")
  362. if is_native_int(exponent):
  363. if exponent < 0:
  364. raise ValueError("Exponent must not be negative")
  365. if exponent < 65536:
  366. _gmp.mpz_powm_ui(self._mpz_p,
  367. self._mpz_p,
  368. c_ulong(exponent),
  369. modulus._mpz_p)
  370. return self
  371. exponent = IntegerGMP(exponent)
  372. elif exponent.is_negative():
  373. raise ValueError("Exponent must not be negative")
  374. _gmp.mpz_powm(self._mpz_p,
  375. self._mpz_p,
  376. exponent._mpz_p,
  377. modulus._mpz_p)
  378. return self
  379. def __pow__(self, exponent, modulus=None):
  380. result = IntegerGMP(self)
  381. return result.inplace_pow(exponent, modulus)
  382. def __abs__(self):
  383. result = IntegerGMP(0)
  384. _gmp.mpz_abs(result._mpz_p, self._mpz_p)
  385. return result
  386. def sqrt(self, modulus=None):
  387. """Return the largest Integer that does not
  388. exceed the square root"""
  389. if modulus is None:
  390. if self < 0:
  391. raise ValueError("Square root of negative value")
  392. result = IntegerGMP(0)
  393. _gmp.mpz_sqrt(result._mpz_p,
  394. self._mpz_p)
  395. else:
  396. if modulus <= 0:
  397. raise ValueError("Modulus must be positive")
  398. modulus = int(modulus)
  399. result = IntegerGMP(self._tonelli_shanks(int(self) % modulus, modulus))
  400. return result
  401. def __iadd__(self, term):
  402. if is_native_int(term):
  403. if 0 <= term < 65536:
  404. _gmp.mpz_add_ui(self._mpz_p,
  405. self._mpz_p,
  406. c_ulong(term))
  407. return self
  408. if -65535 < term < 0:
  409. _gmp.mpz_sub_ui(self._mpz_p,
  410. self._mpz_p,
  411. c_ulong(-term))
  412. return self
  413. term = IntegerGMP(term)
  414. _gmp.mpz_add(self._mpz_p,
  415. self._mpz_p,
  416. term._mpz_p)
  417. return self
  418. def __isub__(self, term):
  419. if is_native_int(term):
  420. if 0 <= term < 65536:
  421. _gmp.mpz_sub_ui(self._mpz_p,
  422. self._mpz_p,
  423. c_ulong(term))
  424. return self
  425. if -65535 < term < 0:
  426. _gmp.mpz_add_ui(self._mpz_p,
  427. self._mpz_p,
  428. c_ulong(-term))
  429. return self
  430. term = IntegerGMP(term)
  431. _gmp.mpz_sub(self._mpz_p,
  432. self._mpz_p,
  433. term._mpz_p)
  434. return self
  435. def __imul__(self, term):
  436. if is_native_int(term):
  437. if 0 <= term < 65536:
  438. _gmp.mpz_mul_ui(self._mpz_p,
  439. self._mpz_p,
  440. c_ulong(term))
  441. return self
  442. if -65535 < term < 0:
  443. _gmp.mpz_mul_ui(self._mpz_p,
  444. self._mpz_p,
  445. c_ulong(-term))
  446. _gmp.mpz_neg(self._mpz_p, self._mpz_p)
  447. return self
  448. term = IntegerGMP(term)
  449. _gmp.mpz_mul(self._mpz_p,
  450. self._mpz_p,
  451. term._mpz_p)
  452. return self
  453. def __imod__(self, divisor):
  454. if not isinstance(divisor, IntegerGMP):
  455. divisor = IntegerGMP(divisor)
  456. comp = _gmp.mpz_cmp(divisor._mpz_p,
  457. divisor._zero_mpz_p)
  458. if comp == 0:
  459. raise ZeroDivisionError("Division by zero")
  460. if comp < 0:
  461. raise ValueError("Modulus must be positive")
  462. _gmp.mpz_mod(self._mpz_p,
  463. self._mpz_p,
  464. divisor._mpz_p)
  465. return self
  466. # Boolean/bit operations
  467. def __and__(self, term):
  468. result = IntegerGMP(0)
  469. if not isinstance(term, IntegerGMP):
  470. term = IntegerGMP(term)
  471. _gmp.mpz_and(result._mpz_p,
  472. self._mpz_p,
  473. term._mpz_p)
  474. return result
  475. def __or__(self, term):
  476. result = IntegerGMP(0)
  477. if not isinstance(term, IntegerGMP):
  478. term = IntegerGMP(term)
  479. _gmp.mpz_ior(result._mpz_p,
  480. self._mpz_p,
  481. term._mpz_p)
  482. return result
  483. def __rshift__(self, pos):
  484. result = IntegerGMP(0)
  485. if pos < 0:
  486. raise ValueError("negative shift count")
  487. if pos > 65536:
  488. if self < 0:
  489. return -1
  490. else:
  491. return 0
  492. _gmp.mpz_tdiv_q_2exp(result._mpz_p,
  493. self._mpz_p,
  494. c_ulong(int(pos)))
  495. return result
  496. def __irshift__(self, pos):
  497. if pos < 0:
  498. raise ValueError("negative shift count")
  499. if pos > 65536:
  500. if self < 0:
  501. return -1
  502. else:
  503. return 0
  504. _gmp.mpz_tdiv_q_2exp(self._mpz_p,
  505. self._mpz_p,
  506. c_ulong(int(pos)))
  507. return self
  508. def __lshift__(self, pos):
  509. result = IntegerGMP(0)
  510. if not 0 <= pos < 65536:
  511. raise ValueError("Incorrect shift count")
  512. _gmp.mpz_mul_2exp(result._mpz_p,
  513. self._mpz_p,
  514. c_ulong(int(pos)))
  515. return result
  516. def __ilshift__(self, pos):
  517. if not 0 <= pos < 65536:
  518. raise ValueError("Incorrect shift count")
  519. _gmp.mpz_mul_2exp(self._mpz_p,
  520. self._mpz_p,
  521. c_ulong(int(pos)))
  522. return self
  523. def get_bit(self, n):
  524. """Return True if the n-th bit is set to 1.
  525. Bit 0 is the least significant."""
  526. if self < 0:
  527. raise ValueError("no bit representation for negative values")
  528. if n < 0:
  529. raise ValueError("negative bit count")
  530. if n > 65536:
  531. return 0
  532. return bool(_gmp.mpz_tstbit(self._mpz_p,
  533. c_ulong(int(n))))
  534. # Extra
  535. def is_odd(self):
  536. return _gmp.mpz_tstbit(self._mpz_p, 0) == 1
  537. def is_even(self):
  538. return _gmp.mpz_tstbit(self._mpz_p, 0) == 0
  539. def size_in_bits(self):
  540. """Return the minimum number of bits that can encode the number."""
  541. if self < 0:
  542. raise ValueError("Conversion only valid for non-negative numbers")
  543. return _gmp.mpz_sizeinbase(self._mpz_p, 2)
  544. def size_in_bytes(self):
  545. """Return the minimum number of bytes that can encode the number."""
  546. return (self.size_in_bits() - 1) // 8 + 1
  547. def is_perfect_square(self):
  548. return _gmp.mpz_perfect_square_p(self._mpz_p) != 0
  549. def fail_if_divisible_by(self, small_prime):
  550. """Raise an exception if the small prime is a divisor."""
  551. if is_native_int(small_prime):
  552. if 0 < small_prime < 65536:
  553. if _gmp.mpz_divisible_ui_p(self._mpz_p,
  554. c_ulong(small_prime)):
  555. raise ValueError("The value is composite")
  556. return
  557. small_prime = IntegerGMP(small_prime)
  558. if _gmp.mpz_divisible_p(self._mpz_p,
  559. small_prime._mpz_p):
  560. raise ValueError("The value is composite")
  561. def multiply_accumulate(self, a, b):
  562. """Increment the number by the product of a and b."""
  563. if not isinstance(a, IntegerGMP):
  564. a = IntegerGMP(a)
  565. if is_native_int(b):
  566. if 0 < b < 65536:
  567. _gmp.mpz_addmul_ui(self._mpz_p,
  568. a._mpz_p,
  569. c_ulong(b))
  570. return self
  571. if -65535 < b < 0:
  572. _gmp.mpz_submul_ui(self._mpz_p,
  573. a._mpz_p,
  574. c_ulong(-b))
  575. return self
  576. b = IntegerGMP(b)
  577. _gmp.mpz_addmul(self._mpz_p,
  578. a._mpz_p,
  579. b._mpz_p)
  580. return self
  581. def set(self, source):
  582. """Set the Integer to have the given value"""
  583. if not isinstance(source, IntegerGMP):
  584. source = IntegerGMP(source)
  585. _gmp.mpz_set(self._mpz_p,
  586. source._mpz_p)
  587. return self
  588. def inplace_inverse(self, modulus):
  589. """Compute the inverse of this number in the ring of
  590. modulo integers.
  591. Raise an exception if no inverse exists.
  592. """
  593. if not isinstance(modulus, IntegerGMP):
  594. modulus = IntegerGMP(modulus)
  595. comp = _gmp.mpz_cmp(modulus._mpz_p,
  596. self._zero_mpz_p)
  597. if comp == 0:
  598. raise ZeroDivisionError("Modulus cannot be zero")
  599. if comp < 0:
  600. raise ValueError("Modulus must be positive")
  601. result = _gmp.mpz_invert(self._mpz_p,
  602. self._mpz_p,
  603. modulus._mpz_p)
  604. if not result:
  605. raise ValueError("No inverse value can be computed")
  606. return self
  607. def inverse(self, modulus):
  608. result = IntegerGMP(self)
  609. result.inplace_inverse(modulus)
  610. return result
  611. def gcd(self, term):
  612. """Compute the greatest common denominator between this
  613. number and another term."""
  614. result = IntegerGMP(0)
  615. if is_native_int(term):
  616. if 0 < term < 65535:
  617. _gmp.mpz_gcd_ui(result._mpz_p,
  618. self._mpz_p,
  619. c_ulong(term))
  620. return result
  621. term = IntegerGMP(term)
  622. _gmp.mpz_gcd(result._mpz_p, self._mpz_p, term._mpz_p)
  623. return result
  624. def lcm(self, term):
  625. """Compute the least common multiplier between this
  626. number and another term."""
  627. result = IntegerGMP(0)
  628. if not isinstance(term, IntegerGMP):
  629. term = IntegerGMP(term)
  630. _gmp.mpz_lcm(result._mpz_p, self._mpz_p, term._mpz_p)
  631. return result
  632. @staticmethod
  633. def jacobi_symbol(a, n):
  634. """Compute the Jacobi symbol"""
  635. if not isinstance(a, IntegerGMP):
  636. a = IntegerGMP(a)
  637. if not isinstance(n, IntegerGMP):
  638. n = IntegerGMP(n)
  639. if n <= 0 or n.is_even():
  640. raise ValueError("n must be positive odd for the Jacobi symbol")
  641. return _gmp.mpz_jacobi(a._mpz_p, n._mpz_p)
  642. @staticmethod
  643. def _mult_modulo_bytes(term1, term2, modulus):
  644. if not isinstance(term1, IntegerGMP):
  645. term1 = IntegerGMP(term1)
  646. if not isinstance(term2, IntegerGMP):
  647. term2 = IntegerGMP(term2)
  648. if not isinstance(modulus, IntegerGMP):
  649. modulus = IntegerGMP(modulus)
  650. if modulus < 0:
  651. raise ValueError("Modulus must be positive")
  652. if modulus == 0:
  653. raise ZeroDivisionError("Modulus cannot be zero")
  654. if (modulus & 1) == 0:
  655. raise ValueError("Odd modulus is required")
  656. numbers_len = len(modulus.to_bytes())
  657. result = ((term1 * term2) % modulus).to_bytes(numbers_len)
  658. return result
  659. # Clean-up
  660. def __del__(self):
  661. try:
  662. if self._mpz_p is not None:
  663. if self._initialized:
  664. _gmp.mpz_clear(self._mpz_p)
  665. self._mpz_p = None
  666. except AttributeError:
  667. pass