some more work on complex addmul and polynomial evaluation

acb.h

@@ -460,21 +460,47 @@ acb_mul_2exp_si(acb_t z, const acb_t x, long e)
 static __inline__ void
 acb_addmul(acb_t z, const acb_t x, const acb_t y, long prec)
 {
-    acb_t t;
-    acb_init(t);
-    acb_mul(t, x, y, prec);
-    acb_add(z, z, t, prec);
-    acb_clear(t);
+    if (arb_is_zero(acb_imagref(y)))
+    {
+        arb_addmul(acb_imagref(z), acb_imagref(x), acb_realref(y), prec);
+        arb_addmul(acb_realref(z), acb_realref(x), acb_realref(y), prec);
+    }
+    else if (arb_is_zero(acb_imagref(x)))
+    {
+        arb_addmul(acb_imagref(z), acb_imagref(y), acb_realref(x), prec);
+        arb_addmul(acb_realref(z), acb_realref(y), acb_realref(x), prec);
+    }
+    else
+    {
+        acb_t t;
+        acb_init(t);
+        acb_mul(t, x, y, prec);
+        acb_add(z, z, t, prec);
+        acb_clear(t);
+    }
 }
 
 static __inline__ void
 acb_submul(acb_t z, const acb_t x, const acb_t y, long prec)
 {
-    acb_t t;
-    acb_init(t);
-    acb_mul(t, x, y, prec);
-    acb_sub(z, z, t, prec);
-    acb_clear(t);
+    if (arb_is_zero(acb_imagref(y)))
+    {
+        arb_submul(acb_imagref(z), acb_imagref(x), acb_realref(y), prec);
+        arb_submul(acb_realref(z), acb_realref(x), acb_realref(y), prec);
+    }
+    else if (arb_is_zero(acb_imagref(x)))
+    {
+        arb_submul(acb_imagref(z), acb_imagref(y), acb_realref(x), prec);
+        arb_submul(acb_realref(z), acb_realref(y), acb_realref(x), prec);
+    }
+    else
+    {
+        acb_t t;
+        acb_init(t);
+        acb_mul(t, x, y, prec);
+        acb_sub(z, z, t, prec);
+        acb_clear(t);
+    }
 }
 
 static __inline__ void
@@ -681,35 +707,17 @@ _acb_vec_sub(acb_ptr res, acb_srcptr vec1, acb_srcptr vec2, long len, long prec)
 static __inline__ void
 _acb_vec_scalar_submul(acb_ptr res, acb_srcptr vec, long len, const acb_t c, long prec)
 {
-    if (len > 0)
-    {
-        long i;
-        acb_t t;
-        acb_init(t);
-        for (i = 0; i < len; i++)
-        {
-            acb_mul(t, vec + i, c, prec);
-            acb_sub(res + i, res + i, t, prec);
-        }
-        acb_clear(t);
-    }
+    long i;
+    for (i = 0; i < len; i++)
+        acb_submul(res + i, vec + i, c, prec);
 }
 
 static __inline__ void
 _acb_vec_scalar_addmul(acb_ptr res, acb_srcptr vec, long len, const acb_t c, long prec)
 {
-    if (len > 0)
-    {
-        long i;
-        acb_t t;
-        acb_init(t);
-        for (i = 0; i < len; i++)
-        {
-            acb_mul(t, vec + i, c, prec);
-            acb_add(res + i, res + i, t, prec);
-        }
-        acb_clear(t);
-    }
+    long i;
+    for (i = 0; i < len; i++)
+        acb_addmul(res + i, vec + i, c, prec);
 }
 
 static __inline__ void

arb_poly.h

@@ -27,6 +27,7 @@
 #define ARB_POLY_H
 
 #include "arb.h"
+#include "acb.h"
 #include "fmpz_poly.h"
 #include "fmpq_poly.h"
 
@@ -526,29 +527,28 @@ void arb_poly_compose_series_brent_kung(arb_poly_t res,
                     const arb_poly_t poly1,
                     const arb_poly_t poly2, long n, long prec);
 
+
+void _arb_poly_evaluate_acb_horner(acb_t res, arb_srcptr f, long len, const acb_t x, long prec);
+void arb_poly_evaluate_acb_horner(acb_t res, const arb_poly_t f, const acb_t a, long prec);
+
+void _arb_poly_evaluate_acb_rectangular(acb_t y, arb_srcptr poly, long len, const acb_t x, long prec);
+void arb_poly_evaluate_acb_rectangular(acb_t res, const arb_poly_t f, const acb_t a, long prec);
+
+void _arb_poly_evaluate_acb(acb_t res, arb_srcptr f, long len, const acb_t x, long prec);
+void arb_poly_evaluate_acb(acb_t res, const arb_poly_t f, const acb_t a, long prec);
+
+void _arb_poly_evaluate2_acb_horner(acb_t y, acb_t z, arb_srcptr f, long len, const acb_t x, long prec);
+void arb_poly_evaluate2_acb_horner(acb_t y, acb_t z, const arb_poly_t f, const acb_t x, long prec);
+
+void _arb_poly_evaluate2_acb_rectangular(acb_t y, acb_t z, arb_srcptr f, long len, const acb_t x, long prec);
+void arb_poly_evaluate2_acb_rectangular(acb_t y, acb_t z, const arb_poly_t f, const acb_t x, long prec);
+
+void _arb_poly_evaluate2_acb(acb_t y, acb_t z, arb_srcptr f, long len, const acb_t x, long prec);
+void arb_poly_evaluate2_acb(acb_t y, acb_t z, const arb_poly_t f, const acb_t x, long prec);
+
 /*
 TBD:
 
-void _arb_poly_evaluate_fmpcb_horner(fmpcb_t res, arb_srcptr f, long len, const fmpcb_t x, long prec);
-void arb_poly_evaluate_fmpcb_horner(fmpcb_t res, const arb_poly_t f, const fmpcb_t a, long prec);
-
-void _arb_poly_evaluate_fmpcb_rectangular(fmpcb_t y, arb_srcptr poly, long len, const fmpcb_t x, long prec);
-void arb_poly_evaluate_fmpcb_rectangular(fmpcb_t res, const arb_poly_t f, const fmpcb_t a, long prec);
-
-void _arb_poly_evaluate_fmpcb(fmpcb_t res, arb_srcptr f, long len, const fmpcb_t x, long prec);
-void arb_poly_evaluate_fmpcb(fmpcb_t res, const arb_poly_t f, const fmpcb_t a, long prec);
-
-void _arb_poly_evaluate2_fmpcb_horner(fmpcb_t y, fmpcb_t z, arb_srcptr f, long len, const fmpcb_t x, long prec);
-void arb_poly_evaluate2_fmpcb_horner(fmpcb_t y, fmpcb_t z, const arb_poly_t f, const fmpcb_t x, long prec);
-
-void _arb_poly_evaluate2_fmpcb_rectangular(fmpcb_t y, fmpcb_t z, arb_srcptr f, long len, const fmpcb_t x, long prec);
-void arb_poly_evaluate2_fmpcb_rectangular(fmpcb_t y, fmpcb_t z, const arb_poly_t f, const fmpcb_t x, long prec);
-
-void _arb_poly_evaluate2_fmpcb(fmpcb_t y, fmpcb_t z, arb_srcptr f, long len, const fmpcb_t x, long prec);
-void arb_poly_evaluate2_fmpcb(fmpcb_t y, fmpcb_t z, const arb_poly_t f, const fmpcb_t x, long prec);
-
-
-
 void _arb_poly_gamma_series(arb_ptr res, arb_srcptr h, long hlen, long len, long prec);
 
 void arb_poly_gamma_series(arb_poly_t res, const arb_poly_t f, long n, long prec);

arb_poly/evaluate2_acb.c

@@ -0,0 +1,39 @@
+/*=============================================================================
+
+    This file is part of ARB.
+
+    ARB is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    ARB is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with ARB; if not, write to the Free Software
+    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+=============================================================================*/
+/******************************************************************************
+
+    Copyright (C) 2013 Fredrik Johansson
+
+******************************************************************************/
+
+#include "arb_poly.h"
+
+void
+_arb_poly_evaluate2_acb(acb_t y, acb_t z, arb_srcptr f, long len, const acb_t x, long prec)
+{
+    _arb_poly_evaluate2_acb_rectangular(y, z, f, len, x, prec);
+}
+
+void
+arb_poly_evaluate2_acb(acb_t r, acb_t s, const arb_poly_t f, const acb_t a, long prec)
+{
+    _arb_poly_evaluate2_acb(r, s, f->coeffs, f->length, a, prec);
+}
+

arb_poly/evaluate2_acb_horner.c

@@ -0,0 +1,87 @@
+/*=============================================================================
+
+    This file is part of ARB.
+
+    ARB is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    ARB is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with ARB; if not, write to the Free Software
+    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+=============================================================================*/
+/******************************************************************************
+
+    Copyright (C) 2013 Fredrik Johansson
+
+******************************************************************************/
+
+#include "arb_poly.h"
+
+void
+_arb_poly_evaluate2_acb_horner(acb_t y, acb_t z,
+    arb_srcptr poly, long len, const acb_t x, long prec)
+{
+    if (len == 0)
+    {
+        acb_zero(y);
+        acb_zero(z);
+    }
+    else if (len == 1)
+    {
+        acb_set_round_arb(y, poly + 0, prec);
+        acb_zero(z);
+    }
+    else if (acb_is_zero(x))
+    {
+        acb_set_round_arb(y, poly + 0, prec);
+        acb_set_round_arb(z, poly + 1, prec);
+    }
+    else if (len == 2)
+    {
+        acb_mul_arb(y, x, poly + 1, prec);
+        acb_add_arb(y, y, poly + 0, prec);
+        acb_set_round_arb(z, poly + 1, prec);
+    }
+    else
+    {
+        acb_t t, u, v;
+        long i;
+
+        acb_init(t);
+        acb_init(u);
+        acb_init(v);
+
+        acb_set_round_arb(u, poly + len - 1, prec);
+        acb_zero(v);
+
+        for (i = len - 2; i >= 0; i--)
+        {
+            acb_mul(t, v, x, prec);
+            acb_add(v, u, t, prec);
+            acb_mul(t, u, x, prec);
+            acb_add_arb(u, t, poly + i, prec);
+        }
+
+        acb_swap(y, u);
+        acb_swap(z, v);
+
+        acb_clear(t);
+        acb_clear(u);
+        acb_clear(v);
+    }
+}
+
+void
+arb_poly_evaluate2_acb_horner(acb_t r, acb_t s, const arb_poly_t f, const acb_t a, long prec)
+{
+    _arb_poly_evaluate2_acb_horner(r, s, f->coeffs, f->length, a, prec);
+}
+

arb_poly/evaluate2_acb_rectangular.c

@@ -0,0 +1,120 @@
+/*=============================================================================
+
+    This file is part of ARB.
+
+    ARB is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    ARB is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with ARB; if not, write to the Free Software
+    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+=============================================================================*/
+/******************************************************************************
+
+    Copyright (C) 2013 Fredrik Johansson
+
+******************************************************************************/
+
+#include "arb_poly.h"
+
+void
+_arb_poly_evaluate2_acb_rectangular(acb_t y, acb_t z,
+    arb_srcptr poly, long len, const acb_t x, long prec)
+{
+    long i, j, m, r;
+    acb_ptr xs;
+    acb_t s, t;
+    arb_t c;
+
+    if (len < 3)
+    {
+        if (len == 0)
+        {
+            acb_zero(y);
+            acb_zero(z);
+        }
+        else if (len == 1)
+        {
+            acb_set_round_arb(y, poly + 0, prec);
+            acb_zero(z);
+        }
+        else if (len == 2)
+        {
+            acb_mul_arb(y, x, poly + 1, prec);
+            acb_add_arb(y, y, poly + 0, prec);
+            acb_set_round_arb(z, poly + 1, prec);
+        }
+        return;
+    }
+
+    m = n_sqrt(len) + 1;
+    m *= 1;
+
+    r = (len + m - 1) / m;
+
+    xs = _acb_vec_init(m + 1);
+    acb_init(s);
+    acb_init(t);
+    arb_init(c);
+
+    _acb_vec_set_powers(xs, x, m + 1, prec);
+
+    acb_set_arb(y, poly + (r - 1) * m);
+    for (j = 1; (r - 1) * m + j < len; j++)
+        acb_addmul_arb(y, xs + j, poly + (r - 1) * m + j, prec);
+
+    for (i = r - 2; i >= 0; i--)
+    {
+        acb_set_arb(s, poly + i * m);
+        for (j = 1; j < m; j++)
+            acb_addmul_arb(s, xs + j, poly + i * m + j, prec);
+
+        acb_mul(y, y, xs + m, prec);
+        acb_add(y, y, s, prec);
+    }
+
+    len -= 1;
+    r = (len + m - 1) / m;
+    arb_mul_ui(acb_realref(z), poly + (r - 1) * m + 1, (r - 1) * m + 1, FMPR_PREC_EXACT);
+    arb_zero(acb_imagref(z));
+    for (j = 1; (r - 1) * m + j < len; j++)
+    {
+        arb_mul_ui(c, poly + (r - 1) * m + j + 1, (r - 1) * m + j + 1, FMPR_PREC_EXACT);
+        acb_addmul_arb(z, xs + j, c, prec);
+    }
+
+    for (i = r - 2; i >= 0; i--)
+    {
+        arb_mul_ui(acb_realref(s), poly + i * m + 1, i * m + 1, FMPR_PREC_EXACT);
+        arb_zero(acb_imagref(s));
+
+        for (j = 1; j < m; j++)
+        {
+            arb_mul_ui(c, poly + i * m + j + 1, i * m + j + 1, FMPR_PREC_EXACT);
+            acb_addmul_arb(s, xs + j, c, prec);
+        }
+
+        acb_mul(z, z, xs + m, prec);
+        acb_add(z, z, s, prec);
+    }
+
+    _acb_vec_clear(xs, m + 1);
+    acb_clear(s);
+    acb_clear(t);
+    arb_clear(c);
+}
+
+void
+arb_poly_evaluate2_acb_rectangular(acb_t r, acb_t s, const arb_poly_t f, const acb_t a, long prec)
+{
+    _arb_poly_evaluate2_acb_rectangular(r, s, f->coeffs, f->length, a, prec);
+}
+

arb_poly/evaluate_acb.c

@@ -0,0 +1,40 @@
+/*=============================================================================
+
+    This file is part of ARB.
+
+    ARB is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    ARB is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with ARB; if not, write to the Free Software
+    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+=============================================================================*/
+/******************************************************************************
+
+    Copyright (C) 2013 Fredrik Johansson
+
+******************************************************************************/
+
+#include "arb_poly.h"
+
+void
+_arb_poly_evaluate_acb(acb_t res, arb_srcptr f, long len,
+                           const acb_t x, long prec)
+{
+    _arb_poly_evaluate_acb_rectangular(res, f, len, x, prec);
+}
+
+void
+arb_poly_evaluate_acb(acb_t res, const arb_poly_t f, const acb_t a, long prec)
+{
+    _arb_poly_evaluate_acb(res, f->coeffs, f->length, a, prec);
+}
+

arb_poly/evaluate_acb_horner.c

@@ -0,0 +1,73 @@
+/*=============================================================================
+
+    This file is part of ARB.
+
+    ARB is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    ARB is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with ARB; if not, write to the Free Software
+    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+=============================================================================*/
+/******************************************************************************
+
+    Copyright (C) 2010 Sebastian Pancratz
+    Copyright (C) 2012 Fredrik Johansson
+
+******************************************************************************/
+
+#include "arb_poly.h"
+
+void
+_arb_poly_evaluate_acb_horner(acb_t y, arb_srcptr f, long len,
+                           const acb_t x, long prec)
+{
+    if (len == 0)
+    {
+        acb_zero(y);
+    }
+    else if (len == 1 || acb_is_zero(x))
+    {
+        acb_set_round_arb(y, f, prec);
+    }
+    else if (len == 2)
+    {
+        acb_mul_arb(y, x, f + 1, prec);
+        acb_add_arb(y, y, f + 0, prec);
+    }
+    else
+    {
+        long i = len - 1;
+        acb_t t, u;
+
+        acb_init(t);
+        acb_init(u);
+        acb_set_arb(u, f + i);
+
+        for (i = len - 2; i >= 0; i--)
+        {
+            acb_mul(t, u, x, prec);
+            acb_add_arb(u, t, f + i, prec);
+        }
+
+        acb_swap(y, u);
+
+        acb_clear(t);
+        acb_clear(u);
+    }
+}
+
+void
+arb_poly_evaluate_acb_horner(acb_t res, const arb_poly_t f, const acb_t a, long prec)
+{
+    _arb_poly_evaluate_acb_horner(res, f->coeffs, f->length, a, prec);
+}
+

arb_poly/evaluate_acb_rectangular.c

@@ -0,0 +1,89 @@
+/*=============================================================================
+
+    This file is part of ARB.
+
+    ARB is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    ARB is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with ARB; if not, write to the Free Software
+    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+=============================================================================*/
+/******************************************************************************
+
+    Copyright (C) 2013 Fredrik Johansson
+
+******************************************************************************/
+
+#include "arb_poly.h"
+
+void
+_arb_poly_evaluate_acb_rectangular(acb_t y, arb_srcptr poly,
+    long len, const acb_t x, long prec)
+{
+    long i, j, m, r;
+    acb_ptr xs;
+    acb_t s, t, c;
+
+    if (len < 3)
+    {
+        if (len == 0)
+        {
+            acb_zero(y);
+        }
+        else if (len == 1)
+        {
+            acb_set_round_arb(y, poly + 0, prec);
+        }
+        else if (len == 2)
+        {
+            acb_mul_arb(y, x, poly + 1, prec);
+            acb_add_arb(y, y, poly + 0, prec);
+        }
+        return;
+    }
+
+    m = n_sqrt(len) + 1;
+    r = (len + m - 1) / m;
+
+    xs = _acb_vec_init(m + 1);
+    acb_init(s);
+    acb_init(t);
+    acb_init(c);
+
+    _acb_vec_set_powers(xs, x, m + 1, prec);
+
+    acb_set_arb(y, poly + (r - 1) * m);
+    for (j = 1; (r - 1) * m + j < len; j++)
+        acb_addmul_arb(y, xs + j, poly + (r - 1) * m + j, prec);
+
+    for (i = r - 2; i >= 0; i--)
+    {
+        acb_set_arb(s, poly + i * m);
+        for (j = 1; j < m; j++)
+            acb_addmul_arb(s, xs + j, poly + i * m + j, prec);
+
+        acb_mul(y, y, xs + m, prec);
+        acb_add(y, y, s, prec);
+    }
+
+    _acb_vec_clear(xs, m + 1);
+    acb_clear(s);
+    acb_clear(t);
+    acb_clear(c);
+}
+
+void
+arb_poly_evaluate_acb_rectangular(acb_t res, const arb_poly_t f, const acb_t a, long prec)
+{
+    _arb_poly_evaluate_acb_rectangular(res, f->coeffs, f->length, a, prec);
+}
+

arb_poly/test/t-evaluate2_acb_rectangular.c

@@ -0,0 +1,87 @@
+/*=============================================================================
+
+    This file is part of ARB.
+
+    ARB is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    ARB is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with ARB; if not, write to the Free Software
+    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+=============================================================================*/
+/******************************************************************************
+
+    Copyright (C) 2012, 2013 Fredrik Johansson
+
+******************************************************************************/
+
+#include "arb_poly.h"
+
+int main()
+{
+    long iter;
+    flint_rand_t state;
+
+    printf("evaluate2_acb_rectangular....");
+    fflush(stdout);
+
+    flint_randinit(state);
+
+    for (iter = 0; iter < 10000; iter++)
+    {
+        arb_poly_t f, g;
+        acb_t x, y1, z1, y2, z2;
+
+        acb_init(x);
+        acb_init(y1);
+        acb_init(z1);
+        acb_init(y2);
+        acb_init(z2);
+        arb_poly_init(f);
+        arb_poly_init(g);
+
+        acb_randtest(x, state, 2 + n_randint(state, 1000), 5);
+        arb_poly_randtest(f, state, 2 + n_randint(state, 100), 2 + n_randint(state, 1000), 5);
+        arb_poly_derivative(g, f, 2 + n_randint(state, 1000));
+
+        arb_poly_evaluate2_acb_rectangular(y1, z1, f, x, 2 + n_randint(state, 1000));
+
+        arb_poly_evaluate_acb_horner(y2, f, x, 2 + n_randint(state, 1000));
+        arb_poly_evaluate_acb_horner(z2, g, x, 2 + n_randint(state, 1000));
+
+        if (!acb_overlaps(y1, y2) || !acb_overlaps(z1, z2))
+        {
+            printf("FAIL\n\n");
+            printf("f = "); arb_poly_printd(f, 15); printf("\n\n");
+            printf("g = "); arb_poly_printd(g, 15); printf("\n\n");
+            printf("x = "); acb_printd(x, 15); printf("\n\n");
+            printf("y1 = "); acb_printd(y1, 15); printf("\n\n");
+            printf("z1 = "); acb_printd(z1, 15); printf("\n\n");
+            printf("y2 = "); acb_printd(y2, 15); printf("\n\n");
+            printf("z2 = "); acb_printd(z2, 15); printf("\n\n");
+            abort();
+        }
+
+        arb_poly_clear(f);
+        arb_poly_clear(g);
+        acb_clear(x);
+        acb_clear(y1);
+        acb_clear(z1);
+        acb_clear(y2);
+        acb_clear(z2);
+    }
+
+    flint_randclear(state);
+    flint_cleanup();
+    printf("PASS\n");
+    return EXIT_SUCCESS;
+}
+

arb_poly/test/t-evaluate_acb_rectangular.c

@@ -0,0 +1,75 @@
+/*=============================================================================
+
+    This file is part of ARB.
+
+    ARB is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    ARB is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with ARB; if not, write to the Free Software
+    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+=============================================================================*/
+/******************************************************************************
+
+    Copyright (C) 2012, 2013 Fredrik Johansson
+
+******************************************************************************/
+
+#include "arb_poly.h"
+
+int main()
+{
+    long iter;
+    flint_rand_t state;
+
+    printf("evaluate_acb_rectangular....");
+    fflush(stdout);
+
+    flint_randinit(state);
+
+    for (iter = 0; iter < 10000; iter++)
+    {
+        arb_poly_t f;
+        acb_t x, y1, y2;
+
+        acb_init(x);
+        acb_init(y1);
+        acb_init(y2);
+        arb_poly_init(f);
+
+        acb_randtest(x, state, 2 + n_randint(state, 1000), 5);
+        arb_poly_randtest(f, state, 2 + n_randint(state, 100), 2 + n_randint(state, 1000), 5);
+
+        arb_poly_evaluate_acb_rectangular(y1, f, x, 2 + n_randint(state, 1000));
+        arb_poly_evaluate_acb_horner(y2, f, x, 2 + n_randint(state, 1000));
+
+        if (!acb_overlaps(y1, y2))
+        {
+            printf("FAIL\n\n");
+            printf("f = "); arb_poly_printd(f, 15); printf("\n\n");
+            printf("x = "); acb_printd(x, 15); printf("\n\n");
+            printf("y1 = "); acb_printd(y1, 15); printf("\n\n");
+            printf("y2 = "); acb_printd(y2, 15); printf("\n\n");
+            abort();
+        }
+
+        arb_poly_clear(f);
+        acb_clear(x);
+        acb_clear(y1);
+        acb_clear(y2);
+    }
+
+    flint_randclear(state);
+    flint_cleanup();
+    printf("PASS\n");
+    return EXIT_SUCCESS;
+}
+