1 /**
2  * @file ropt_str.c
3  * Basic structs used in ropt.
4  */
5 
6 
7 #include "cado.h" // IWYU pragma: keep
8 #include <stdint.h>
9 #include <stdio.h>
10 #include <stdlib.h> // free malloc exit
11 #include <float.h> // DBL_MAX
12 #include <math.h> // fmin
13 #include <gmp.h>
14 #include "cado_poly.h"
15 #include "mpz_poly.h"
16 #include "ropt_param.h"        // NUM_SUBLATTICE_PRIMES
17 #include "ropt_str.h"
18 #include "ropt_arith.h" // ab2uv
19 #include "auxiliary.h"  // ALG_SIDE RAT_SIDE
20 #include "area.h"
21 #include "auxiliary.h"  // L2_skew_lognorm // IWYU pragma: keep
22 #include "size_optimization.h"
23 
24 
25 /* -----------------*/
26 /* Static Functions */
27 /* -----------------*/
28 
29 
30 /**
31  * Find bound V for constant rotation.
32  */
33 static inline double
rotate_bounds_V_mpz(mpz_t * f,mpz_t * g,int d,ropt_bound_t bound)34 rotate_bounds_V_mpz ( mpz_t *f,
35                       mpz_t *g,
36                       int d,
37                       ropt_bound_t bound )
38 {
39   mpz_t V;
40   mpz_poly F0, G0, F, G;
41   double lognorm, exp_E = DBL_MAX, min_exp_E = DBL_MAX;
42   mpz_poly_init (F0, d);
43   mpz_poly_init (F, d);
44   mpz_poly_init (G0, 1);
45   mpz_poly_init (G, 1);
46   mpz_poly_setcoeffs (F0, f, d);
47   mpz_poly_setcoeffs (G0, g, 1);
48 
49   mpz_poly_ptr fptr, gptr;
50   cado_poly poly;
51   cado_poly_init (poly);
52   fptr = poly->pols[ALG_SIDE];
53   fptr->deg = d;
54   gptr = poly->pols[RAT_SIDE];
55   gptr->deg = 1;
56   for (int i = 0; i < (d+1); i++)
57     mpz_set(poly->pols[ALG_SIDE]->coeff[i], f[i]);
58   for (int i = 0; i < 2; i++)
59     mpz_set(poly->pols[RAT_SIDE]->coeff[i], g[i]);
60 
61 
62   /* look for positive V: 2, 4, 8, ... */
63   mpz_init_set_ui (V, 1);
64   for (unsigned int i = 0; i < 150; i++, mpz_mul_2exp (V, V, 1) )
65   {
66     /* F = F0 + V*G0 */
67     mpz_poly_rotation (F, F0, G0, V, 0);
68 
69     /* translation-optimize the rotated polynomial and compute exp_E */
70     sopt_local_descent (F, G, F, G0, 1, -1, SOPT_DEFAULT_MAX_STEPS, 0);
71     lognorm = L2_skew_lognorm (F, SKEWNESS_DEFAULT_PREC);
72     exp_E = lognorm + expected_rotation_gain (F, G);
73     min_exp_E = fmin(min_exp_E, exp_E);
74 #if 0
75     gmp_printf (" i: %u, translation: %Zd --", i, V);
76     printf (" lognorm: %f, exp_apha: %f, exp_E: %f\n", lognorm,
77             expected_rotation_gain (F, G),  exp_E);
78 #endif
79     if (exp_E > bound->bound_E) break;
80     if (lognorm > bound->bound_lognorm) break;
81 
82   }
83   mpz_set (bound->global_v_boundr, V);
84 
85   /* look for negative V: -2, -4, -8, ... */
86   mpz_set_si (V, -1);
87   for (unsigned int i = 0; i < 150; i++, mpz_mul_2exp (V, V, 1))
88   {
89     /* F = F0 + V*G0 */
90     mpz_poly_rotation (F, F0, G0, V, 0);
91 
92     /* translation-optimize the rotated polynomial */
93     sopt_local_descent (F, G, F, G0, 1, -1, SOPT_DEFAULT_MAX_STEPS, 0);
94     lognorm = L2_skew_lognorm (F, SKEWNESS_DEFAULT_PREC);
95     exp_E = lognorm + expected_rotation_gain (F, G);
96     min_exp_E = fmin(min_exp_E, exp_E);
97 #if 0
98     gmp_printf (" i: %u, translation: %Zd --", i, V);
99     printf (" lognorm: %f, exp_apha: %f, exp_E: %f\n", lognorm,
100             expected_rotation_gain (F, G),  exp_E);
101 #endif
102     if (exp_E > bound->bound_E) break;
103     if (lognorm > bound->bound_lognorm) break;
104   }
105   mpz_set (bound->global_v_boundl, V);
106 
107   cado_poly_clear (poly);
108   mpz_poly_clear (F0);
109   mpz_poly_clear (F);
110   mpz_poly_clear (G0);
111   mpz_poly_clear (G);
112   mpz_clear (V);
113   return min_exp_E;
114 }
115 
116 
117 /**
118  * Find bound U for linear rotation.
119  */
120 static inline double
rotate_bounds_U_lu(mpz_t * f,mpz_t * g,int d,ropt_bound_t bound)121 rotate_bounds_U_lu ( mpz_t *f,
122                      mpz_t *g,
123                      int d,
124                      ropt_bound_t bound )
125 {
126   mpz_poly F0, G0, F, G;
127   double lognorm, exp_E = DBL_MAX, min_exp_E = DBL_MAX;
128   mpz_poly_init (F0, d);
129   mpz_poly_init (F, d);
130   mpz_poly_init (G0, 1);
131   mpz_poly_init (G, 1);
132   mpz_poly_setcoeffs (F0, f, d);
133   mpz_poly_setcoeffs (G0, g, 1);
134 
135   /* Look for positive u: 1, 2, 4, 8, ... */
136   int64_t u = 1;
137   for (unsigned int i = 0; i < 62; i++, u *= 2)
138   {
139     /* F = F0 + u*x*G0 */
140     mpz_poly_rotation_int64 (F, F0, G0, u, 1);
141 
142     /* translation-optimize the rotated polynomial */
143     sopt_local_descent (F, G, F, G0, 1, -1, SOPT_DEFAULT_MAX_STEPS, 0);
144     lognorm = L2_skew_lognorm (F, SKEWNESS_DEFAULT_PREC);
145     exp_E = lognorm + expected_rotation_gain (F, G);
146     min_exp_E = fmin(min_exp_E, exp_E);
147 #if 0
148     printf (" i: %u, translation: %ld --", i, u);
149     printf (" lognorm: %f, exp_apha: %f, exp_E: %f\n", lognorm,
150             expected_rotation_gain (F, G),  exp_E);
151 #endif
152     if (exp_E > bound->bound_E) break;
153     if (lognorm > bound->bound_lognorm) break;
154   }
155   bound->global_u_boundr = u;
156 
157   /* Look for negative u: -1, -2, -4, -8, ... */
158   u = -1;
159   for (unsigned int i = 0; i < 62; i++, u *= 2)
160   {
161     /* F = F0 + u*x*G0 */
162     mpz_poly_rotation_int64 (F, F0, G0, u, 1);
163 
164 
165     /* translation-optimize the rotated polynomial */
166     sopt_local_descent (F, G, F, G0, 1, -1, SOPT_DEFAULT_MAX_STEPS, 0);
167     lognorm = L2_skew_lognorm (F, SKEWNESS_DEFAULT_PREC);
168     exp_E = lognorm + expected_rotation_gain (F, G);
169     min_exp_E = fmin(min_exp_E, exp_E);
170 #if 0
171     printf (" i: %u, translation: %ld --", i, u);
172     printf (" lognorm: %f, exp_apha: %f, exp_E: %f\n", lognorm,
173             expected_rotation_gain (F, G),  exp_E);
174 #endif
175     if (exp_E > bound->bound_E) break;
176     if (lognorm > bound->bound_lognorm) break;
177   }
178   bound->global_u_boundl = u;
179 
180   mpz_poly_clear (F0);
181   mpz_poly_clear (F);
182   mpz_poly_clear (G0);
183   mpz_poly_clear (G);
184   return exp_E;
185 }
186 
187 
188 /**
189  * Find bound W for quadratic rotation.
190  */
191 static inline void
rotate_bounds_W_lu(ropt_poly_t poly,ropt_bound_t bound)192 rotate_bounds_W_lu ( ropt_poly_t poly,
193                      ropt_bound_t bound )
194 {
195   mpz_poly F0, G0, F, G;
196 
197   mpz_poly_init (F0, poly->d);
198   mpz_poly_init (F, poly->d);
199   mpz_poly_init (G0, 1);
200   mpz_poly_init (G, 1);
201   mpz_poly_setcoeffs (F0, poly->f, poly->d);
202   mpz_poly_setcoeffs (G0, poly->g, 1);
203 
204   /* look for positive w: 0, 1, 2, ... */
205   int64_t w = 0;
206   for (int64_t w = 0; w < 4096; w++)
207   {
208     /* F = F0 + w*x^2*G0 */
209     mpz_poly_rotation_int64 (F, F0, G0, w, 2);
210 
211     /* translation-optimize the rotated polynomial */
212     sopt_local_descent (F, G, F, G0, 1, -1, SOPT_DEFAULT_MAX_STEPS, 0);
213     double lognorm = L2_skew_lognorm (F, SKEWNESS_DEFAULT_PREC);
214 
215     if (lognorm > bound->bound_lognorm)
216       break;
217   }
218   bound->global_w_boundr = w;
219 
220   for (int64_t w = 0; w > -4096; w--)
221   {
222     /* F = F0 + w*x^2*G0 */
223     mpz_poly_rotation_int64 (F, F0, G0, w, 2);
224 
225     /* translation-optimize the rotated polynomial */
226     sopt_local_descent (F, G, F, G0, 1, -1, SOPT_DEFAULT_MAX_STEPS, 0);
227     double lognorm = L2_skew_lognorm (F, SKEWNESS_DEFAULT_PREC);
228 
229     if (lognorm > bound->bound_lognorm)
230       break;
231   }
232   bound->global_w_boundl = w;
233 
234   mpz_poly_clear (F0);
235   mpz_poly_clear (F);
236   mpz_poly_clear (G0);
237   mpz_poly_clear (G);
238 }
239 
240 
241 /* -----------------*/
242 /* Public Functions */
243 /* -----------------*/
244 
245 
246 /**
247  * Init ropt_poly_t.
248  */
249 void
ropt_poly_init(ropt_poly_t poly)250 ropt_poly_init ( ropt_poly_t poly )
251 {
252   unsigned int i;
253   mpz_init (poly->n);
254   mpz_init (poly->m);
255 
256   /* fx, gx holds pre-computed values f(r), g(r) where 0 <= r < p. */
257   poly->f = (mpz_t*) malloc ((MAX_DEGREE + 1) * sizeof (mpz_t));
258   poly->g = (mpz_t*) malloc ((MAX_DEGREE + 1) * sizeof (mpz_t));
259   (poly->fx) = (mpz_t *) malloc ((primes[ROPT_NPRIMES-1]+1) * sizeof (mpz_t));
260   (poly->gx) = (mpz_t *) malloc ((primes[ROPT_NPRIMES-1]+1) * sizeof (mpz_t));
261   (poly->numerator) = (mpz_t *) malloc ((primes[ROPT_NPRIMES-1]+1) * sizeof (mpz_t));
262 
263   if ( (poly->f == NULL) || (poly->g == NULL) ||
264        ((poly->fx) == NULL) || ((poly->gx) == NULL) ||
265        ((poly->numerator) == NULL) ) {
266     fprintf (stderr, "Error, cannot allocate memory for polynomial"
267              " coefficients in ropt_poly_init().\n");
268     exit(1);
269   }
270 
271   for (i = 0; i <= MAX_DEGREE; i++) {
272     mpz_init (poly->f[i]);
273     mpz_init (poly->g[i]);
274   }
275 
276   for (i = 0; i <= primes[ROPT_NPRIMES-1]; i++) {
277     mpz_init (poly->fx[i]);
278     mpz_init (poly->gx[i]);
279     mpz_init (poly->numerator[i]);
280   }
281 }
282 
283 
284 /**
285  * clean coefficients
286  */
287 void
ropt_poly_refresh(ropt_poly_t poly)288 ropt_poly_refresh ( ropt_poly_t poly )
289 {
290   unsigned int i;
291   mpz_set_ui (poly->n, 0);
292   mpz_set_ui (poly->m, 0);
293   for (i = 0; i <= MAX_DEGREE; i++) {
294     mpz_set_ui (poly->f[i], 0);
295     mpz_set_ui (poly->g[i], 0);
296   }
297   for (i = 0; i <= primes[ROPT_NPRIMES-1]; i++) {
298     mpz_set_ui (poly->fx[i], 0);
299     mpz_set_ui (poly->gx[i], 0);
300     mpz_set_ui (poly->numerator[i], 0);
301   }
302 }
303 
304 
305 /**
306  * Evaluation polynomials at many points.
307  */
308 static inline void
ropt_poly_setup_eval(mpz_t * fr,mpz_t * gr,mpz_t * numerator,mpz_poly_srcptr f,mpz_poly_srcptr g,const unsigned int * p)309 ropt_poly_setup_eval (mpz_t *fr, mpz_t *gr, mpz_t *numerator, mpz_poly_srcptr f,
310                       mpz_poly_srcptr g, const unsigned int *p)
311 {
312   unsigned long i;
313   mpz_t tmp;
314   mpz_init (tmp);
315 
316   for ( i = 0; i <= p[ROPT_NPRIMES-1]; i ++ ) {
317     mpz_set_ui (fr[i], 0);
318     mpz_set_ui (gr[i], 0);
319     mpz_poly_eval_ui (fr[i], f, i);
320     mpz_poly_eval_ui (gr[i], g, i);
321     mpz_poly_eval_diff_ui (numerator[i], f, i);
322     mpz_mul (numerator[i], numerator[i], gr[i]);
323     mpz_neg (numerator[i], numerator[i]);
324     mpz_mul (tmp, fr[i], g->coeff[1]);
325     mpz_add (numerator[i], tmp, numerator[i]);
326   }
327 
328   mpz_clear (tmp);
329 }
330 
331 
332 /**
333  * check if g[0]/g[1] is a root of f (mod n)
334  */
335 bool
ropt_poly_setup_check(ropt_poly_t poly)336 ropt_poly_setup_check ( ropt_poly_t poly )
337 {
338   int i;
339   mpz_t t;
340   mpz_init (t);
341 
342   /* degree */
343   for ( (poly->d) = MAX_DEGREE;
344         (poly->d) > 0 && mpz_cmp_ui ((poly->f[poly->d]), 0) == 0;
345         poly->d -- );
346 
347   /* m = -Y0/Y1 mod n */
348   mpz_invert (poly->m, poly->g[1], poly->n);
349   mpz_neg (poly->m, poly->m);
350   mpz_mul (poly->m, poly->m, poly->g[0]);
351   mpz_mod (poly->m, poly->m, poly->n);
352 
353   /* check if m is a root of f mod n */
354   mpz_set (t, poly->f[poly->d]);
355   for (i = poly->d - 1; i >= 0; i --) {
356     mpz_mul (t, t, poly->m);
357     mpz_add (t, t, poly->f[i]);
358     mpz_mod (t, t, poly->n);
359   }
360 
361   if (mpz_cmp_ui (t, 0) != 0) {
362     mpz_clear (t);
363     return 0;
364   }
365   mpz_clear (t);
366   return 1;
367 }
368 
369 
370 /**
371  * Precompute fx, gx and numerator in ropt_poly_t. Note: poly->f,
372  * poly->g, poly->d, poly->n must be set in prior.
373  * This function can be called to reset poly after rotation.
374  */
375 void
ropt_poly_setup(ropt_poly_t poly)376 ropt_poly_setup ( ropt_poly_t poly )
377 {
378   if (!ropt_poly_setup_check (poly)) {
379     fprintf (stderr, "ERROR: The following polynomial have no common"
380              " root. \n");
381     print_cadopoly_fg (stderr, poly->f, poly->d, poly->g, 1, poly->n);
382     exit (1);
383   }
384 
385   mpz_poly F, G;
386   F->coeff = poly->f;
387   F->deg = poly->d;
388   G->coeff = poly->g;
389   G->deg = 1;
390 
391   /* pre-compute f(r) for all r < B */
392   ropt_poly_setup_eval (poly->fx, poly->gx, poly->numerator,
393                         F, G, primes);
394   /* projective alpha */
395   poly->alpha_proj = get_alpha_projective (F, get_alpha_bound ());
396 }
397 
398 
399 /**
400  * Free ropt_poly_t.
401  */
402 void
ropt_poly_free(ropt_poly_t poly)403 ropt_poly_free ( ropt_poly_t poly )
404 {
405   unsigned int i;
406 
407   for (i = 0; i <= primes[ROPT_NPRIMES-1]; i ++) {
408     mpz_clear(poly->fx[i]);
409     mpz_clear(poly->gx[i]);
410     mpz_clear(poly->numerator[i]);
411   }
412 
413   for (i = 0; i <= MAX_DEGREE; i++) {
414     mpz_clear (poly->f[i]);
415     mpz_clear (poly->g[i]);
416   }
417 
418   mpz_clear (poly->n);
419   mpz_clear (poly->m);
420 
421   free (poly->fx);
422   free (poly->gx);
423   free (poly->numerator);
424   free (poly->f);
425   free (poly->g);
426 }
427 
428 
429 /**
430  * Init ropt_bound_t.
431  */
432 void
ropt_bound_init(ropt_bound_t bound)433 ropt_bound_init ( ropt_bound_t bound )
434 {
435   bound->global_w_boundl = 0;
436   bound->global_w_boundr = 0;
437   bound->global_u_boundl = 0;
438   bound->global_u_boundr = 0;
439   mpz_init_set_ui (bound->global_v_boundl, 0UL);
440   mpz_init_set_ui (bound->global_v_boundr, 0UL);
441   bound->init_lognorm = 0.0;
442   bound->bound_lognorm = DBL_MAX;
443   bound->bound_E = 0.0;
444   bound->exp_min_alpha = 0.0;
445 }
446 
447 
448 /**
449  * Subroutine for ropt_bound_setup().
450  */
451 static inline void
ropt_bound_setup_normbound(ropt_poly_t poly,ropt_bound_t bound,ropt_param_t param,double incr)452 ropt_bound_setup_normbound ( ropt_poly_t poly,
453                              ropt_bound_t bound,
454                              ropt_param_t param,
455                              double incr )
456 {
457   mpz_poly F, G;
458   double skewness;
459   F->coeff = poly->f;
460   F->deg = poly->d;
461   G->coeff = poly->g;
462   G->deg = 1;
463   skewness = L2_skewness (F, SKEWNESS_DEFAULT_PREC);
464   bound->init_lognorm = L2_lognorm (F, skewness);
465   bound->bound_E = (bound->init_lognorm + expected_rotation_gain (F, G))
466     * incr;
467   /* setup lognorm bound, either from input or by default
468      By default, we do not use bound_lognorm to bound the
469      rotation; instead we use (lognorm + exp_E)          */
470   if (param->bound_lognorm > 0)
471     bound->bound_lognorm = param->bound_lognorm;
472 }
473 
474 
475 #if 0
476 double
477 ropt_bound_expected_E (mpz_poly F, mpz_poly G)
478 {
479   mpz_t *f = F->coeff;
480   mpz_t *g = G->coeff;
481   unsigned int d = F->deg;
482   ropt_bound_t bound;
483   double exp_E;
484 
485   ASSERT_ALWAYS(G->deg == 1);
486   ropt_bound_init (bound);
487   F->coeff = f;
488   F->deg = d;
489   bound->init_lognorm = L2_skew_lognorm (F, SKEWNESS_DEFAULT_PREC);
490   bound->bound_lognorm = bound->init_lognorm * BOUND_LOGNORM_INCR_MAX;
491   exp_E = rotate_bounds_U_lu (f, g, d, bound);
492   exp_E = fmin(exp_E, rotate_bounds_V_mpz (f, g, d, bound));
493   ropt_bound_free (bound);
494   return exp_E;
495 }
496 #endif
497 
498 
499 /**
500  * Subroutine for ropt_bound_setup().
501  */
502 static inline void
ropt_bound_setup_globalbound(ropt_poly_t poly,ropt_bound_t bound,ropt_param_t param)503 ropt_bound_setup_globalbound ( ropt_poly_t poly,
504                                ropt_bound_t bound,
505                                ropt_param_t param )
506 {
507   /* w bound */
508   if (poly->d == 6) {
509     if (param->w_length > 0) {
510       bound->global_w_boundl = param->w_left_bound;
511       bound->global_w_boundr = param->w_left_bound + param->w_length - 1;
512     }
513     else
514       rotate_bounds_W_lu (poly, bound);
515   }
516   else {
517     bound->global_w_boundr = 0;
518     bound->global_w_boundl = 0;
519     param->w_left_bound = 0;
520     param->w_length = 1;
521   }
522 
523   /* u bound */
524   rotate_bounds_U_lu (poly->f, poly->g, poly->d, bound);
525 
526   /* v bound */
527   rotate_bounds_V_mpz (poly->f, poly->g, poly->d, bound);
528 }
529 
530 
531 /**
532  * Subroutine for ropt_bound_setup().
533  */
534 static inline void
ropt_bound_setup_others(ropt_poly_t poly,ropt_bound_t bound)535 ropt_bound_setup_others ( ropt_poly_t poly,
536                           ropt_bound_t bound )
537 {
538   mpz_poly F, G;
539   mpz_poly_init (F, poly->d);
540   mpz_poly_init (G, 1);
541   mpz_poly_setcoeffs (F, poly->f, poly->d);
542   mpz_poly_setcoeffs (G, poly->g, 1);
543   bound->exp_min_alpha = expected_rotation_gain (F, G);
544   mpz_poly_clear (F);
545   mpz_poly_clear (G);
546 }
547 
548 
549 /**
550  * Setup ropt_bound_t (independent of manually-input param).
551  * Note, this function should be called in the very beginning
552  * before doing any rotation, since the init_lognorm parameter
553  * will be set to decide the rotation range in the rest.
554  */
555 void
ropt_bound_setup(ropt_poly_t poly,ropt_bound_t bound,ropt_param_t param,double incr)556 ropt_bound_setup ( ropt_poly_t poly,
557                    ropt_bound_t bound,
558                    ropt_param_t param,
559                    double incr )
560 {
561   /* set bound->bound_lognorm */
562   ropt_bound_setup_normbound (poly, bound, param, incr);
563 
564   /* set w, u, v bounds */
565   ropt_bound_setup_globalbound (poly, bound, param);
566 
567   /* set exp_min_alpha */
568   ropt_bound_setup_others (poly, bound);
569 
570   if (param->verbose >= 1) {
571     gmp_fprintf ( stderr, "# Info: global bounds (%d:%d, %ld:%ld, %Zd:%Zd)"
572                   " gives:\n",
573                   bound->global_w_boundl,
574                   bound->global_w_boundr,
575                   bound->global_u_boundl,
576                   bound->global_u_boundr,
577                   bound->global_v_boundl,
578                   bound->global_v_boundr );
579     if (bound->bound_lognorm==DBL_MAX) {
580       gmp_fprintf ( stderr, "# Info: exp_alpha: %.3f, norm bound: DBL_MAX, "
581                     "E bound: %.3f (init. norm: %.3f)\n",
582                     bound->exp_min_alpha,
583                     bound->bound_E,
584                     bound->init_lognorm );
585     }
586     else {
587       gmp_fprintf ( stderr, "# Info: exp_alpha: %.3f, norm bound: %.3f, "
588                     "E bound: %.3f (init. norm: %.3f)\n",
589                     bound->exp_min_alpha,
590                     bound->bound_lognorm,
591                     bound->bound_E,
592                     bound->init_lognorm );
593     }
594   }
595 }
596 
597 
598 /**
599  * This is for tuning function
600  */
601 void
ropt_bound_setup_incr(ropt_poly_t poly,ropt_bound_t bound,ropt_param_t param,double incr)602 ropt_bound_setup_incr ( ropt_poly_t poly,
603                         ropt_bound_t bound,
604                         ropt_param_t param,
605                         double incr )
606 {
607   ropt_bound_setup_normbound (poly, bound, param, incr);
608   ropt_bound_setup_globalbound (poly, bound, param);
609   ropt_bound_setup_others (poly, bound);
610   if (param->verbose >= 2) {
611     gmp_fprintf ( stderr, "# Info: tune (%d:%d, %ld:%ld, %Zd:%Zd)\n",
612                   bound->global_w_boundl,
613                   bound->global_w_boundr,
614                   bound->global_u_boundl,
615                   bound->global_u_boundr,
616                   bound->global_v_boundl,
617                   bound->global_v_boundr );
618   }
619 }
620 
621 
622 /**
623  * For existing rsparam and when rs is changed,
624  */
625 void
ropt_bound_reset(ropt_poly_t poly,ropt_bound_t bound,ropt_param_t param)626 ropt_bound_reset ( ropt_poly_t poly,
627                    ropt_bound_t bound,
628                    ropt_param_t param )
629 {
630 
631   /* u bound */
632   rotate_bounds_U_lu (poly->f, poly->g, poly->d, bound);
633 
634   /* v bound */
635   rotate_bounds_V_mpz (poly->f, poly->g, poly->d, bound);
636 
637   /* set exp_min_alpha */
638   ropt_bound_setup_others (poly, bound);
639 
640   if (param->verbose >= 2) {
641     gmp_fprintf ( stderr, "# Info: reset bounds (%d:%d, %ld:%ld, %Zd:%Zd)"
642                   " gives:\n",
643                   bound->global_w_boundl,
644                   bound->global_w_boundr,
645                   bound->global_u_boundl,
646                   bound->global_u_boundr,
647                   bound->global_v_boundl,
648                   bound->global_v_boundr );
649     gmp_fprintf ( stderr, "# Info: exp_alpha: %.3f, norm bound: %.3f, "
650                   "E bound: %.3f (init norm: %.3f)\n",
651                   bound->exp_min_alpha,
652                   bound->bound_lognorm,
653                   bound->bound_E,
654                   bound->init_lognorm );
655   }
656 }
657 
658 
659 /**
660  * Free ropt_bound_t.
661  */
662 void
ropt_bound_free(ropt_bound_t bound)663 ropt_bound_free ( ropt_bound_t bound )
664 {
665   mpz_clear (bound->global_v_boundl);
666   mpz_clear (bound->global_v_boundr);
667 }
668 
669 
670 /**
671  * Init stage 1 parameters. The default customisation/parameters
672  * happens in ropt_s1param_setup().
673  */
674 void
ropt_s1param_init(ropt_s1param_t s1param)675 ropt_s1param_init ( ropt_s1param_t s1param )
676 {
677   int i;
678 
679   /* will be set either from param (stdin) or by default */
680   s1param->len_e_sl = 0;
681   s1param->tlen_e_sl = 0;
682   s1param->nbest_sl = 0;
683   s1param->nbest_sl_tune = 1;
684   s1param->nbest_sieve = 0;
685   s1param->modbound = 0;
686 
687   /* set to 1 for using quicker, smaller nbest_sl for tuning
688      sublattices */
689   s1param->nbest_sl_tunemode = 0;
690 
691   mpz_init_set_ui (s1param->modulus, 1UL);
692 
693   s1param->e_sl = (unsigned int*)
694     malloc ( NUM_SUBLATTICE_PRIMES * sizeof (unsigned int) );
695 
696   s1param->individual_nbest_sl = (unsigned int*)
697     malloc ( NUM_SUBLATTICE_PRIMES * sizeof (unsigned int) );
698 
699   if (s1param->e_sl == NULL || s1param->individual_nbest_sl == NULL) {
700     fprintf (stderr, "Error, cannot allocate memory in "
701              "ropt_s1param_init().\n");
702     exit (1);
703   }
704 
705   for (i = 0; i < NUM_SUBLATTICE_PRIMES; i ++) {
706     s1param->e_sl[i] = 0;
707     s1param->individual_nbest_sl[i] = 0;
708   }
709 }
710 
711 
712 /**
713  * Helper function to setup "e_sl[]".
714  * Note the factor 16 in the following function seems
715  * important. It has the effect to reduce the bounds
716  * which is better in practice.
717  */
718 static inline void
ropt_s1param_setup_e_sl(ropt_poly_t poly,ropt_s1param_t s1param,ropt_bound_t bound,ropt_param_t param)719 ropt_s1param_setup_e_sl ( ropt_poly_t poly,
720                           ropt_s1param_t s1param,
721                           ropt_bound_t bound,
722                           ropt_param_t param )
723 {
724   unsigned int i, j;
725   unsigned long modbound;
726   mpz_t bound_by_v;
727   mpz_init (bound_by_v);
728 
729   /* check bound_by_u */
730   unsigned long bound_by_u = (unsigned long)
731     (bound->global_u_boundr < (-bound->global_u_boundl)) ?
732     bound->global_u_boundr : (-bound->global_u_boundl);
733 
734   /* check bound_by_v */
735   mpz_fdiv_q_ui (bound_by_v, bound->global_v_boundr,
736                  SIZE_SIEVEARRAY_V_MIN);
737 
738   /* compare two bounds */
739   if ( mpz_cmp_ui (bound_by_v, bound_by_u) < 0 ) {
740     modbound = mpz_get_ui (bound_by_v);
741     if (modbound > bound_by_u / 8)
742       modbound = bound_by_u / 8;
743   }
744   else
745     modbound = bound_by_u / 8;
746 
747   /* adjust for small skewness but large bound (not sure if this is good) */
748   mpz_poly F;
749   F->coeff = poly->f;
750   F->deg = poly->d;
751   double skew = L2_skewness (F, SKEWNESS_DEFAULT_PREC);
752   if ((double) modbound > skew)
753     modbound = (unsigned long) skew;
754 
755   /* find i */
756   for (i = 0; i < NUM_DEFAULT_SUBLATTICE; i ++)
757     if (default_sublattice_prod[i] > modbound)
758       break;
759 
760   /* check if i > limit and fix if so */
761   i = (i >= NUM_DEFAULT_SUBLATTICE) ?
762     (NUM_DEFAULT_SUBLATTICE - 1) : i;
763 
764   /* set e_sl[] from default array */
765   for (j = 0; j < NUM_SUBLATTICE_PRIMES; j++) {
766     s1param->e_sl[j] = default_sublattice_pe[i][j];
767   }
768 
769   /* get mod bound */
770   s1param->modbound = default_sublattice_prod[i];
771 
772   /* overwrite e_sl[] from from stdin, if needed */
773   if (param->s1_num_e_sl != 0) {
774     s1param->modbound = 1;
775     for (i = 0; i < NUM_SUBLATTICE_PRIMES; i++) {
776       s1param->e_sl[i] = param->s1_e_sl[i];
777       s1param->modbound *= param->s1_e_sl[i];
778     }
779   }
780 
781   mpz_clear (bound_by_v);
782 }
783 
784 
785 /**
786  *  Function to setup "individual_nbest_sl[]".
787  */
788 void
ropt_s1param_setup_individual_nbest_sl(ropt_s1param_t s1param)789 ropt_s1param_setup_individual_nbest_sl (ropt_s1param_t s1param)
790 {
791   unsigned int i;
792   for (i = 0; i < s1param->len_e_sl; i ++)
793     s1param->individual_nbest_sl[i] =
794         s1_size_each_sublattice[s1param->tlen_e_sl - 1][i];
795 }
796 
797 
798 /**
799  *  Function to setup shorter "individual_nbest_sl[]".
800  */
801 void
ropt_s1param_setup_individual_nbest_sl_tune(ropt_s1param_t s1param)802 ropt_s1param_setup_individual_nbest_sl_tune (ropt_s1param_t s1param)
803 {
804   unsigned int i;
805   for (i = 0; i < s1param->len_e_sl; i ++)
806     s1param->individual_nbest_sl[i]
807         = s1_size_each_sublattice_tune[i];
808 }
809 
810 
811 /**
812  * Setup s1param by default parameters and/or param (stdin).
813  */
814 void
ropt_s1param_setup(ropt_poly_t poly,ropt_s1param_t s1param,ropt_bound_t bound,ropt_param_t param)815 ropt_s1param_setup ( ropt_poly_t poly,
816                      ropt_s1param_t s1param,
817                      ropt_bound_t bound,
818                      ropt_param_t param )
819 {
820   unsigned int i, j;
821 
822   /* Set 1: "len_e_sl" and "tlen_e_sl" */
823   s1param->len_e_sl = NUM_SUBLATTICE_PRIMES;
824   s1param->tlen_e_sl = s1param->len_e_sl;
825 
826   /* Set 2: "nbest_sl", depending on the size of n */
827   j = mpz_sizeinbase (poly->n, 10);
828   for (i = 0; i < NUM_DEFAULT_DIGITS-1; i ++)
829     if (size_total_sublattices[i][0] > j)
830       break;
831 
832   s1param->nbest_sl = (unsigned int) ((double) size_total_sublattices[i][1]);
833 
834   s1param->nbest_sieve = (unsigned int) ((double) size_total_sublattices[i][2]);
835 
836   s1param->nbest_sl_tune = (unsigned int) ((double) size_total_sublattices[i][3]);
837 
838 
839   if (s1param->nbest_sl < 4)
840     s1param->nbest_sl = 4;
841   if (param->verbose >= 1) {
842     printf ("# Info: s1param->nbest_sl: %u (= size_total_sublattices*ropteffort)\n", s1param->nbest_sl);
843     printf ("# Info: s1param->nbest_sieve: %u (for sieving)\n", s1param->nbest_sieve);
844   }
845 
846   /* Set 3: set "e_sl[]" */
847   ropt_s1param_setup_e_sl (poly, s1param, bound, param);
848 
849   /* Set 4: set "individual_nbest_sl[]" */
850   ropt_s1param_setup_individual_nbest_sl (s1param);
851 
852 }
853 
854 
855 /**
856  * Setup s1param by default parameters and/or param (stdin).
857  */
858 void
ropt_s1param_resetup(ropt_poly_t poly,ropt_s1param_t s1param,ropt_bound_t bound,ropt_param_t param,unsigned int nbest)859 ropt_s1param_resetup ( ropt_poly_t poly,
860                        ropt_s1param_t s1param,
861                        ropt_bound_t bound,
862                        ropt_param_t param,
863                        unsigned int nbest )
864 {
865   unsigned int i, j;
866   s1param->len_e_sl = NUM_SUBLATTICE_PRIMES;
867   s1param->tlen_e_sl = s1param->len_e_sl;
868   j = mpz_sizeinbase (poly->n, 10);
869   for (i = 0; i < NUM_DEFAULT_DIGITS; i ++)
870     if (size_total_sublattices[i][0] > j)
871       break;
872   s1param->nbest_sl = (unsigned long) nbest;
873   s1param->nbest_sieve = (unsigned long) ((double) size_total_sublattices[i][2]);
874   if (s1param->nbest_sl < 4)
875     s1param->nbest_sl = 4;
876 
877   /* Set 3: set "e_sl[]" */
878   ropt_s1param_setup_e_sl (poly, s1param, bound, param);
879 
880   /* Set 4: set "individual_nbest_sl[]" */
881   ropt_s1param_setup_individual_nbest_sl (s1param);
882 
883 }
884 
885 
886 /**
887  * Setup s1param with modbound as extra input
888  */
889 void
ropt_s1param_resetup_modbound(ropt_poly_t poly,ropt_s1param_t s1param,ropt_bound_t bound,ropt_param_t param,unsigned int nbest,unsigned long modbound)890 ropt_s1param_resetup_modbound ( ropt_poly_t poly,
891                                 ropt_s1param_t s1param,
892                                 ropt_bound_t bound,
893                                 ropt_param_t param,
894                                 unsigned int nbest,
895                                 unsigned long modbound)
896 {
897   ropt_s1param_resetup (poly, s1param, bound, param, nbest);
898 
899   /* find i */
900   int i, j;
901   for (i = 0; i < NUM_DEFAULT_SUBLATTICE; i ++)
902     if (default_sublattice_prod[i] > modbound)
903       break;
904 
905   /* check if i > limit and fix if so */
906   i = (i >= NUM_DEFAULT_SUBLATTICE) ?
907     (NUM_DEFAULT_SUBLATTICE - 1) : i;
908 
909   /* set e_sl[] from default array */
910   for (j = 0; j < NUM_SUBLATTICE_PRIMES; j++) {
911     s1param->e_sl[j] = default_sublattice_pe[i][j];
912   }
913 
914   /* get mod bound */
915   s1param->modbound = default_sublattice_prod[i];
916 
917 }
918 
919 
920 /**
921  * Free s1param.
922  */
923 void
ropt_s1param_free(ropt_s1param_t s1param)924 ropt_s1param_free ( ropt_s1param_t s1param )
925 {
926   free(s1param->e_sl);
927   free(s1param->individual_nbest_sl);
928   mpz_clear (s1param->modulus);
929 }
930 
931 
932 /**
933  * Init ropt_s2param.
934  */
935 void
ropt_s2param_init(ropt_poly_t poly,ropt_s2param_t s2param)936 ropt_s2param_init ( ropt_poly_t poly,
937                     ropt_s2param_t s2param )
938 {
939   int i;
940   s2param->len_p_rs = ROPT_NPRIMES - 1;
941   s2param->Amax = s2param->Amin = 0;
942   s2param->Bmax = s2param->Bmin = 0;
943 
944   mpz_init_set_ui (s2param->Umax, 0UL);
945   mpz_init_set_ui (s2param->Umin, 0UL);
946   mpz_init_set_ui (s2param->Vmax, 0UL);
947   mpz_init_set_ui (s2param->Vmin, 0UL);
948   mpz_init_set_ui (s2param->A, 0UL);
949   mpz_init_set_ui (s2param->B, 0UL);
950   mpz_init_set_ui (s2param->MOD, 0UL);
951 
952   s2param->f = (mpz_t*) malloc ((poly->d + 1) * sizeof (mpz_t));
953   s2param->g = (mpz_t*) malloc (2 * sizeof (mpz_t));
954   if (s2param->f == NULL || s2param->g == NULL) {
955     fprintf ( stderr, "Error, cannot allocate memory in "
956               "ropt_s2param_init().\n" );
957     exit (1);
958   }
959 
960   for (i = 0; i <= poly->d; i++)
961     mpz_init (s2param->f[i]);
962   mpz_init (s2param->g[0]);
963   mpz_init (s2param->g[1]);
964 }
965 
966 
967 /**
968  * Free ropt_s2param.
969  */
970 void
ropt_s2param_free(ropt_poly_t poly,ropt_s2param_t s2param)971 ropt_s2param_free ( ropt_poly_t poly,
972                     ropt_s2param_t s2param )
973 {
974   int i;
975 
976   mpz_clear (s2param->Umax);
977   mpz_clear (s2param->Umin);
978   mpz_clear (s2param->Vmax);
979   mpz_clear (s2param->Vmin);
980   mpz_clear (s2param->A);
981   mpz_clear (s2param->B);
982   mpz_clear (s2param->MOD);
983   for (i = 0; i <= poly->d; i++)
984     mpz_clear (s2param->f[i]);
985   mpz_clear (s2param->g[0]);
986   mpz_clear (s2param->g[1]);
987   free (s2param->f);
988   free (s2param->g);
989 }
990 
991 
992 /**
993  * Setup sublattice (u, v), mod and Umax, Umin, Vmax, Vmin.
994  */
995 static inline void
ropt_s2param_setup_sublattice(ropt_s2param_t s2param,mpz_t A,mpz_t B,mpz_t MOD)996 ropt_s2param_setup_sublattice ( ropt_s2param_t s2param,
997                                 mpz_t A,
998                                 mpz_t B,
999                                 mpz_t MOD )
1000 {
1001   mpz_set (s2param->A, A);
1002   mpz_set (s2param->B, B);
1003   /* the s1param->modulus might be not true since rsparam might be
1004      changed for different quadratic rotations. Instead, the true mod
1005      is recorded in the priority queue, now in 'MOD'. */
1006   mpz_set (s2param->MOD, MOD);
1007   ab2uv (s2param->A, s2param->MOD, s2param->Amax, s2param->Umax);
1008   ab2uv (s2param->A, s2param->MOD, s2param->Amin, s2param->Umin);
1009   ab2uv (s2param->B, s2param->MOD, s2param->Bmax, s2param->Vmax);
1010   ab2uv (s2param->B, s2param->MOD, s2param->Bmin, s2param->Vmin);
1011 }
1012 
1013 
1014 /**
1015  * Set sieving region for s2param -> Amax, Amin, Amax, Amin.
1016  */
1017 static inline void
ropt_s2param_setup_range(ropt_bound_t bound,ropt_s2param_t s2param,ropt_param_t param,mpz_t mod)1018 ropt_s2param_setup_range ( ropt_bound_t bound,
1019                            ropt_s2param_t s2param,
1020                            ropt_param_t param,
1021                            mpz_t mod )
1022 {
1023   /* read sieve length from param (stdin) */
1024   if (param->s2_Amax >= 0 && param->s2_Bmax > 0) {
1025     s2param->Amax = param->s2_Amax;
1026     s2param->Bmax = param->s2_Bmax;
1027   }
1028   /* or compute sieve V length */
1029   else {
1030     s2param->Amax = 0;
1031     unsigned long len;
1032     mpz_t q;
1033     mpz_init (q);
1034     mpz_fdiv_q (q, bound->global_v_boundr, mod);
1035     len =  mpz_get_ui (q);
1036 
1037     /* upper bound */
1038     s2param->Bmax = ( (len > SIZE_SIEVEARRAY_V_MAX) ?
1039                       SIZE_SIEVEARRAY_V_MAX : (long) len );
1040 
1041     /* fix if len is too small */
1042     if (s2param->Bmax == 0)
1043       s2param->Bmax = 128;
1044 
1045     mpz_clear (q);
1046   }
1047 
1048   s2param->Bmin = -s2param->Bmax;
1049   s2param->Amin = -s2param->Amax;
1050 }
1051 
1052 
1053 /**
1054  * Set up sieving length and sublattice for s2param.
1055  */
1056 void
ropt_s2param_setup(ropt_bound_t bound,ropt_s1param_t s1param,ropt_s2param_t s2param,ropt_param_t param,mpz_t A,mpz_t B,mpz_t MOD)1057 ropt_s2param_setup ( ropt_bound_t bound,
1058                      ropt_s1param_t s1param,
1059                      ropt_s2param_t s2param,
1060                      ropt_param_t param,
1061                      mpz_t A,
1062                      mpz_t B,
1063                      mpz_t MOD )
1064 {
1065   /* normally we should have more primes than those in the sublattice */
1066   if (s2param->len_p_rs < s1param->len_e_sl) {
1067     fprintf ( stderr, "# Warning: number of primes considered "
1068               "in stage 2 is smaller than that in "
1069               "stage 1. See ropt_s2param_setup().\n" );
1070   }
1071 
1072   /* set sieving length */
1073   ropt_s2param_setup_range (bound, s2param, param, MOD);
1074 
1075   /* set sublattice */
1076   ropt_s2param_setup_sublattice (s2param, A, B, MOD);
1077 }
1078 
1079 
1080 /**
1081  * Set up s2param without s1param.
1082  */
1083 void
ropt_s2param_setup_stage2_only(ropt_bound_t bound,ropt_s2param_t s2param,ropt_param_t param,mpz_t A,mpz_t B,mpz_t MOD)1084 ropt_s2param_setup_stage2_only ( ropt_bound_t bound,
1085                                  ropt_s2param_t s2param,
1086                                  ropt_param_t param,
1087                                  mpz_t A,
1088                                  mpz_t B,
1089                                  mpz_t MOD )
1090 {
1091   s2param->len_p_rs = ROPT_NPRIMES - 1;
1092 
1093   /* set sieving length */
1094   ropt_s2param_setup_range (bound, s2param, param, MOD);
1095 
1096   /* set sublattice */
1097   ropt_s2param_setup_sublattice (s2param, A, B, MOD);
1098 }
1099 
1100 
1101 /**
1102  * Set tune (short) sieving region for s2param -> Amax, Amin, Amax, Amin.
1103  */
1104 static inline void
ropt_s2param_setup_tune_range(ropt_s2param_t s2param,unsigned long Amax,unsigned long Bmax)1105 ropt_s2param_setup_tune_range  ( ropt_s2param_t s2param,
1106                                  unsigned long Amax,
1107                                  unsigned long Bmax )
1108 {
1109   s2param->Amax = (long) Amax;
1110   s2param->Bmax = (long) Bmax;
1111   s2param->Bmin = -s2param->Bmax;
1112   s2param->Amin = -s2param->Amax;
1113 }
1114 
1115 
1116 /**
1117  * Set up tune sieving length and sublattice for s2param.
1118  */
1119 void
ropt_s2param_setup_tune(ropt_s1param_t s1param,ropt_s2param_t s2param,mpz_t A,mpz_t B,mpz_t MOD,unsigned long Amax,unsigned long Bmax,unsigned int len_p_rs)1120 ropt_s2param_setup_tune ( ropt_s1param_t s1param,
1121                           ropt_s2param_t s2param,
1122                           mpz_t A,
1123                           mpz_t B,
1124                           mpz_t MOD,
1125                           unsigned long Amax,
1126                           unsigned long Bmax,
1127                           unsigned int len_p_rs )
1128 {
1129   /* setup s2param->len_p_rs */
1130   if (len_p_rs < s1param->tlen_e_sl) {
1131     fprintf ( stderr, "# Warning: number of primes considered "
1132               "in stage 2 (%d) is smaller than that (%d) in "
1133               "stage 1. See ropt_s2param_setup().\n",
1134               len_p_rs, s1param->tlen_e_sl );
1135   }
1136   s2param->len_p_rs = len_p_rs;
1137 
1138   /* set tune sieving length */
1139   ropt_s2param_setup_tune_range (s2param, Amax, Bmax);
1140 
1141   /* set sublattice */
1142   ropt_s2param_setup_sublattice (s2param, A, B, MOD);
1143 }
1144 
1145 
1146 /**
1147  * Print s2param.
1148  */
1149 void
ropt_s2param_print(ropt_s2param_t s2param)1150 ropt_s2param_print ( ropt_s2param_t s2param )
1151 {
1152   fprintf ( stderr, "# Info: sieving matrix: "
1153             "[%ld, %ld] x [%ld, %ld], len. prime = %d.\n",
1154             s2param->Amin, s2param->Amax,
1155             s2param->Bmin, s2param->Bmax,
1156             s2param->len_p_rs );
1157 
1158   gmp_fprintf ( stderr,
1159                 "# Info: (u, v) = (%Zd + i * %Zd, %Zd + j * %Zd)\n",
1160                 s2param->A, s2param->MOD, s2param->B,
1161                 s2param->MOD );
1162 
1163   gmp_fprintf ( stderr,
1164                 "# Info: (Amin: %4ld, Amax: %4ld) -> "
1165                 "(Umin: %6Zd, Umax: %6Zd)\n",
1166                 s2param->Amin, s2param->Amax,
1167                 s2param->Umin, s2param->Umax );
1168 
1169   gmp_fprintf ( stderr,
1170                 "# Info: (Bmin: %4ld, Bmax: %4ld) -> "
1171                 "(Vmin: %6Zd, Vmax: %6Zd)\n",
1172                 s2param->Bmin, s2param->Bmax,
1173                 s2param->Vmin, s2param->Vmax );
1174 }
1175 
1176 
1177 /**
1178  * Init bestpoly.
1179  */
1180 void
ropt_bestpoly_init(ropt_bestpoly_t bestpoly,int d)1181 ropt_bestpoly_init ( ropt_bestpoly_t bestpoly,
1182                      int d )
1183 {
1184   int i;
1185   bestpoly->f = (mpz_t*) malloc ((d + 1) * sizeof (mpz_t));
1186   bestpoly->g = (mpz_t*) malloc (2 * sizeof (mpz_t));
1187 
1188   if ((bestpoly->f == NULL) || (bestpoly->g == NULL)) {
1189     fprintf ( stderr, "Error, cannot allocate memory for"
1190               " ropt_bestpoly_init()\n" );
1191     exit (1);
1192   }
1193 
1194   for (i = 0; i <= d; i++)
1195     mpz_init (bestpoly->f[i]);
1196   mpz_init (bestpoly->g[0]);
1197   mpz_init (bestpoly->g[1]);
1198 }
1199 
1200 
1201 /**
1202  * Setup bestpoly.
1203  */
1204 void
ropt_bestpoly_setup(ropt_bestpoly_t bestpoly,mpz_t * f,mpz_t * g,int d)1205 ropt_bestpoly_setup ( ropt_bestpoly_t bestpoly,
1206                       mpz_t *f,
1207                       mpz_t *g,
1208                       int d )
1209 {
1210   int i;
1211   for (i = 0; i <= d; i++)
1212     mpz_set (bestpoly->f[i], f[i]);
1213   mpz_set (bestpoly->g[0], g[0]);
1214   mpz_set (bestpoly->g[1], g[1]);
1215 }
1216 
1217 
1218 /**
1219  * Free bestpoly.
1220  */
1221 void
ropt_bestpoly_free(ropt_bestpoly_t bestpoly,int d)1222 ropt_bestpoly_free ( ropt_bestpoly_t bestpoly,
1223                      int d )
1224 {
1225   int i;
1226   for (i = 0; i <= d; i++)
1227     mpz_clear (bestpoly->f[i]);
1228 
1229   mpz_clear (bestpoly->g[0]);
1230   mpz_clear (bestpoly->g[1]);
1231   free (bestpoly->f);
1232   free (bestpoly->g);
1233 }
1234 
1235 
1236 /**
1237  * Init param.
1238  */
1239 void
ropt_param_init(ropt_param_t param)1240 ropt_param_init ( ropt_param_t param )
1241 {
1242   mpz_init (param->s2_u);
1243   mpz_init (param->s2_v);
1244   mpz_init (param->s2_mod);
1245   mpz_init (param->n);
1246   mpz_set_ui (param->s2_u, 0);
1247   mpz_set_ui (param->s2_v, 0);
1248   mpz_set_ui (param->s2_mod, 0);
1249   mpz_set_ui (param->n, 0);
1250   param->w_left_bound = 0;
1251   param->w_length = 0;
1252   param->s1_num_e_sl = 0;
1253   param->s2_Amax = 0;
1254   param->s2_Bmax = 0;
1255   param->s2_w = 0;
1256   param->bound_lognorm = -1;
1257   param->s1_e_sl = (unsigned short*)
1258     malloc ( NUM_SUBLATTICE_PRIMES * sizeof (unsigned short) );
1259   int i;
1260   for (i = 0; i < NUM_SUBLATTICE_PRIMES; i ++)
1261     param->s1_e_sl[i] = 0;
1262   param->d = 0;
1263   param->verbose = 0;
1264   param->effort = DEFAULT_ROPTEFFORT;
1265   param->skip_ropt = 0;
1266   param->gen_raw = 0;
1267   param->sopt = 0;
1268 }
1269 
1270 
1271 /**
1272  * Free param.
1273  */
1274 void
ropt_param_free(ropt_param_t param)1275 ropt_param_free ( ropt_param_t param )
1276 {
1277   mpz_clear (param->s2_u);
1278   mpz_clear (param->s2_v);
1279   mpz_clear (param->s2_mod);
1280   mpz_clear (param->n);
1281   free (param->s1_e_sl);
1282 }
1283 
1284 
1285 /**
1286  * Init info.
1287  */
1288 void
ropt_info_init(ropt_info_t info)1289 ropt_info_init ( ropt_info_t info )
1290 {
1291   info->ave_MurphyE = 0.0;
1292   info->best_MurphyE = 0.0;
1293   info->mode = 0;
1294   info->w = 0;
1295   info->ropt_time_stage1 = 0.0;
1296   info->ropt_time_tuning = 0.0;
1297   info->ropt_time_stage2 = 0.0;
1298 }
1299 
1300 
1301 /**
1302  * Free info.
1303  */
1304 void
ropt_info_free(ropt_info_t info)1305 ropt_info_free ( ropt_info_t info )
1306 {
1307   info->ave_MurphyE = 0.0;
1308   info->best_MurphyE = 0.0;
1309   info->mode = 0;
1310   info->w = 0;
1311 }
1312