1 #include "cado.h" // IWYU pragma: keep
2 #include <cstdio>
3 #include <cstdlib>
4 #include <cstdint>
5 #include <fcntl.h>   /* for _O_BINARY */
6 #include <fstream>
7 #include "cado_poly.h"  // cado_poly
8 #include "filter_io.h"  // earlyparsed_relation_ptr
9 #include "gzip.h"       // set_antebuffer_path
10 #include "misc.h"       // filelist
11 #include "params.h"     // param_list
12 #include "renumber.hpp" // renumber_t
13 #include "rootfinder.h" // mpz_poly_roots_ulong
14 #include "typedefs.h"   // p_r_values_t
15 #include "macros.h"
16 
17 /* This program creates an "appendix" to a renumber table. It's never
18  * used in cado, it seems, so that the code is _most likely_ buggy.
19  */
20 
21 char *argv0; /* = argv[0] */
22 
23 uint64_t *very_large_primes = NULL;
24 uint64_t nb_very_large_primes = 0;
25 uint64_t very_large_primes_alloc = 0;
26 
27 
28 void *
look_for_very_large_primes(void * context_data,earlyparsed_relation_ptr rel)29 look_for_very_large_primes (void * context_data, earlyparsed_relation_ptr rel)
30 {
31 
32   uint64_t *lpb = (uint64_t *) context_data;
33   for (weight_t i = 0; i < rel->nb; i++)
34   {
35     int side = (int) rel->primes[i].side;
36     p_r_values_t p = rel->primes[i].p;
37     if (p > lpb[side])
38     {
39       if (nb_very_large_primes == very_large_primes_alloc)
40       {
41         very_large_primes_alloc = very_large_primes_alloc << 1;
42         very_large_primes = (uint64_t *) realloc
43               (very_large_primes, very_large_primes_alloc * sizeof (uint64_t));
44       }
45       very_large_primes[nb_very_large_primes++] = p;
46     }
47   }
48 
49   return NULL;
50 }
51 
declare_usage(param_list pl)52 static void declare_usage(param_list pl)
53 {
54   param_list_decl_usage(pl, "out", "output file for extra renumbering table");
55   param_list_decl_usage(pl, "poly", "input polynomial file");
56   param_list_decl_usage(pl, "lpb0", "large prime bound on side 0");
57   param_list_decl_usage(pl, "lpb1", "large prime bound on side 1");
58   param_list_decl_usage(pl, "filelist", "file containing a list of input files");
59   param_list_decl_usage(pl, "basepath", "path added to all file in filelist");
60   param_list_decl_usage(pl, "force-posix-threads", "force the use of posix threads, do not rely on platform memory semantics");
61   param_list_decl_usage(pl, "path_antebuffer", "path to antebuffer program");
62 }
63 
64 static void
usage(param_list pl,char * argv0)65 usage (param_list pl, char *argv0)
66 {
67     param_list_print_usage(pl, argv0, stderr);
68     exit(EXIT_FAILURE);
69 }
70 
71 uint64_t
cmp64(void * a,void * b)72 cmp64 (void * a, void *b)
73 {
74   uint64_t u = *((uint64_t *) a);
75   uint64_t v = *((uint64_t *) b);
76   if (u < v) return -1;
77   else if (u > v) return 1;
78   else return 0;
79 }
80 
81 int
main(int argc,char * argv[])82 main (int argc, char *argv[])
83 {
84   argv0 = argv[0];
85   std::vector<unsigned int> lpb_arg(2,0);
86   uint64_t lpb[2] = {0, 0};
87   cado_poly poly;
88 
89   param_list pl;
90   param_list_init(pl);
91   declare_usage(pl);
92   argv++,argc--;
93 
94   param_list_configure_switch(pl, "force-posix-threads", &filter_rels_force_posix_threads);
95 
96 #ifdef HAVE_MINGW
97   _fmode = _O_BINARY;     /* Binary open for all files */
98 #endif
99 
100   if (argc == 0)
101     usage (pl, argv0);
102 
103   for( ; argc ; ) {
104       if (param_list_update_cmdline(pl, &argc, &argv)) { continue; }
105       /* Since we accept file names freeform, we decide to never abort
106        * on unrecognized options */
107       break;
108       // fprintf (stderr, "Unknown option: %s\n", argv[0]);
109       // abort();
110   }
111   /* print command-line arguments */
112   param_list_print_command_line (stdout, pl);
113   fflush(stdout);
114 
115   param_list_parse_uint(pl, "lpb0", &lpb_arg[0]);
116   param_list_parse_uint(pl, "lpb1", &lpb_arg[1]);
117 
118   const char * filelist = param_list_lookup_string(pl, "filelist");
119   const char * basepath = param_list_lookup_string(pl, "basepath");
120   const char * polyfilename = param_list_lookup_string(pl, "poly");
121   const char * outname = param_list_lookup_string(pl, "out");
122   const char * path_antebuffer = param_list_lookup_string(pl, "path_antebuffer");
123 
124   if (param_list_warn_unused(pl))
125   {
126     fprintf(stderr, "Error, unused parameters are given\n");
127     usage(pl, argv0);
128   }
129 
130   if (lpb_arg[0] == 0 || lpb_arg[1] == 0)
131   {
132     fprintf (stderr, "Error, missing -lpb0 or -lpb1 command line argument\n");
133     usage (pl, argv0);
134   }
135   if (outname == NULL)
136   {
137     fprintf (stderr, "Error, missing -out command line argument\n");
138     usage(pl, argv0);
139   }
140   if (polyfilename == NULL)
141   {
142     fprintf (stderr, "Error, missing -poly command line argument\n");
143     usage(pl, argv0);
144   }
145   if (basepath && !filelist)
146   {
147     fprintf(stderr, "Error, -basepath only valid with -filelist\n");
148     usage(pl, argv0);
149   }
150   if ((filelist != NULL) + (argc != 0) != 1) {
151     fprintf(stderr, "Error, provide either -filelist or freeform file names\n");
152     usage(pl, argv0);
153   }
154 
155   set_antebuffer_path (argv0, path_antebuffer);
156 
157   cado_poly_init(poly);
158   if (!cado_poly_read (poly, polyfilename))
159   {
160     fprintf (stderr, "Error reading polynomial file\n");
161     exit (EXIT_FAILURE);
162   }
163 
164   lpb[0] = 1 << lpb_arg[0];
165   lpb[1] = 1 << lpb_arg[1];
166 
167   very_large_primes_alloc = 1024;
168   very_large_primes = (uint64_t *)
169                           malloc (very_large_primes_alloc * sizeof (uint64_t));
170   ASSERT_ALWAYS (very_large_primes != NULL);
171 
172   renumber_t renumber_table(poly);
173   renumber_table.set_lpb(lpb_arg);
174 
175   char ** files = filelist ? filelist_from_file (basepath, filelist, 0) : argv;
176 
177   /* compute the list of very large primes appearing in the relations. */
178   filter_rels(files, (filter_rels_callback_t) &look_for_very_large_primes,
179               (void *) lpb, EARLYPARSE_NEED_AB_DECIMAL | EARLYPARSE_NEED_PRIMES,
180               NULL, NULL);
181 
182   qsort (very_large_primes, nb_very_large_primes, sizeof (uint64_t),
183                                                         (__compar_fn_t) cmp64);
184 
185   std::ofstream r_out(outname);
186   renumber_table.write_header(r_out);
187   /* true, we didn't compute the bad ideals. But we need to stick in
188    * there the proper format info saying "no bad ideals there"
189    */
190   renumber_table.write_bad_ideals(r_out);
191 
192   unsigned long prev = 0;
193   std::vector<std::vector<unsigned long>> roots;
194 
195   for (uint64_t i = 0; i < nb_very_large_primes; i++)
196   {
197       unsigned long p = very_large_primes[i];
198       if (p == prev) continue;
199       printf ("%lu\n", p);
200       for (int side = 0; side < 2; side++) {
201           mpz_poly_srcptr f = poly->pols[side];
202           roots.push_back(mpz_poly_roots(f, p));
203           // Check for a projective root
204           if (mpz_divisible_ui_p (f->coeff[f->deg], p))
205               roots.back().push_back(p);
206       }
207 
208       auto cook = renumber_table.cook(p, roots);
209       renumber_table.use_cooked(cook);
210       r_out << cook.text;
211       prev = p;
212   }
213 
214 
215   if (filelist)
216     filelist_clear(files);
217 
218   cado_poly_clear (poly);
219   param_list_clear(pl);
220   free (very_large_primes);
221   return 0;
222 }
223