Merge pull request #2312 from mezzarobba/shiftless

Shift equivalence, dispersion, and shiftless decomposition of generic polynomials

Author: fredrik-johansson

doc/source/gr.rst

@@ -230,6 +230,14 @@ that the value is not computable:
 
 .. function:: void truth_println(truth_t x)
 
+.. function:: truth_t gr_in_domain(int status)
+
+   Returns ``T_TRUE`` when ``status`` is ``GR_SUCCESS``, ``T_FALSE`` when ``GR_DOMAIN`` is set but ``GR_UNABLE`` is not, and ``T_UNKNOWN`` otherwise.
+
+.. function:: int gr_check(truth_t t)
+
+   Maps ``T_TRUE`` to ``GR_SUCCESS``, ``T_FALSE`` to ``GR_DOMAIN``, and ``T_UNKNOWN`` to ``GR_UNABLE``.
+
 Context operations
 -------------------------------------------------------------------------------
 

doc/source/gr_poly.rst

@@ -780,6 +780,82 @@ TODO: currently only fields of characteristic 0 are supported.
 
     Sets *res* to the squarefreepart of *poly*.
 
+Shift equivalence
+-------------------------------------------------------------------------------
+
+.. function:: truth_t gr_poly_shift_equivalent(fmpz_t shift, const gr_poly_t p, const gr_poly_t q, gr_ctx_t ctx)
+
+    Returns whether there exists an integer *n* such that `p(x + n) = q(x)`. If
+    the result is ``T_TRUE`` and *shift* is not ``NULL``, *shift* is set to
+    such an *n*.
+
+.. function:: int gr_poly_leading_taylor_shift(gr_ptr shift, const gr_poly_t p, const gr_poly_t q, gr_ctx_t ctx)
+
+    Computes (if possible) *s* such that `p(x+s) = q(x)(1+O(x^2))`.
+
+.. function:: int gr_poly_dispersion_resultant(fmpz_t disp, gr_vec_t disp_set, const gr_poly_t f, const gr_poly_t g, gr_ctx_t ctx);
+              int gr_poly_dispersion_factor(fmpz_t disp, gr_vec_t disp_set, const gr_poly_t f, const gr_poly_t g, gr_ctx_t ctx);
+              int gr_poly_dispersion(fmpz_t disp, gr_vec_t disp_set, const gr_poly_t f, const gr_poly_t g, gr_ctx_t ctx);
+
+    Computes the dispersion and/or the dispersion set of *f* and *g*.
+
+    The dispersion set of two polynomials *f* and *g* (over a unique
+    factorization domain of characteristic zero) is the set of nonnegative
+    integers *n* such that `f(x + n)` and `g(x)` have a nonconstant common
+    factor. The dispersion is the largest element of the dispersion set.
+
+    The output variables *disp* and/or *disp_set* can be ``NULL``, in which case
+    the corresponding result is not stored.
+    When the dispersion set is empty, *disp* is left unchanged.
+    The elements of *disp_set* are sorted in increasing order.
+
+    The *factor* version uses the algorithm described in [ManWright1994]_.
+    The *resultant* version computes the integer roots of a bivariate resultant
+    and is mainly intended for testing.
+
+.. function:: int gr_poly_dispersion_from_factors(fmpz_t disp, gr_vec_t disp_set, const gr_vec_t ffac, const gr_vec_t gfac, gr_ctx_t ctx);
+
+    Same as :func:`gr_poly_dispersion_factor` for nonzero *f* and *g* but takes
+    as input their nonconstant irreducible factors (without multiplicities)
+    instead of the polynomials themselves.
+
+.. function:: int gr_poly_shiftless_decomposition_factor(gr_ptr c, gr_vec_t slfac, gr_vec_t slshifts, gr_vec_t slmult, const gr_poly_t f, gr_ctx_t ctx)
+              int gr_poly_shiftless_decomposition(gr_ptr c, gr_vec_t slfac, gr_vec_t slshifts, gr_vec_t slmult, const gr_poly_t f, gr_ctx_t ctx)
+
+
+    Computes a decomposition of *f* of the form
+
+        .. math:: c \prod_i \prod_j g_i(x + h_{i,j})^{e_{i,j}}
+
+    where
+
+    * `c` is a constant,
+    * the `g_i` are squarefree polynomials of degree at least one,
+    * `g_i(x)` and `g_j(x + h)` (with `i \neq j`) are coprime for all
+      `h \in \mathbb Z`,
+    * `g_i(x)` and `g_i(x + h)` are coprime for all nonzero `h \in \mathbb Z`,
+    * `e_{i,j}` and `h_{i,j}` are integers with `e_{i,j} \geq 1`
+      and `0 = h_{i,1} < h_{i,2} < \cdots`.
+
+    The output variable *slfac* must be initialized to a vector of polynomials
+    of the same type as *f*. The other two output vectors *slshift* and
+    *slmult* must be initialized to vectors *of vectors* with entries of type
+    *fmpz*.
+
+    The *factor* version computes an irreducible factorization and sorts the
+    factors into shift-equivalence classes.
+
+    No algorithm avoiding a full irreducible factorization is currently
+    implemented.
+
+.. function:: int _gr_poly_shiftless_decomposition_from_factors(gr_vec_t slfac, gr_vec_t slshifts, gr_vec_t slmult, const gr_vec_t fac, const gr_vec_t mult, gr_ctx_t ctx)
+              int gr_poly_shiftless_decomposition_from_factors(gr_vec_t slfac, gr_vec_t slshifts, gr_vec_t slmult, const gr_vec_t fac, const gr_vec_t mult, gr_ctx_t ctx)
+
+    Same as :func:`gr_poly_shiftless_decomposition_factor` but takes as input
+    an irreducible factorization (*fac*, *mult*) of *f* (without the
+    prefactor *c*). The underscore method does not support aliasing of *slfac*
+    with *fac*.
+
 Roots
 -------------------------------------------------------------------------------
 

doc/source/references.rst

@@ -223,6 +223,8 @@ References
 
 .. [Lüb2004] \F. Lübeck, "Conway polynomials for finite fields", RTWH Aachen, https://www.math.rwth-aachen.de/~Frank.Luebeck/data/ConwayPol/index.html, (accessed 2024-01-12)
 
+.. [ManWright1994] \Y.-K. Man and F. J. Wright. "Fast polynomial dispersion computation and its application to indefinite summation". Proceedings of the International Symposium on Symbolic and Algebraic Computation (1994), 175-180. https://doi.org/10.1145/190347.190413
+
 .. [MN2019] \P. Molin and C. Neurohr, "Computing period matrices and the Abel--Jacobi map of superelliptic curves", Math. Comp. 88:316 (2019), 847--888.
 
 .. [MP2006] \M. Monagan and R. Pearce. "Rational simplification modulo a polynomial ideal". Proceedings of the 2006 international symposium on Symbolic and algebraic computation - ISSAC '06. https://doi.org/10.1145/1145768.1145809

src/gr.h

@@ -51,6 +51,28 @@ GR_INLINE truth_t truth_not(truth_t x)
     return T_UNKNOWN;
 }
 
+GR_INLINE truth_t gr_in_domain(int status)
+{
+    if (status == GR_SUCCESS)
+	return T_TRUE;
+    else if (status & GR_UNABLE)
+	return T_UNKNOWN;
+    else if (status & GR_DOMAIN)
+	return T_FALSE;
+    else
+	return T_UNKNOWN;
+}
+
+GR_INLINE int gr_check(truth_t t)
+{
+    if (t == T_TRUE)
+	return GR_SUCCESS;
+    else if (t == T_FALSE)
+	return GR_DOMAIN;
+    else
+	return GR_UNABLE;
+}
+
 GR_INLINE void truth_println(truth_t x)
 {
     if (x == T_TRUE) flint_printf("T_TRUE\n");

src/gr_poly.h

@@ -429,6 +429,22 @@ WARN_UNUSED_RESULT int gr_poly_evaluate_vec_iter(gr_vec_t ys, const gr_poly_t po
 int gr_poly_factor_squarefree(gr_ptr c, gr_vec_t fac, gr_vec_t exp, const gr_poly_t F, gr_ctx_t ctx);
 int gr_poly_squarefree_part(gr_poly_t res, const gr_poly_t poly, gr_ctx_t ctx);
 
+/* Shift factorization */
+
+WARN_UNUSED_RESULT int gr_poly_leading_taylor_shift(gr_ptr shift, const gr_poly_t p, const gr_poly_t q, gr_ctx_t ctx);
+
+truth_t gr_poly_shift_equivalent(fmpz_t shift, const gr_poly_t p, const gr_poly_t q, gr_ctx_t ctx);
+
+WARN_UNUSED_RESULT int gr_poly_dispersion_resultant(fmpz_t disp, gr_vec_t disp_set, const gr_poly_t f, const gr_poly_t g, gr_ctx_t ctx);
+WARN_UNUSED_RESULT int gr_poly_dispersion_from_factors(fmpz_t disp, gr_vec_t disp_set, const gr_vec_t ffac, const gr_vec_t gfac, gr_ctx_t ctx);
+WARN_UNUSED_RESULT int gr_poly_dispersion_factor(fmpz_t disp, gr_vec_t disp_set, const gr_poly_t f, const gr_poly_t g, gr_ctx_t ctx);
+WARN_UNUSED_RESULT int gr_poly_dispersion(fmpz_t disp, gr_vec_t disp_set, const gr_poly_t f, const gr_poly_t g, gr_ctx_t ctx);
+
+WARN_UNUSED_RESULT int _gr_poly_shiftless_decomposition_from_factors(gr_vec_t slfac, gr_vec_t slshifts, gr_vec_t slmult, const gr_vec_t fac, const gr_vec_t mult, gr_ctx_t ctx);
+WARN_UNUSED_RESULT int gr_poly_shiftless_decomposition_from_factors(gr_vec_t slfac, gr_vec_t slshifts, gr_vec_t slmult, const gr_vec_t fac, const gr_vec_t mult, gr_ctx_t ctx);
+WARN_UNUSED_RESULT int gr_poly_shiftless_decomposition_factor(gr_ptr c, gr_vec_t slfac, gr_vec_t slshifts, gr_vec_t slmult, const gr_poly_t pol, gr_ctx_t ctx);
+WARN_UNUSED_RESULT int gr_poly_shiftless_decomposition(gr_ptr c, gr_vec_t slfac, gr_vec_t slshifts, gr_vec_t slmult, const gr_poly_t pol, gr_ctx_t ctx);
+
 /* Roots */
 
 int _gr_poly_refine_roots_aberth(gr_ptr w, gr_srcptr f, gr_srcptr f_prime, slong deg, gr_srcptr z, int progressive, gr_ctx_t ctx);