diff --git a/Makefile.in b/Makefile.in
index 6050b08189..9a2a413f04 100644
--- a/Makefile.in
+++ b/Makefile.in
@@ -10,7 +10,7 @@ BUILD_DIRS = aprcl ulong_extras long_extras perm fmpz fmpz_vec fmpz_poly \
    fmpq_poly fmpz_mat fmpz_lll mpfr_vec mpfr_mat mpf_vec mpf_mat nmod_vec nmod_poly \
    nmod_poly_factor arith mpn_extras nmod_mat fmpq fmpq_vec fmpq_mat padic \
    fmpz_poly_q fmpz_poly_mat nmod_poly_mat fmpz_mod_poly \
-   fmpz_mod_poly_factor fmpz_factor fmpz_poly_factor fft qsieve \
+   fmpz_mod_poly_factor fmpz_factor fmpz_poly_factor fmpz_poly_roots fft qsieve \
    double_extras d_vec d_mat padic_poly padic_mat qadic  \
    fq fq_vec fq_mat fq_poly fq_poly_factor fq_embed \
    fq_nmod fq_nmod_vec fq_nmod_mat fq_nmod_poly fq_nmod_poly_factor fq_nmod_embed \
diff --git a/doc/source/fmpz_poly_roots.rst b/doc/source/fmpz_poly_roots.rst
new file mode 100644
index 0000000000..b8abb2617c
--- /dev/null
+++ b/doc/source/fmpz_poly_roots.rst
@@ -0,0 +1,182 @@
+.. _fmpz-poly-roots:
+
+**fmpz_poly_roots.h** -- roots of univariate polynomials over the integers
+==========================================================================
+
+Find roots of univariate polynomials over the integers over finite
+fields, p-adic or q-adic extensions.
+
+-----------------------------
+
+Roots over finite fields
+------------------------
+
+Types, macros and constants
+___________________________
+
+Type containing the information of roots over a finite field.
+
+.. type:: fmpz_polynomial_roots_fq_t
+	  
+
+Memory management
+_________________
+
+.. function:: void fmpz_poly_roots_fq_init2 (fmpz_poly_roots_fq_t roots, slong n, fq_ctx_t fctx)
+
+	Initializes ``roots`` for use, with context ``fctx``
+	to contain at most ``n`` roots. A corresponding call to
+	:func:`fmpz_poly_roots_fq_clear` must be made to free the memory
+	used bye the roots.
+
+.. function:: void fmpz_poly_roots_fq_clear (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+
+	Clears the given roots, releasing any memory used. It must be
+	reinitialized in order to be used again.
+
+Calculation of roots
+____________________
+
+.. function:: void fmpz_poly_roots_fq (fmpz_poly_roots_fq_t roots, fmpz_poly_t poly,  fq_ctx_t fctx)
+	      
+	 Takes as input a polynomial ``poly`` and a context of of a
+	 finite field ``fctx`` and initialized ``roots``.
+
+Input and output
+________________
+
+.. function:: char * fmpz_poly_roots_fq_get_str (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+	      
+	 Returns  a representation of ``roots``
+	 as a zero terminated string.
+
+.. function:: int fmpz_poly_roots_fq_print (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+
+	 Prints a representation of ``roots`` to ``stdout``.
+	 In case of success, returns a positive value.  In case of
+	 failure, returns a non-positive value.
+
+.. function:: char * fmpz_poly_roots_fq_get_str_pretty (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+
+	 Returns  a pretty representation of ``roots``
+	 as a zero terminated string.
+
+.. function:: int fmpz_poly_roots_fq_fprint_pretty (FILE * file, fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+
+	 Prints a pretty representation of ``roots`` to ``file``.
+	 In case of success, returns a positive value.  In case of
+	 failure, returns a non-positive value.
+	 
+.. function:: int fmpz_poly_roots_fq_print_pretty (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+
+	 Prints a pretty representation of ``roots`` to ``stdout``.
+	 In case of success, returns a positive value.  In case of
+	 failure, returns a non-positive value.
+
+
+Roots over padic fields
+-----------------------
+
+Types, macros and constants
+___________________________
+
+Type containing the information of roots over a p-adic field.
+
+.. type:: fmpz_polynomial_roots_padic_t
+	  
+Memory management
+_________________
+
+.. function:: void fmpz_poly_roots_padic_init2 (fmpz_poly_roots_padic_t roots, slong n, fq_ctx_t fctx)
+
+	Initializes ``roots`` for use, with context ``fctx``
+	to contain at most ``n`` roots. A corresponding call to :func:`fmpz_poly_roots_padic_clear` must be made to free the memory
+	used bye the roots.
+
+.. function:: void fmpz_poly_roots_padic_clear (fmpz_poly_roots_padic_t roots, fq_ctx_t fctx)
+
+	Clears the given roots, releasing any memory used. It must be
+	reinitialized in order to be used again.
+
+Calculation of roots
+____________________
+
+.. function:: void fmpz_poly_roots_padic_fmpz_poly (fmpz_poly_roots_padic_t roots, fmpz_poly_t poly,  fq_ctx_t fctx)
+	      
+	 Takes as input a polynomial ``poly`` and a context of of a
+	 finite field ``fctx`` and initialized ``roots``.
+
+Input and output
+________________
+
+
+.. function:: char * fmpz_poly_roots_padic_get_str (fmpz_poly_roots_padic_t roots, fq_ctx_t fctx)
+
+	 Returns  a representation of ``roots``
+	 as a zero terminated string.
+	 
+.. function:: int fmpz_poly_roots_padic_fprint_pretty (FILE * file, fmpz_poly_roots_padic_t roots, fq_ctx_t fctx)
+
+	 Prints a representation of ``roots`` to ``file``.
+	 In case of success, returns a positive value.  In case of
+	 failure, returns a non-positive value.
+	 
+.. function:: int fmpz_poly_roots_padic_print_pretty (fmpz_poly_roots_padic_t roots, fq_ctx_t fctx)
+
+	 Prints a representation of ``roots`` to ``stdout``.
+	 In case of success, returns a positive value.  In case of
+	 failure, returns a non-positive value.
+
+Roots over q-adic fields
+------------------------
+
+Types, macros and constants
+___________________________
+
+Type containing the information of roots over an unramified
+extension over p-adic field.
+
+.. type:: fmpz_polynomial_roots_qadic_t
+ 
+Memory management
+_________________
+
+.. function:: void fmpz_poly_roots_qadic_init2 (fmpz_poly_roots_qadic_t roots, slong n, fq_ctx_t fctx)
+
+	Initializes ``roots`` for use, with context ``fctx``
+	to contain at most ``n`` roots. A corresponding call to
+	:func:`fmpz_poly_roots_qadic_clear` must be made to free the memory
+	used bye the roots.
+
+.. function:: void fmpz_poly_roots_qadic_clear (fmpz_poly_roots_qadic_t roots, fq_ctx_t fctx)
+
+	Clears the given roots, releasing any memory used. It must be
+	reinitialized in order to be used again.
+
+Calculation of roots
+________________________
+
+.. function:: void fmpz_poly_roots_qadic (fmpz_poly_roots_qadic_t roots, fmpz_poly_t poly,  fq_ctx_t fctx)
+	      
+	 Takes as input a polynomial ``poly`` and a context of of a
+	 finite field ``fctx`` and initialized ``roots``.
+
+Input and output
+________________
+
+.. function:: char * fmpz_poly_roots_qadic_get_str_pretty (fmpz_poly_roots_qadic_t roots, fq_ctx_t fctx)
+
+	 Returns a pretty representation of ``roots`` as zero
+	 terminated string.
+
+.. function:: int fmpz_poly_roots_qadic_fprint_pretty (FILE *file, fmpz_poly_roots_qadic_t roots, fq_ctx_t fctx)
+
+	 Prints a pretty representation of ``roots`` to ``file``.
+	 In case of success, returns a positive value.  In case of
+	 failure, returns a non-positive value.
+
+.. function:: int fmpz_poly_roots_qadic_print_pretty (fmpz_poly_roots_qadic_t roots, fq_ctx_t fctx)
+
+	 Prints a pretty representation of ``roots`` to ``stdout``.
+	 In case of success, returns a positive value.  In case of
+	 failure, returns a non-positive value.
diff --git a/doc/source/index.rst b/doc/source/index.rst
index f6d5dcb633..d0a6486a09 100644
--- a/doc/source/index.rst
+++ b/doc/source/index.rst
@@ -59,6 +59,7 @@ Integers
    fmpz_poly.rst
    fmpz_poly_mat.rst
    fmpz_poly_factor.rst
+   fmpz_poly_roots.rst
    fmpz_mpoly.rst
    fmpz_mpoly_factor.rst
    long_extras.rst
diff --git a/fmpz_poly_roots.h b/fmpz_poly_roots.h
new file mode 100644
index 0000000000..d98ff71f59
--- /dev/null
+++ b/fmpz_poly_roots.h
@@ -0,0 +1,64 @@
+#ifndef FMPZ_POLY_ROOTS_H_
+#define FMPZ_POLY_ROOTS_H_
+
+#include "fmpz_poly.h"
+
+#include "fq.h"
+#include "padic.h"
+#include "qadic.h"
+
+typedef struct
+{
+  fq_struct *x0;
+  slong *multiplicity;
+  slong num;
+  slong alloc;
+} fmpz_poly_roots_fq_struct;
+
+typedef fmpz_poly_roots_fq_struct fmpz_poly_roots_fq_t[1];
+
+void fmpz_poly_roots_fq_init2 (fmpz_poly_roots_fq_t roots, slong n, fq_ctx_t fctx);
+void fmpz_poly_roots_fq_clear (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx);
+void fmpz_poly_roots_fq_print (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx);
+char* fmpz_poly_roots_fq_get_str_pretty (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx);
+int fmpz_poly_roots_fq_fprint_pretty (FILE *file, fmpz_poly_roots_fq_t roots, fq_ctx_t fctx);
+int fmpz_poly_roots_fq_print_pretty (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx);
+void fmpz_poly_roots_fq(fmpz_poly_roots_fq_t roots, fmpz_poly_t poly, fq_ctx_t fctx);
+
+typedef struct
+{
+  padic_struct *x0;
+  slong *multiplicity;
+  slong num;
+  slong alloc;
+} fmpz_poly_roots_padic_struct;
+
+typedef fmpz_poly_roots_padic_struct fmpz_poly_roots_padic_t[1];
+
+void fmpz_poly_roots_padic_init2 (fmpz_poly_roots_padic_t roots, slong n);
+void fmpz_poly_roots_padic_clear (fmpz_poly_roots_padic_t roots);
+char* fmpz_poly_roots_padic_get_str (fmpz_poly_roots_padic_t roots, padic_ctx_t fctx);
+int fmpz_poly_roots_padic_fprint (FILE* file, fmpz_poly_roots_padic_t roots, padic_ctx_t fctx);
+int fmpz_poly_roots_padic_print (fmpz_poly_roots_padic_t roots, padic_ctx_t fctx);
+void fmpz_poly_roots_padic (fmpz_poly_roots_padic_t roots, fmpz_poly_t poly,
+			    padic_ctx_t fctx);
+
+typedef struct
+{
+  qadic_struct *x0;
+  slong *multiplicity;
+  slong num;
+  slong alloc;
+} fmpz_poly_roots_qadic_struct;
+
+typedef fmpz_poly_roots_qadic_struct fmpz_poly_roots_qadic_t[1];
+
+void fmpz_poly_roots_qadic_init2 (fmpz_poly_roots_qadic_t roots, slong n);
+void fmpz_poly_roots_qadic_clear (fmpz_poly_roots_qadic_t roots);
+char* fmpz_poly_roots_qadic_get_str_pretty (fmpz_poly_roots_qadic_t roots, qadic_ctx_t fctx);
+int fmpz_poly_roots_qadic_fprint_pretty (FILE * file, fmpz_poly_roots_qadic_t roots, qadic_ctx_t fctx);
+int fmpz_poly_roots_qadic_print_pretty (fmpz_poly_roots_qadic_t roots, qadic_ctx_t fctx);
+void fmpz_poly_roots_qadic (fmpz_poly_roots_qadic_t roots, fmpz_poly_t poly,
+			    qadic_ctx_t fctx);
+
+#endif
diff --git a/fmpz_poly_roots/fq_roots.c b/fmpz_poly_roots/fq_roots.c
new file mode 100644
index 0000000000..315aeac89d
--- /dev/null
+++ b/fmpz_poly_roots/fq_roots.c
@@ -0,0 +1,132 @@
+#include "flint.h"
+#include "fq_vec.h"
+#include "fq_poly.h"
+#include "fmpz_poly.h"
+#include "fmpz_mod_poly.h"
+#include "fmpz_poly_roots.h"
+
+void
+fmpz_poly_roots_fq_init2 (fmpz_poly_roots_fq_t roots, slong n, fq_ctx_t fctx)
+{
+  roots->x0 = _fq_vec_init (n, fctx);
+  roots->multiplicity = flint_malloc (sizeof (slong) * n);
+  roots->num = n;
+  roots->alloc = n;
+}
+
+void
+fmpz_poly_roots_fq_clear (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+{
+  _fq_vec_clear (roots->x0, roots->alloc, fctx);
+  flint_free (roots->multiplicity);
+}
+
+
+char *
+fmpz_poly_roots_fq_get_str (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+{
+  char * buffer = NULL;
+  size_t buffer_size = 0;
+  FILE * out = open_memstream(&buffer, &buffer_size);
+  
+  _fq_vec_fprint (out, roots->x0, roots->num, fctx);
+  fclose(out);
+
+  return buffer;
+}
+
+void
+fmpz_poly_roots_fq_print (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+{
+  _fq_vec_print (roots->x0, roots->num, fctx);
+}
+
+char *
+fmpz_poly_roots_fq_get_str_pretty (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+{
+  char * buffer = NULL;
+  size_t buffer_size = 0;
+  FILE * out = open_memstream(&buffer, &buffer_size);
+  slong j;
+  
+  for (j = 0; j < roots->num; j++)
+    {
+      fq_fprint_pretty (out, roots->x0 + j, fctx);
+      flint_fprintf (out, " %wd\n", roots->multiplicity[j]);
+    }
+
+  fclose(out);
+  
+  return buffer;
+}
+
+int
+fmpz_poly_roots_fq_fprint_pretty (FILE* file, fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+{
+  slong j;
+  
+  for (j = 0; j < roots->num; j++)
+    {
+      fq_fprint_pretty (file, roots->x0 + j, fctx);
+      flint_fprintf (file, " %wd\n", roots->multiplicity[j]);
+    }
+
+  return 1;
+}
+
+int
+fmpz_poly_roots_fq_print_pretty (fmpz_poly_roots_fq_t roots, fq_ctx_t fctx)
+{
+  fmpz_poly_roots_fq_fprint_pretty(stdout, roots, fctx);
+  return 1;
+}
+
+void
+fmpz_poly_roots_fq (fmpz_poly_roots_fq_t roots, fmpz_poly_t poly, fq_ctx_t fctx)
+{
+  slong j, k, num;
+
+  fmpz_mod_ctx_t fmctx;
+  fmpz_mod_poly_t mpoly;
+  fq_poly_factor_t f;
+  fq_t lead;
+  fq_poly_t fpoly;
+  
+  fmpz_mod_ctx_init(fmctx, fq_ctx_prime(fctx));
+  fmpz_mod_poly_init (mpoly, fmctx);
+  fmpz_mod_poly_set_fmpz_poly (mpoly, poly, fmctx);
+
+  fq_init (lead, fctx);
+  fq_poly_factor_init (f, fctx);
+  fq_poly_init (fpoly, fctx);
+  fq_poly_set_fmpz_mod_poly (fpoly, mpoly, fctx);
+  fq_poly_factor(f, lead, fpoly, fctx);
+  fq_clear(lead, fctx);
+  
+  fmpz_mod_poly_clear (mpoly, fmctx);
+  fmpz_mod_ctx_clear(fmctx);
+  fq_poly_clear (fpoly, fctx);
+ 
+  num = 0;
+  
+  for (j = 0; j < f->num; j++)
+    {
+      if (fq_poly_degree (f->poly + j, fctx) == 1)
+	num++;
+    }
+  
+  fmpz_poly_roots_fq_init2 (roots, num, fctx);
+  
+  k = 0;
+  
+  for (j = 0; j < f->num; j++)
+    {
+      if (fq_poly_degree (f->poly + j, fctx) == 1)
+	{
+	  fq_poly_get_coeff (roots->x0 + k, f->poly + j, 0, fctx);
+	  fq_neg(roots->x0 + k, roots->x0 + k, fctx);
+	  roots->multiplicity[k] = f->exp[j];
+	  k++;
+	}
+    }
+}
diff --git a/fmpz_poly_roots/padic_roots.c b/fmpz_poly_roots/padic_roots.c
new file mode 100644
index 0000000000..7abc8e9f90
--- /dev/null
+++ b/fmpz_poly_roots/padic_roots.c
@@ -0,0 +1,148 @@
+#include <stdio.h>
+
+#include "padic.h"
+#include "padic_poly.h"
+#include "fmpz_poly_roots.h"
+
+
+void
+fmpz_poly_roots_padic_init2 (fmpz_poly_roots_padic_t roots, slong n)
+{
+  slong j;
+  
+  roots->x0 = flint_malloc (sizeof (padic_t) * n);
+  
+  for (j = 0; j < n; j++)
+    {
+      padic_init (roots->x0 + j);
+    }
+  
+  roots->multiplicity = flint_malloc (sizeof (slong) * n);
+  roots->num = n;
+  roots->alloc = n;
+}
+
+void
+fmpz_poly_roots_padic_clear (fmpz_poly_roots_padic_t roots)
+{
+  slong j;
+  
+  for (j = 0; j < roots->alloc; j++)
+    {
+      padic_clear (roots->x0 + j);
+    }
+  
+  flint_free (roots->x0);
+  flint_free (roots->multiplicity);
+}
+
+char*
+fmpz_poly_roots_padic_get_str (fmpz_poly_roots_padic_t roots, padic_ctx_t pctx)
+{
+  char * buffer = NULL;
+  size_t buffer_size = 0;
+  FILE *out = open_memstream(&buffer, &buffer_size);
+  
+  slong j;
+  
+  for (j = 0; j < roots->num; j++)
+    {
+      padic_fprint(out, roots->x0 + j, pctx);
+      flint_fprintf (out, " %wd\n", roots->multiplicity[j]);
+    }
+  
+  fclose(out);
+  
+  return buffer;
+}
+
+int
+fmpz_poly_roots_padic_fprint (FILE *file, fmpz_poly_roots_padic_t roots, padic_ctx_t pctx)
+{
+  slong j;
+  
+  for (j = 0; j < roots->num; j++)
+    {
+      padic_fprint(file, roots->x0 + j, pctx);
+      flint_fprintf (file, " %wd\n", roots->multiplicity[j]);
+    }
+
+  return 1;
+}
+
+int
+fmpz_poly_roots_padic_print (fmpz_poly_roots_padic_t roots, padic_ctx_t pctx)
+{
+  fmpz_poly_roots_padic_fprint(stdout, roots, pctx);
+  return 1;
+}
+
+static void
+padic_hensel_iteration (fmpz_poly_t poly, padic_t x, padic_ctx_t ctx)
+{
+  slong j;
+  padic_t tmp, y0, y1;
+  
+  padic_init (tmp);
+  padic_init (y0);
+  padic_init (y1);
+  
+  do
+    {
+      /* Horner evaluation of poly and poly' at x */
+      padic_set_fmpz (y0, poly->coeffs + poly->length - 1, ctx);
+      padic_zero (y1);
+      for (j = poly->length - 2; j >= 0; j--)
+	{
+	  padic_mul (y1, y1, x, ctx);
+	  padic_add (y1, y1, y0, ctx);
+	  padic_mul (y0, y0, x, ctx);
+	  padic_set_fmpz (tmp, poly->coeffs + j, ctx);
+	  padic_add (y0, y0, tmp, ctx);
+	}
+      /* Newton step: x -> x - poly / poly' */
+      padic_inv (y1, y1, ctx);
+      padic_mul (y1, y1, y0, ctx);
+      padic_sub (x, x, y1, ctx);
+    }
+  while (padic_val (y0));
+  
+  padic_clear (tmp);
+  padic_clear (y0);
+  padic_clear (y1);
+}
+
+void
+fmpz_poly_roots_padic (fmpz_poly_roots_padic_t roots, fmpz_poly_t poly,
+		       padic_ctx_t pctx)
+{
+  slong j, k, multiplicity;
+  fq_ctx_t fctx;
+  fmpz_poly_roots_fq_t froots;
+  fmpz_poly_t tmp_poly;
+  
+  fq_ctx_init (fctx, pctx->p, 1, "a");
+  
+  fmpz_poly_roots_fq (froots, poly, fctx);
+  fmpz_poly_roots_padic_init2 (roots, froots->num);
+
+  fmpz_poly_init(tmp_poly);
+  
+  for (j = 0; j < roots->num; j++)
+    {
+      multiplicity = *(froots->multiplicity + j);
+      *(roots->multiplicity + j) = multiplicity;
+      fmpz_poly_set(tmp_poly, poly);
+      if( multiplicity > 1 ) {
+	for(k = 1; k < multiplicity; k++) {
+	  fmpz_poly_derivative(tmp_poly, tmp_poly);
+	}
+      }
+      padic_set_fmpz (roots->x0 + j, (froots->x0 + j)->coeffs, pctx);
+      padic_hensel_iteration (tmp_poly, roots->x0 + j, pctx);
+    }
+
+  fmpz_poly_clear(tmp_poly);
+  fq_ctx_clear (fctx);
+  fmpz_poly_roots_fq_clear (froots, fctx);
+}
diff --git a/fmpz_poly_roots/qadic_roots.c b/fmpz_poly_roots/qadic_roots.c
new file mode 100644
index 0000000000..2730c934e5
--- /dev/null
+++ b/fmpz_poly_roots/qadic_roots.c
@@ -0,0 +1,162 @@
+#include "qadic.h"
+#include "fmpz_poly_roots.h"
+
+void
+fmpz_poly_roots_qadic_init2 (fmpz_poly_roots_qadic_t roots, slong n)
+{
+  slong j;
+  
+  roots->x0 = flint_malloc (sizeof (qadic_t) * n);
+  
+  for (j = 0; j < n; j++)
+    {
+      qadic_init (roots->x0 + j);
+    }
+  roots->multiplicity = flint_malloc (sizeof (slong) * n);
+  roots->num = n;
+  roots->alloc = n;
+}
+
+void
+fmpz_poly_roots_qadic_clear (fmpz_poly_roots_qadic_t roots)
+{
+  slong j;
+  
+  for (j = 0; j < roots->alloc; j++)
+    {
+      qadic_clear (roots->x0 + j);
+    }
+  
+  flint_free (roots->x0);
+  flint_free (roots->multiplicity);
+}
+
+char *
+fmpz_poly_roots_qadic_get_str_pretty (fmpz_poly_roots_qadic_t roots, qadic_ctx_t qctx)
+{
+  char * buffer = NULL;
+  size_t buffer_size = 0;
+  FILE *out = open_memstream(&buffer, &buffer_size);
+  slong j;
+  
+  for (j = 0; j < roots->num; j++)
+    {
+      qadic_fprint_pretty (out, roots->x0 + j, qctx);
+      flint_fprintf (out, " %wd\n", roots->multiplicity[j]);
+    }
+
+  fclose(out);
+
+  return buffer;
+}
+
+
+int
+fmpz_poly_roots_qadic_fprint_pretty (FILE *file, fmpz_poly_roots_qadic_t roots, qadic_ctx_t qctx)
+{
+  slong j;
+  
+  for (j = 0; j < roots->num; j++)
+    {
+      qadic_fprint_pretty (file, roots->x0 + j, qctx);
+      flint_fprintf (file, " %wd\n", roots->multiplicity[j]);
+    }
+
+  return 1;
+}
+
+int
+fmpz_poly_roots_qadic_print_pretty (fmpz_poly_roots_qadic_t roots, qadic_ctx_t qctx)
+{
+  fmpz_poly_roots_qadic_fprint_pretty (stdout, roots, qctx);
+  
+  return 1;
+}
+
+static void
+qadic_hensel_iteration (fmpz_poly_t poly, qadic_t x, qadic_ctx_t ctx)
+{
+  slong j;
+  
+  qadic_t tmp, y0, y1;
+  
+  qadic_init (tmp);
+  qadic_init (y0);
+  qadic_init (y1);
+  
+  do
+    {
+      /* Horner evaluation of poly and poly' at x */
+      padic_poly_set_fmpz (y0, poly->coeffs + poly->length - 1, &(ctx->pctx));
+      qadic_zero (y1);
+      for (j = poly->length - 2; j >= 0; j--)
+	{
+	  qadic_mul (y1, y1, x, ctx);
+	  qadic_add (y1, y1, y0, ctx);
+	  qadic_mul (y0, y0, x, ctx);
+	  padic_poly_set_fmpz (tmp, poly->coeffs + j, &(ctx->pctx));
+	  qadic_add (y0, y0, tmp, ctx);
+	}
+      /* Newton step: x -> x - poly / poly' */
+      qadic_inv (y1, y1, ctx);
+      qadic_mul (y1, y1, y0, ctx);
+      qadic_sub (x, x, y1, ctx);
+    }
+  
+  while (qadic_val (y0));
+  
+  qadic_clear (tmp);
+  qadic_clear (y0);
+  qadic_clear (y1);
+}
+
+void
+fmpz_poly_roots_qadic (fmpz_poly_roots_qadic_t roots, fmpz_poly_t poly,
+		       qadic_ctx_t qctx)
+{
+  slong j, k, multiplicity;
+  fq_ctx_t fctx;
+  fmpz_poly_roots_fq_t froots;
+  fmpz_poly_t tmp_poly;
+  qadic_t q, a;
+  
+  fq_ctx_init (fctx, (qctx->pctx).p, qadic_ctx_degree (qctx), "a");
+  fmpz_poly_roots_fq (froots, poly, fctx);
+  fmpz_poly_roots_qadic_init2 (roots, froots->num);
+  
+  fmpz_poly_init(tmp_poly);
+  
+  qadic_init (q);
+  qadic_init (a);
+  qadic_gen (q, qctx);
+  
+  for (j = 0; j < roots->num; j++)
+    {
+      multiplicity = *(froots->multiplicity + j);
+      *(roots->multiplicity + j) = multiplicity;
+      padic_poly_set_fmpz (roots->x0 + j,
+			   (froots->x0 + j)->coeffs + (froots->x0 +
+						       j)->length - 1,
+			   &(qctx->pctx));
+      for (k = (froots->x0 + j)->length - 2; k >= 0; k--)
+	{
+	  qadic_mul (roots->x0 + j, roots->x0 + j, q, qctx);
+	  padic_poly_set_fmpz (a, (froots->x0 + j)->coeffs + k,
+			       &(qctx->pctx));
+	  qadic_add (roots->x0 + j, roots->x0 + j, a, qctx);
+	}
+      fmpz_poly_set(tmp_poly, poly);
+      if( multiplicity > 1 ) {
+	for(k = 1; k < multiplicity; k++) {
+	  fmpz_poly_derivative(tmp_poly, tmp_poly);
+	}
+      }
+      qadic_hensel_iteration (tmp_poly, roots->x0 + j, qctx);
+    }
+  
+  qadic_clear (a);
+  qadic_clear (q);
+  
+  fq_ctx_clear (fctx);
+  fmpz_poly_roots_fq_clear (froots, fctx);
+}