diff --git a/lib/pari.c b/lib/pari.c index f4032226c78aa9daf8e086241daa33f2c66d6eaa..680510e2e2f93f3bc8c282f718b77945970ad824 100644 --- a/lib/pari.c +++ b/lib/pari.c @@ -3,7 +3,7 @@ pari.c - functions using pari; for factoring polynomials and for computing generators of class groups -Copyright (C) 2010, 2015, 2018, 2021, 2022, 2023 Andreas Enge +Copyright (C) 2010, 2015, 2018, 2021, 2022, 2023, 2024 Andreas Enge This file is part of CM. @@ -30,6 +30,9 @@ static GEN icl_get_Z (int_cl_t z); static int_cl_t Z_get_icl (GEN x); static GEN mpzx_get_FpX (mpzx_srcptr f, mpz_srcptr p); static void FpX_get_mpzx (mpzx_ptr f, GEN x); +static GEN cert_get_V (mpz_t **cert, int depth); +static bool primecertentryisvalid (GEN c, int i, bool type, bool verbose); +static int primecertpartialisvalid (GEN c, bool type, bool verbose); /*****************************************************************************/ /* */ @@ -149,6 +152,28 @@ static void FpX_get_mpzx (mpzx_ptr f, GEN x) Z_get_mpz (f->coeff [i], gel (x, i+2)); } +/*****************************************************************************/ + +static GEN cert_get_V (mpz_t **cert, int depth) + /* Given a (partial or complete) ECPP certificate of the given depth + in cert, return a corresponding PARI object. */ +{ + GEN c, ci; + int i, j; + + c = cgetg (depth + 1, t_VEC); + for (i = 0; i < depth; i++) { + ci = cgetg (6, t_VEC); + gel (c, i+1) = ci; + for (j = 0; j < 4; j++) + gel (ci, j+1) = mpz_get_Z (cert [i][j]); + gel (ci, 5) = mkvec2 (mpz_get_Z (cert [i][4]), + mpz_get_Z (cert [i][5])); + } + + return c; +} + /*****************************************************************************/ /* */ /* Various simple functions. */ @@ -1149,28 +1174,173 @@ bool cm_pari_cornacchia (mpz_ptr t, mpz_ptr v, mpz_srcptr p, /*****************************************************************************/ +static bool primecertentryisvalid (GEN c, int i, bool type, bool verbose) + /* Check whether line i of the partial ECPP certificate c is valid + (where numbering starts with 1 as conventional in PARI/GP). + c is supposed to be of type t_VEC. + If type is true, assume a CM style certificate, that is, the points + are of prime order; + if type is false, assume a PARI/GP style certificate, that is, + the points must first be multiplied by the cofactor to reach prime + order. */ +{ + GEN ci, N, t, s, a, P, x, y, cip1, N1, m, q, PP, Q; + int depth; + bool res; + pari_sp av; + + av = avma; + + /* Extract the line and check all types in passing. */ + depth = glength (c); + if (i > depth) + return false; + ci = gel (c, i); + if (glength (ci) != 5) + return false; + N = gel (ci, 1); + t = gel (ci, 2); + s = gel (ci, 3); + a = gel (ci, 4); + P = gel (ci, 5); + if (typ (N) != t_INT || typ (t) != t_INT || typ (s) != t_INT + || typ (a) != t_INT || typ (P) != t_VEC || glength (P) != 2) + return false; + x = gel (P, 1); + y = gel (P, 2); + if ( typ (x) != t_INT || typ (y) != t_INT) + return false; + + /* Do the same for the first entry of the next line if applicable. */ + if (i < depth) { + cip1 = gel (c, i+1); + if (typ (cip1) != t_VEC || glength (cip1) != 5) + return false; + N1 = gel (cip1, 1); + if (typ (N1) != t_INT) + return false; + } + + /* Check for numerical boundaries. + From here on, a simple return is not possible any more, we need + garbage collection. */ + res = true; + if (cmpui (2, N) > 0 + || cmpii (sqri (t), mului (4, N)) > 0 + || cmpui (0, s) >= 0 + || cmpui (0, a) > 0 || cmpii (a, N) >= 0 + || cmpui (0, x) > 0 || cmpii (x, N) >= 0 + || cmpui (0, y) > 0 || cmpii (y, N) >= 0) + res = false; + else { + m = subii (addui (1, N), t); + if (!equalui (0, modii (m, s))) + res = false; + else { + /* Check for consistency with next line if applicable. */ + q = diviiexact (m, s); + if (i < depth && !equalii (q, N1)) + res = false; + else { + /* Check that q > (\sqrt [4]{N} + 1)^2; the following + equivalent formula has been obtained by repeated sorting + and squaring. */ + if (cmpii (sqri (subii (sqri (subiu (q, 1)), N)), + mului (16, mulii (q, N))) < 0) + res = false; + else { + PP = cgetg (4, t_VEC); + gel (PP, 1) = x; + gel (PP, 2) = y; + gel (PP, 3) = gen_1; + /* If type is false, replace PP by s*PP and check that it + is different from the point at infinity modulo all prime + divisors of N, otherwise said that its third projective + coordinate is coprime to N. */ + if (!type) { + PP = FpJ_mul (PP, s, a, N); + if (!equalui (1, gcdii (gel (PP, 3), N))) + res = false; + } + /* Check that q*PP = O. */ + if (res) { + Q = FpJ_mul (PP, q, a, N); + res = (equalui (0, modii (gel (Q, 3), N)) != 0); + } + } + } + } + } + + if (verbose) + printf (" ECPP line %5i: %i\n", i, res); + + avma = av; + + return res; +} + +/*****************************************************************************/ + +static int primecertpartialisvalid (GEN c, bool type, bool verbose) + /* Check whether c is a valid partial ECPP certificate. + type indicates whether the certificate is of type created by CM + or not. + The return value is 0, 1 or 2. + 0 means the certificate is invalid. + 2 means the certificate is valid, but only partial. + 1 means the certificate is valid and its first entry is a + proven prime. */ +{ + int depth, val, res, i; + GEN ci, N, t, s, q; + pari_sp av; + + av = avma; + + if (typ (c) != t_VEC) + return 0; + + depth = glength (c); + val = 0; + for (i = 1; i <= depth; i++) + val += (primecertentryisvalid (c, i, type, verbose) ? 1 : 0); + + if (val != depth) + return 0; + + ci = gel (c, depth); + N = gel (ci, 1); + t = gel (ci, 2); + s = gel (ci, 3); + q = diviiexact (subii (addui (1, N), t), s); + if (expi (q) >= 64) + res = 2; + else if (equalui (1, gispseudoprime (N, 0))) + res = 1; + else + res = 0; + + avma = av; + + return res; +} + +/*****************************************************************************/ + bool cm_pari_ecpp_check (mpz_t **cert, int depth) /* Given a complete ECPP certificate (after step 2) of length depth in cert, use PARI to check it and return its validity. */ { pari_sp av; bool res; - GEN c, ci; - int i, j; + GEN c; av = avma; - c = cgetg (depth + 1, t_VEC); - for (i = 0; i < depth; i++) { - ci = cgetg (6, t_VEC); - gel (c, i+1) = ci; - for (j = 0; j < 4; j++) - gel (ci, j+1) = mpz_get_Z (cert [i][j]); - gel (ci, 5) = mkvec2 (mpz_get_Z (cert [i][4]), - mpz_get_Z (cert [i][5])); - } + c = cert_get_V (cert, depth); - res = (primecertisvalid (c) == 1); + res = (primecertpartialisvalid (c, true, true) == 1); avma = av;