1 /*
2     Copyright (C) 2016 Pascal Molin
3 
4     This file is part of Arb.
5 
6     Arb is free software: you can redistribute it and/or modify it under
7     the terms of the GNU Lesser General Public License (LGPL) as published
8     by the Free Software Foundation; either version 2.1 of the License, or
9     (at your option) any later version.  See <http://www.gnu.org/licenses/>.
10 */
11 
12 #include <string.h>
13 #include "dirichlet.h"
14 #include "flint/profiler.h"
15 
16 #define LOG 0
17 #define CSV 1
18 #define JSON 2
19 
20 typedef ulong (*do_f) (ulong q1, ulong q2);
21 
22 static ulong
do_gcd(ulong q1,ulong q2)23 do_gcd(ulong q1, ulong q2)
24 {
25     ulong n, q, k;
26 
27     for (n = 0, q = q1; q <= q2; q++)
28         for (k = 1; k < q; k++)
29             n += (n_gcd(k, q) == 1);
30 
31     return n;
32 }
33 
34 static ulong
do_char(ulong q1,ulong q2)35 do_char(ulong q1, ulong q2)
36 {
37     ulong n, q;
38 
39     for (n = 0, q = q1; q <= q2; q++)
40     {
41         dirichlet_group_t G;
42         dirichlet_char_t x;
43 
44         dirichlet_group_init(G, q);
45         dirichlet_char_init(x, G);
46 
47         dirichlet_char_one(x, G);
48         n++;
49         for (; dirichlet_char_next(x, G) >= 0; n++);
50 
51         dirichlet_char_clear(x);
52         dirichlet_group_clear(G);
53     }
54 
55     return n;
56 }
57 
58 static ulong
do_gcdpluscond(ulong q1,ulong q2)59 do_gcdpluscond(ulong q1, ulong q2)
60 {
61     ulong n, q, k;
62 
63     for (n = 0, q = q1; q <= q2; q++)
64     {
65         dirichlet_group_t G;
66         dirichlet_group_init(G, q);
67 
68         for (k = 1; k < q; k++)
69         {
70             /* known factors -> faster gcd */
71             slong i;
72             if (G->q_even > 1 && k % 2 == 0)
73                 continue;
74             for (i = G->neven; i < G->num; i++)
75                 if (k % G->P[i].p == 0)
76                     break;
77 
78             if (i == G->num)
79                 dirichlet_conductor_ui(G, k);
80         }
81 
82         n += G->phi_q;
83         dirichlet_group_clear(G);
84 
85     }
86 
87     return n;
88 }
89 
90 static ulong
do_charpluscond(ulong q1,ulong q2)91 do_charpluscond(ulong q1, ulong q2)
92 {
93     ulong n, q, k;
94 
95     for (n = 0, q = q1; q <= q2; q++)
96     {
97         dirichlet_group_t G;
98         dirichlet_char_t x;
99 
100         dirichlet_group_init(G, q);
101         dirichlet_char_init(x, G);
102 
103         dirichlet_char_one(x, G);
104         n++;
105 
106         for (; dirichlet_char_next(x, G) >= 0; n++)
107             dirichlet_conductor_char(G, x);
108 
109         dirichlet_char_clear(x);
110         dirichlet_group_clear(G);
111 
112     }
113 
114     return n;
115 }
116 
117 static ulong
do_charplusorder(ulong q1,ulong q2)118 do_charplusorder(ulong q1, ulong q2)
119 {
120     ulong n, q, k;
121 
122     for (n = 0, q = q1; q <= q2; q++)
123     {
124         dirichlet_group_t G;
125         dirichlet_char_t x;
126 
127         dirichlet_group_init(G, q);
128         dirichlet_char_init(x, G);
129 
130         dirichlet_char_one(x, G);
131         n++;
132 
133         for (; dirichlet_char_next(x, G) >= 0; n++)
134             dirichlet_order_char(G, x);
135 
136         dirichlet_char_clear(x);
137         dirichlet_group_clear(G);
138 
139     }
140 
141     return n;
142 }
143 
main(int argc,char * argv[])144 int main(int argc, char *argv[])
145 {
146     int out;
147     ulong n, nref, maxq = 5000;
148 
149     int l, nf = 5;
150     do_f func[5] = {
151         do_gcd,
152         do_char,
153         do_gcdpluscond,
154         do_charpluscond,
155         do_charplusorder
156     };
157     char * name[5] = {
158         "gcd",
159         "char",
160         "gcd + cond",
161         "char + cond",
162         "char + order"
163     };
164 
165     int i, ni = 5;
166     ulong qmin[5] = { 2,   1000, 10000, 100000, 1000000 };
167     ulong qmax[5] = { 500, 3000, 11000, 100100, 1000010 };
168 
169     if (argc < 2)
170         out = LOG;
171     else if (!strcmp(argv[1], "json"))
172         out = JSON;
173     else if (!strcmp(argv[1], "csv"))
174         out = CSV;
175     else if (!strcmp(argv[1], "log"))
176         out = LOG;
177     else
178     {
179         printf("usage: %s [log|csv|json]\n", argv[0]);
180         flint_abort();
181     }
182 
183     if (out == CSV)
184         flint_printf("# %-12s, %7s, %7s, %7s\n","name", "qmin", "qmax", "time");
185 
186     for (i = 0; i < ni; i++)
187     {
188 
189         if (out == LOG)
190         {
191             flint_printf("loop over all (Z/q)* for %wu<=q<=%wu\n", qmin[i], qmax[i]);
192         }
193 
194         for (l = 0; l < nf; l++)
195         {
196 
197             if (out == LOG)
198                 flint_printf("%-14s ...  ",name[l]);
199             else if (out == CSV)
200                 flint_printf("%-12s, %7d, %7d,   ",name[l],qmin[i],qmax[i]);
201             else if (out == JSON)
202                 flint_printf("{ \"name\": \"%s\", \"qmin\": %d, \"qmax\": %d, \"time\": ",
203                         name[l],qmin[i],qmax[i]);
204 
205             TIMEIT_ONCE_START
206                 (func[l])(qmin[i], qmax[i]);
207             TIMEIT_ONCE_STOP
208 
209             if (l == 0)
210                 nref = n;
211             else if (n != nref)
212             {
213                 flint_printf("FAIL: wrong number of elements %wu != %wu\n\n",n, nref);
214                 flint_abort();
215             }
216 
217             if (out == JSON)
218                 flint_printf("}\n");
219             else
220                 flint_printf("\n");
221         }
222 
223     }
224 
225     flint_cleanup();
226     return EXIT_SUCCESS;
227 }
228