1 /*
2     Copyright (C) 2020 Daniel Schultz
3     This file is part of FLINT.
4     FLINT is free software: you can redistribute it and/or modify it under
5     the terms of the GNU Lesser General Public License (LGPL) as published
6     by the Free Software Foundation; either version 2.1 of the License, or
7     (at your option) any later version.  See <https://www.gnu.org/licenses/>.
8 */
9 
10 #include <stdio.h>
11 #include <stdlib.h>
12 #include "flint.h"
13 #include "fmpz_poly.h"
14 #include "profiler.h"
15 
16 
count_omega(const fmpz_poly_factor_t fac)17 static slong count_omega(const fmpz_poly_factor_t fac)
18 {
19     slong i, omega = 0;
20     for (i = 0; i < fac->num; i++)
21         omega += fac->exp[i];
22     return omega;
23 }
24 
25 
26 int
main(void)27 main(void)
28 {
29     flint_bitcnt_t max_bits;
30     slong it, i, j, omega;
31     ulong tab[11][20];
32     timeit_t timer;
33     FLINT_TEST_INIT(state);
34 
35     printf("        bits : ");
36 
37     omega = 0;
38     for (max_bits = 10, it = 0; max_bits <= 3000; max_bits *= 2, it++)
39     {
40         fmpz_t d;
41         fmpz_poly_t g;
42         fmpz_poly_struct * f;
43         fmpz_poly_factor_t fac;
44         slong reps = 1000000/(max_bits + 100);
45         slong outer_reps = 10;
46 
47         f = flint_malloc(reps * sizeof(fmpz_poly_struct));
48         for (i = 0; i < reps; i++)
49             fmpz_poly_init(f + i);
50         fmpz_poly_init(g);
51         fmpz_poly_factor_init(fac);
52         fmpz_init(d);
53 
54         tab[0][it] = max_bits;
55         flint_printf(" %06wu", tab[0][it]);
56         fflush(stdout);
57 
58         /* (degree 1)^2 */
59         for (i = 0; i < reps; i++)
60         {
61             do {
62                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
63             } while (g->length != 2);
64             fmpz_poly_pow(f + i, g, 2);
65         }
66         timeit_start(timer);
67         for (j = 0; j < outer_reps; j++)
68         for (i = 0; i < reps; i++)
69             fmpz_poly_factor(fac, f + i);
70         timeit_stop(timer);
71         tab[1][it] = timer->cpu;
72 
73         /* (degree 1)*(degree 1) */
74         for (i = 0; i < reps; i++)
75         {
76             do {
77                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
78             } while (g->length != 2);
79             fmpz_poly_set(f + i, g);
80             do {
81                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
82             } while (g->length != 2);
83             fmpz_poly_mul(f + i, f + i, g);
84         }
85         timeit_start(timer);
86         for (i = 0; i < reps; i++)
87         {
88             fmpz_poly_factor(fac, f + i);
89             omega += count_omega(fac);
90         }
91         for (j = 0; j < outer_reps; j++)
92         for (i = 0; i < reps; i++)
93             fmpz_poly_factor(fac, f + i);
94         timeit_stop(timer);
95         tab[2][it] = timer->cpu;
96 
97         /* (degree 2 +) */
98         for (i = 0; i < reps; i++)
99         {
100             do {
101                fmpz_poly_randtest(g, state, 3, n_randint(state, max_bits) + 2);
102             } while (g->length != 3 || (fmpz_poly_discriminant(d, g), fmpz_sgn(d) < 0));
103             fmpz_poly_set(f + i, g);
104         }
105         timeit_start(timer);
106         for (i = 0; i < reps; i++)
107         {
108             fmpz_poly_factor(fac, f + i);
109             omega += count_omega(fac);
110         }
111         for (j = 0; j < outer_reps; j++)
112         for (i = 0; i < reps; i++)
113             fmpz_poly_factor(fac, f + i);
114         timeit_stop(timer);
115         tab[3][it] = timer->cpu;
116 
117         /* (degree 2 -) */
118         for (i = 0; i < reps; i++)
119         {
120             do {
121                fmpz_poly_randtest(g, state, 3, n_randint(state, max_bits) + 2);
122             } while (g->length != 3 || (fmpz_poly_discriminant(d, g), fmpz_sgn(d) > 0));
123             fmpz_poly_set(f + i, g);
124         }
125         timeit_start(timer);
126         for (i = 0; i < reps; i++)
127         {
128             fmpz_poly_factor(fac, f + i);
129             omega += count_omega(fac);
130         }
131         for (j = 0; j < outer_reps; j++)
132         for (i = 0; i < reps; i++)
133             fmpz_poly_factor(fac, f + i);
134         timeit_stop(timer);
135         tab[4][it] = timer->cpu;
136 
137         /* (degree 1)^3 */
138         for (i = 0; i< reps; i++)
139         {
140             do {
141                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
142             } while (g->length != 2);
143             fmpz_poly_pow(f + i, g, 3);
144         }
145         timeit_start(timer);
146         for (i = 0; i < reps; i++)
147         {
148             fmpz_poly_factor(fac, f + i);
149             omega += count_omega(fac);
150         }
151         for (j = 0; j < outer_reps; j++)
152         for (i = 0; i < reps; i++)
153             fmpz_poly_factor(fac, f + i);
154         timeit_stop(timer);
155         tab[5][it] = timer->cpu;
156 
157         /* (degree 1)*(degree 1)^2 */
158         for (i = 0; i < reps; i++)
159         {
160             do {
161                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
162             } while (g->length != 2);
163             fmpz_poly_pow(f + i, g, 2);
164             do {
165                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
166             } while (g->length != 2);
167             fmpz_poly_mul(f + i, f + i, g);
168         }
169         timeit_start(timer);
170         for (i = 0; i < reps; i++)
171         {
172             fmpz_poly_factor(fac, f + i);
173             omega += count_omega(fac);
174         }
175         for (j = 0; j < outer_reps; j++)
176         for (i = 0; i < reps; i++)
177             fmpz_poly_factor(fac, f + i);
178         timeit_stop(timer);
179         tab[6][it] = timer->cpu;
180 
181         /* (degree 1)*(degree 1)*(degree 1) */
182         for (i = 0; i < reps; i++)
183         {
184             do {
185                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
186             } while (g->length != 2);
187             fmpz_poly_set(f + i, g);
188             do {
189                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
190             } while (g->length != 2);
191             fmpz_poly_mul(f + i, f + i, g);
192             do {
193                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
194             } while (g->length != 2);
195             fmpz_poly_mul(f + i, f + i, g);
196         }
197         timeit_start(timer);
198         for (i = 0; i < reps; i++)
199         {
200             fmpz_poly_factor(fac, f + i);
201             omega += count_omega(fac);
202         }
203         for (j = 0; j < outer_reps; j++)
204         for (i = 0; i < reps; i++)
205             fmpz_poly_factor(fac, f + i);
206         timeit_stop(timer);
207         tab[7][it] = timer->cpu;
208 
209         /* (degree 1)*(degree 2 +) */
210         for (i = 0; i < reps; i++)
211         {
212             do {
213                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
214             } while (g->length != 2);
215             fmpz_poly_set(f + i, g);
216             do {
217                 fmpz_poly_randtest(g, state, 3, n_randint(state, max_bits) + 2);
218             } while (g->length != 3 || (fmpz_poly_discriminant(d, g), fmpz_sgn(d) < 0));
219             fmpz_poly_mul(f + i, f + i, g);
220         }
221         timeit_start(timer);
222         for (i = 0; i < reps; i++)
223         {
224             fmpz_poly_factor(fac, f + i);
225             omega += count_omega(fac);
226         }
227         for (j = 0; j < outer_reps; j++)
228         for (i = 0; i < reps; i++)
229             fmpz_poly_factor(fac, f + i);
230         timeit_stop(timer);
231         tab[8][it] = timer->cpu;
232 
233         /* (degree 1)*(degree 2 -) */
234         for (i = 0; i < reps; i++)
235         {
236             do {
237                fmpz_poly_randtest(g, state, 2, n_randint(state, max_bits) + 2);
238             } while (g->length != 2);
239             fmpz_poly_set(f + i, g);
240             do {
241                 fmpz_poly_randtest(g, state, 3, n_randint(state, max_bits) + 2);
242             } while (g->length != 3 || (fmpz_poly_discriminant(d, g), fmpz_sgn(d) > 0));
243             fmpz_poly_mul(f + i, f + i, g);
244         }
245         timeit_start(timer);
246         for (i = 0; i < reps; i++)
247         {
248             fmpz_poly_factor(fac, f + i);
249             omega += count_omega(fac);
250         }
251         for (j = 0; j < outer_reps; j++)
252         for (i = 0; i < reps; i++)
253             fmpz_poly_factor(fac, f + i);
254         timeit_stop(timer);
255         tab[9][it] = timer->cpu;
256 
257         /* (degree 3) */
258         for (i = 0; i < reps; i++)
259         {
260             do {
261                fmpz_poly_randtest(g, state, 4, n_randint(state, max_bits) + 2);
262             } while (g->length != 4);
263             fmpz_poly_set(f + i, g);
264         }
265         timeit_start(timer);
266         for (i = 0; i < reps; i++)
267         {
268             fmpz_poly_factor(fac, f + i);
269             omega += count_omega(fac);
270         }
271         for (j = 0; j < outer_reps; j++)
272         for (i = 0; i < reps; i++)
273             fmpz_poly_factor(fac, f + i);
274         timeit_stop(timer);
275         tab[10][it] = timer->cpu;
276 
277         for (i = 0; i < reps; i++)
278             fmpz_poly_clear(f + i);
279         flint_free(f);
280         fmpz_poly_clear(g);
281         fmpz_poly_factor_clear(fac);
282         fmpz_clear(d);
283     }
284 
285     flint_printf("\n", max_bits);
286 
287     printf("---------------");
288     for (i = 0; i < it; i++)
289         flint_printf("-------");
290     printf("\n");
291 
292     printf("      (d1)^2 : ");
293     for (i = 0; i < it; i++)
294         flint_printf(" %06wu", tab[1][i]);
295     printf("\n");
296 
297     printf("    (d1)(d1) : ");
298     for (i = 0; i < it; i++)
299         flint_printf(" %06wu", tab[2][i]);
300     printf("\n");
301 
302     printf("       (d2+) : ");
303     for (i = 0; i < it; i++)
304         flint_printf(" %06wu", tab[3][i]);
305     printf("\n");
306 
307     printf("       (d2-) : ");
308     for (i = 0; i < it; i++)
309         flint_printf(" %06wu", tab[4][i]);
310     printf("\n");
311 
312     printf("      (d1)^3 : ");
313     for (i = 0; i < it; i++)
314         flint_printf(" %06wu", tab[5][i]);
315     printf("\n");
316 
317     printf("  (d1)(d1)^2 : ");
318     for (i = 0; i < it; i++)
319         flint_printf(" %06wu", tab[6][i]);
320     printf("\n");
321 
322     printf("(d1)(d1)(d1) : ");
323     for (i = 0; i < it; i++)
324         flint_printf(" %06wu", tab[7][i]);
325     printf("\n");
326 
327     printf("   (d1)(d2+) : ");
328     for (i = 0; i < it; i++)
329         flint_printf(" %06wu", tab[8][i]);
330     printf("\n");
331 
332     printf("   (d1)(d2-) : ");
333     for (i = 0; i < it; i++)
334         flint_printf(" %06wu", tab[9][i]);
335     printf("\n");
336 
337     printf("        (d3) : ");
338     for (i = 0; i < it; i++)
339         flint_printf(" %06wu", tab[10][i]);
340     printf("\n");
341 
342     FLINT_TEST_CLEANUP(state);
343 
344     flint_printf("factors found: %wd\n", omega);
345     return 0;
346 }
347