1 /*
2 Copyright (C) 2011 Fredrik Johansson
3 Copyright (C) 2013, 2014 Martin Lee
4
5 This file is part of FLINT.
6
7 FLINT 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 #undef ulong
14 #define ulong ulongxx/* interferes with system includes */
15
16 #include <stdlib.h>
17 #include <stdio.h>
18
19 #undef ulong
20
21 #include <gmp.h>
22
23 #define ulong mp_limb_t
24
25 #include "flint.h"
26 #include "nmod_poly.h"
27 #include "ulong_extras.h"
28
29 int
main(void)30 main(void)
31 {
32 int i;
33 FLINT_TEST_INIT(state);
34
35 flint_printf("compose_mod_brent_kung_vec_preinv....");
36 fflush(stdout);
37
38 for (i = 0; i < 100 * flint_test_multiplier(); i++)
39 {
40 nmod_poly_t a, ainv, b, c;
41 mp_limb_t m = n_randtest_prime(state, 0);
42 slong j, k, l;
43 nmod_poly_struct * pow, * res;
44
45 nmod_poly_init(a, m);
46 nmod_poly_init(b, m);
47 nmod_poly_init(c, m);
48 nmod_poly_init(ainv, m);
49
50 nmod_poly_randtest(b, state, 1+n_randint(state, 20));
51 nmod_poly_randtest_not_zero(a, state, 1+n_randint(state, 20));
52 l= n_randint(state, 20) + 1;
53 k= n_randint(state, l ) + 1;
54
55 nmod_poly_rem(b, b, a);
56 nmod_poly_reverse(ainv, a, a->length);
57 nmod_poly_inv_series(ainv, ainv, a->length);
58 pow = (nmod_poly_struct *) flint_malloc((l + k)*sizeof(nmod_poly_struct));
59 res = pow + l;
60
61 for (j = 0; j < l; j++)
62 {
63 nmod_poly_init(pow + j, m);
64 nmod_poly_randtest(pow + j, state, n_randint(state, 20) + 1);
65 nmod_poly_rem(pow + j, pow + j, a);
66 }
67
68 for (j = 0; j < k; j++)
69 nmod_poly_init(res + j, m);
70
71 nmod_poly_compose_mod_brent_kung_vec_preinv(res, pow, l, k, b, a, ainv);
72
73 for (j = 0; j < k; j++)
74 {
75 nmod_poly_compose_mod(c, pow + j, b, a);
76 if (!nmod_poly_equal(res + j, c))
77 {
78 flint_printf("FAIL (composition):\n");
79 flint_printf("a:\n"); nmod_poly_print(a); flint_printf("\n");
80 flint_printf("res:\n"); nmod_poly_print(res + j); flint_printf("\n");
81 flint_printf("pow:\n"); nmod_poly_print(pow + j); flint_printf("\n");
82 flint_printf("b:\n"); nmod_poly_print(b); flint_printf("\n");
83 flint_printf("c:\n"); nmod_poly_print(c); flint_printf("\n");
84 flint_printf("j: %wd\n", j);
85 abort();
86 }
87 }
88
89 nmod_poly_clear(a);
90 nmod_poly_clear(ainv);
91 nmod_poly_clear(b);
92 nmod_poly_clear(c);
93 for (j = 0; j < l; j++)
94 nmod_poly_clear(pow + j);
95 for (j = 0; j < k; j++)
96 nmod_poly_clear(res + j);
97 flint_free(pow);
98 }
99
100
101 FLINT_TEST_CLEANUP(state);
102
103 flint_printf("PASS\n");
104 return 0;
105 }
106