Merge pull request #2702 from fredrik-johansson/divlowsmodp

Base fmpz_poly_div[low/high]_smodp on _fmpz_mod_poly_div_series

Author: fredrik-johansson

doc/source/fmpz_mod_vec.rst

@@ -11,6 +11,13 @@ Conversions
     Set the `fmpz_mod_vec` `(A, len)` to the `fmpz_vec` `(B, len)` after
     reduction of each entry modulo the modulus..
 
+.. function:: void _fmpz_mod_vec_get_fmpz_vec_smod(fmpz * res, const fmpz * vec, slong len, const fmpz_mod_ctx_t ctx)
+
+    Given `(vec, len)` of residues in `[0, n)`, set `(res, len)` to the symmetric
+    residues in `(-n/2, n/2]`. This is equivalent to applying
+    :func:`_fmpz_vec_scalar_smod_fmpz` but optimized for the case where
+    the input vector is known to be reduced mod `n` ahead of time.
+
 Arithmetic
 --------------------------------------------------------------------------------
 

src/fmpz_mod_vec.h

@@ -28,6 +28,9 @@ extern "C" {
 void _fmpz_mod_vec_set_fmpz_vec(fmpz * A, const fmpz * B, slong len,
                                                      const fmpz_mod_ctx_t ctx);
 
+void _fmpz_mod_vec_get_fmpz_vec_smod(fmpz * res, const fmpz * vec, slong len,
+    const fmpz_mod_ctx_t ctx);
+
 void _fmpz_mod_vec_neg(fmpz * A, const fmpz * B, slong len,
                                                      const fmpz_mod_ctx_t ctx);
 

src/fmpz_mod_vec/get_fmpz_vec_smod.c

@@ -0,0 +1,33 @@
+/*
+    Copyright (C) 2026 Fredrik Johansson
+
+    This file is part of FLINT.
+
+    FLINT is free software: you can redistribute it and/or modify it under
+    the terms of the GNU Lesser General Public License (LGPL) as published
+    by the Free Software Foundation; either version 3 of the License, or
+    (at your option) any later version.  See <https://www.gnu.org/licenses/>.
+*/
+#include "fmpz.h"
+#include "fmpz_mod.h"
+#include "fmpz_mod_vec.h"
+
+void _fmpz_mod_vec_get_fmpz_vec_smod(fmpz * res, const fmpz * vec, slong len, const fmpz_mod_ctx_t ctx)
+{
+    const fmpz * p = fmpz_mod_ctx_modulus(ctx);
+    fmpz_t phalf;
+    slong i;
+
+    fmpz_init(phalf);
+    fmpz_fdiv_q_2exp(phalf, p, 1);
+
+    for (i = 0; i < len; i++)
+    {
+        if (fmpz_cmp(vec + i, phalf) > 0)
+            fmpz_sub(res + i, vec + i, p);
+        else
+            fmpz_set(res + i, vec + i);
+    }
+
+    fmpz_clear(phalf);
+}

src/fmpz_poly/divhigh_smodp.c

@@ -1,5 +1,6 @@
 /*
     Copyright (C) 2016 William Hart
+    Copyright (C) 2026 Fredrik Johansson
 
     This file is part of FLINT.
 
@@ -12,43 +13,43 @@
 #include "fmpz.h"
 #include "fmpz_vec.h"
 #include "fmpz_poly.h"
+#include "fmpz_mod.h"
+#include "fmpz_mod_vec.h"
+#include "fmpz_mod_poly.h"
 
 void fmpz_poly_divhigh_smodp(fmpz * res, const fmpz_poly_t f,
-                                  const fmpz_poly_t g, const fmpz_t p, slong n)
+                              const fmpz_poly_t g, const fmpz_t p, slong n)
 {
-   fmpz_t d, cinv;
-   slong i = 0, k, start = 0, len_g = g->length;
-   fmpz_poly_t tf;
-
-   fmpz_init(d);
-   fmpz_init(cinv);
-
-   fmpz_poly_init2(tf, f->length);
-
-   fmpz_poly_set(tf, f);
-
-   fmpz_gcdinv(d, cinv, g->coeffs + len_g - 1, p);
-   if (!fmpz_is_one(d))
-   {
-      flint_throw(FLINT_ERROR, "Exception (fmpz_poly_divhigh_smodp). Impossible inverse.\n");
-   }
-
-   for (k = n - 1, i = f->length - len_g; k >= 0; i--, k--)
-   {
-      if (i < f->length - n)
-         start++;
-
-      fmpz_mul(res + k, tf->coeffs + i + len_g - 1, cinv);
-
-      fmpz_smod(res + k, res + k, p);
-
-      _fmpz_vec_scalar_submul_fmpz(tf->coeffs + i + start,
-                                    g->coeffs + start, len_g - start, res + k);
-      _fmpz_vec_scalar_smod_fmpz(tf->coeffs + i + start,
-                                     tf->coeffs + i + start, len_g - start, p);
-   }
-
-   fmpz_poly_clear(tf);
-   fmpz_clear(cinv);
-   fmpz_clear(d);
+    fmpz_mod_ctx_t ctx;
+    fmpz *tA, *tB;
+    slong glen, flen, Alen, Blen;
+
+    glen = g->length;
+    flen = f->length;
+    Blen = FLINT_MIN(glen, n);
+    Alen = FLINT_MIN(flen, n);
+
+    if (Alen == 0)
+    {
+        _fmpz_vec_zero(res, n);
+        return;
+    }
+
+    fmpz_mod_ctx_init(ctx, p);
+    tA = _fmpz_vec_init(Alen);
+    tB = _fmpz_vec_init(Blen);
+
+    _fmpz_mod_vec_set_fmpz_vec(tA, f->coeffs + (flen - Alen), Alen, ctx);
+    _fmpz_poly_reverse(tA, tA, Alen, Alen);
+    _fmpz_mod_vec_set_fmpz_vec(tB, g->coeffs + (glen - Blen), Blen, ctx);
+    _fmpz_poly_reverse(tB, tB, Blen, Blen);
+
+    _fmpz_mod_poly_div_series(res, tA, Alen, tB, Blen, n, ctx);
+    _fmpz_poly_reverse(res, res, n, n);
+    _fmpz_mod_vec_get_fmpz_vec_smod(res, res, n, ctx);
+
+    _fmpz_vec_clear(tA, Alen);
+    _fmpz_vec_clear(tB, Blen);
+    fmpz_mod_ctx_clear(ctx);
 }
+

src/fmpz_poly/divlow_smodp.c

@@ -1,5 +1,6 @@
 /*
     Copyright (C) 2016 William Hart
+    Copyright (C) 2026 Fredrik Johansson
 
     This file is part of FLINT.
 
@@ -12,58 +13,43 @@
 #include "fmpz.h"
 #include "fmpz_vec.h"
 #include "fmpz_poly.h"
+#include "fmpz_mod.h"
+#include "fmpz_mod_vec.h"
+#include "fmpz_mod_poly.h"
 
 void fmpz_poly_divlow_smodp(fmpz * res, const fmpz_poly_t f,
-                                  const fmpz_poly_t g, const fmpz_t p, slong n)
+                             const fmpz_poly_t g, const fmpz_t p, slong n)
 {
-   fmpz_t d, cinv;
-   slong i = 0, k, zeroes;
-   fmpz_poly_t tf;
-
-   fmpz_init(d);
-   fmpz_init(cinv);
-
-   while (fmpz_is_zero(g->coeffs + i))
-      i++;
-
-   zeroes = i;
-
-   fmpz_poly_init2(tf, n + zeroes);
-
-   fmpz_poly_set(tf, f);
-
-   if (fmpz_sgn(g->coeffs + zeroes) >= 0)
-      fmpz_gcdinv(d, cinv, g->coeffs + zeroes, p);
-   else
-   {
-      fmpz_t temp;
-
-      fmpz_init(temp);
-
-      fmpz_add(temp, g->coeffs + zeroes, p);
-      fmpz_gcdinv(d, cinv, temp, p);
-
-      fmpz_clear(temp);
-   }
-
-   if (!fmpz_is_one(d))
-   {
-      flint_throw(FLINT_ERROR, "Exception (fmpz_poly_divlow_smodp). Impossible inverse.\n");
-   }
-
-   for (k = 0; k < n; i++, k++)
-   {
-      fmpz_mul(res + k, tf->coeffs + i, cinv);
-
-      fmpz_smod(res + k, res + k, p);
-
-      _fmpz_vec_scalar_submul_fmpz(tf->coeffs + i,  g->coeffs + zeroes,
-                                FLINT_MIN(g->length - zeroes, n - k), res + k);
-      _fmpz_vec_scalar_smod_fmpz(tf->coeffs + i, tf->coeffs + i,
-                                      FLINT_MIN(g->length - zeroes, n - k), p);
-   }
-
-   fmpz_poly_clear(tf);
-   fmpz_clear(cinv);
-   fmpz_clear(d);
+    fmpz_mod_ctx_t ctx;
+    fmpz *tA, *tB;
+    slong zeroes, glen, flen, Alen, Blen;
+
+    zeroes = 0;
+    while (fmpz_is_zero(g->coeffs + zeroes))
+        zeroes++;
+
+    flen = f->length - zeroes;
+    glen = g->length - zeroes;
+    Alen = FLINT_MAX(0, FLINT_MIN(flen, n));
+    Blen = FLINT_MIN(glen, n);
+
+    if (Alen == 0)
+    {
+        _fmpz_vec_zero(res, n);
+        return;
+    }
+
+    fmpz_mod_ctx_init(ctx, p);
+    tA = _fmpz_vec_init(Alen);
+    tB = _fmpz_vec_init(Blen);
+
+    _fmpz_mod_vec_set_fmpz_vec(tA, f->coeffs + zeroes, Alen, ctx);
+    _fmpz_mod_vec_set_fmpz_vec(tB, g->coeffs + zeroes, Blen, ctx);
+    _fmpz_mod_poly_div_series(res, tA, Alen, tB, Blen, n, ctx);
+    _fmpz_mod_vec_get_fmpz_vec_smod(res, res, n, ctx);
+
+    _fmpz_vec_clear(tA, Alen);
+    _fmpz_vec_clear(tB, Blen);
+    fmpz_mod_ctx_clear(ctx);
 }
+

src/fmpz_poly/test/t-divlow_smodp.c

@@ -63,6 +63,9 @@ TEST_FUNCTION_START(fmpz_poly_divlow_smodp, state)
            }
         } while (b->length < 2 || !fmpz_is_one(d1) || !fmpz_is_one(d2));
 
+        /* Test with some leading zeros */
+        fmpz_poly_shift_left(b, b, n_randint(state, 3));
+
         fmpz_poly_mul(c, a, b);
 
         fmpz_poly_scalar_mod_fmpz(d, b, P);