1 /*
2 Copyright (C) 2011 William Hart
3 Copyright (C) 2012 Fredrik Johansson
4
5 This file is part of Arb.
6
7 Arb is free software: you can redistribute it and/or modify it under
8 the terms of the GNU Lesser General Public License (LGPL) as published
9 by the Free Software Foundation; either version 2.1 of the License, or
10 (at your option) any later version. See <http://www.gnu.org/licenses/>.
11 */
12
13 #include "acb_poly.h"
14
15 void
_acb_poly_product_roots(acb_ptr poly,acb_srcptr xs,slong n,slong prec)16 _acb_poly_product_roots(acb_ptr poly, acb_srcptr xs, slong n, slong prec)
17 {
18 if (n == 0)
19 {
20 acb_one(poly);
21 }
22 else if (n == 1)
23 {
24 acb_neg(poly, xs);
25 acb_one(poly + 1);
26 }
27 else if (n == 2)
28 {
29 acb_mul(poly, xs + 0, xs + 1, prec);
30 acb_add(poly + 1, xs + 0, xs + 1, prec);
31 acb_neg(poly + 1, poly + 1);
32 acb_one(poly + 2);
33 }
34 else if (n == 3)
35 {
36 acb_mul(poly + 1, xs, xs + 1, prec);
37 acb_mul(poly, poly + 1, xs + 2, prec);
38 acb_neg(poly, poly);
39 acb_add(poly + 2, xs, xs + 1, prec);
40 acb_addmul(poly + 1, poly + 2, xs + 2, prec);
41 acb_add(poly + 2, poly + 2, xs + 2, prec);
42 acb_neg(poly + 2, poly + 2);
43 acb_one(poly + 3);
44 }
45 else
46 {
47 const slong m = (n + 1) / 2;
48 acb_ptr tmp;
49
50 tmp = _acb_vec_init(n + 2);
51
52 _acb_poly_product_roots(tmp, xs, m, prec);
53 _acb_poly_product_roots(tmp + m + 1, xs + m, n - m, prec);
54 _acb_poly_mul_monic(poly, tmp, m + 1, tmp + m + 1, n - m + 1, prec);
55
56 _acb_vec_clear(tmp, n + 2);
57 }
58 }
59
60 void
acb_poly_product_roots(acb_poly_t poly,acb_srcptr xs,slong n,slong prec)61 acb_poly_product_roots(acb_poly_t poly, acb_srcptr xs, slong n, slong prec)
62 {
63 acb_poly_fit_length(poly, n + 1);
64 _acb_poly_product_roots(poly->coeffs, xs, n, prec);
65 _acb_poly_set_length(poly, n + 1);
66 }
67