From d5bd34e2df1f75c7c857ee1423bd659ae0586ea6 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Tue, 9 Jul 2024 11:43:49 +0200 Subject: [PATCH 01/29] Implement algorithm 2.2 from BS2007 for EllipticCurve_with_prime_order() --- src/doc/en/reference/references/index.rst | 3 + src/sage/schemes/elliptic_curves/all.py | 2 +- .../elliptic_curves/ell_finite_field.py | 114 ++++++++++++++++++ 3 files changed, 118 insertions(+), 1 deletion(-) diff --git a/src/doc/en/reference/references/index.rst b/src/doc/en/reference/references/index.rst index d7fe85422ba..35bc023c2db 100644 --- a/src/doc/en/reference/references/index.rst +++ b/src/doc/en/reference/references/index.rst @@ -1326,6 +1326,9 @@ REFERENCES: no. 1 (2003): 97-111, http://www.moi.math.bas.bg/moiuser/~iliya/pdf_site/gf5srev.pdf. +.. [BS2007] \R. Bröker and P. Stevenhagen. *Constructing elliptic curves of + prime order*. [math.NT] (2007), :arXiv:`0712.2022`. + .. [BS2010] \P. Baseilhac and K. Shigechi. *A new current algebra and the reflection equation*. Lett. Math. Phys. **92** (2010), pp. 47-65. :arxiv:`0906.1482`. diff --git a/src/sage/schemes/elliptic_curves/all.py b/src/sage/schemes/elliptic_curves/all.py index 84f7b0d5a50..82df6f45a96 100644 --- a/src/sage/schemes/elliptic_curves/all.py +++ b/src/sage/schemes/elliptic_curves/all.py @@ -27,7 +27,7 @@ lazy_import('sage.schemes.elliptic_curves.jacobian', 'Jacobian') lazy_import('sage.schemes.elliptic_curves.ell_finite_field', 'special_supersingular_curve') - +lazy_import('sage.schemes.elliptic_curves.ell_finite_field', 'EllipticCurve_with_prime_order') lazy_import('sage.schemes.elliptic_curves.ell_rational_field', ['cremona_curves', 'cremona_optimal_curves']) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index fcd23e7243f..e81011494b7 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -15,6 +15,8 @@ - Lorenz Panny, John Cremona (2023-02): ``.twists()`` - Lorenz Panny (2023): ``special_supersingular_curve()`` + +- Martin Grenouilloux (2024): ``EllipticCurve_with_prime_order()`` """ # **************************************************************************** @@ -2786,3 +2788,115 @@ def find_q(m, m4_fac, D): yield Et except ValueError: pass + +def EllipticCurve_with_prime_order(N): + r""" + Given a prime number ``N``, find another prime number `p` and construct an + elliptic curve `E` defined over `\mathbb F_p` such that + `\#E(\mathbb F_p) = N`. + + INPUT: + + - ``N`` -- integer; the order for which we seek an elliptic curve. Must be a + prime number. + + OUTPUT: an elliptic curve `E/\mathbb F_p` of order ``N`` + + EXAMPLES:: + + sage: N = next_prime(int(b'sagemath'.hex(), 16)) + sage: E = EllipticCurve_with_prime_order(N) + sage: E + Elliptic Curve defined by y^2 = x^3 + 4757897140353078952*x + + 1841350074072114366 over Finite Field of size 8314040074357871443 + sage: E.order() == N + True + + The execution time largely depends on the input:: + + sage: N = 125577861263605878504082476745517446213 + sage: E = EllipticCurve_with_prime_order(N) # Takes ~1 second. + sage: E.order() == N + True + + TESTS:: + + sage: N = next_prime(123456789) + sage: E = EllipticCurve_with_prime_order(N) + sage: E.order() == N + True + + sage: for N in prime_range(3, 100): + ....: E = EllipticCurve_with_prime_order(N) + ....: assert E.order() == N + + sage: N = 123456789 + sage: E = EllipticCurve_with_prime_order(N) + Traceback (most recent call last): + ... + ValueError: input order is not prime + + sage: E = EllipticCurve_with_prime_order(0) + Traceback (most recent call last): + ... + ValueError: input order is not prime + + ..NOTE:: + + Depending on the input, this function may run forever. + + ALGORITHM: [BS2007]_, Algorithm 2.2 + """ + from sage.arith.misc import is_prime, next_prime + from sage.combinat.subset import powerset + from sage.functions.other import floor + from sage.misc.functional import symbolic_prod as product + from sage.quadratic_forms.binary_qf import BinaryQF + from sage.schemes.elliptic_curves.cm import hilbert_class_polynomial + + if not is_prime(N): + raise ValueError("input order is not prime") + + # The algorithm consists of multiple ‘search rounds’ for a suitable + # discriminant `D`, `r` defines the number of rounds. We expect this + # algorithm to terminate after a number of rounds that is polynomial in + # loglog N. + r = 0 + S = [] + # First prime. Iterating over the primes by chunks of size log(`N`). + p = next_prime(0) + + while True: + while p < (r + 1) * floor(N.log()): + S.append(p) + p = next_prime(p) + r += 1 + + # Every possible products of distinct elements of `S`. + # There probably is a more optimal way to compute all possible products + # of elements of S than using a powerset. Here many multiplications are + # done multiple times. + for e in powerset(S): + D = product(e) + if -D % 8 != 5: + continue + + Q = BinaryQF([1, 0, D]) + sol = Q.solve_integer(4 * N) + if sol is None: + continue + + x, _ = sol + p1 = N + 1 - x + p2 = N + 1 + x + for p_i in [p1, p2]: + if is_prime(p_i): + H = hilbert_class_polynomial(-D) + for j0, _ in H.roots(ring=GF(p_i)): + E = EllipticCurve(j=j0) + if E.order() == N: + return E + else: + for Et in E.twists(): + if Et.order() == N: + return Et From 1767c5da1abfb71f35a988afcc04adae0a5681e0 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Tue, 9 Jul 2024 11:48:10 +0200 Subject: [PATCH 02/29] Changed $ sign for backticks in EllipticCurve_with_order() docstring --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index e81011494b7..624efe45ecf 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2725,8 +2725,8 @@ def EllipticCurve_with_order(m, *, D=None): sage: all(E.order() == 21 for E in Es) True - Indeed, we can verify that this is correct. Hasse's bounds tell us that $p \leq 50$ - (approximately), and the rest can be checked via bruteforce:: + Indeed, we can verify that this is correct. Hasse's bounds tell us that + `p \leq 50` (approximately), and the rest can be checked via bruteforce:: sage: for p in prime_range(50): ....: for j in range(p): From f26730ba9986b3cd77f9b385056be8d46ded8ba6 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Tue, 9 Jul 2024 12:26:07 +0200 Subject: [PATCH 03/29] Docstring "NOTE::" block format --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 624efe45ecf..94005b9bcb2 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2841,7 +2841,7 @@ def EllipticCurve_with_prime_order(N): ... ValueError: input order is not prime - ..NOTE:: + .. NOTE:: Depending on the input, this function may run forever. From 1d8a783fbcc5f88be265a16af21706c3c5ea8888 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Wed, 10 Jul 2024 11:12:06 +0200 Subject: [PATCH 04/29] Fixed to odd primes iteration and optimized their generation using prime_range() --- .../elliptic_curves/ell_finite_field.py | 36 ++++++++++++------- 1 file changed, 23 insertions(+), 13 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 94005b9bcb2..650209d2f4f 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2821,14 +2821,19 @@ def EllipticCurve_with_prime_order(N): TESTS:: - sage: N = next_prime(123456789) + sage: for N in prime_range(3, 100): + ....: E = EllipticCurve_with_prime_order(N) + ....: assert E.order() == N + + sage: N = 15175980689839334471 sage: E = EllipticCurve_with_prime_order(N) sage: E.order() == N True - sage: for N in prime_range(3, 100): - ....: E = EllipticCurve_with_prime_order(N) - ....: assert E.order() == N + sage: N = next_prime(123456789) + sage: E = EllipticCurve_with_prime_order(N) + sage: E.order() == N + True sage: N = 123456789 sage: E = EllipticCurve_with_prime_order(N) @@ -2841,17 +2846,23 @@ def EllipticCurve_with_prime_order(N): ... ValueError: input order is not prime + sage: E = EllipticCurve_with_prime_order(-7) + Traceback (most recent call last): + ... + ValueError: input order is not prime + .. NOTE:: Depending on the input, this function may run forever. ALGORITHM: [BS2007]_, Algorithm 2.2 """ - from sage.arith.misc import is_prime, next_prime + from sage.arith.misc import is_prime from sage.combinat.subset import powerset - from sage.functions.other import floor + from sage.functions.other import ceil from sage.misc.functional import symbolic_prod as product from sage.quadratic_forms.binary_qf import BinaryQF + from sage.rings.fast_arith import prime_range from sage.schemes.elliptic_curves.cm import hilbert_class_polynomial if not is_prime(N): @@ -2862,15 +2873,10 @@ def EllipticCurve_with_prime_order(N): # algorithm to terminate after a number of rounds that is polynomial in # loglog N. r = 0 - S = [] - # First prime. Iterating over the primes by chunks of size log(`N`). - p = next_prime(0) while True: - while p < (r + 1) * floor(N.log()): - S.append(p) - p = next_prime(p) - r += 1 + # Iterating over the odd primes by chunks of size log(`N`). + S = prime_range(3, ceil((r + 1) * N.log())) # Every possible products of distinct elements of `S`. # There probably is a more optimal way to compute all possible products @@ -2900,3 +2906,7 @@ def EllipticCurve_with_prime_order(N): for Et in E.twists(): if Et.order() == N: return Et + + # At this point, no discriminant has been found, moving to next round + # and extending the prime list. + r += 1 From dd08f6c384925ce5dbbe907568502314c2f559f1 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Wed, 10 Jul 2024 21:02:37 +0200 Subject: [PATCH 05/29] Better wording for NOTE block in docstring --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 650209d2f4f..5778d90c108 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2853,7 +2853,11 @@ def EllipticCurve_with_prime_order(N): .. NOTE:: - Depending on the input, this function may run forever. + Depending on the input, this function may run for a *very* long time. + This algorithm consists of multiple "search rounds" for a suitable + discriminant `D`. We expect this algorithm to terminate after a number + of rounds that is polynomial in `loglog N`. In practice (cf. Section 5), + this number is usually 1. ALGORITHM: [BS2007]_, Algorithm 2.2 """ From 95be4d80f5d45fa56847cf8ec1df268d5c7f9091 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Thu, 11 Jul 2024 10:29:37 +0200 Subject: [PATCH 06/29] Dynamic and precise bounds for the prime table at each round --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 5778d90c108..c7d6ce1c9b0 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2877,10 +2877,13 @@ def EllipticCurve_with_prime_order(N): # algorithm to terminate after a number of rounds that is polynomial in # loglog N. r = 0 + prime_start = 3 + prime_end = ceil((r + 1) * N.log()) + S = [] while True: # Iterating over the odd primes by chunks of size log(`N`). - S = prime_range(3, ceil((r + 1) * N.log())) + S.extend(prime_range(prime_start, prime_end)) # Every possible products of distinct elements of `S`. # There probably is a more optimal way to compute all possible products @@ -2914,3 +2917,6 @@ def EllipticCurve_with_prime_order(N): # At this point, no discriminant has been found, moving to next round # and extending the prime list. r += 1 + + prime_start = prime_end + prime_end = ceil((r + 1) * N.log()) From c84d6e661c9c7e11f3ba55d23af9befd8d003713 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Thu, 11 Jul 2024 17:48:48 +0200 Subject: [PATCH 07/29] Support python native int type --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index c7d6ce1c9b0..28cbd302631 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2864,7 +2864,7 @@ def EllipticCurve_with_prime_order(N): from sage.arith.misc import is_prime from sage.combinat.subset import powerset from sage.functions.other import ceil - from sage.misc.functional import symbolic_prod as product + from sage.misc.functional import symbolic_prod as product, log from sage.quadratic_forms.binary_qf import BinaryQF from sage.rings.fast_arith import prime_range from sage.schemes.elliptic_curves.cm import hilbert_class_polynomial @@ -2878,7 +2878,7 @@ def EllipticCurve_with_prime_order(N): # loglog N. r = 0 prime_start = 3 - prime_end = ceil((r + 1) * N.log()) + prime_end = ceil((r + 1) * log(N)) S = [] while True: @@ -2919,4 +2919,4 @@ def EllipticCurve_with_prime_order(N): r += 1 prime_start = prime_end - prime_end = ceil((r + 1) * N.log()) + prime_end = ceil((r + 1) * log(N)) From 4e78e2cacaaa59f866a3f8e8e0e1577f3aac1cf3 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Thu, 11 Jul 2024 20:09:36 +0200 Subject: [PATCH 08/29] Improved latex formatting in NOTE docstring --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 28cbd302631..7e31389ce80 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2856,8 +2856,8 @@ def EllipticCurve_with_prime_order(N): Depending on the input, this function may run for a *very* long time. This algorithm consists of multiple "search rounds" for a suitable discriminant `D`. We expect this algorithm to terminate after a number - of rounds that is polynomial in `loglog N`. In practice (cf. Section 5), - this number is usually 1. + of rounds that is polynomial in `\log\log N`. In practice (cf. Section + 5), this number is usually 1. ALGORITHM: [BS2007]_, Algorithm 2.2 """ From f2b9164d8a89ab5c808b0f6f85d7cd056a6c12ed Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Thu, 15 Aug 2024 17:01:57 +0200 Subject: [PATCH 09/29] Blazingly fast optimization by filtering primes with kronecker symbol --- .../elliptic_curves/ell_finite_field.py | 18 +++++++++++++----- 1 file changed, 13 insertions(+), 5 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 7e31389ce80..776b56ffa22 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2812,10 +2812,17 @@ def EllipticCurve_with_prime_order(N): sage: E.order() == N True - The execution time largely depends on the input:: + It works for large primes:: - sage: N = 125577861263605878504082476745517446213 - sage: E = EllipticCurve_with_prime_order(N) # Takes ~1 second. + sage: N = 0x6cbc824032974516623e732462f4b74b56c4ffbd984380d9 + sage: E = EllipticCurve_with_prime_order(N) + sage: E.order() == N + True + + But the execution time largely depends on the input:: + + sage: N = 200396817641911230625970463749415493753 + sage: E = EllipticCurve_with_prime_order(N) sage: E.order() == N True @@ -2861,7 +2868,7 @@ def EllipticCurve_with_prime_order(N): ALGORITHM: [BS2007]_, Algorithm 2.2 """ - from sage.arith.misc import is_prime + from sage.arith.misc import is_prime, kronecker from sage.combinat.subset import powerset from sage.functions.other import ceil from sage.misc.functional import symbolic_prod as product, log @@ -2883,7 +2890,8 @@ def EllipticCurve_with_prime_order(N): while True: # Iterating over the odd primes by chunks of size log(`N`). - S.extend(prime_range(prime_start, prime_end)) + S.extend(p for p in prime_range(prime_start, prime_end) + if kronecker(N, p) == 1) # Every possible products of distinct elements of `S`. # There probably is a more optimal way to compute all possible products From 2e37e99e43ccd18ee6101f3a41ed3bc2a7930474 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Wed, 4 Sep 2024 15:22:04 +0200 Subject: [PATCH 10/29] Perfect compliance with algorithm from paper --- .../elliptic_curves/ell_finite_field.py | 39 ++++++++++--------- 1 file changed, 20 insertions(+), 19 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 776b56ffa22..8afb1ddc840 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2868,7 +2868,7 @@ def EllipticCurve_with_prime_order(N): ALGORITHM: [BS2007]_, Algorithm 2.2 """ - from sage.arith.misc import is_prime, kronecker + from sage.arith.misc import is_prime, legendre_symbol from sage.combinat.subset import powerset from sage.functions.other import ceil from sage.misc.functional import symbolic_prod as product, log @@ -2879,7 +2879,7 @@ def EllipticCurve_with_prime_order(N): if not is_prime(N): raise ValueError("input order is not prime") - # The algorithm consists of multiple ‘search rounds’ for a suitable + # The algorithm consists of multiple "search rounds" for a suitable # discriminant `D`, `r` defines the number of rounds. We expect this # algorithm to terminate after a number of rounds that is polynomial in # loglog N. @@ -2890,8 +2890,10 @@ def EllipticCurve_with_prime_order(N): while True: # Iterating over the odd primes by chunks of size log(`N`). - S.extend(p for p in prime_range(prime_start, prime_end) - if kronecker(N, p) == 1) + for p in prime_range(prime_start, prime_end): + if legendre_symbol(N, p) == 1: + # Equivalent to p* = (-1)^((p - 1) / 2) * p in [BS2007]_ page 5. + S.append(-p if p >> 1 & 1 else p) # Every possible products of distinct elements of `S`. # There probably is a more optimal way to compute all possible products @@ -2899,32 +2901,31 @@ def EllipticCurve_with_prime_order(N): # done multiple times. for e in powerset(S): D = product(e) - if -D % 8 != 5: + if D % 8 != 5 or D >= 0 or D >= prime_start^2: continue - Q = BinaryQF([1, 0, D]) - sol = Q.solve_integer(4 * N) + Q = BinaryQF([1, 0, -D]) + sol = Q.solve_integer(4 * N, algorithm='cornacchia') if sol is None: continue x, _ = sol - p1 = N + 1 - x - p2 = N + 1 + x - for p_i in [p1, p2]: + for p_i in [N + 1 - x, N + 1 + x]: if is_prime(p_i): - H = hilbert_class_polynomial(-D) - for j0, _ in H.roots(ring=GF(p_i)): + H = hilbert_class_polynomial(D) + for j0 in H.roots(ring=GF(p_i), multiplicities=False): E = EllipticCurve(j=j0) - if E.order() == N: - return E - else: - for Et in E.twists(): - if Et.order() == N: - return Et + # `E.twists()` also contains E. + for Et in E.twists(): + if Et.order() == N: + return Et # At this point, no discriminant has been found, moving to next round # and extending the prime list. r += 1 - prime_start = prime_end + # For small `N`, the value `(r+1)log(N)` can be less than 3. When that's + # the case, we don't need to worry about prime duplicates in `S` since + # it will be empty. + prime_start = max(3, prime_end) prime_end = ceil((r + 1) * log(N)) From 9085cc91a0fbdae6e57cdcd0f4e9ab4fbc9df131 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Fri, 6 Sep 2024 13:39:12 +0200 Subject: [PATCH 11/29] Changing to iterator output --- .../elliptic_curves/ell_finite_field.py | 35 ++++++++++++------- 1 file changed, 23 insertions(+), 12 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 8afb1ddc840..b85336ff843 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2800,60 +2800,71 @@ def EllipticCurve_with_prime_order(N): - ``N`` -- integer; the order for which we seek an elliptic curve. Must be a prime number. - OUTPUT: an elliptic curve `E/\mathbb F_p` of order ``N`` + OUTPUT: an iterator of elliptic curves `E/\mathbb F_p` of order ``N`` EXAMPLES:: sage: N = next_prime(int(b'sagemath'.hex(), 16)) - sage: E = EllipticCurve_with_prime_order(N) + sage: E = next(EllipticCurve_with_prime_order(N)) sage: E Elliptic Curve defined by y^2 = x^3 + 4757897140353078952*x + 1841350074072114366 over Finite Field of size 8314040074357871443 sage: E.order() == N True + You can directly iterate over the function call; here only on the first 10 + curves:: + + sage: N = 54675917 + sage: for _, E in zip(range(10), EllipticCurve_with_prime_order(N)): + ....: assert E.order() == N + It works for large primes:: sage: N = 0x6cbc824032974516623e732462f4b74b56c4ffbd984380d9 - sage: E = EllipticCurve_with_prime_order(N) + sage: E = next(EllipticCurve_with_prime_order(N)) sage: E.order() == N True But the execution time largely depends on the input:: sage: N = 200396817641911230625970463749415493753 - sage: E = EllipticCurve_with_prime_order(N) + sage: E = next(EllipticCurve_with_prime_order(N)) sage: E.order() == N True TESTS:: sage: for N in prime_range(3, 100): - ....: E = EllipticCurve_with_prime_order(N) + ....: E = next(EllipticCurve_with_prime_order(N)) + ....: assert E.order() == N + + sage: N = 113 + sage: for _, E in zip(range(30), EllipticCurve_with_prime_order(N)): ....: assert E.order() == N sage: N = 15175980689839334471 - sage: E = EllipticCurve_with_prime_order(N) + sage: E = next(EllipticCurve_with_prime_order(N)) sage: E.order() == N True sage: N = next_prime(123456789) - sage: E = EllipticCurve_with_prime_order(N) + sage: E = next(EllipticCurve_with_prime_order(N)) sage: E.order() == N True sage: N = 123456789 - sage: E = EllipticCurve_with_prime_order(N) + sage: E = next(EllipticCurve_with_prime_order(N)) Traceback (most recent call last): ... ValueError: input order is not prime - sage: E = EllipticCurve_with_prime_order(0) + sage: E = next(EllipticCurve_with_prime_order(0)) Traceback (most recent call last): ... ValueError: input order is not prime - sage: E = EllipticCurve_with_prime_order(-7) + sage: E = next(EllipticCurve_with_prime_order(-7)) Traceback (most recent call last): ... ValueError: input order is not prime @@ -2901,7 +2912,7 @@ def EllipticCurve_with_prime_order(N): # done multiple times. for e in powerset(S): D = product(e) - if D % 8 != 5 or D >= 0 or D >= prime_start^2: + if D % 8 != 5 or D >= 0: continue Q = BinaryQF([1, 0, -D]) @@ -2918,7 +2929,7 @@ def EllipticCurve_with_prime_order(N): # `E.twists()` also contains E. for Et in E.twists(): if Et.order() == N: - return Et + yield Et # At this point, no discriminant has been found, moving to next round # and extending the prime list. From 36cb2af4441200df4d12110b7eee76f1e5a583c9 Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Mon, 9 Sep 2024 20:32:05 +0100 Subject: [PATCH 12/29] factor out `has_order` --- src/sage/schemes/elliptic_curves/all.py | 2 +- .../elliptic_curves/ell_finite_field.py | 117 ++++++++++++++---- 2 files changed, 92 insertions(+), 27 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/all.py b/src/sage/schemes/elliptic_curves/all.py index 82df6f45a96..df527315527 100644 --- a/src/sage/schemes/elliptic_curves/all.py +++ b/src/sage/schemes/elliptic_curves/all.py @@ -27,10 +27,10 @@ lazy_import('sage.schemes.elliptic_curves.jacobian', 'Jacobian') lazy_import('sage.schemes.elliptic_curves.ell_finite_field', 'special_supersingular_curve') -lazy_import('sage.schemes.elliptic_curves.ell_finite_field', 'EllipticCurve_with_prime_order') lazy_import('sage.schemes.elliptic_curves.ell_rational_field', ['cremona_curves', 'cremona_optimal_curves']) +from sage.schemes.elliptic_curves.ell_finite_field import EllipticCurve_with_prime_order from sage.schemes.elliptic_curves.cm import (cm_orders, cm_j_invariants, cm_j_invariants_and_orders, diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index b85336ff843..4d4e01a0bff 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -1343,6 +1343,87 @@ def is_ordinary(self, proof=True): """ return not is_j_supersingular(self.j_invariant(), proof=proof) + def has_order(self, value, num_checks=8): + r""" + Return ``True`` if the curve has order ``value``. + + INPUT: + + - ``value`` -- integer in the Hasse-Weil range for this curve + + - ``num_checks``-- integer (default: `8`); the number of times to check + whether ``value`` times a random point on this curve equals the + identity + + .. NOTE:: + + Since the method is probabilistic, there is a possibility for the + method to yield false positives (i.e. returning ``True`` even when + the result is ``False``). Even worse, it is possible for this to + happen even when ``num_checks`` is increased arbitrarily. See below + for an example and :issue:`38617` for an open discussion. + + EXAMPLES: + + For curves over small finite fields, the order is computed and compared + directly:: + + sage: E = EllipticCurve(GF(7), [0, 1]) + sage: E.order() + 12 + sage: E.has_order(12, num_checks=0) + True + sage: E.has_order(11, num_checks=0) + False + sage: E.has_order(13, num_checks=0) + False + + This tests the method on a random curve:: + + sage: # long time (10s) + sage: p = random_prime(2**128, lbound=2**127) + sage: K = GF((p, 2), name="a") + sage: E = EllipticCurve(K, [K.random_element() for _ in range(2)]) + sage: N = E.order() + sage: E.has_order(N, num_checks=20) + True + sage: E.has_order(N + 1) + False + + This demonstrates the bug mentioned in the NOTE above. The return value + should be ``False`` after :issue:`38617` is fixed:: + + sage: E = EllipticCurve(GF(127), [0, 1]) + sage: E.order() + 108 + sage: E.has_order(126) + True + + AUTHORS: + + - Mariah Lenox (2011-02-16): Initial implementation + + - Gareth Ma (2024-01-21): Fix bug for small curves + """ + # This method does *not* use the cached value of `_order` even when available + q = self.base_field().order() + a, b = Hasse_bounds(q, 1) + if not a <= value <= b: + return False + + # For really small values, the random tests are too weak to detect wrong orders + # So we go with computing directly instead. + if q <= 100: + return self.order() == value + + # Is value * random == identity? + for _ in range(num_checks): + G = self.random_point() + if not (value * G).is_zero(): + return False + + return True + def set_order(self, value, *, check=True, num_checks=8): r""" Set the value of ``self._order`` to ``value``. @@ -1357,7 +1438,7 @@ def set_order(self, value, *, check=True, num_checks=8): - ``check``-- boolean (default: ``True``); whether or not to run sanity checks on the input - - ``num_checks``-- integer (default: 8); if ``check`` is + - ``num_checks``-- integer (default: `8`); if ``check`` is ``True``, the number of times to check whether ``value`` times a random point on this curve equals the identity @@ -1403,11 +1484,11 @@ def set_order(self, value, *, check=True, num_checks=8): sage: E.set_order(0) Traceback (most recent call last): ... - ValueError: Value 0 illegal (not an integer in the Hasse range) + ValueError: Elliptic Curve defined by y^2 = x^3 + 1 over Finite Field of size 7 does not have order 0 sage: E.set_order(1000) Traceback (most recent call last): ... - ValueError: Value 1000 illegal (not an integer in the Hasse range) + ValueError: Elliptic Curve defined by y^2 = x^3 + 1 over Finite Field of size 7 does not have order 1000 It is also very likely an error to pass a value which is not the actual order of this curve. How unlikely is determined by @@ -1418,7 +1499,7 @@ def set_order(self, value, *, check=True, num_checks=8): sage: E.set_order(947) Traceback (most recent call last): ... - ValueError: Value 947 illegal (multiple of random point not the identity) + ValueError: Elliptic Curve defined by y^2 = x^3 + 1 over Finite Field of size 1009 does not have order 947 For curves over small finite fields, the order is cheap to compute, so it is computed directly and compared:: @@ -1427,7 +1508,7 @@ def set_order(self, value, *, check=True, num_checks=8): sage: E.set_order(11) Traceback (most recent call last): ... - ValueError: Value 11 illegal (correct order is 12) + ValueError: Elliptic Curve defined by y^2 = x^3 + 1 over Finite Field of size 7 does not have order 11 TESTS: @@ -1438,7 +1519,7 @@ def set_order(self, value, *, check=True, num_checks=8): sage: E.set_order(3) Traceback (most recent call last): ... - ValueError: Value 3 illegal (correct order is 1) + ValueError: Elliptic Curve defined by y^2 + y = x^3 + x + 1 over Finite Field of size 2 does not have order 3 :: @@ -1446,7 +1527,7 @@ def set_order(self, value, *, check=True, num_checks=8): sage: E.set_order(4, num_checks=0) Traceback (most recent call last): ... - ValueError: Value 4 illegal (correct order is 12) + ValueError: Elliptic Curve defined by y^2 = x^3 + 1 over Finite Field of size 7 does not have order 4 sage: E.order() 12 @@ -1461,24 +1542,8 @@ def set_order(self, value, *, check=True, num_checks=8): """ value = Integer(value) - if check: - # Is value in the Hasse range? - q = self.base_field().order() - a,b = Hasse_bounds(q,1) - if not a <= value <= b: - raise ValueError(f"Value {value} illegal (not an integer in the Hasse range)") - - # For really small values, the random tests are too weak to detect wrong orders - # So we go with computing directly instead. - if q <= 100: - if self.order() != value: - raise ValueError(f"Value {value} illegal (correct order is {self.order()})") - - # Is value*random == identity? - for _ in range(num_checks): - G = self.random_point() - if value * G != self(0): - raise ValueError(f"Value {value} illegal (multiple of random point not the identity)") + if check and not self.has_order(value, num_checks=num_checks): + raise ValueError(f"{self} does not have order {value}") # TODO: It might help some of PARI's algorithms if we # could copy this over to the .pari_curve() as well. @@ -2678,7 +2743,7 @@ def EllipticCurve_with_order(m, *, D=None): Return an iterator for elliptic curves over finite fields with the given order. The curves are computed using the Complex Multiplication (CM) method. - A `:sage:`~sage.structure.factorization.Factorization` can be passed for ``m``, in which case + A :sage:`~sage.structure.factorization.Factorization` can be passed for ``m``, in which case the algorithm is more efficient. If ``D`` is specified, it is used as the discriminant. From 57cf88b019371c355d920f8118dac322f057f48e Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Mon, 9 Sep 2024 21:50:34 +0100 Subject: [PATCH 13/29] add extensive docs and gray code optimisation --- .../elliptic_curves/ell_finite_field.py | 138 ++++++++++++------ 1 file changed, 97 insertions(+), 41 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 4d4e01a0bff..8095c755041 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -1416,9 +1416,17 @@ def has_order(self, value, num_checks=8): if q <= 100: return self.order() == value + # This might be slow + # if value.is_prime(): + # num_checks = 1 + # Is value * random == identity? for _ in range(num_checks): - G = self.random_point() + while True: + G = self.random_point() + if not G.is_zero(): + break + if not (value * G).is_zero(): return False @@ -2887,17 +2895,50 @@ def EllipticCurve_with_prime_order(N): It works for large primes:: sage: N = 0x6cbc824032974516623e732462f4b74b56c4ffbd984380d9 - sage: E = next(EllipticCurve_with_prime_order(N)) + sage: E = next(EllipticCurve_with_prime_order(N)); E + Elliptic Curve defined by y^2 = x^3 + 2666207849820848272386538889427721639173508298483739490459*x + + 77986137112576 over Finite Field of size 2666207849820848272386538889427721639173508298487130585243 sage: E.order() == N True - But the execution time largely depends on the input:: + The execution time largely depends on the input, specifically the smallest + discriminant ``D`` for which we can apply CM method on. Here it takes + slightly longer, though still within `1` second:: sage: N = 200396817641911230625970463749415493753 - sage: E = next(EllipticCurve_with_prime_order(N)) + sage: E = next(EllipticCurve_with_prime_order(N)); E sage: E.order() == N True + Note that the iterator does *not* return all curves with the given order:: + + sage: any(E.base_ring() is GF(7) for E in EllipticCurve_with_prime_order(7)) + False + sage: EllipticCurve(GF(7), [0, 5]).order() + 7 + + However, experimentally it returns many of them. Here it returns all of them:: + + sage: N = 23 + sage: curves = list(EllipticCurve_with_prime_order(N)); curves + [Elliptic Curve defined by y^2 = x^3 + 20*x + 20 over Finite Field of size 23, + Elliptic Curve defined by y^2 = x^3 + 10*x + 16 over Finite Field of size 23, + Elliptic Curve defined by y^2 = x^3 + 14*x + 14 over Finite Field of size 17, + Elliptic Curve defined by y^2 = x^3 + 16*x + 25 over Finite Field of size 31, + Elliptic Curve defined by y^2 = x^3 + 13*x + 6 over Finite Field of size 19, + Elliptic Curve defined by y^2 = x^3 + 24*x + 3 over Finite Field of size 29] + sage: import itertools + sage: # These are the only primes, by the Weil-Hasse bound + sage: for q in prime_range(17, 35): + ....: K = GF(q) + ....: for u in itertools.product(range(q), repeat=2): + ....: try: E = EllipticCurve(GF(q), u) + ....: except ArithmeticError: continue + ....: if E.order() == N: + ....: assert any(E.is_isomorphic(E_) for E_ in curves) + + + TESTS:: sage: for N in prime_range(3, 100): @@ -2944,42 +2985,63 @@ def EllipticCurve_with_prime_order(N): ALGORITHM: [BS2007]_, Algorithm 2.2 """ + import itertools + from sage.misc.verbose import verbose + from sage.combinat import gray_codes + from sage.sets.primes import Primes from sage.arith.misc import is_prime, legendre_symbol - from sage.combinat.subset import powerset - from sage.functions.other import ceil - from sage.misc.functional import symbolic_prod as product, log from sage.quadratic_forms.binary_qf import BinaryQF - from sage.rings.fast_arith import prime_range from sage.schemes.elliptic_curves.cm import hilbert_class_polynomial if not is_prime(N): raise ValueError("input order is not prime") - # The algorithm consists of multiple "search rounds" for a suitable - # discriminant `D`, `r` defines the number of rounds. We expect this - # algorithm to terminate after a number of rounds that is polynomial in - # loglog N. - r = 0 - prime_start = 3 - prime_end = ceil((r + 1) * log(N)) + if N < 1000: + # Log how long this algorithm will take + from sage.rings.fast_arith import prime_range + num = sum(1 for p in prime_range(3, 4 * N) if legendre_symbol(N, p) == 1) + verbose(f"Total work to enumerate curves with order {N}: 2^{num}", level=2) + + # The algorithm considers smooth discriminants `D`, sorted by their largest + # prime factor. We expect this algorithm to terminate after a number of + # rounds that is polynomial in loglog N. S = [] - while True: - # Iterating over the odd primes by chunks of size log(`N`). - for p in prime_range(prime_start, prime_end): - if legendre_symbol(N, p) == 1: - # Equivalent to p* = (-1)^((p - 1) / 2) * p in [BS2007]_ page 5. - S.append(-p if p >> 1 & 1 else p) - - # Every possible products of distinct elements of `S`. - # There probably is a more optimal way to compute all possible products - # of elements of S than using a powerset. Here many multiplications are - # done multiple times. - for e in powerset(S): - D = product(e) - if D % 8 != 5 or D >= 0: + for p in Primes(): + if p == 2: + continue + + if legendre_symbol(N, p) != 1: + continue + + verbose(f"Considering {len(S) + 1}th valid prime {p}", level=3) + + # Equivalent to p* = (-1)^((p - 1) / 2) * p in [BS2007]_ page 5. + p_star = -p if p >> 1 & 1 else p + + # later we need x^2 + (-D)y^2 = 4N, and since y = 0 has no solution, we need + # p = |p_star| <= |-D| <= 4N + if p > 4 * N: + break + + # Although this is O(2^n) instead of O(n2^n), it misses out on the optimisation where + # when D becomes too large (see -D > 4 * N below), one can cut the branch entirely + # so perhaps this is slower + D = p_star + for idx, dx in itertools.chain([(None, 1)], gray_codes.product([2] * len(S))): + assert abs(dx) == 1, f"hmm, gray codes shouldn't do that. file a bug report" + + if idx is not None: + if dx == 1: + D *= S[idx] + else: + D //= S[idx] + + if D % 8 != 5 or D >= 0 or -D > 4 * N: continue + verbose(f"Testing {D=}", level=3) + Q = BinaryQF([1, 0, -D]) sol = Q.solve_integer(4 * N, algorithm='cornacchia') if sol is None: @@ -2989,19 +3051,13 @@ def EllipticCurve_with_prime_order(N): for p_i in [N + 1 - x, N + 1 + x]: if is_prime(p_i): H = hilbert_class_polynomial(D) - for j0 in H.roots(ring=GF(p_i), multiplicities=False): - E = EllipticCurve(j=j0) + K = GF(p_i) + for j0 in H.roots(ring=K, multiplicities=False): + E = EllipticCurve(K, j=j0) # `E.twists()` also contains E. for Et in E.twists(): - if Et.order() == N: + # num_checks=1 is sufficient for prime order + if Et.has_order(N, num_checks=1): yield Et - # At this point, no discriminant has been found, moving to next round - # and extending the prime list. - r += 1 - - # For small `N`, the value `(r+1)log(N)` can be less than 3. When that's - # the case, we don't need to worry about prime duplicates in `S` since - # it will be empty. - prime_start = max(3, prime_end) - prime_end = ceil((r + 1) * log(N)) + S.append(p_star) From bae2c8488b1ea6daf24cbbf5a79d2ad8376797a5 Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Mon, 9 Sep 2024 22:49:12 +0100 Subject: [PATCH 14/29] use recursive implementation instead of gray code --- .../elliptic_curves/ell_finite_field.py | 113 +++++++++++++----- 1 file changed, 83 insertions(+), 30 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 8095c755041..212093cdca0 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2862,7 +2862,7 @@ def find_q(m, m4_fac, D): except ValueError: pass -def EllipticCurve_with_prime_order(N): +def EllipticCurve_with_prime_order(N, max_D=10**8): r""" Given a prime number ``N``, find another prime number `p` and construct an elliptic curve `E` defined over `\mathbb F_p` such that @@ -2877,15 +2877,24 @@ def EllipticCurve_with_prime_order(N): EXAMPLES:: - sage: N = next_prime(int(b'sagemath'.hex(), 16)) + sage: N = next_prime(int.from_bytes(b'sagemath', 'big')) sage: E = next(EllipticCurve_with_prime_order(N)) sage: E - Elliptic Curve defined by y^2 = x^3 + 4757897140353078952*x + - 1841350074072114366 over Finite Field of size 8314040074357871443 + Elliptic Curve defined by y^2 = x^3 + 4757897140353078952*x + 1841350074072114366 + over Finite Field of size 8314040074357871443 sage: E.order() == N True - You can directly iterate over the function call; here only on the first 10 + The returned curves are sometimes random because + :meth:`~sage.schemes.elliptic_curves.ell_finite_field.EllipticCurve_finite_field.twists` is + not deterministic. However, it's always isomorphic:: + + sage: E = next(EllipticCurve_with_prime_order(23)); E # random + Elliptic Curve defined by y^2 = x^3 + 12*x + 11 over Finite Field of size 17 + sage: E.is_isomorphic(EllipticCurve(GF(17), [3, 5])) + True + + You can directly iterate over the iterator; here only on the first 10 curves:: sage: N = 54675917 @@ -2920,13 +2929,12 @@ def EllipticCurve_with_prime_order(N): However, experimentally it returns many of them. Here it returns all of them:: sage: N = 23 + sage: set_random_seed(1337) # as the function returns random twists of curves sage: curves = list(EllipticCurve_with_prime_order(N)); curves - [Elliptic Curve defined by y^2 = x^3 + 20*x + 20 over Finite Field of size 23, - Elliptic Curve defined by y^2 = x^3 + 10*x + 16 over Finite Field of size 23, - Elliptic Curve defined by y^2 = x^3 + 14*x + 14 over Finite Field of size 17, - Elliptic Curve defined by y^2 = x^3 + 16*x + 25 over Finite Field of size 31, - Elliptic Curve defined by y^2 = x^3 + 13*x + 6 over Finite Field of size 19, - Elliptic Curve defined by y^2 = x^3 + 24*x + 3 over Finite Field of size 29] + [Elliptic Curve defined by y^2 = x^3 + 3*x + 5 over Finite Field of size 17, + Elliptic Curve defined by y^2 = x^3 + 19*x + 14 over Finite Field of size 31, + Elliptic Curve defined by y^2 = x^3 + 2*x + 9 over Finite Field of size 19, + Elliptic Curve defined by y^2 = x^3 + 7*x + 18 over Finite Field of size 29] sage: import itertools sage: # These are the only primes, by the Weil-Hasse bound sage: for q in prime_range(17, 35): @@ -2937,7 +2945,43 @@ def EllipticCurve_with_prime_order(N): ....: if E.order() == N: ....: assert any(E.is_isomorphic(E_) for E_ in curves) - + The algorithm is efficient for small ``N`` due to the low number of + suitable discriminants (see the ``recur`` internal function of the code for + details). The following runs in less than a second:: + + sage: len(list(EllipticCurve_with_prime_order(next_prime(5000)))) + 402 + + There is different verbose data for level `2` and `3`. Note that if you try + to compute the 5th curve in the iterator, the algorithm takes a long time + to run, as it is trying to compute a large Hilbert class polynomial, as + noted in the verbose log:: + + sage: set_random_seed(1337) # as the function returns random twists of curves + sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): + ....: print(E) + Elliptic Curve defined by y^2 = x^3 + 703734957*x + 232615553 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 278352808*x + 442354703 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 22138181*x + 343211325 over Finite Field of size 999969307 + sage: from sage.misc.verbose import set_verbose + sage: set_verbose(2) + sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): + ....: print(E) + verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 + Elliptic Curve defined by y^2 = x^3 + 276281332*x + 798205466 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 515293837*x + 588191363 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 684893162*x + 679980762 over Finite Field of size 999969307 + sage: set_verbose(3) + sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): + ....: print(E) + verbose 3 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Considering 1th valid prime 19 + verbose 3 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Considering 2th valid prime 23 + verbose 3 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Considering 3th valid prime 29 + verbose 3 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-667 + verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 + Elliptic Curve defined by y^2 = x^3 + 620755650*x + 816246309 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 651416562*x + 875284469 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 149283707*x + 72132991 over Finite Field of size 999969307 TESTS:: @@ -2986,6 +3030,7 @@ def EllipticCurve_with_prime_order(N): ALGORITHM: [BS2007]_, Algorithm 2.2 """ import itertools + from sage.misc.misc_c import prod from sage.misc.verbose import verbose from sage.combinat import gray_codes from sage.sets.primes import Primes @@ -2996,17 +3041,33 @@ def EllipticCurve_with_prime_order(N): if not is_prime(N): raise ValueError("input order is not prime") - if N < 1000: - # Log how long this algorithm will take - from sage.rings.fast_arith import prime_range - num = sum(1 for p in prime_range(3, 4 * N) if legendre_symbol(N, p) == 1) - verbose(f"Total work to enumerate curves with order {N}: 2^{num}", level=2) + # if N < 1000: + # # Log how long this algorithm will take + # from sage.rings.fast_arith import prime_range + # num = sum(1 for p in prime_range(3, 4 * N) if legendre_symbol(N, p) == 1) + # verbose(f"Total work to enumerate curves with order {N}: 2^{num}", level=2) # The algorithm considers smooth discriminants `D`, sorted by their largest # prime factor. We expect this algorithm to terminate after a number of # rounds that is polynomial in loglog N. S = [] + def recur(idx, bound): + """ + This function returns all subsets of S[idx:] (at the time of running) + with product of absolute value <= bound. It's extremely quick since S + is basically prime_range(N) while bound is basically 4 * N. + """ + if idx >= len(S): + return + yield [] + for nxt in range(idx + 1, len(S)): + if abs(S[nxt]) <= bound: + for r in recur(nxt, bound // abs(S[nxt])): + yield [S[nxt]] + r + else: + break + for p in Primes(): if p == 2: continue @@ -3024,20 +3085,11 @@ def EllipticCurve_with_prime_order(N): if p > 4 * N: break - # Although this is O(2^n) instead of O(n2^n), it misses out on the optimisation where - # when D becomes too large (see -D > 4 * N below), one can cut the branch entirely - # so perhaps this is slower - D = p_star - for idx, dx in itertools.chain([(None, 1)], gray_codes.product([2] * len(S))): - assert abs(dx) == 1, f"hmm, gray codes shouldn't do that. file a bug report" - - if idx is not None: - if dx == 1: - D *= S[idx] - else: - D //= S[idx] + for e in recur(0, 4 * N // p): + D = p_star * prod(e) + assert abs(D) <= 4 * N - if D % 8 != 5 or D >= 0 or -D > 4 * N: + if D % 8 != 5 or D >= 0: continue verbose(f"Testing {D=}", level=3) @@ -3050,6 +3102,7 @@ def EllipticCurve_with_prime_order(N): x, _ = sol for p_i in [N + 1 - x, N + 1 + x]: if is_prime(p_i): + verbose(f"Computing the Hilbert class polynomial H_{D}", level=2) H = hilbert_class_polynomial(D) K = GF(p_i) for j0 in H.roots(ring=K, multiplicities=False): From bf706c5d390ca52d1720e099fde78e0a17dd9695 Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Tue, 10 Sep 2024 00:55:16 +0100 Subject: [PATCH 15/29] fix docs --- .../elliptic_curves/ell_finite_field.py | 175 ++++++++++-------- 1 file changed, 95 insertions(+), 80 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 212093cdca0..9d2666ac989 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2862,7 +2862,7 @@ def find_q(m, m4_fac, D): except ValueError: pass -def EllipticCurve_with_prime_order(N, max_D=10**8): +def EllipticCurve_with_prime_order(N): r""" Given a prime number ``N``, find another prime number `p` and construct an elliptic curve `E` defined over `\mathbb F_p` such that @@ -2875,6 +2875,16 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): OUTPUT: an iterator of elliptic curves `E/\mathbb F_p` of order ``N`` + .. NOTE:: + + Depending on the input, this function may run for a *very* long time. + This algorithm consists of multiple "search rounds" for a suitable + discriminant `D`. We expect this algorithm to terminate after a number + of rounds that is polynomial in `\log\log N`. In practice (cf. Section + 5), this number is usually 1. + + ALGORITHM: [BS2007]_, Algorithm 2.2 + EXAMPLES:: sage: N = next_prime(int.from_bytes(b'sagemath', 'big')) @@ -2882,7 +2892,7 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): sage: E Elliptic Curve defined by y^2 = x^3 + 4757897140353078952*x + 1841350074072114366 over Finite Field of size 8314040074357871443 - sage: E.order() == N + sage: E.has_order(N) True The returned curves are sometimes random because @@ -2890,7 +2900,7 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): not deterministic. However, it's always isomorphic:: sage: E = next(EllipticCurve_with_prime_order(23)); E # random - Elliptic Curve defined by y^2 = x^3 + 12*x + 11 over Finite Field of size 17 + Elliptic Curve defined by y^2 = x^3 + 12*x + 6 over Finite Field of size 17 sage: E.is_isomorphic(EllipticCurve(GF(17), [3, 5])) True @@ -2899,7 +2909,7 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): sage: N = 54675917 sage: for _, E in zip(range(10), EllipticCurve_with_prime_order(N)): - ....: assert E.order() == N + ....: assert E.has_order(N) It works for large primes:: @@ -2907,16 +2917,19 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): sage: E = next(EllipticCurve_with_prime_order(N)); E Elliptic Curve defined by y^2 = x^3 + 2666207849820848272386538889427721639173508298483739490459*x + 77986137112576 over Finite Field of size 2666207849820848272386538889427721639173508298487130585243 - sage: E.order() == N + sage: E.has_order(N) True - The execution time largely depends on the input, specifically the smallest - discriminant ``D`` for which we can apply CM method on. Here it takes - slightly longer, though still within `1` second:: + :: - sage: N = 200396817641911230625970463749415493753 - sage: E = next(EllipticCurve_with_prime_order(N)); E - sage: E.order() == N + sage: N = next_prime(2^256) + sage: E = next(EllipticCurve_with_prime_order(N)); E # random + Elliptic Curve defined by y^2 = x^3 + 6056521267553273205988520276135607487700943205131813669424576873701361709521*x + + 86942739955486781674010637133214195706465136689012129911736706024465988573567 over Finite Field of size + 115792089237316195423570985008687907853847329310253429036565151476471048389761 + sage: E.j_invariant() + 111836223967433630316209796253554285080540088646141285337487360944738698436350 + sage: E.has_order(N) True Note that the iterator does *not* return all curves with the given order:: @@ -2930,11 +2943,13 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): sage: N = 23 sage: set_random_seed(1337) # as the function returns random twists of curves - sage: curves = list(EllipticCurve_with_prime_order(N)); curves + sage: curves = list(EllipticCurve_with_prime_order(N)); curves # random [Elliptic Curve defined by y^2 = x^3 + 3*x + 5 over Finite Field of size 17, Elliptic Curve defined by y^2 = x^3 + 19*x + 14 over Finite Field of size 31, Elliptic Curve defined by y^2 = x^3 + 2*x + 9 over Finite Field of size 19, - Elliptic Curve defined by y^2 = x^3 + 7*x + 18 over Finite Field of size 29] + Elliptic Curve defined by y^2 = x^3 + 7*x + 18 over Finite Field of size 29, + Elliptic Curve defined by y^2 = x^3 + 20*x + 20 over Finite Field of size 23, + Elliptic Curve defined by y^2 = x^3 + 10*x + 16 over Finite Field of size 23] sage: import itertools sage: # These are the only primes, by the Weil-Hasse bound sage: for q in prime_range(17, 35): @@ -2942,7 +2957,7 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): ....: for u in itertools.product(range(q), repeat=2): ....: try: E = EllipticCurve(GF(q), u) ....: except ArithmeticError: continue - ....: if E.order() == N: + ....: if E.has_order(N): ....: assert any(E.is_isomorphic(E_) for E_ in curves) The algorithm is efficient for small ``N`` due to the low number of @@ -2950,57 +2965,63 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): details). The following runs in less than a second:: sage: len(list(EllipticCurve_with_prime_order(next_prime(5000)))) - 402 + 534 - There is different verbose data for level `2` and `3`. Note that if you try - to compute the 5th curve in the iterator, the algorithm takes a long time - to run, as it is trying to compute a large Hilbert class polynomial, as - noted in the verbose log:: + There is different verbose data for level `2` to `4`, though level `3` + rarely logs anything (it logs when a new prime `p` is added to the + smoothness bound):: + sage: from sage.misc.verbose import set_verbose sage: set_random_seed(1337) # as the function returns random twists of curves sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): ....: print(E) - Elliptic Curve defined by y^2 = x^3 + 703734957*x + 232615553 over Finite Field of size 999969307 - Elliptic Curve defined by y^2 = x^3 + 278352808*x + 442354703 over Finite Field of size 999969307 - Elliptic Curve defined by y^2 = x^3 + 22138181*x + 343211325 over Finite Field of size 999969307 - sage: from sage.misc.verbose import set_verbose + Elliptic Curve defined by y^2 = x^3 + 265977778*x + 120868502 over Finite Field of size 1000041437 + Elliptic Curve defined by y^2 = x^3 + 665393686*x + 948152000 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 572311614*x + 178583984 over Finite Field of size 999969307 sage: set_verbose(2) + sage: set_random_seed(1337) sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): ....: print(E) + verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 + Elliptic Curve defined by y^2 = x^3 + 265977778*x + 120868502 over Finite Field of size 1000041437 verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 - Elliptic Curve defined by y^2 = x^3 + 276281332*x + 798205466 over Finite Field of size 999969307 - Elliptic Curve defined by y^2 = x^3 + 515293837*x + 588191363 over Finite Field of size 999969307 - Elliptic Curve defined by y^2 = x^3 + 684893162*x + 679980762 over Finite Field of size 999969307 - sage: set_verbose(3) + Elliptic Curve defined by y^2 = x^3 + 665393686*x + 948152000 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 572311614*x + 178583984 over Finite Field of size 999969307 + sage: set_verbose(4) + sage: set_random_seed(1337) sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): ....: print(E) - verbose 3 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Considering 1th valid prime 19 - verbose 3 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Considering 2th valid prime 23 - verbose 3 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Considering 3th valid prime 29 - verbose 3 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-667 + verbose 4 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-19 + ... + verbose 4 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-163 + verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 + Elliptic Curve defined by y^2 = x^3 + 265977778*x + 120868502 over Finite Field of size 1000041437 + verbose 4 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-179 + ... + verbose 4 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-667 verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 - Elliptic Curve defined by y^2 = x^3 + 620755650*x + 816246309 over Finite Field of size 999969307 - Elliptic Curve defined by y^2 = x^3 + 651416562*x + 875284469 over Finite Field of size 999969307 - Elliptic Curve defined by y^2 = x^3 + 149283707*x + 72132991 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 665393686*x + 948152000 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 572311614*x + 178583984 over Finite Field of size 999969307 TESTS:: + sage: set_verbose(0) sage: for N in prime_range(3, 100): ....: E = next(EllipticCurve_with_prime_order(N)) - ....: assert E.order() == N + ....: assert E.has_order(N) sage: N = 113 sage: for _, E in zip(range(30), EllipticCurve_with_prime_order(N)): - ....: assert E.order() == N + ....: assert E.has_order(N) sage: N = 15175980689839334471 sage: E = next(EllipticCurve_with_prime_order(N)) - sage: E.order() == N + sage: E.has_order(N) True sage: N = next_prime(123456789) sage: E = next(EllipticCurve_with_prime_order(N)) - sage: E.order() == N + sage: E.has_order(N) True sage: N = 123456789 @@ -3018,25 +3039,14 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): Traceback (most recent call last): ... ValueError: input order is not prime - - .. NOTE:: - - Depending on the input, this function may run for a *very* long time. - This algorithm consists of multiple "search rounds" for a suitable - discriminant `D`. We expect this algorithm to terminate after a number - of rounds that is polynomial in `\log\log N`. In practice (cf. Section - 5), this number is usually 1. - - ALGORITHM: [BS2007]_, Algorithm 2.2 """ import itertools - from sage.misc.misc_c import prod - from sage.misc.verbose import verbose - from sage.combinat import gray_codes - from sage.sets.primes import Primes from sage.arith.misc import is_prime, legendre_symbol + from sage.misc.verbose import verbose from sage.quadratic_forms.binary_qf import BinaryQF + from sage.rings.fast_arith import prime_range from sage.schemes.elliptic_curves.cm import hilbert_class_polynomial + from sage.sets.primes import Primes if not is_prime(N): raise ValueError("input order is not prime") @@ -3050,49 +3060,54 @@ def EllipticCurve_with_prime_order(N, max_D=10**8): # The algorithm considers smooth discriminants `D`, sorted by their largest # prime factor. We expect this algorithm to terminate after a number of # rounds that is polynomial in loglog N. - S = [] + # We start with small primes directly to accelerate the search + S = [(-p if p >> 1 & 1 else p) for p in prime_range(3, min(1000, 4 * N)) + if legendre_symbol(N, p) == 1] - def recur(idx, bound): + def recur(bound): """ - This function returns all subsets of S[idx:] (at the time of running) - with product of absolute value <= bound. It's extremely quick since S - is basically prime_range(N) while bound is basically 4 * N. + This function returns an iterator of all numbers with absolute value + not exceeding ``bound`` expressable as product of distinct elements in + ``S`` in ascending order. """ - if idx >= len(S): - return - yield [] - for nxt in range(idx + 1, len(S)): - if abs(S[nxt]) <= bound: - for r in recur(nxt, bound // abs(S[nxt])): - yield [S[nxt]] + r - else: - break + import heapq + hq = [(1, 1, -1)] + while len(hq): + abs_n, n, idx = heapq.heappop(hq) + yield n + for nxt in range(idx + 1, len(S)): + if abs_n * abs(S[nxt]) <= bound: + heapq.heappush(hq, (abs_n * abs(S[nxt]), n * S[nxt], nxt)) + else: + break - for p in Primes(): - if p == 2: - continue + for p in itertools.chain([1], Primes()): + # We add p = 1 to process the small primes + if p != 1: + if p < 1000: + continue - if legendre_symbol(N, p) != 1: - continue + if legendre_symbol(N, p) != 1: + continue - verbose(f"Considering {len(S) + 1}th valid prime {p}", level=3) + # later we need x^2 + (-D)y^2 = 4N, and since y = 0 has no solution, we need + # p = |p_star| <= |-D| <= 4N + if p > 4 * N: + break + + verbose(f"Considering {len(S) + 1}th valid prime {p}", level=3) # Equivalent to p* = (-1)^((p - 1) / 2) * p in [BS2007]_ page 5. p_star = -p if p >> 1 & 1 else p - # later we need x^2 + (-D)y^2 = 4N, and since y = 0 has no solution, we need - # p = |p_star| <= |-D| <= 4N - if p > 4 * N: - break - - for e in recur(0, 4 * N // p): - D = p_star * prod(e) + for e in recur(4 * N // p): + D = p_star * e assert abs(D) <= 4 * N if D % 8 != 5 or D >= 0: continue - verbose(f"Testing {D=}", level=3) + verbose(f"Testing {D=}", level=4) Q = BinaryQF([1, 0, -D]) sol = Q.solve_integer(4 * N, algorithm='cornacchia') From 94c4866d83e7035cb3352b2eeece13e46932afbf Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Tue, 10 Sep 2024 11:59:06 +0100 Subject: [PATCH 16/29] remove bithack --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 9d2666ac989..49779804c1d 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -3084,7 +3084,7 @@ def recur(bound): for p in itertools.chain([1], Primes()): # We add p = 1 to process the small primes if p != 1: - if p < 1000: + if p < S[-1]: continue if legendre_symbol(N, p) != 1: @@ -3098,7 +3098,7 @@ def recur(bound): verbose(f"Considering {len(S) + 1}th valid prime {p}", level=3) # Equivalent to p* = (-1)^((p - 1) / 2) * p in [BS2007]_ page 5. - p_star = -p if p >> 1 & 1 else p + p_star = -p if p % 4 == 3 else p for e in recur(4 * N // p): D = p_star * e From 459433c7a92eab7fe6420bccfccfe154553d974f Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Tue, 10 Sep 2024 12:10:28 +0100 Subject: [PATCH 17/29] add note about magic constant 1000 --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 49779804c1d..889ff1aa6e1 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -3061,6 +3061,8 @@ def EllipticCurve_with_prime_order(N): # prime factor. We expect this algorithm to terminate after a number of # rounds that is polynomial in loglog N. # We start with small primes directly to accelerate the search + # Note: 1000 is a magic constant, it's just fast enough to compute without + # sacrificing much speed S = [(-p if p >> 1 & 1 else p) for p in prime_range(3, min(1000, 4 * N)) if legendre_symbol(N, p) == 1] From 5fb4bc4ef2bb5090448056ec32b3df3518adeec0 Mon Sep 17 00:00:00 2001 From: Martin Grenouilloux Date: Tue, 10 Sep 2024 16:20:09 +0200 Subject: [PATCH 18/29] Coding style, naming and skip p=2 --- .../elliptic_curves/ell_finite_field.py | 96 ++++++++++--------- 1 file changed, 50 insertions(+), 46 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 889ff1aa6e1..e1d960d9fc1 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -1405,14 +1405,15 @@ def has_order(self, value, num_checks=8): - Gareth Ma (2024-01-21): Fix bug for small curves """ - # This method does *not* use the cached value of `_order` even when available + # This method does *not* use the cached value of `_order` even when + # available. q = self.base_field().order() a, b = Hasse_bounds(q, 1) if not a <= value <= b: return False - # For really small values, the random tests are too weak to detect wrong orders - # So we go with computing directly instead. + # For really small values, the random tests are too weak to detect wrong + # orders So we go with computing directly instead. if q <= 100: return self.order() == value @@ -1498,10 +1499,9 @@ def set_order(self, value, *, check=True, num_checks=8): ... ValueError: Elliptic Curve defined by y^2 = x^3 + 1 over Finite Field of size 7 does not have order 1000 - It is also very likely an error to pass a value which is not - the actual order of this curve. How unlikely is determined by - ``num_checks``, the factorization of the actual order, and the - actual group structure:: + It is also very likely an error to pass a value which is not the actual + order of this curve. How unlikely is determined by ``num_checks``, the + factorization of the actual order, and the actual group structure:: sage: E = EllipticCurve(GF(1009), [0, 1]) # This curve has order 948 sage: E.set_order(947) @@ -1509,8 +1509,8 @@ def set_order(self, value, *, check=True, num_checks=8): ... ValueError: Elliptic Curve defined by y^2 = x^3 + 1 over Finite Field of size 1009 does not have order 947 - For curves over small finite fields, the order is cheap to compute, so it is computed - directly and compared:: + For curves over small finite fields, the order is cheap to compute, so + it is computed directly and compared:: sage: E = EllipticCurve(GF(7), [0, 1]) # This curve has order 12 sage: E.set_order(11) @@ -1520,8 +1520,8 @@ def set_order(self, value, *, check=True, num_checks=8): TESTS: - The previous version's random tests are not strong enough. In particular, the following used - to work:: + The previous version's random tests are not strong enough. In particular, + the following used to work:: sage: E = EllipticCurve(GF(2), [0, 0, 1, 1, 1]) # This curve has order 1 sage: E.set_order(3) @@ -1539,8 +1539,8 @@ def set_order(self, value, *, check=True, num_checks=8): sage: E.order() 12 - .. TODO:: Add provable correctness check by computing the abelian group structure and - comparing. + .. TODO:: Add provable correctness check by computing the abelian group + structure and comparing. AUTHORS: @@ -2896,10 +2896,10 @@ def EllipticCurve_with_prime_order(N): True The returned curves are sometimes random because - :meth:`~sage.schemes.elliptic_curves.ell_finite_field.EllipticCurve_finite_field.twists` is - not deterministic. However, it's always isomorphic:: + :meth:`~sage.schemes.elliptic_curves.ell_finite_field.EllipticCurve_finite_field.twists` + is not deterministic. However, it's always isomorphic:: - sage: E = next(EllipticCurve_with_prime_order(23)); E # random + sage: E = next(EllipticCurve_with_prime_order(23)); E # random Elliptic Curve defined by y^2 = x^3 + 12*x + 6 over Finite Field of size 17 sage: E.is_isomorphic(EllipticCurve(GF(17), [3, 5])) True @@ -2939,7 +2939,8 @@ def EllipticCurve_with_prime_order(N): sage: EllipticCurve(GF(7), [0, 5]).order() 7 - However, experimentally it returns many of them. Here it returns all of them:: + However, experimentally it returns many of them. Here it returns all of + them:: sage: N = 23 sage: set_random_seed(1337) # as the function returns random twists of curves @@ -2960,12 +2961,12 @@ def EllipticCurve_with_prime_order(N): ....: if E.has_order(N): ....: assert any(E.is_isomorphic(E_) for E_ in curves) - The algorithm is efficient for small ``N`` due to the low number of - suitable discriminants (see the ``recur`` internal function of the code for - details). The following runs in less than a second:: + The algorithm is efficient for small ``N`` due to the low number of suitable + discriminants (see the ``abs_products_under`` internal function of the code + for details). The following runs in less than a second:: sage: len(list(EllipticCurve_with_prime_order(next_prime(5000)))) - 534 + 945 There is different verbose data for level `2` to `4`, though level `3` rarely logs anything (it logs when a new prime `p` is added to the @@ -3024,21 +3025,26 @@ def EllipticCurve_with_prime_order(N): sage: E.has_order(N) True + sage: E = next(EllipticCurve_with_prime_order(2)) + Traceback (most recent call last): + ... + ValueError: input order is not an odd prime + sage: N = 123456789 sage: E = next(EllipticCurve_with_prime_order(N)) Traceback (most recent call last): ... - ValueError: input order is not prime + ValueError: input order is not an odd prime sage: E = next(EllipticCurve_with_prime_order(0)) Traceback (most recent call last): ... - ValueError: input order is not prime + ValueError: input order is not an odd prime sage: E = next(EllipticCurve_with_prime_order(-7)) Traceback (most recent call last): ... - ValueError: input order is not prime + ValueError: input order is not an odd prime """ import itertools from sage.arith.misc import is_prime, legendre_symbol @@ -3048,29 +3054,23 @@ def EllipticCurve_with_prime_order(N): from sage.schemes.elliptic_curves.cm import hilbert_class_polynomial from sage.sets.primes import Primes - if not is_prime(N): - raise ValueError("input order is not prime") - - # if N < 1000: - # # Log how long this algorithm will take - # from sage.rings.fast_arith import prime_range - # num = sum(1 for p in prime_range(3, 4 * N) if legendre_symbol(N, p) == 1) - # verbose(f"Total work to enumerate curves with order {N}: 2^{num}", level=2) + if not is_prime(N) or N < 3: + raise ValueError("input order is not an odd prime") # The algorithm considers smooth discriminants `D`, sorted by their largest # prime factor. We expect this algorithm to terminate after a number of # rounds that is polynomial in loglog N. - # We start with small primes directly to accelerate the search + # We start with small primes directly to accelerate the search. # Note: 1000 is a magic constant, it's just fast enough to compute without - # sacrificing much speed + # sacrificing much speed. S = [(-p if p >> 1 & 1 else p) for p in prime_range(3, min(1000, 4 * N)) if legendre_symbol(N, p) == 1] - def recur(bound): + def abs_products_under(bound): """ - This function returns an iterator of all numbers with absolute value - not exceeding ``bound`` expressable as product of distinct elements in - ``S`` in ascending order. + This function returns an iterator of all numbers with absolute value not + exceeding ``bound`` expressable as product of distinct elements in ``S`` + in ascending order. """ import heapq hq = [(1, 1, -1)] @@ -3083,8 +3083,9 @@ def recur(bound): else: break + # We add p = 1 to process the small primes. for p in itertools.chain([1], Primes()): - # We add p = 1 to process the small primes + if p == 2: continue if p != 1: if p < S[-1]: continue @@ -3092,17 +3093,18 @@ def recur(bound): if legendre_symbol(N, p) != 1: continue - # later we need x^2 + (-D)y^2 = 4N, and since y = 0 has no solution, we need - # p = |p_star| <= |-D| <= 4N + # Later we need x^2 + (-D)y^2 = 4N, and since y = 0 has no solution, + # we need p = |p_star| <= |-D| <= 4N. This is a stopping condition + # for the algorithm. if p > 4 * N: break verbose(f"Considering {len(S) + 1}th valid prime {p}", level=3) - # Equivalent to p* = (-1)^((p - 1) / 2) * p in [BS2007]_ page 5. - p_star = -p if p % 4 == 3 else p + # Equivalent to p_star = (-1)^((p - 1) / 2) * p in [BS2007]_ page 5. + p_star = -p if p >> 1 & 1 else p - for e in recur(4 * N // p): + for e in abs_products_under(4 * N // p): D = p_star * e assert abs(D) <= 4 * N @@ -3119,15 +3121,17 @@ def recur(bound): x, _ = sol for p_i in [N + 1 - x, N + 1 + x]: if is_prime(p_i): - verbose(f"Computing the Hilbert class polynomial H_{D}", level=2) + verbose(f"Computing the Hilbert class polynomial H_{D}", + level=2) H = hilbert_class_polynomial(D) K = GF(p_i) for j0 in H.roots(ring=K, multiplicities=False): E = EllipticCurve(K, j=j0) # `E.twists()` also contains E. for Et in E.twists(): - # num_checks=1 is sufficient for prime order + # `num_checks=1` is sufficient for prime order. if Et.has_order(N, num_checks=1): yield Et + # Extending our prime list and continuing onto the next round. S.append(p_star) From d7d681adc3e8f24cf042f234739f6472d9c67e1d Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Tue, 10 Sep 2024 16:07:49 +0100 Subject: [PATCH 19/29] fix things --- .../elliptic_curves/ell_finite_field.py | 69 +++++++++++-------- 1 file changed, 41 insertions(+), 28 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index e1d960d9fc1..59835b784d9 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -16,7 +16,7 @@ - Lorenz Panny (2023): ``special_supersingular_curve()`` -- Martin Grenouilloux (2024): ``EllipticCurve_with_prime_order()`` +- Martin Grenouilloux, Gareth Ma (2024-09): ``EllipticCurve_with_prime_order()`` """ # **************************************************************************** @@ -2923,7 +2923,7 @@ def EllipticCurve_with_prime_order(N): :: sage: N = next_prime(2^256) - sage: E = next(EllipticCurve_with_prime_order(N)); E # random + sage: E = next(EllipticCurve_with_prime_order(N)); E # random Elliptic Curve defined by y^2 = x^3 + 6056521267553273205988520276135607487700943205131813669424576873701361709521*x + 86942739955486781674010637133214195706465136689012129911736706024465988573567 over Finite Field of size 115792089237316195423570985008687907853847329310253429036565151476471048389761 @@ -2966,7 +2966,9 @@ def EllipticCurve_with_prime_order(N): for details). The following runs in less than a second:: sage: len(list(EllipticCurve_with_prime_order(next_prime(5000)))) - 945 + 534 + sage: len(list(EllipticCurve_with_prime_order(next_prime(50000)))) # long time (6s) + 3841 There is different verbose data for level `2` to `4`, though level `3` rarely logs anything (it logs when a new prime `p` is added to the @@ -3006,6 +3008,11 @@ def EllipticCurve_with_prime_order(N): TESTS:: + sage: list(EllipticCurve_with_prime_order(2)) + [Elliptic Curve defined by y^2 + x*y + y = x^3 + 1 over Finite Field of size 2, + Elliptic Curve defined by y^2 = x^3 + 2*x^2 + 2 over Finite Field of size 3, + Elliptic Curve defined by y^2 = x^3 + 2*x over Finite Field of size 5] + sage: set_verbose(0) sage: for N in prime_range(3, 100): ....: E = next(EllipticCurve_with_prime_order(N)) @@ -3025,26 +3032,25 @@ def EllipticCurve_with_prime_order(N): sage: E.has_order(N) True - sage: E = next(EllipticCurve_with_prime_order(2)) - Traceback (most recent call last): - ... - ValueError: input order is not an odd prime - sage: N = 123456789 sage: E = next(EllipticCurve_with_prime_order(N)) Traceback (most recent call last): ... - ValueError: input order is not an odd prime + ValueError: input order is not a prime sage: E = next(EllipticCurve_with_prime_order(0)) Traceback (most recent call last): ... - ValueError: input order is not an odd prime + ValueError: input order is not a prime sage: E = next(EllipticCurve_with_prime_order(-7)) Traceback (most recent call last): ... - ValueError: input order is not an odd prime + ValueError: input order is not a prime + + AUTHORS: + + - Martin Grenouilloux, Gareth Ma (2024-09): initial implementation """ import itertools from sage.arith.misc import is_prime, legendre_symbol @@ -3054,16 +3060,22 @@ def EllipticCurve_with_prime_order(N): from sage.schemes.elliptic_curves.cm import hilbert_class_polynomial from sage.sets.primes import Primes - if not is_prime(N) or N < 3: - raise ValueError("input order is not an odd prime") + if not is_prime(N): + raise ValueError("input order is not a prime") - # The algorithm considers smooth discriminants `D`, sorted by their largest - # prime factor. We expect this algorithm to terminate after a number of - # rounds that is polynomial in loglog N. - # We start with small primes directly to accelerate the search. - # Note: 1000 is a magic constant, it's just fast enough to compute without + if N == 2: + yield from [ + EllipticCurve(GF(2), [1, 0, 1, 0, 1]), + EllipticCurve(GF(3), [0, 2, 0, 0, 2]), + EllipticCurve(GF(5), [2, 0]) + ] + return + + # We start with small primes directly to accelerate the search. Note that + # 1000 is a magic constant, it's just fast enough to compute without # sacrificing much speed. - S = [(-p if p >> 1 & 1 else p) for p in prime_range(3, min(1000, 4 * N)) + # The if-then-else term is (-1)^((p - 1) / 2) * p in [BS2007]_ page 5. + S = [(-p if p % 4 == 3 else p) for p in prime_range(3, min(1000, 4 * N)) if legendre_symbol(N, p) == 1] def abs_products_under(bound): @@ -3085,26 +3097,26 @@ def abs_products_under(bound): # We add p = 1 to process the small primes. for p in itertools.chain([1], Primes()): - if p == 2: continue if p != 1: - if p < S[-1]: + if p < abs(S[-1]): continue if legendre_symbol(N, p) != 1: continue - # Later we need x^2 + (-D)y^2 = 4N, and since y = 0 has no solution, - # we need p = |p_star| <= |-D| <= 4N. This is a stopping condition - # for the algorithm. + # Later we need x^2 + (-D)y^2 = 4N, and since y = 0 has no + # solution, we need p = |p_star| <= |-D| <= 4N. This is a stopping + # condition for the algorithm. if p > 4 * N: break verbose(f"Considering {len(S) + 1}th valid prime {p}", level=3) - # Equivalent to p_star = (-1)^((p - 1) / 2) * p in [BS2007]_ page 5. - p_star = -p if p >> 1 & 1 else p + p_star = -p if p % 4 == 3 else p for e in abs_products_under(4 * N // p): + # According to the paper, the expected minimum D to work is + # O(log(N)^2) D = p_star * e assert abs(D) <= 4 * N @@ -3133,5 +3145,6 @@ def abs_products_under(bound): if Et.has_order(N, num_checks=1): yield Et - # Extending our prime list and continuing onto the next round. - S.append(p_star) + if p != 1: + # Extending our prime list and continuing onto the next round. + S.append(p_star) From 2854e94fcdc447f22b44d4ccc27e0deba906cc0d Mon Sep 17 00:00:00 2001 From: grhkm21 <83517584+grhkm21@users.noreply.github.com> Date: Wed, 11 Sep 2024 18:48:16 +0100 Subject: [PATCH 20/29] Apply suggestions from code review Co-authored-by: Vincent Macri --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 59835b784d9..b249a62885f 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2873,7 +2873,7 @@ def EllipticCurve_with_prime_order(N): - ``N`` -- integer; the order for which we seek an elliptic curve. Must be a prime number. - OUTPUT: an iterator of elliptic curves `E/\mathbb F_p` of order ``N`` + OUTPUT: an iterator of (some) elliptic curves `E/\mathbb F_p` of order ``N`` .. NOTE:: @@ -2883,11 +2883,11 @@ def EllipticCurve_with_prime_order(N): of rounds that is polynomial in `\log\log N`. In practice (cf. Section 5), this number is usually 1. - ALGORITHM: [BS2007]_, Algorithm 2.2 + ALGORITHM: Based on [BS2007]_, Algorithm 2.2 EXAMPLES:: - sage: N = next_prime(int.from_bytes(b'sagemath', 'big')) + sage: N = 8314040072427107567 sage: E = next(EllipticCurve_with_prime_order(N)) sage: E Elliptic Curve defined by y^2 = x^3 + 4757897140353078952*x + 1841350074072114366 @@ -2913,7 +2913,7 @@ def EllipticCurve_with_prime_order(N): It works for large primes:: - sage: N = 0x6cbc824032974516623e732462f4b74b56c4ffbd984380d9 + sage: N = 2666207849820848272386538889527600954292544013630953455833 sage: E = next(EllipticCurve_with_prime_order(N)); E Elliptic Curve defined by y^2 = x^3 + 2666207849820848272386538889427721639173508298483739490459*x + 77986137112576 over Finite Field of size 2666207849820848272386538889427721639173508298487130585243 From 4b9e197ec179539e2fdf17b73eda45d120d474ff Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Wed, 11 Sep 2024 21:23:18 +0100 Subject: [PATCH 21/29] apply review changes (explain ALGORITHM changes) --- .../elliptic_curves/ell_finite_field.py | 42 ++++++++++++++----- 1 file changed, 32 insertions(+), 10 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 59835b784d9..43366a793d7 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -1414,7 +1414,9 @@ def has_order(self, value, num_checks=8): # For really small values, the random tests are too weak to detect wrong # orders So we go with computing directly instead. - if q <= 100: + # In #38341, the bound has been increased to a large value (2^64), but + # it should be decreased (to ~100) after bug #38617 is fixed. + if q <= 2**64: return self.order() == value # This might be slow @@ -2875,15 +2877,35 @@ def EllipticCurve_with_prime_order(N): OUTPUT: an iterator of elliptic curves `E/\mathbb F_p` of order ``N`` - .. NOTE:: - - Depending on the input, this function may run for a *very* long time. - This algorithm consists of multiple "search rounds" for a suitable - discriminant `D`. We expect this algorithm to terminate after a number - of rounds that is polynomial in `\log\log N`. In practice (cf. Section - 5), this number is usually 1. + ALGORITHM: - ALGORITHM: [BS2007]_, Algorithm 2.2 + Our algorithm is based on [BS2007]_, Algorithm 2.2, but we deviate for + several key steps. Firstly, the authors in the paper perform the search for + a suitable `D` *incrementally*, by enlarging the table `S` by `log(N)`-size + interval of primes `p` and testing all products of distinct primes `p` (or + rather `p^*`). We find this difficult to implement without testing + duplicate `D`s, so we instead enlarge the table one prime at a time + (effectively replacing `[r\log(N), (r + 1)\log(N)]` in the paper by `[r, + r]`). To compensate for the speed loss, we begin the algorithm by + prefilling `S` with the primes below `1000` (satisfying quadratic + reciprocity properties). The constant `1000` is determined experimentally + to be fast for many purposes, and for most `N` we tested we are able to + find a suitable small `D` without increasing the size of `S`. + + The paper also doesn't specify how to enumerate such `D`s, which recall + should be product of distinct values in the table `S`. We implement this + with a priority queue (min heap), which also allows us to search for the + suitable `D`s in increasing (absolute value) order. This is suitable for + the algorithm because smaller `D` means the Hilbert class polynomial is + computed quicker. + + Finally, to avoid repeatedly testing the same `D`s, we require the latest + prime to be added to the table to be included as a factor of `D` (see code + for more explanation). As we need to find integers `x, y` such that `x^2 + + (-D)y^2 = 4N` with `D < 0` and `N` prime, we actually need `|D| \leq 4N`, + so we terminate the algorithm when the primes in the table are larger than + that bound. This makes the iterator return all curves it can find in finite + time :) EXAMPLES:: @@ -2963,7 +2985,7 @@ def EllipticCurve_with_prime_order(N): The algorithm is efficient for small ``N`` due to the low number of suitable discriminants (see the ``abs_products_under`` internal function of the code - for details). The following runs in less than a second:: + for details):: sage: len(list(EllipticCurve_with_prime_order(next_prime(5000)))) 534 From 7941d8b458b674eadf0f173cf8f226977f8bd00c Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Wed, 11 Sep 2024 21:27:59 +0100 Subject: [PATCH 22/29] linter --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 43366a793d7..e2cba07647a 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2753,7 +2753,7 @@ def EllipticCurve_with_order(m, *, D=None): Return an iterator for elliptic curves over finite fields with the given order. The curves are computed using the Complex Multiplication (CM) method. - A :sage:`~sage.structure.factorization.Factorization` can be passed for ``m``, in which case + A :class:`~sage.structure.factorization.Factorization` can be passed for ``m``, in which case the algorithm is more efficient. If ``D`` is specified, it is used as the discriminant. From 4840e9a2f09048b482e84c2f20deb2be49126cd5 Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Wed, 11 Sep 2024 21:58:09 +0100 Subject: [PATCH 23/29] update test to use larger valued curves --- .../elliptic_curves/ell_finite_field.py | 30 +++++++++---------- 1 file changed, 15 insertions(+), 15 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 8ddbdfdb8c2..47c630650a8 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -1390,13 +1390,16 @@ def has_order(self, value, num_checks=8): sage: E.has_order(N + 1) False - This demonstrates the bug mentioned in the NOTE above. The return value - should be ``False`` after :issue:`38617` is fixed:: + This demonstrates the bug mentioned in the NOTE above. The last return + value should be ``False`` after :issue:`38617` is fixed:: - sage: E = EllipticCurve(GF(127), [0, 1]) - sage: E.order() - 108 - sage: E.has_order(126) + sage: E = EllipticCurve(GF(5443568676570036275321323), [0, 13]) + sage: N = 2333145661241 + sage: E.order() == N^2 + True + sage: E.has_order(N^2) + True + sage: E.has_order(N^2 + N) True AUTHORS: @@ -2837,7 +2840,7 @@ def find_q(m, m4_fac, D): m_val = m if D is None: - Ds = (D for D in range(-4 * m_val, 0) if D % 4 in [0, 1]) + Ds = (D for D in range(-1, -4 * m_val - 1, -1) if D % 4 in [0, 1]) else: assert D < 0 and D % 4 in [0, 1] Ds = [D] @@ -2849,20 +2852,16 @@ def find_q(m, m4_fac, D): continue H = hilbert_class_polynomial(D) - K = GF(q) - roots = H.roots(ring=K) - for j0, _ in roots: + for j0 in H.roots(ring=GF(q), multiplicities=False): E = EllipticCurve(j=j0) for Et in E.twists(): if any(Et.is_isomorphic(E) for E in seen): continue - try: - # This tests whether the curve has given order - Et.set_order(m_val) + # This tests whether the curve has given order + if Et.has_order(m_val): + Et.set_order(m_val, check=False) seen.add(Et) yield Et - except ValueError: - pass def EllipticCurve_with_prime_order(N): r""" @@ -3168,6 +3167,7 @@ def abs_products_under(bound): for Et in E.twists(): # `num_checks=1` is sufficient for prime order. if Et.has_order(N, num_checks=1): + Et.set_order(N, check=False) yield Et if p != 1: From e258c18eac18b7231549e0a3a99c054f93c0dda8 Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Wed, 11 Sep 2024 22:20:47 +0100 Subject: [PATCH 24/29] fix tests --- .../elliptic_curves/ell_finite_field.py | 28 +++++++++---------- 1 file changed, 14 insertions(+), 14 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 47c630650a8..841f90f9931 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -3003,32 +3003,32 @@ def EllipticCurve_with_prime_order(N): sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): ....: print(E) Elliptic Curve defined by y^2 = x^3 + 265977778*x + 120868502 over Finite Field of size 1000041437 - Elliptic Curve defined by y^2 = x^3 + 665393686*x + 948152000 over Finite Field of size 999969307 - Elliptic Curve defined by y^2 = x^3 + 572311614*x + 178583984 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 689795416*x + 188156157 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 999178436*x + 900579394 over Finite Field of size 999969307 sage: set_verbose(2) sage: set_random_seed(1337) sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): ....: print(E) - verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 + verbose 2 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 Elliptic Curve defined by y^2 = x^3 + 265977778*x + 120868502 over Finite Field of size 1000041437 - verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 - Elliptic Curve defined by y^2 = x^3 + 665393686*x + 948152000 over Finite Field of size 999969307 - Elliptic Curve defined by y^2 = x^3 + 572311614*x + 178583984 over Finite Field of size 999969307 + verbose 2 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 + Elliptic Curve defined by y^2 = x^3 + 689795416*x + 188156157 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 999178436*x + 900579394 over Finite Field of size 999969307 sage: set_verbose(4) sage: set_random_seed(1337) sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): ....: print(E) - verbose 4 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-19 + verbose 4 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-19 ... - verbose 4 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-163 - verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 + verbose 4 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-163 + verbose 2 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 Elliptic Curve defined by y^2 = x^3 + 265977778*x + 120868502 over Finite Field of size 1000041437 - verbose 4 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-179 + verbose 4 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-179 ... - verbose 4 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-667 - verbose 2 (2865: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 - Elliptic Curve defined by y^2 = x^3 + 665393686*x + 948152000 over Finite Field of size 999969307 - Elliptic Curve defined by y^2 = x^3 + 572311614*x + 178583984 over Finite Field of size 999969307 + verbose 4 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-667 + verbose 2 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 + Elliptic Curve defined by y^2 = x^3 + 689795416*x + 188156157 over Finite Field of size 999969307 + Elliptic Curve defined by y^2 = x^3 + 999178436*x + 900579394 over Finite Field of size 999969307 TESTS:: From c9e33ea7d9a3959ee3066044ebf334e8c4f36805 Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Wed, 11 Sep 2024 22:29:07 +0100 Subject: [PATCH 25/29] apply review changes --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 841f90f9931..50b0fb4727c 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -1436,6 +1436,7 @@ def has_order(self, value, num_checks=8): if not (value * G).is_zero(): return False + self._order = value return True def set_order(self, value, *, check=True, num_checks=8): @@ -2859,7 +2860,6 @@ def find_q(m, m4_fac, D): continue # This tests whether the curve has given order if Et.has_order(m_val): - Et.set_order(m_val, check=False) seen.add(Et) yield Et @@ -2891,7 +2891,6 @@ def EllipticCurve_with_prime_order(N): to be fast for many purposes, and for most `N` we tested we are able to find a suitable small `D` without increasing the size of `S`. -<<<<<<< HEAD The paper also doesn't specify how to enumerate such `D`s, which recall should be product of distinct values in the table `S`. We implement this with a priority queue (min heap), which also allows us to search for the @@ -2905,7 +2904,7 @@ def EllipticCurve_with_prime_order(N): (-D)y^2 = 4N` with `D < 0` and `N` prime, we actually need `|D| \leq 4N`, so we terminate the algorithm when the primes in the table are larger than that bound. This makes the iterator return all curves it can find in finite - time :) + time. ALGORITHM: Based on [BS2007]_, Algorithm 2.2 @@ -3167,7 +3166,6 @@ def abs_products_under(bound): for Et in E.twists(): # `num_checks=1` is sufficient for prime order. if Et.has_order(N, num_checks=1): - Et.set_order(N, check=False) yield Et if p != 1: From 41cb6e912e3c0b668aaabc41b39d43e18bcc5ba3 Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Wed, 11 Sep 2024 22:48:18 +0100 Subject: [PATCH 26/29] use _order in has_order when possible --- .../elliptic_curves/ell_finite_field.py | 19 ++++++++++--------- 1 file changed, 10 insertions(+), 9 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 50b0fb4727c..87fa0ec67b8 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -1399,6 +1399,7 @@ def has_order(self, value, num_checks=8): True sage: E.has_order(N^2) True + sage: del E._order sage: E.has_order(N^2 + N) True @@ -1419,7 +1420,7 @@ def has_order(self, value, num_checks=8): # orders So we go with computing directly instead. # In #38341, the bound has been increased to a large value (2^64), but # it should be decreased (to ~100) after bug #38617 is fixed. - if q <= 2**64: + if q <= 2**64 or hasattr(self, "_order"): return self.order() == value # This might be slow @@ -3008,24 +3009,24 @@ def EllipticCurve_with_prime_order(N): sage: set_random_seed(1337) sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): ....: print(E) - verbose 2 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 + verbose 2 (...: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 Elliptic Curve defined by y^2 = x^3 + 265977778*x + 120868502 over Finite Field of size 1000041437 - verbose 2 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 + verbose 2 (...: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 Elliptic Curve defined by y^2 = x^3 + 689795416*x + 188156157 over Finite Field of size 999969307 Elliptic Curve defined by y^2 = x^3 + 999178436*x + 900579394 over Finite Field of size 999969307 sage: set_verbose(4) sage: set_random_seed(1337) sage: for _, E in zip(range(3), EllipticCurve_with_prime_order(10^9 + 7)): ....: print(E) - verbose 4 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-19 + verbose 4 (...: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-19 ... - verbose 4 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-163 - verbose 2 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 + verbose 4 (...: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-163 + verbose 2 (...: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-163 Elliptic Curve defined by y^2 = x^3 + 265977778*x + 120868502 over Finite Field of size 1000041437 - verbose 4 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-179 + verbose 4 (...: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-179 ... - verbose 4 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-667 - verbose 2 (2866: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 + verbose 4 (...: ell_finite_field.py, EllipticCurve_with_prime_order) Testing D=-667 + verbose 2 (...: ell_finite_field.py, EllipticCurve_with_prime_order) Computing the Hilbert class polynomial H_-667 Elliptic Curve defined by y^2 = x^3 + 689795416*x + 188156157 over Finite Field of size 999969307 Elliptic Curve defined by y^2 = x^3 + 999178436*x + 900579394 over Finite Field of size 999969307 From 228a2bf301601baef9cbec59650f1deb96593c65 Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Wed, 11 Sep 2024 22:50:18 +0100 Subject: [PATCH 27/29] revert --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index 87fa0ec67b8..b877f605c77 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -1437,7 +1437,9 @@ def has_order(self, value, num_checks=8): if not (value * G).is_zero(): return False - self._order = value + # TODO: uncomment this and remove the line in `set_order` after 38617 is fixed. + # self._order = value + return True def set_order(self, value, *, check=True, num_checks=8): From 163191fdb7f63e62d19bb70ca56a5662c476a26d Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Wed, 11 Sep 2024 22:54:10 +0100 Subject: [PATCH 28/29] apply review changes --- src/sage/schemes/elliptic_curves/ell_finite_field.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index b877f605c77..b1f77da41e8 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -2863,6 +2863,8 @@ def find_q(m, m4_fac, D): continue # This tests whether the curve has given order if Et.has_order(m_val): + # TODO: remove after 38617 + Et.set_order(m_val, check=False) seen.add(Et) yield Et @@ -3167,8 +3169,10 @@ def abs_products_under(bound): E = EllipticCurve(K, j=j0) # `E.twists()` also contains E. for Et in E.twists(): - # `num_checks=1` is sufficient for prime order. + # `num_checks=1` is sufficient for prime order if Et.has_order(N, num_checks=1): + # TODO: remove after 38617 + Et.set_order(N, check=False) yield Et if p != 1: From b60e81f323703f3d641cf5d78c519fd8701f1449 Mon Sep 17 00:00:00 2001 From: Gareth Ma Date: Thu, 12 Sep 2024 18:32:04 +0100 Subject: [PATCH 29/29] apply review changes --- .../schemes/elliptic_curves/ell_finite_field.py | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) diff --git a/src/sage/schemes/elliptic_curves/ell_finite_field.py b/src/sage/schemes/elliptic_curves/ell_finite_field.py index b1f77da41e8..65b49cfe1bd 100755 --- a/src/sage/schemes/elliptic_curves/ell_finite_field.py +++ b/src/sage/schemes/elliptic_curves/ell_finite_field.py @@ -1353,7 +1353,9 @@ def has_order(self, value, num_checks=8): - ``num_checks``-- integer (default: `8`); the number of times to check whether ``value`` times a random point on this curve equals the - identity + identity. If ``value`` is a prime and the curve is over a field of + order at least `5`, it is sufficient to pass in ``num_checks=1`` - + see the examples below. .. NOTE:: @@ -1403,14 +1405,20 @@ def has_order(self, value, num_checks=8): sage: E.has_order(N^2 + N) True + Due to the nature of the algorithm (testing multiple of points) and the Hasse-Weil bound, we see that for testing prime orders, ``num_checks=1`` is sufficient:: + + sage: p = random_prime(1000) + sage: E = EllipticCurve(GF(p), j=randrange(p)) + sage: q = random_prime(p + 20, lbound=p - 20) + sage: E.has_order(q, num_checks=20) == E.has_order(q, num_checks=1) + True + AUTHORS: - Mariah Lenox (2011-02-16): Initial implementation - Gareth Ma (2024-01-21): Fix bug for small curves """ - # This method does *not* use the cached value of `_order` even when - # available. q = self.base_field().order() a, b = Hasse_bounds(q, 1) if not a <= value <= b: @@ -2980,7 +2988,7 @@ def EllipticCurve_with_prime_order(N): Elliptic Curve defined by y^2 = x^3 + 20*x + 20 over Finite Field of size 23, Elliptic Curve defined by y^2 = x^3 + 10*x + 16 over Finite Field of size 23] sage: import itertools - sage: # These are the only primes, by the Weil-Hasse bound + sage: # These are the only primes, by the Hasse-Weil bound sage: for q in prime_range(17, 35): ....: K = GF(q) ....: for u in itertools.product(range(q), repeat=2):