_Numbers_gmp.py 26 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720
  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. from Crypto.Util.py3compat import tobytes, b, bchr
  31. from Crypto.Util._raw_api import (backend, load_lib,
  32. get_raw_buffer, get_c_string,
  33. null_pointer, create_string_buffer,
  34. c_ulong, c_ulonglong, c_size_t)
  35. # GMP uses unsigned longs in several functions prototypes.
  36. # On a UNIX 64 bit platform that type takes 64 bits but in Windows 64
  37. # it is still 32 bits.
  38. # The intention of the MPIR developers is to maintain binary compatibility
  39. # so they probably assumed that that GMP would compile on Windows 64
  40. # by treating it as a UNIX platform.
  41. gmp_defs_common = """
  42. typedef struct { int a; int b; void *c; } MPZ;
  43. typedef MPZ mpz_t[1];
  44. typedef UNIX_ULONG mp_bitcnt_t;
  45. void __gmpz_init (mpz_t x);
  46. void __gmpz_init_set (mpz_t rop, const mpz_t op);
  47. void __gmpz_init_set_ui (mpz_t rop, UNIX_ULONG op);
  48. int __gmp_sscanf (const char *s, const char *fmt, ...);
  49. void __gmpz_set (mpz_t rop, const mpz_t op);
  50. int __gmp_snprintf (uint8_t *buf, size_t size, const char *fmt, ...);
  51. void __gmpz_add (mpz_t rop, const mpz_t op1, const mpz_t op2);
  52. void __gmpz_add_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  53. void __gmpz_sub_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  54. void __gmpz_addmul (mpz_t rop, const mpz_t op1, const mpz_t op2);
  55. void __gmpz_addmul_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  56. void __gmpz_submul_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  57. void __gmpz_import (mpz_t rop, size_t count, int order, size_t size,
  58. int endian, size_t nails, const void *op);
  59. void * __gmpz_export (void *rop, size_t *countp, int order,
  60. size_t size,
  61. int endian, size_t nails, const mpz_t op);
  62. size_t __gmpz_sizeinbase (const mpz_t op, int base);
  63. void __gmpz_sub (mpz_t rop, const mpz_t op1, const mpz_t op2);
  64. void __gmpz_mul (mpz_t rop, const mpz_t op1, const mpz_t op2);
  65. void __gmpz_mul_ui (mpz_t rop, const mpz_t op1, UNIX_ULONG op2);
  66. int __gmpz_cmp (const mpz_t op1, const mpz_t op2);
  67. void __gmpz_powm (mpz_t rop, const mpz_t base, const mpz_t exp, const
  68. mpz_t mod);
  69. void __gmpz_powm_ui (mpz_t rop, const mpz_t base, UNIX_ULONG exp,
  70. const mpz_t mod);
  71. void __gmpz_pow_ui (mpz_t rop, const mpz_t base, UNIX_ULONG exp);
  72. void __gmpz_sqrt(mpz_t rop, const mpz_t op);
  73. void __gmpz_mod (mpz_t r, const mpz_t n, const mpz_t d);
  74. void __gmpz_neg (mpz_t rop, const mpz_t op);
  75. void __gmpz_abs (mpz_t rop, const mpz_t op);
  76. void __gmpz_and (mpz_t rop, const mpz_t op1, const mpz_t op2);
  77. void __gmpz_ior (mpz_t rop, const mpz_t op1, const mpz_t op2);
  78. void __gmpz_clear (mpz_t x);
  79. void __gmpz_tdiv_q_2exp (mpz_t q, const mpz_t n, mp_bitcnt_t b);
  80. void __gmpz_fdiv_q (mpz_t q, const mpz_t n, const mpz_t d);
  81. void __gmpz_mul_2exp (mpz_t rop, const mpz_t op1, mp_bitcnt_t op2);
  82. int __gmpz_tstbit (const mpz_t op, mp_bitcnt_t bit_index);
  83. int __gmpz_perfect_square_p (const mpz_t op);
  84. int __gmpz_jacobi (const mpz_t a, const mpz_t b);
  85. void __gmpz_gcd (mpz_t rop, const mpz_t op1, const mpz_t op2);
  86. UNIX_ULONG __gmpz_gcd_ui (mpz_t rop, const mpz_t op1,
  87. UNIX_ULONG op2);
  88. void __gmpz_lcm (mpz_t rop, const mpz_t op1, const mpz_t op2);
  89. int __gmpz_invert (mpz_t rop, const mpz_t op1, const mpz_t op2);
  90. int __gmpz_divisible_p (const mpz_t n, const mpz_t d);
  91. int __gmpz_divisible_ui_p (const mpz_t n, UNIX_ULONG d);
  92. """
  93. try:
  94. gmp_defs = "typedef unsigned long UNIX_ULONG;" + gmp_defs_common
  95. lib = load_lib("gmp", gmp_defs)
  96. implementation = { "library":"gmp", "api":backend }
  97. except OSError:
  98. import platform
  99. bits, linkage = platform.architecture()
  100. if bits.startswith("64") and linkage.startswith("Win"):
  101. # MPIR uses unsigned long long where GMP uses unsigned long
  102. # (LLP64 vs LP64)
  103. gmp_defs = "typedef unsigned long long UNIX_ULONG;" + gmp_defs_common
  104. c_ulong = c_ulonglong
  105. # Try to load private MPIR lib first (wheel)
  106. try:
  107. from Crypto.Util._file_system import pycryptodome_filename
  108. mpir_dll = pycryptodome_filename(("Crypto", "Math"), "mpir.dll")
  109. lib = load_lib(mpir_dll, gmp_defs)
  110. except OSError:
  111. lib = load_lib("mpir", gmp_defs)
  112. implementation = { "library":"mpir", "api":backend }
  113. # In order to create a function that returns a pointer to
  114. # a new MPZ structure, we need to break the abstraction
  115. # and know exactly what ffi backend we have
  116. if implementation["api"] == "ctypes":
  117. from ctypes import Structure, c_int, c_void_p, byref
  118. class _MPZ(Structure):
  119. _fields_ = [('_mp_alloc', c_int),
  120. ('_mp_size', c_int),
  121. ('_mp_d', c_void_p)]
  122. def new_mpz():
  123. return byref(_MPZ())
  124. else:
  125. # We are using CFFI
  126. from Crypto.Util._raw_api import ffi
  127. def new_mpz():
  128. return ffi.new("MPZ*")
  129. # Unfortunately, all symbols exported by the GMP library start with "__"
  130. # and have no trailing underscore.
  131. # You cannot directly refer to them as members of the ctypes' library
  132. # object from within any class because Python will replace the double
  133. # underscore with "_classname_".
  134. class _GMP(object):
  135. pass
  136. _gmp = _GMP()
  137. _gmp = _GMP()
  138. _gmp.mpz_init = lib.__gmpz_init
  139. _gmp.mpz_init_set = lib.__gmpz_init_set
  140. _gmp.mpz_init_set_ui = lib.__gmpz_init_set_ui
  141. _gmp.mpz_set = lib.__gmpz_set
  142. _gmp.gmp_snprintf = lib.__gmp_snprintf
  143. _gmp.gmp_sscanf = lib.__gmp_sscanf
  144. _gmp.mpz_add = lib.__gmpz_add
  145. _gmp.mpz_add_ui = lib.__gmpz_add_ui
  146. _gmp.mpz_sub_ui = lib.__gmpz_sub_ui
  147. _gmp.mpz_addmul = lib.__gmpz_addmul
  148. _gmp.mpz_addmul_ui = lib.__gmpz_addmul_ui
  149. _gmp.mpz_submul_ui = lib.__gmpz_submul_ui
  150. _gmp.mpz_import = lib.__gmpz_import
  151. _gmp.mpz_export = lib.__gmpz_export
  152. _gmp.mpz_sizeinbase = lib.__gmpz_sizeinbase
  153. _gmp.mpz_sub = lib.__gmpz_sub
  154. _gmp.mpz_mul = lib.__gmpz_mul
  155. _gmp.mpz_mul_ui = lib.__gmpz_mul_ui
  156. _gmp.mpz_cmp = lib.__gmpz_cmp
  157. _gmp.mpz_powm = lib.__gmpz_powm
  158. _gmp.mpz_powm_ui = lib.__gmpz_powm_ui
  159. _gmp.mpz_pow_ui = lib.__gmpz_pow_ui
  160. _gmp.mpz_sqrt = lib.__gmpz_sqrt
  161. _gmp.mpz_mod = lib.__gmpz_mod
  162. _gmp.mpz_neg = lib.__gmpz_neg
  163. _gmp.mpz_abs = lib.__gmpz_abs
  164. _gmp.mpz_and = lib.__gmpz_and
  165. _gmp.mpz_ior = lib.__gmpz_ior
  166. _gmp.mpz_clear = lib.__gmpz_clear
  167. _gmp.mpz_tdiv_q_2exp = lib.__gmpz_tdiv_q_2exp
  168. _gmp.mpz_fdiv_q = lib.__gmpz_fdiv_q
  169. _gmp.mpz_mul_2exp = lib.__gmpz_mul_2exp
  170. _gmp.mpz_tstbit = lib.__gmpz_tstbit
  171. _gmp.mpz_perfect_square_p = lib.__gmpz_perfect_square_p
  172. _gmp.mpz_jacobi = lib.__gmpz_jacobi
  173. _gmp.mpz_gcd = lib.__gmpz_gcd
  174. _gmp.mpz_gcd_ui = lib.__gmpz_gcd_ui
  175. _gmp.mpz_lcm = lib.__gmpz_lcm
  176. _gmp.mpz_invert = lib.__gmpz_invert
  177. _gmp.mpz_divisible_p = lib.__gmpz_divisible_p
  178. _gmp.mpz_divisible_ui_p = lib.__gmpz_divisible_ui_p
  179. class Integer(object):
  180. """A fast, arbitrary precision integer"""
  181. _zero_mpz_p = new_mpz()
  182. _gmp.mpz_init_set_ui(_zero_mpz_p, c_ulong(0))
  183. def __init__(self, value):
  184. """Initialize the integer to the given value."""
  185. self._mpz_p = new_mpz()
  186. self._initialized = False
  187. if isinstance(value, float):
  188. raise ValueError("A floating point type is not a natural number")
  189. self._initialized = True
  190. if isinstance(value, (int, long)):
  191. _gmp.mpz_init(self._mpz_p)
  192. result = _gmp.gmp_sscanf(tobytes(str(value)), b("%Zd"), self._mpz_p)
  193. if result != 1:
  194. raise ValueError("Error converting '%d'" % value)
  195. else:
  196. _gmp.mpz_init_set(self._mpz_p, value._mpz_p)
  197. # Conversions
  198. def __int__(self):
  199. # buf will contain the integer encoded in decimal plus the trailing
  200. # zero, and possibly the negative sign.
  201. # dig10(x) < log10(x) + 1 = log2(x)/log2(10) + 1 < log2(x)/3 + 1
  202. buf_len = _gmp.mpz_sizeinbase(self._mpz_p, 2) // 3 + 3
  203. buf = create_string_buffer(buf_len)
  204. _gmp.gmp_snprintf(buf, c_size_t(buf_len), b("%Zd"), self._mpz_p)
  205. return int(get_c_string(buf))
  206. def __str__(self):
  207. return str(int(self))
  208. def __repr__(self):
  209. return "Integer(%s)" % str(self)
  210. def to_bytes(self, block_size=0):
  211. """Convert the number into a byte string.
  212. This method encodes the number in network order and prepends
  213. as many zero bytes as required. It only works for non-negative
  214. values.
  215. :Parameters:
  216. block_size : integer
  217. The exact size the output byte string must have.
  218. If zero, the string has the minimal length.
  219. :Returns:
  220. A byte string.
  221. :Raise ValueError:
  222. If the value is negative or if ``block_size`` is
  223. provided and the length of the byte string would exceed it.
  224. """
  225. if self < 0:
  226. raise ValueError("Conversion only valid for non-negative numbers")
  227. buf_len = (_gmp.mpz_sizeinbase(self._mpz_p, 2) + 7) // 8
  228. if buf_len > block_size > 0:
  229. raise ValueError("Number is too big to convert to byte string"
  230. "of prescribed length")
  231. buf = create_string_buffer(buf_len)
  232. _gmp.mpz_export(
  233. buf,
  234. null_pointer, # Ignore countp
  235. 1, # Big endian
  236. c_size_t(1), # Each word is 1 byte long
  237. 0, # Endianess within a word - not relevant
  238. c_size_t(0), # No nails
  239. self._mpz_p)
  240. return bchr(0) * max(0, block_size - buf_len) + get_raw_buffer(buf)
  241. @staticmethod
  242. def from_bytes(byte_string):
  243. """Convert a byte string into a number.
  244. :Parameters:
  245. byte_string : byte string
  246. The input number, encoded in network order.
  247. It can only be non-negative.
  248. :Return:
  249. The ``Integer`` object carrying the same value as the input.
  250. """
  251. result = Integer(0)
  252. _gmp.mpz_import(
  253. result._mpz_p,
  254. c_size_t(len(byte_string)), # Amount of words to read
  255. 1, # Big endian
  256. c_size_t(1), # Each word is 1 byte long
  257. 0, # Endianess within a word - not relevant
  258. c_size_t(0), # No nails
  259. byte_string)
  260. return result
  261. # Relations
  262. def _apply_and_return(self, func, term):
  263. if not isinstance(term, Integer):
  264. term = Integer(term)
  265. return func(self._mpz_p, term._mpz_p)
  266. def __eq__(self, term):
  267. if not isinstance(term, (Integer, int, long)):
  268. return False
  269. return self._apply_and_return(_gmp.mpz_cmp, term) == 0
  270. def __ne__(self, term):
  271. if not isinstance(term, (Integer, int, long)):
  272. return True
  273. return self._apply_and_return(_gmp.mpz_cmp, term) != 0
  274. def __lt__(self, term):
  275. return self._apply_and_return(_gmp.mpz_cmp, term) < 0
  276. def __le__(self, term):
  277. return self._apply_and_return(_gmp.mpz_cmp, term) <= 0
  278. def __gt__(self, term):
  279. return self._apply_and_return(_gmp.mpz_cmp, term) > 0
  280. def __ge__(self, term):
  281. return self._apply_and_return(_gmp.mpz_cmp, term) >= 0
  282. def __nonzero__(self):
  283. return _gmp.mpz_cmp(self._mpz_p, self._zero_mpz_p) != 0
  284. def is_negative(self):
  285. return _gmp.mpz_cmp(self._mpz_p, self._zero_mpz_p) < 0
  286. # Arithmetic operations
  287. def __add__(self, term):
  288. result = Integer(0)
  289. if not isinstance(term, Integer):
  290. term = Integer(term)
  291. _gmp.mpz_add(result._mpz_p,
  292. self._mpz_p,
  293. term._mpz_p)
  294. return result
  295. def __sub__(self, term):
  296. result = Integer(0)
  297. if not isinstance(term, Integer):
  298. term = Integer(term)
  299. _gmp.mpz_sub(result._mpz_p,
  300. self._mpz_p,
  301. term._mpz_p)
  302. return result
  303. def __mul__(self, term):
  304. result = Integer(0)
  305. if not isinstance(term, Integer):
  306. term = Integer(term)
  307. _gmp.mpz_mul(result._mpz_p,
  308. self._mpz_p,
  309. term._mpz_p)
  310. return result
  311. def __floordiv__(self, divisor):
  312. if not isinstance(divisor, Integer):
  313. divisor = Integer(divisor)
  314. if _gmp.mpz_cmp(divisor._mpz_p,
  315. self._zero_mpz_p) == 0:
  316. raise ZeroDivisionError("Division by zero")
  317. result = Integer(0)
  318. _gmp.mpz_fdiv_q(result._mpz_p,
  319. self._mpz_p,
  320. divisor._mpz_p)
  321. return result
  322. def __mod__(self, divisor):
  323. if not isinstance(divisor, Integer):
  324. divisor = Integer(divisor)
  325. comp = _gmp.mpz_cmp(divisor._mpz_p,
  326. self._zero_mpz_p)
  327. if comp == 0:
  328. raise ZeroDivisionError("Division by zero")
  329. if comp < 0:
  330. raise ValueError("Modulus must be positive")
  331. result = Integer(0)
  332. _gmp.mpz_mod(result._mpz_p,
  333. self._mpz_p,
  334. divisor._mpz_p)
  335. return result
  336. def inplace_pow(self, exponent, modulus=None):
  337. if modulus is None:
  338. if exponent < 0:
  339. raise ValueError("Exponent must not be negative")
  340. # Normal exponentiation
  341. if exponent > 256:
  342. raise ValueError("Exponent is too big")
  343. _gmp.mpz_pow_ui(self._mpz_p,
  344. self._mpz_p, # Base
  345. c_ulong(int(exponent))
  346. )
  347. else:
  348. # Modular exponentiation
  349. if not isinstance(modulus, Integer):
  350. modulus = Integer(modulus)
  351. if not modulus:
  352. raise ZeroDivisionError("Division by zero")
  353. if modulus.is_negative():
  354. raise ValueError("Modulus must be positive")
  355. if isinstance(exponent, (int, long)):
  356. if exponent < 0:
  357. raise ValueError("Exponent must not be negative")
  358. if exponent < 65536:
  359. _gmp.mpz_powm_ui(self._mpz_p,
  360. self._mpz_p,
  361. c_ulong(exponent),
  362. modulus._mpz_p)
  363. return self
  364. exponent = Integer(exponent)
  365. elif exponent.is_negative():
  366. raise ValueError("Exponent must not be negative")
  367. _gmp.mpz_powm(self._mpz_p,
  368. self._mpz_p,
  369. exponent._mpz_p,
  370. modulus._mpz_p)
  371. return self
  372. def __pow__(self, exponent, modulus=None):
  373. result = Integer(self)
  374. return result.inplace_pow(exponent, modulus)
  375. def __abs__(self):
  376. result = Integer(0)
  377. _gmp.mpz_abs(result._mpz_p, self._mpz_p)
  378. return result
  379. def sqrt(self):
  380. """Return the largest Integer that does not
  381. exceed the square root"""
  382. if self < 0:
  383. raise ValueError("Square root of negative value")
  384. result = Integer(0)
  385. _gmp.mpz_sqrt(result._mpz_p,
  386. self._mpz_p)
  387. return result
  388. def __iadd__(self, term):
  389. if isinstance(term, (int, long)):
  390. if 0 <= term < 65536:
  391. _gmp.mpz_add_ui(self._mpz_p,
  392. self._mpz_p,
  393. c_ulong(term))
  394. return self
  395. if -65535 < term < 0:
  396. _gmp.mpz_sub_ui(self._mpz_p,
  397. self._mpz_p,
  398. c_ulong(-term))
  399. return self
  400. term = Integer(term)
  401. _gmp.mpz_add(self._mpz_p,
  402. self._mpz_p,
  403. term._mpz_p)
  404. return self
  405. def __isub__(self, term):
  406. if isinstance(term, (int, long)):
  407. if 0 <= term < 65536:
  408. _gmp.mpz_sub_ui(self._mpz_p,
  409. self._mpz_p,
  410. c_ulong(term))
  411. return self
  412. if -65535 < term < 0:
  413. _gmp.mpz_add_ui(self._mpz_p,
  414. self._mpz_p,
  415. c_ulong(-term))
  416. return self
  417. term = Integer(term)
  418. _gmp.mpz_sub(self._mpz_p,
  419. self._mpz_p,
  420. term._mpz_p)
  421. return self
  422. def __imul__(self, term):
  423. if isinstance(term, (int, long)):
  424. if 0 <= term < 65536:
  425. _gmp.mpz_mul_ui(self._mpz_p,
  426. self._mpz_p,
  427. c_ulong(term))
  428. return self
  429. if -65535 < term < 0:
  430. _gmp.mpz_mul_ui(self._mpz_p,
  431. self._mpz_p,
  432. c_ulong(-term))
  433. _gmp.mpz_neg(self._mpz_p, self._mpz_p)
  434. return self
  435. term = Integer(term)
  436. _gmp.mpz_mul(self._mpz_p,
  437. self._mpz_p,
  438. term._mpz_p)
  439. return self
  440. def __imod__(self, divisor):
  441. if not isinstance(divisor, Integer):
  442. divisor = Integer(divisor)
  443. comp = _gmp.mpz_cmp(divisor._mpz_p,
  444. divisor._zero_mpz_p)
  445. if comp == 0:
  446. raise ZeroDivisionError("Division by zero")
  447. if comp < 0:
  448. raise ValueError("Modulus must be positive")
  449. _gmp.mpz_mod(self._mpz_p,
  450. self._mpz_p,
  451. divisor._mpz_p)
  452. return self
  453. # Boolean/bit operations
  454. def __and__(self, term):
  455. result = Integer(0)
  456. if not isinstance(term, Integer):
  457. term = Integer(term)
  458. _gmp.mpz_and(result._mpz_p,
  459. self._mpz_p,
  460. term._mpz_p)
  461. return result
  462. def __or__(self, term):
  463. result = Integer(0)
  464. if not isinstance(term, Integer):
  465. term = Integer(term)
  466. _gmp.mpz_ior(result._mpz_p,
  467. self._mpz_p,
  468. term._mpz_p)
  469. return result
  470. def __rshift__(self, pos):
  471. result = Integer(0)
  472. if not 0 <= pos < 65536:
  473. raise ValueError("Incorrect shift count")
  474. _gmp.mpz_tdiv_q_2exp(result._mpz_p,
  475. self._mpz_p,
  476. c_ulong(int(pos)))
  477. return result
  478. def __irshift__(self, pos):
  479. if not 0 <= pos < 65536:
  480. raise ValueError("Incorrect shift count")
  481. _gmp.mpz_tdiv_q_2exp(self._mpz_p,
  482. self._mpz_p,
  483. c_ulong(int(pos)))
  484. return self
  485. def __lshift__(self, pos):
  486. result = Integer(0)
  487. if not 0 <= pos < 65536:
  488. raise ValueError("Incorrect shift count")
  489. _gmp.mpz_mul_2exp(result._mpz_p,
  490. self._mpz_p,
  491. c_ulong(int(pos)))
  492. return result
  493. def __ilshift__(self, pos):
  494. if not 0 <= pos < 65536:
  495. raise ValueError("Incorrect shift count")
  496. _gmp.mpz_mul_2exp(self._mpz_p,
  497. self._mpz_p,
  498. c_ulong(int(pos)))
  499. return self
  500. def get_bit(self, n):
  501. """Return True if the n-th bit is set to 1.
  502. Bit 0 is the least significant."""
  503. if not 0 <= n < 65536:
  504. raise ValueError("Incorrect bit position")
  505. return bool(_gmp.mpz_tstbit(self._mpz_p,
  506. c_ulong(int(n))))
  507. # Extra
  508. def is_odd(self):
  509. return _gmp.mpz_tstbit(self._mpz_p, 0) == 1
  510. def is_even(self):
  511. return _gmp.mpz_tstbit(self._mpz_p, 0) == 0
  512. def size_in_bits(self):
  513. """Return the minimum number of bits that can encode the number."""
  514. if self < 0:
  515. raise ValueError("Conversion only valid for non-negative numbers")
  516. return _gmp.mpz_sizeinbase(self._mpz_p, 2)
  517. def size_in_bytes(self):
  518. """Return the minimum number of bytes that can encode the number."""
  519. return (self.size_in_bits() - 1) // 8 + 1
  520. def is_perfect_square(self):
  521. return _gmp.mpz_perfect_square_p(self._mpz_p) != 0
  522. def fail_if_divisible_by(self, small_prime):
  523. """Raise an exception if the small prime is a divisor."""
  524. if isinstance(small_prime, (int, long)):
  525. if 0 < small_prime < 65536:
  526. if _gmp.mpz_divisible_ui_p(self._mpz_p,
  527. c_ulong(small_prime)):
  528. raise ValueError("The value is composite")
  529. return
  530. small_prime = Integer(small_prime)
  531. if _gmp.mpz_divisible_p(self._mpz_p,
  532. small_prime._mpz_p):
  533. raise ValueError("The value is composite")
  534. def multiply_accumulate(self, a, b):
  535. """Increment the number by the product of a and b."""
  536. if not isinstance(a, Integer):
  537. a = Integer(a)
  538. if isinstance(b, (int, long)):
  539. if 0 < b < 65536:
  540. _gmp.mpz_addmul_ui(self._mpz_p,
  541. a._mpz_p,
  542. c_ulong(b))
  543. return self
  544. if -65535 < b < 0:
  545. _gmp.mpz_submul_ui(self._mpz_p,
  546. a._mpz_p,
  547. c_ulong(-b))
  548. return self
  549. b = Integer(b)
  550. _gmp.mpz_addmul(self._mpz_p,
  551. a._mpz_p,
  552. b._mpz_p)
  553. return self
  554. def set(self, source):
  555. """Set the Integer to have the given value"""
  556. if not isinstance(source, Integer):
  557. source = Integer(source)
  558. _gmp.mpz_set(self._mpz_p,
  559. source._mpz_p)
  560. return self
  561. def inplace_inverse(self, modulus):
  562. """Compute the inverse of this number in the ring of
  563. modulo integers.
  564. Raise an exception if no inverse exists.
  565. """
  566. if not isinstance(modulus, Integer):
  567. modulus = Integer(modulus)
  568. comp = _gmp.mpz_cmp(modulus._mpz_p,
  569. self._zero_mpz_p)
  570. if comp == 0:
  571. raise ZeroDivisionError("Modulus cannot be zero")
  572. if comp < 0:
  573. raise ValueError("Modulus must be positive")
  574. result = _gmp.mpz_invert(self._mpz_p,
  575. self._mpz_p,
  576. modulus._mpz_p)
  577. if not result:
  578. raise ValueError("No inverse value can be computed")
  579. return self
  580. def inverse(self, modulus):
  581. result = Integer(self)
  582. result.inplace_inverse(modulus)
  583. return result
  584. def gcd(self, term):
  585. """Compute the greatest common denominator between this
  586. number and another term."""
  587. result = Integer(0)
  588. if isinstance(term, (int, long)):
  589. if 0 < term < 65535:
  590. _gmp.mpz_gcd_ui(result._mpz_p,
  591. self._mpz_p,
  592. c_ulong(term))
  593. return result
  594. term = Integer(term)
  595. _gmp.mpz_gcd(result._mpz_p, self._mpz_p, term._mpz_p)
  596. return result
  597. def lcm(self, term):
  598. """Compute the least common multiplier between this
  599. number and another term."""
  600. result = Integer(0)
  601. if not isinstance(term, Integer):
  602. term = Integer(term)
  603. _gmp.mpz_lcm(result._mpz_p, self._mpz_p, term._mpz_p)
  604. return result
  605. @staticmethod
  606. def jacobi_symbol(a, n):
  607. """Compute the Jacobi symbol"""
  608. if not isinstance(a, Integer):
  609. a = Integer(a)
  610. if not isinstance(n, Integer):
  611. n = Integer(n)
  612. if n <= 0 or n.is_even():
  613. raise ValueError("n must be positive even for the Jacobi symbol")
  614. return _gmp.mpz_jacobi(a._mpz_p, n._mpz_p)
  615. # Clean-up
  616. def __del__(self):
  617. try:
  618. if self._mpz_p is not None:
  619. if self._initialized:
  620. _gmp.mpz_clear(self._mpz_p)
  621. self._mpz_p = None
  622. except AttributeError:
  623. pass