Une MAJ de sécurité est nécessaire sur notre version actuelle. Elle sera effectuée lundi 02/08 entre 12h30 et 13h. L'interruption de service devrait durer quelques minutes (probablement moins de 5 minutes).

Commit 7c10b68f authored by Emmanuel Thomé's avatar Emmanuel Thomé
Browse files

scale algebraic square by $`{\hat f}'(\hat\alpha)^2`$

parent cfc03e7e
Pipeline #179642 passed with stage
in 13 minutes and 11 seconds
......@@ -70,6 +70,7 @@ struct cxx_mpz_polymod_scaled {
unsigned long v;
cxx_mpz_polymod_scaled(int deg) : p(deg), v(0) {}
cxx_mpz_polymod_scaled() = default;
cxx_mpz_polymod_scaled(cxx_mpz_poly const & p, unsigned long v = 0) : p(p), v(v) {}
cxx_mpz_polymod_scaled(cxx_mpz_polymod_scaled const &) = delete;
cxx_mpz_polymod_scaled(cxx_mpz_polymod_scaled &&) = default;
cxx_mpz_polymod_scaled& operator=(cxx_mpz_polymod_scaled const &) = delete;
......@@ -1005,6 +1006,16 @@ cxx_mpz_polymod_scaled_sqrt (cxx_mpz_polymod_scaled & res, cxx_mpz_polymod_scale
the coefficients. */
target_size = mpz_poly_sizeinbase (AA.p, 2);
target_size = target_size / 2;
/* note that we scale by the derivative of f_hat, which might increase
* the coefficient a little bit. I suppose that the bitsize of the
* discriminant is a safe bet. Better take into account */
{
cxx_mpz disc;
mpz_poly_discriminant(disc, F);
target_size += mpz_sizeinbase(disc, 2);
}
target_size += target_size / 10;
#pragma omp critical
{
......@@ -1301,8 +1312,23 @@ calculateSqrtAlg (const char *prefix, int numdep,
// Init F to be the corresponding polynomial
cxx_mpz_poly F(pol->pols[side]);
cxx_mpz_poly F_hat(F->deg);
cxx_mpz_polymod_scaled prod;
/* {{{ create F_hat, the minimal polynomial of alpha_hat = lc(F) *
* alpha */
{
cxx_mpz tmp;
mpz_init_set_ui(tmp, 1);
mpz_set_ui(F_hat->coeff[F->deg], 1);
for(int i = F->deg - 1 ; i >= 0 ; i--) {
mpz_mul(F_hat->coeff[i], tmp, F->coeff[i]);
mpz_mul(tmp, tmp, F->coeff[F->deg]);
}
mpz_poly_cleandeg(F_hat, F->deg);
}
/* }}} */
// Accumulate product with a subproduct tree
{
cxx_mpz_polymod_scaled_functions M(F, pinf);
......@@ -1312,6 +1338,19 @@ calculateSqrtAlg (const char *prefix, int numdep,
std::vector<cxx_mpz_polymod_scaled> prd = read_ab_pairs_from_depfile(depname, M, message, nab, nfree);
free(depname);
{
cxx_mpz_poly scale, alpha_hat;
mpz_poly_set_xi(alpha_hat, 1);
mpz_poly_mul_mpz(alpha_hat, alpha_hat, mpz_poly_lc(F));
mpz_poly_eval_diff_poly(scale, F_hat, alpha_hat);
/* Add these two to the subproduct tree. ok, they're slightly
* larger than the other elements, but it won't make a big
* difference in the long run, and is certainly not worse than
* doing the full unbalanced mul after the accumulation */
prd.emplace_back(scale,0);
prd.emplace_back(scale,0);
}
ASSERT_ALWAYS ((nab & 1) == 0);
ASSERT_ALWAYS ((nfree & 1) == 0);
/* nfree being even is forced by a specific character column added
......@@ -1411,7 +1450,20 @@ calculateSqrtAlg (const char *prefix, int numdep,
}
} while (!ret);
mpz_poly_eval_mod_mpz(algsqrt, prod.p, m, Np);
/* now divide by f'(alpha_hat), but mod N. alpha_hat mod N is lc*m
*/
{
cxx_mpz scale_modN, alpha_hat_modN;
mpz_mul(alpha_hat_modN, mpz_poly_lc(F), m);
mpz_poly_eval_diff(scale_modN, F_hat, alpha_hat_modN);
mpz_invert(scale_modN, scale_modN, Np);
mpz_mul(algsqrt, algsqrt, scale_modN);
mpz_mod(algsqrt, algsqrt, Np);
}
mpz_clear(m);
mpz_invert(aux, F->coeff[F->deg], Np); // 1/fd mod n
mpz_powm_ui(aux, aux, prod.v, Np); // 1/fd^v mod n
mpz_mul(algsqrt, algsqrt, aux);
......
......@@ -1258,6 +1258,31 @@ mpz_poly_sub_ui (mpz_poly_ptr g, mpz_poly_srcptr f, unsigned long a)
}
}
void mpz_poly_add_mpz(mpz_poly_ptr f, mpz_poly_srcptr g, mpz_srcptr a)
{
mpz_poly_set(f, g);
if (f->deg == -1) {
mpz_poly_set_mpz(f, a);
} else {
mpz_add(f->coeff[0], f->coeff[0], a);
if (f->deg == 0)
mpz_poly_cleandeg(f, 0);
}
}
void mpz_poly_sub_mpz(mpz_poly_ptr f, mpz_poly_srcptr g, mpz_srcptr a)
{
mpz_poly_set(f, g);
if (f->deg == -1) {
mpz_poly_set_mpz(f, a);
mpz_neg(f->coeff[0], f->coeff[0]);
} else {
mpz_sub(f->coeff[0], f->coeff[0], a);
if (f->deg == 0)
mpz_poly_cleandeg(f, 0);
}
}
/* Set f=g-h (mod m). */
void
mpz_poly_sub_mod_mpz (mpz_poly_ptr f, mpz_poly_srcptr g, mpz_poly_srcptr h, mpz_srcptr m)
......@@ -1794,6 +1819,7 @@ void mpz_poly_parallel_interface<inf>::mpz_poly_div_2_mod_mpz (mpz_poly_ptr f, m
/* Set res=f(x). Assumes res and x are different variables. */
void mpz_poly_eval(mpz_ptr res, mpz_poly_srcptr f, mpz_srcptr x) {
int i, d;
ASSERT_ALWAYS(res != x);
d = f->deg;
if (d == -1) {
mpz_set_ui(res, 0);
......@@ -1806,6 +1832,24 @@ void mpz_poly_eval(mpz_ptr res, mpz_poly_srcptr f, mpz_srcptr x) {
}
}
/* Set res=f(x). Assumes res and x are different variables. */
void mpz_poly_eval_poly(mpz_poly_ptr res, mpz_poly_srcptr f, mpz_poly_srcptr x)
{
int i, d;
ASSERT_ALWAYS(res != x);
ASSERT_ALWAYS(res != f);
d = f->deg;
if (d == -1) {
mpz_poly_set_zero(res);
return;
}
mpz_poly_set_mpz(res, f->coeff[d]);
for (i = d-1; i>=0; --i) {
mpz_poly_mul(res, res, x);
mpz_poly_add_mpz(res, res, f->coeff[i]);
}
}
/* Set res=f(x) where x is an unsigned long. */
void mpz_poly_eval_ui (mpz_ptr res, mpz_poly_srcptr f, unsigned long x)
{
......@@ -1833,6 +1877,41 @@ mpz_poly_eval_diff_ui (mpz_ptr res, mpz_poly_srcptr f, unsigned long x)
}
}
void
mpz_poly_eval_diff (mpz_ptr res, mpz_poly_srcptr f, mpz_srcptr x)
{
int d = f->deg;
ASSERT_ALWAYS(res != x);
mpz_mul_ui (res, f->coeff[d], d);
for (int i = d - 1; i >= 1; i--)
{
mpz_mul (res, res, x);
mpz_addmul_ui (res, f->coeff[i], i); /* res <- res + i*f[i] */
}
}
/* Set res=f'(x), where x is a polynomial */
void mpz_poly_eval_diff_poly (mpz_poly_ptr res, mpz_poly_srcptr f, mpz_poly_srcptr x)
{
ASSERT_ALWAYS(res != x);
ASSERT_ALWAYS(res != f);
int d = f->deg;
mpz_poly_realloc(res, f->deg * x->deg + 1);
mpz_t t;
mpz_init(t);
mpz_mul_ui (t, f->coeff[d], d);
mpz_poly_add_mpz(res, res, t);
for (int i = d - 1; i >= 1; i--)
{
mpz_poly_mul (res, res, x);
mpz_mul_ui (t, f->coeff[i], i); /* res <- res + i*f[i] */
mpz_poly_add_mpz(res, res, t);
}
mpz_clear(t);
}
/* Return 1 if poly(root) % modulus == 0, return 0 otherwise */
/* Coefficients of f(x) need not be reduced mod m */
......
......@@ -153,6 +153,8 @@ void mpz_poly_add(mpz_poly_ptr f, mpz_poly_srcptr g, mpz_poly_srcptr h);
void mpz_poly_sub(mpz_poly_ptr f, mpz_poly_srcptr g, mpz_poly_srcptr h);
void mpz_poly_add_ui(mpz_poly_ptr g, mpz_poly_srcptr f, unsigned long a);
void mpz_poly_sub_ui(mpz_poly_ptr g, mpz_poly_srcptr f, unsigned long a);
void mpz_poly_add_mpz(mpz_poly_ptr f, mpz_poly_srcptr g, mpz_srcptr a);
void mpz_poly_sub_mpz(mpz_poly_ptr f, mpz_poly_srcptr g, mpz_srcptr a);
void mpz_poly_sub_mod_mpz(mpz_poly_ptr f, mpz_poly_srcptr g, mpz_poly_srcptr h, mpz_srcptr m);
void mpz_poly_mul(mpz_poly_ptr f, mpz_poly_srcptr g, mpz_poly_srcptr h);
void mpz_poly_mul_mpz(mpz_poly_ptr Q, mpz_poly_srcptr P, mpz_srcptr a);
......@@ -186,6 +188,9 @@ void mpz_poly_mul_xplusa(mpz_poly_ptr g, mpz_poly_srcptr f, mpz_srcptr a);
void mpz_poly_eval(mpz_ptr res, mpz_poly_srcptr f, mpz_srcptr x);
void mpz_poly_eval_ui (mpz_ptr res, mpz_poly_srcptr f, unsigned long x);
void mpz_poly_eval_diff_ui (mpz_ptr res, mpz_poly_srcptr f, unsigned long x);
void mpz_poly_eval_diff (mpz_ptr res, mpz_poly_srcptr f, mpz_srcptr x);
void mpz_poly_eval_poly(mpz_poly_ptr res, mpz_poly_srcptr f, mpz_poly_srcptr x);
void mpz_poly_eval_diff_poly (mpz_poly_ptr res, mpz_poly_srcptr f, mpz_poly_srcptr x);
void mpz_poly_eval_mod_mpz(mpz_t res, mpz_poly_srcptr f, mpz_srcptr x, mpz_srcptr m);
int mpz_poly_is_root(mpz_poly_srcptr poly, mpz_srcptr root, mpz_srcptr modulus);
void mpz_poly_eval_several_mod_mpz(mpz_ptr * res, mpz_poly_srcptr * f, int k, mpz_srcptr x, mpz_srcptr m);
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment