_IntegerGMP.py 25 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689
  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 Cryptodome.Util.py3compat import tobytes, is_native_int
  32. from Cryptodome.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)
  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. int __gmp_sscanf (const char *s, const char *fmt, ...);
  45. void __gmpz_set (mpz_t rop, const mpz_t op);
  46. int __gmp_snprintf (uint8_t *buf, size_t size, const char *fmt, ...);
  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. lib = load_lib("gmp", gmp_defs)
  90. implementation = {"library": "gmp", "api": backend}
  91. if hasattr(lib, "__mpir_version"):
  92. raise ImportError("MPIR library detected")
  93. if sys.platform == "win32":
  94. raise ImportError("Not using GMP on Windows")
  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 Cryptodome.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. self._initialized = True
  135. if is_native_int(value):
  136. _gmp.mpz_init(self._mpz_p)
  137. result = _gmp.gmp_sscanf(tobytes(str(value)), b"%Zd", self._mpz_p)
  138. if result != 1:
  139. raise ValueError("Error converting '%d'" % value)
  140. else:
  141. _gmp.mpz_init_set(self._mpz_p, value._mpz_p)
  142. # Conversions
  143. def __int__(self):
  144. # buf will contain the integer encoded in decimal plus the trailing
  145. # zero, and possibly the negative sign.
  146. # dig10(x) < log10(x) + 1 = log2(x)/log2(10) + 1 < log2(x)/3 + 1
  147. buf_len = _gmp.mpz_sizeinbase(self._mpz_p, 2) // 3 + 3
  148. buf = create_string_buffer(buf_len)
  149. _gmp.gmp_snprintf(buf, c_size_t(buf_len), b"%Zd", self._mpz_p)
  150. return int(get_c_string(buf))
  151. def __str__(self):
  152. return str(int(self))
  153. def __repr__(self):
  154. return "Integer(%s)" % str(self)
  155. def to_bytes(self, block_size=0):
  156. """Convert the number into a byte string.
  157. This method encodes the number in network order and prepends
  158. as many zero bytes as required. It only works for non-negative
  159. values.
  160. :Parameters:
  161. block_size : integer
  162. The exact size the output byte string must have.
  163. If zero, the string has the minimal length.
  164. :Returns:
  165. A byte string.
  166. :Raise ValueError:
  167. If the value is negative or if ``block_size`` is
  168. provided and the length of the byte string would exceed it.
  169. """
  170. if self < 0:
  171. raise ValueError("Conversion only valid for non-negative numbers")
  172. buf_len = (_gmp.mpz_sizeinbase(self._mpz_p, 2) + 7) // 8
  173. if buf_len > block_size > 0:
  174. raise ValueError("Number is too big to convert to byte string"
  175. "of prescribed length")
  176. buf = create_string_buffer(buf_len)
  177. _gmp.mpz_export(
  178. buf,
  179. null_pointer, # Ignore countp
  180. 1, # Big endian
  181. c_size_t(1), # Each word is 1 byte long
  182. 0, # Endianess within a word - not relevant
  183. c_size_t(0), # No nails
  184. self._mpz_p)
  185. return b'\x00' * max(0, block_size - buf_len) + get_raw_buffer(buf)
  186. @staticmethod
  187. def from_bytes(byte_string):
  188. """Convert a byte string into a number.
  189. :Parameters:
  190. byte_string : byte string
  191. The input number, encoded in network order.
  192. It can only be non-negative.
  193. :Return:
  194. The ``Integer`` object carrying the same value as the input.
  195. """
  196. result = IntegerGMP(0)
  197. _gmp.mpz_import(
  198. result._mpz_p,
  199. c_size_t(len(byte_string)), # Amount of words to read
  200. 1, # Big endian
  201. c_size_t(1), # Each word is 1 byte long
  202. 0, # Endianess within a word - not relevant
  203. c_size_t(0), # No nails
  204. byte_string)
  205. return result
  206. # Relations
  207. def _apply_and_return(self, func, term):
  208. if not isinstance(term, IntegerGMP):
  209. term = IntegerGMP(term)
  210. return func(self._mpz_p, term._mpz_p)
  211. def __eq__(self, term):
  212. if not (isinstance(term, IntegerGMP) or is_native_int(term)):
  213. return False
  214. return self._apply_and_return(_gmp.mpz_cmp, term) == 0
  215. def __ne__(self, term):
  216. if not (isinstance(term, IntegerGMP) or is_native_int(term)):
  217. return True
  218. return self._apply_and_return(_gmp.mpz_cmp, term) != 0
  219. def __lt__(self, term):
  220. return self._apply_and_return(_gmp.mpz_cmp, term) < 0
  221. def __le__(self, term):
  222. return self._apply_and_return(_gmp.mpz_cmp, term) <= 0
  223. def __gt__(self, term):
  224. return self._apply_and_return(_gmp.mpz_cmp, term) > 0
  225. def __ge__(self, term):
  226. return self._apply_and_return(_gmp.mpz_cmp, term) >= 0
  227. def __nonzero__(self):
  228. return _gmp.mpz_cmp(self._mpz_p, self._zero_mpz_p) != 0
  229. __bool__ = __nonzero__
  230. def is_negative(self):
  231. return _gmp.mpz_cmp(self._mpz_p, self._zero_mpz_p) < 0
  232. # Arithmetic operations
  233. def __add__(self, term):
  234. result = IntegerGMP(0)
  235. if not isinstance(term, IntegerGMP):
  236. term = IntegerGMP(term)
  237. _gmp.mpz_add(result._mpz_p,
  238. self._mpz_p,
  239. term._mpz_p)
  240. return result
  241. def __sub__(self, term):
  242. result = IntegerGMP(0)
  243. if not isinstance(term, IntegerGMP):
  244. term = IntegerGMP(term)
  245. _gmp.mpz_sub(result._mpz_p,
  246. self._mpz_p,
  247. term._mpz_p)
  248. return result
  249. def __mul__(self, term):
  250. result = IntegerGMP(0)
  251. if not isinstance(term, IntegerGMP):
  252. term = IntegerGMP(term)
  253. _gmp.mpz_mul(result._mpz_p,
  254. self._mpz_p,
  255. term._mpz_p)
  256. return result
  257. def __floordiv__(self, divisor):
  258. if not isinstance(divisor, IntegerGMP):
  259. divisor = IntegerGMP(divisor)
  260. if _gmp.mpz_cmp(divisor._mpz_p,
  261. self._zero_mpz_p) == 0:
  262. raise ZeroDivisionError("Division by zero")
  263. result = IntegerGMP(0)
  264. _gmp.mpz_fdiv_q(result._mpz_p,
  265. self._mpz_p,
  266. divisor._mpz_p)
  267. return result
  268. def __mod__(self, divisor):
  269. if not isinstance(divisor, IntegerGMP):
  270. divisor = IntegerGMP(divisor)
  271. comp = _gmp.mpz_cmp(divisor._mpz_p,
  272. self._zero_mpz_p)
  273. if comp == 0:
  274. raise ZeroDivisionError("Division by zero")
  275. if comp < 0:
  276. raise ValueError("Modulus must be positive")
  277. result = IntegerGMP(0)
  278. _gmp.mpz_mod(result._mpz_p,
  279. self._mpz_p,
  280. divisor._mpz_p)
  281. return result
  282. def inplace_pow(self, exponent, modulus=None):
  283. if modulus is None:
  284. if exponent < 0:
  285. raise ValueError("Exponent must not be negative")
  286. # Normal exponentiation
  287. if exponent > 256:
  288. raise ValueError("Exponent is too big")
  289. _gmp.mpz_pow_ui(self._mpz_p,
  290. self._mpz_p, # Base
  291. c_ulong(int(exponent))
  292. )
  293. else:
  294. # Modular exponentiation
  295. if not isinstance(modulus, IntegerGMP):
  296. modulus = IntegerGMP(modulus)
  297. if not modulus:
  298. raise ZeroDivisionError("Division by zero")
  299. if modulus.is_negative():
  300. raise ValueError("Modulus must be positive")
  301. if is_native_int(exponent):
  302. if exponent < 0:
  303. raise ValueError("Exponent must not be negative")
  304. if exponent < 65536:
  305. _gmp.mpz_powm_ui(self._mpz_p,
  306. self._mpz_p,
  307. c_ulong(exponent),
  308. modulus._mpz_p)
  309. return self
  310. exponent = IntegerGMP(exponent)
  311. elif exponent.is_negative():
  312. raise ValueError("Exponent must not be negative")
  313. _gmp.mpz_powm(self._mpz_p,
  314. self._mpz_p,
  315. exponent._mpz_p,
  316. modulus._mpz_p)
  317. return self
  318. def __pow__(self, exponent, modulus=None):
  319. result = IntegerGMP(self)
  320. return result.inplace_pow(exponent, modulus)
  321. def __abs__(self):
  322. result = IntegerGMP(0)
  323. _gmp.mpz_abs(result._mpz_p, self._mpz_p)
  324. return result
  325. def sqrt(self, modulus=None):
  326. """Return the largest Integer that does not
  327. exceed the square root"""
  328. if modulus is None:
  329. if self < 0:
  330. raise ValueError("Square root of negative value")
  331. result = IntegerGMP(0)
  332. _gmp.mpz_sqrt(result._mpz_p,
  333. self._mpz_p)
  334. else:
  335. if modulus <= 0:
  336. raise ValueError("Modulus must be positive")
  337. modulus = int(modulus)
  338. result = IntegerGMP(self._tonelli_shanks(int(self) % modulus, modulus))
  339. return result
  340. def __iadd__(self, term):
  341. if is_native_int(term):
  342. if 0 <= term < 65536:
  343. _gmp.mpz_add_ui(self._mpz_p,
  344. self._mpz_p,
  345. c_ulong(term))
  346. return self
  347. if -65535 < term < 0:
  348. _gmp.mpz_sub_ui(self._mpz_p,
  349. self._mpz_p,
  350. c_ulong(-term))
  351. return self
  352. term = IntegerGMP(term)
  353. _gmp.mpz_add(self._mpz_p,
  354. self._mpz_p,
  355. term._mpz_p)
  356. return self
  357. def __isub__(self, term):
  358. if is_native_int(term):
  359. if 0 <= term < 65536:
  360. _gmp.mpz_sub_ui(self._mpz_p,
  361. self._mpz_p,
  362. c_ulong(term))
  363. return self
  364. if -65535 < term < 0:
  365. _gmp.mpz_add_ui(self._mpz_p,
  366. self._mpz_p,
  367. c_ulong(-term))
  368. return self
  369. term = IntegerGMP(term)
  370. _gmp.mpz_sub(self._mpz_p,
  371. self._mpz_p,
  372. term._mpz_p)
  373. return self
  374. def __imul__(self, term):
  375. if is_native_int(term):
  376. if 0 <= term < 65536:
  377. _gmp.mpz_mul_ui(self._mpz_p,
  378. self._mpz_p,
  379. c_ulong(term))
  380. return self
  381. if -65535 < term < 0:
  382. _gmp.mpz_mul_ui(self._mpz_p,
  383. self._mpz_p,
  384. c_ulong(-term))
  385. _gmp.mpz_neg(self._mpz_p, self._mpz_p)
  386. return self
  387. term = IntegerGMP(term)
  388. _gmp.mpz_mul(self._mpz_p,
  389. self._mpz_p,
  390. term._mpz_p)
  391. return self
  392. def __imod__(self, divisor):
  393. if not isinstance(divisor, IntegerGMP):
  394. divisor = IntegerGMP(divisor)
  395. comp = _gmp.mpz_cmp(divisor._mpz_p,
  396. divisor._zero_mpz_p)
  397. if comp == 0:
  398. raise ZeroDivisionError("Division by zero")
  399. if comp < 0:
  400. raise ValueError("Modulus must be positive")
  401. _gmp.mpz_mod(self._mpz_p,
  402. self._mpz_p,
  403. divisor._mpz_p)
  404. return self
  405. # Boolean/bit operations
  406. def __and__(self, term):
  407. result = IntegerGMP(0)
  408. if not isinstance(term, IntegerGMP):
  409. term = IntegerGMP(term)
  410. _gmp.mpz_and(result._mpz_p,
  411. self._mpz_p,
  412. term._mpz_p)
  413. return result
  414. def __or__(self, term):
  415. result = IntegerGMP(0)
  416. if not isinstance(term, IntegerGMP):
  417. term = IntegerGMP(term)
  418. _gmp.mpz_ior(result._mpz_p,
  419. self._mpz_p,
  420. term._mpz_p)
  421. return result
  422. def __rshift__(self, pos):
  423. result = IntegerGMP(0)
  424. if pos < 0:
  425. raise ValueError("negative shift count")
  426. if pos > 65536:
  427. if self < 0:
  428. return -1
  429. else:
  430. return 0
  431. _gmp.mpz_tdiv_q_2exp(result._mpz_p,
  432. self._mpz_p,
  433. c_ulong(int(pos)))
  434. return result
  435. def __irshift__(self, pos):
  436. if pos < 0:
  437. raise ValueError("negative shift count")
  438. if pos > 65536:
  439. if self < 0:
  440. return -1
  441. else:
  442. return 0
  443. _gmp.mpz_tdiv_q_2exp(self._mpz_p,
  444. self._mpz_p,
  445. c_ulong(int(pos)))
  446. return self
  447. def __lshift__(self, pos):
  448. result = IntegerGMP(0)
  449. if not 0 <= pos < 65536:
  450. raise ValueError("Incorrect shift count")
  451. _gmp.mpz_mul_2exp(result._mpz_p,
  452. self._mpz_p,
  453. c_ulong(int(pos)))
  454. return result
  455. def __ilshift__(self, pos):
  456. if not 0 <= pos < 65536:
  457. raise ValueError("Incorrect shift count")
  458. _gmp.mpz_mul_2exp(self._mpz_p,
  459. self._mpz_p,
  460. c_ulong(int(pos)))
  461. return self
  462. def get_bit(self, n):
  463. """Return True if the n-th bit is set to 1.
  464. Bit 0 is the least significant."""
  465. if self < 0:
  466. raise ValueError("no bit representation for negative values")
  467. if n < 0:
  468. raise ValueError("negative bit count")
  469. if n > 65536:
  470. return 0
  471. return bool(_gmp.mpz_tstbit(self._mpz_p,
  472. c_ulong(int(n))))
  473. # Extra
  474. def is_odd(self):
  475. return _gmp.mpz_tstbit(self._mpz_p, 0) == 1
  476. def is_even(self):
  477. return _gmp.mpz_tstbit(self._mpz_p, 0) == 0
  478. def size_in_bits(self):
  479. """Return the minimum number of bits that can encode the number."""
  480. if self < 0:
  481. raise ValueError("Conversion only valid for non-negative numbers")
  482. return _gmp.mpz_sizeinbase(self._mpz_p, 2)
  483. def size_in_bytes(self):
  484. """Return the minimum number of bytes that can encode the number."""
  485. return (self.size_in_bits() - 1) // 8 + 1
  486. def is_perfect_square(self):
  487. return _gmp.mpz_perfect_square_p(self._mpz_p) != 0
  488. def fail_if_divisible_by(self, small_prime):
  489. """Raise an exception if the small prime is a divisor."""
  490. if is_native_int(small_prime):
  491. if 0 < small_prime < 65536:
  492. if _gmp.mpz_divisible_ui_p(self._mpz_p,
  493. c_ulong(small_prime)):
  494. raise ValueError("The value is composite")
  495. return
  496. small_prime = IntegerGMP(small_prime)
  497. if _gmp.mpz_divisible_p(self._mpz_p,
  498. small_prime._mpz_p):
  499. raise ValueError("The value is composite")
  500. def multiply_accumulate(self, a, b):
  501. """Increment the number by the product of a and b."""
  502. if not isinstance(a, IntegerGMP):
  503. a = IntegerGMP(a)
  504. if is_native_int(b):
  505. if 0 < b < 65536:
  506. _gmp.mpz_addmul_ui(self._mpz_p,
  507. a._mpz_p,
  508. c_ulong(b))
  509. return self
  510. if -65535 < b < 0:
  511. _gmp.mpz_submul_ui(self._mpz_p,
  512. a._mpz_p,
  513. c_ulong(-b))
  514. return self
  515. b = IntegerGMP(b)
  516. _gmp.mpz_addmul(self._mpz_p,
  517. a._mpz_p,
  518. b._mpz_p)
  519. return self
  520. def set(self, source):
  521. """Set the Integer to have the given value"""
  522. if not isinstance(source, IntegerGMP):
  523. source = IntegerGMP(source)
  524. _gmp.mpz_set(self._mpz_p,
  525. source._mpz_p)
  526. return self
  527. def inplace_inverse(self, modulus):
  528. """Compute the inverse of this number in the ring of
  529. modulo integers.
  530. Raise an exception if no inverse exists.
  531. """
  532. if not isinstance(modulus, IntegerGMP):
  533. modulus = IntegerGMP(modulus)
  534. comp = _gmp.mpz_cmp(modulus._mpz_p,
  535. self._zero_mpz_p)
  536. if comp == 0:
  537. raise ZeroDivisionError("Modulus cannot be zero")
  538. if comp < 0:
  539. raise ValueError("Modulus must be positive")
  540. result = _gmp.mpz_invert(self._mpz_p,
  541. self._mpz_p,
  542. modulus._mpz_p)
  543. if not result:
  544. raise ValueError("No inverse value can be computed")
  545. return self
  546. def inverse(self, modulus):
  547. result = IntegerGMP(self)
  548. result.inplace_inverse(modulus)
  549. return result
  550. def gcd(self, term):
  551. """Compute the greatest common denominator between this
  552. number and another term."""
  553. result = IntegerGMP(0)
  554. if is_native_int(term):
  555. if 0 < term < 65535:
  556. _gmp.mpz_gcd_ui(result._mpz_p,
  557. self._mpz_p,
  558. c_ulong(term))
  559. return result
  560. term = IntegerGMP(term)
  561. _gmp.mpz_gcd(result._mpz_p, self._mpz_p, term._mpz_p)
  562. return result
  563. def lcm(self, term):
  564. """Compute the least common multiplier between this
  565. number and another term."""
  566. result = IntegerGMP(0)
  567. if not isinstance(term, IntegerGMP):
  568. term = IntegerGMP(term)
  569. _gmp.mpz_lcm(result._mpz_p, self._mpz_p, term._mpz_p)
  570. return result
  571. @staticmethod
  572. def jacobi_symbol(a, n):
  573. """Compute the Jacobi symbol"""
  574. if not isinstance(a, IntegerGMP):
  575. a = IntegerGMP(a)
  576. if not isinstance(n, IntegerGMP):
  577. n = IntegerGMP(n)
  578. if n <= 0 or n.is_even():
  579. raise ValueError("n must be positive even for the Jacobi symbol")
  580. return _gmp.mpz_jacobi(a._mpz_p, n._mpz_p)
  581. # Clean-up
  582. def __del__(self):
  583. try:
  584. if self._mpz_p is not None:
  585. if self._initialized:
  586. _gmp.mpz_clear(self._mpz_p)
  587. self._mpz_p = None
  588. except AttributeError:
  589. pass