test_qs.py 3.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110
  1. from __future__ import annotations
  2. import math
  3. from sympy.core.random import _randint
  4. from sympy.ntheory import qs, qs_factor
  5. from sympy.ntheory.qs import SievePolynomial, _generate_factor_base, \
  6. _generate_polynomial, \
  7. _gen_sieve_array, _check_smoothness, _trial_division_stage, _find_factor
  8. from sympy.testing.pytest import slow
  9. @slow
  10. def test_qs_1():
  11. assert qs(10009202107, 100, 10000) == {100043, 100049}
  12. assert qs(211107295182713951054568361, 1000, 10000) == \
  13. {13791315212531, 15307263442931}
  14. assert qs(980835832582657*990377764891511, 2000, 10000) == \
  15. {980835832582657, 990377764891511}
  16. assert qs(18640889198609*20991129234731, 1000, 50000) == \
  17. {18640889198609, 20991129234731}
  18. def test_qs_2() -> None:
  19. n = 10009202107
  20. M = 50
  21. sieve_poly = SievePolynomial(10, 80, n)
  22. assert sieve_poly.eval_v(10) == sieve_poly.eval_u(10)**2 - n == -10009169707
  23. assert sieve_poly.eval_v(5) == sieve_poly.eval_u(5)**2 - n == -10009185207
  24. idx_1000, idx_5000, factor_base = _generate_factor_base(2000, n)
  25. assert idx_1000 == 82
  26. assert [factor_base[i].prime for i in range(15)] == \
  27. [2, 3, 7, 11, 17, 19, 29, 31, 43, 59, 61, 67, 71, 73, 79]
  28. assert [factor_base[i].tmem_p for i in range(15)] == \
  29. [1, 1, 3, 5, 3, 6, 6, 14, 1, 16, 24, 22, 18, 22, 15]
  30. assert [factor_base[i].log_p for i in range(5)] == \
  31. [710, 1125, 1993, 2455, 2901]
  32. it = _generate_polynomial(
  33. n, M, factor_base, idx_1000, idx_5000, _randint(0))
  34. g = next(it)
  35. assert g.a == 1133107
  36. assert g.b == 682543
  37. assert [factor_base[i].soln1 for i in range(15)] == \
  38. [0, 0, 3, 7, 13, 0, 8, 19, 9, 43, 27, 25, 63, 29, 19]
  39. assert [factor_base[i].soln2 for i in range(15)] == \
  40. [0, 1, 1, 3, 12, 16, 15, 6, 15, 1, 56, 55, 61, 58, 16]
  41. assert [factor_base[i].b_ainv for i in range(5)] == \
  42. [[0, 0], [0, 2], [3, 0], [3, 9], [13, 13]]
  43. g_1 = next(it)
  44. assert g_1.a == 1133107
  45. assert g_1.b == 136765
  46. sieve_array = _gen_sieve_array(M, factor_base)
  47. assert sieve_array[0:5] == [8424, 13603, 1835, 5335, 710]
  48. assert _check_smoothness(9645, factor_base) == (36028797018963972, 5)
  49. assert _check_smoothness(210313, factor_base) == (20992, 1)
  50. partial_relations: dict[int, tuple[int, int]] = {}
  51. smooth_relation, proper_factor = _trial_division_stage(
  52. n, M, factor_base, sieve_array, sieve_poly, partial_relations,
  53. ERROR_TERM=25*2**10)
  54. assert partial_relations == {
  55. 8699: (440, -10009008507, 75557863761098695507973),
  56. 166741: (490, -10008962007, 524341),
  57. 131449: (530, -10008921207, 664613997892457936451903530140172325),
  58. 6653: (550, -10008899607, 19342813113834066795307021)
  59. }
  60. assert [smooth_relation[i][0] for i in range(5)] == [
  61. -250, 1064469, 72819, 231957, 44167]
  62. assert [smooth_relation[i][1] for i in range(5)] == [
  63. -10009139607, 1133094251961, 5302606761, 53804049849, 1950723889]
  64. assert smooth_relation[0][2] == 89213869829863962596973701078031812362502145
  65. assert proper_factor == set()
  66. def test_qs_3():
  67. N = 1817
  68. smooth_relations = [
  69. (2455024, 637, 8),
  70. (-27993000, 81536, 10),
  71. (11461840, 12544, 0),
  72. (149, 20384, 10),
  73. (-31138074, 19208, 2)
  74. ]
  75. assert next(_find_factor(N, smooth_relations, 4)) == 23
  76. def test_qs_4():
  77. N = 10007**2 * 10009 * 10037**3 * 10039
  78. for factor in qs(N, 1000, 2000):
  79. assert N % factor == 0
  80. N //= factor
  81. def test_qs_factor():
  82. assert qs_factor(1009 * 100003, 2000, 10000) == {1009: 1, 100003: 1}
  83. n = 1009**2 * 2003**2*30011*400009
  84. factors = qs_factor(n, 2000, 10000)
  85. assert len(factors) > 1
  86. assert math.prod(p**e for p, e in factors.items()) == n
  87. def test_issue_27616():
  88. #https://github.com/sympy/sympy/issues/27616
  89. N = 9804659461513846513 + 1
  90. assert qs(N, 5000, 20000) is not None