Line data Source code
1 : /* Auxiliary routines for the ecm library.
2 :
3 : Copyright 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2011 Paul Zimmermann,
4 : Alexander Kruppa.
5 :
6 : This file is part of the ECM Library.
7 :
8 : The ECM Library is free software; you can redistribute it and/or modify
9 : it under the terms of the GNU Lesser General Public License as published by
10 : the Free Software Foundation; either version 3 of the License, or (at your
11 : option) any later version.
12 :
13 : The ECM Library is distributed in the hope that it will be useful, but
14 : WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15 : or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16 : License for more details.
17 :
18 : You should have received a copy of the GNU Lesser General Public License
19 : along with the ECM Library; see the file COPYING.LIB. If not, see
20 : http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
21 : 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
22 :
23 : /* need stdio.h and stdarg.h for gmp.h to declare gmp_vfprintf */
24 : #include <stdio.h>
25 : #include <stdarg.h>
26 : #include <gmp.h>
27 : #include "ecm-impl.h"
28 :
29 : #if TIME_WITH_SYS_TIME
30 : # include <sys/time.h>
31 : # include <time.h>
32 : #else
33 : # if HAVE_SYS_TIME_H
34 : # include <sys/time.h>
35 : # else
36 : # include <time.h>
37 : # endif
38 : #endif
39 :
40 : #ifdef HAVE_LIMITS_H
41 : # include <limits.h>
42 : #else
43 : # ifndef ULONG_MAX
44 : # define LONG_MAX (__GMP_ULONG_MAX / 2)
45 : # endif
46 : #endif
47 :
48 : #ifdef HAVE_STDINT
49 : #include <stdint.h>
50 : #else
51 : /* size_t is an unsigned integer so this ought to work */
52 : #ifndef SIZE_MAX
53 : #define SIZE_MAX (~((size_t) 0))
54 : #endif
55 : #endif
56 :
57 : #define VERBOSE __ECM(verbose)
58 : static int VERBOSE = OUTPUT_NORMAL;
59 :
60 : void
61 2034 : mpz_add_si (mpz_t r, mpz_t s, long i)
62 : {
63 2034 : if (i >= 0)
64 1336 : mpz_add_ui (r, s, (unsigned long) i);
65 : else
66 698 : mpz_sub_ui (r, s, (unsigned long) (-i));
67 2034 : }
68 :
69 : void
70 30752 : mpz_sub_si (mpz_t r, mpz_t s, long i)
71 : {
72 30752 : if (i >= 0)
73 400 : mpz_sub_ui (r, s, (unsigned long) i);
74 : else
75 30352 : mpz_add_ui (r, s, (unsigned long) (-i));
76 30752 : }
77 :
78 : /* Divide RS by 3 */
79 : void
80 24318463 : mpz_divby3_1op (mpz_t RS)
81 : {
82 24318463 : mp_size_t abssize = mpz_size (RS);
83 :
84 24318463 : if (abssize == 0)
85 0 : return;
86 :
87 24318463 : mpn_divexact_by3 (RS->_mp_d, RS->_mp_d, abssize);
88 :
89 24318463 : if (RS->_mp_d[abssize - 1] == 0)
90 502202 : RS->_mp_size -= mpz_sgn (RS);
91 : }
92 :
93 : /* Convert a double d to a size_t.
94 : Assumes d >= 0. If d > SIZE_MAX, returns SIZE_MAX. */
95 : size_t
96 58 : double_to_size (double d)
97 : {
98 : ASSERT(d >= 0.0);
99 58 : return (d > (double) SIZE_MAX) ? SIZE_MAX : (size_t) d;
100 : }
101 :
102 : /* cputime () gives the elapsed time in milliseconds */
103 :
104 : #if defined (_WIN32)
105 : /* First case - GetProcessTimes () is the only known way of getting process
106 : * time (as opposed to calendar time) under mingw32 */
107 :
108 : #include <windows.h>
109 :
110 : long
111 : cputime ()
112 : {
113 : FILETIME lpCreationTime, lpExitTime, lpKernelTime, lpUserTime;
114 : ULARGE_INTEGER n;
115 :
116 : HANDLE hProcess = GetCurrentProcess();
117 :
118 : GetProcessTimes (hProcess, &lpCreationTime, &lpExitTime, &lpKernelTime,
119 : &lpUserTime);
120 :
121 : /* copy FILETIME to a ULARGE_INTEGER as recommended by MSDN docs */
122 : n.u.LowPart = lpUserTime.dwLowDateTime;
123 : n.u.HighPart = lpUserTime.dwHighDateTime;
124 :
125 : /* lpUserTime is in units of 100 ns. Return time in milliseconds */
126 : return (long) (n.QuadPart / 10000);
127 : }
128 :
129 : #elif defined (HAVE_GETRUSAGE)
130 : /* Next case: getrusage () has higher resolution than clock () and so is
131 : preferred. */
132 :
133 : #ifdef HAVE_SYS_TYPES_H
134 : # include <sys/types.h>
135 : #endif
136 : #ifdef HAVE_SYS_RESOURCE_H
137 : # include <sys/resource.h>
138 : #endif
139 :
140 : long
141 463639 : cputime ()
142 : {
143 : struct rusage rus;
144 :
145 463639 : getrusage (RUSAGE_SELF, &rus);
146 : /* This overflows a 32 bit signed int after 2147483s = 24.85 days */
147 463639 : return rus.ru_utime.tv_sec * 1000L + rus.ru_utime.tv_usec / 1000L;
148 : }
149 :
150 : #else
151 : /* Resort to clock (), which on some systems may return calendar time. */
152 :
153 : long
154 : cputime ()
155 : {
156 : /* Return time in milliseconds */
157 : return (long) (clock () * (1000. / (double) CLOCKS_PER_SEC));
158 : }
159 :
160 : #endif /* defining cputime () */
161 :
162 : /* ellapsed time (in milliseconds) between st0 and st1 (values of cputime) */
163 : long
164 422338 : elltime (long st0, long st1)
165 : {
166 422338 : return st1 - st0; /* assumes no wrap around */
167 : }
168 :
169 : /* Get real (wall-clock) time in milliseconds */
170 : long
171 4406 : realtime ()
172 : {
173 : #ifdef HAVE_GETTIMEOFDAY
174 : struct timeval tv;
175 4406 : int ret = gettimeofday (&tv, NULL);
176 4406 : ASSERT_ALWAYS(ret == 0); /* if HAVE_GETTIMEOFDAY, it should be functional */
177 4406 : return (long) tv.tv_sec * 1000L + (long) tv.tv_usec / 1000L;
178 : #else
179 : return 0L;
180 : #endif
181 : }
182 :
183 : /* Tests if loglevel gets printed with the current verbose setting */
184 :
185 : int
186 6730584 : test_verbose (int loglevel)
187 : {
188 6730584 : return (loglevel <= VERBOSE);
189 : }
190 :
191 : void
192 2129 : set_verbose (int v)
193 : {
194 2129 : VERBOSE = v;
195 2129 : }
196 :
197 : int
198 29953752 : outputf (int loglevel, const char *format, ...)
199 : {
200 : va_list ap;
201 29953752 : int n = 0;
202 :
203 29953752 : va_start (ap, format);
204 :
205 29953752 : if (loglevel != OUTPUT_ERROR && loglevel <= VERBOSE)
206 : {
207 499553 : n = gmp_vfprintf (ECM_STDOUT, format, ap);
208 499553 : fflush (ECM_STDOUT);
209 : }
210 29454199 : else if (loglevel == OUTPUT_ERROR)
211 72 : n = gmp_vfprintf (ECM_STDERR, format, ap);
212 :
213 29953752 : va_end (ap);
214 :
215 29953752 : return n;
216 : }
217 :
218 : /* for P-1 and P+1 we have A = y = z = NULL */
219 : void
220 20 : writechkfile (char *chkfilename, int method, double p, mpmod_t modulus,
221 : mpres_t A, mpres_t x, mpres_t y, mpres_t z)
222 : {
223 : FILE *chkfile;
224 : char *methodname;
225 : mpz_t t;
226 :
227 20 : outputf (OUTPUT_VERBOSE, "Writing checkpoint to %s at p = %.0f\n",
228 : chkfilename, p);
229 :
230 20 : switch (method)
231 : {
232 10 : case ECM_ECM : methodname = "ECM"; break;
233 5 : case ECM_PM1 : methodname = "P-1"; break;
234 5 : case ECM_PP1 : methodname = "P+1"; break;
235 0 : default:
236 0 : outputf (OUTPUT_ERROR, "writechkfile: Invalid method\n");
237 0 : return;
238 : }
239 :
240 20 : chkfile = fopen (chkfilename, "w");
241 20 : ASSERT_ALWAYS(chkfile != NULL);
242 :
243 20 : mpz_init (t);
244 :
245 20 : gmp_fprintf (chkfile, "METHOD=%s; B1=%.0f; N=%Zd;",
246 20 : methodname, p, modulus->orig_modulus);
247 20 : mpres_get_z (t, x, modulus);
248 20 : gmp_fprintf (chkfile, " X=0x%Zx;", t);
249 20 : if (method == ECM_ECM)
250 : {
251 10 : if (y != NULL) /* this should mean Weierstrass form */
252 : {
253 : /* actually, we want to print (x:y:1) */
254 0 : mpres_get_z (t, y, modulus);
255 0 : gmp_fprintf (chkfile, " Y=0x%Zx;", t);
256 0 : fprintf (chkfile, " Z=0x1;");
257 : }
258 : else /* one day, we could have some homogeneous form to deal with */
259 : {
260 10 : mpres_get_z (t, z, modulus);
261 10 : gmp_fprintf (chkfile, " Z=0x%Zx;", t);
262 : }
263 10 : mpres_get_z (t, A, modulus);
264 10 : gmp_fprintf (chkfile, " A=0x%Zx;", t);
265 : }
266 20 : fprintf (chkfile, "\n");
267 20 : mpz_clear (t);
268 20 : fflush (chkfile);
269 20 : fclose (chkfile);
270 : }
271 :
272 : #if 0 /* currently unused (only used in listz_handle.c, currently inactive) */
273 : int
274 : aux_fseek64(FILE *f, const int64_t offset, const int whence)
275 : {
276 : #ifdef HAVE__FSEEKI64
277 : return _fseeki64(f, offset, whence);
278 : #endif
279 : #if LONG_MAX == INT64_MAX
280 : return fseek (f, (long) offset, whence);
281 : #endif
282 : ASSERT_ALWAYS (offset <= LONG_MAX);
283 : return fseek (f, (long) offset, whence);
284 : }
285 : #endif
286 :
287 : int
288 21151796 : ecm_tstbit (mpz_srcptr u, ecm_uint bit_index)
289 : {
290 21151796 : mp_srcptr u_ptr = PTR(u);
291 21151796 : ecm_int size = SIZ(u);
292 21151796 : ecm_uint abs_size = ABS(size);
293 21151796 : ecm_uint limb_index = bit_index / GMP_NUMB_BITS;
294 21151796 : mp_srcptr p = u_ptr + limb_index;
295 : mp_limb_t limb;
296 :
297 21151796 : if (limb_index >= abs_size)
298 216328 : return (size < 0);
299 :
300 20935468 : limb = *p;
301 20935468 : if (size < 0)
302 : {
303 0 : limb = -limb; /* twos complement */
304 :
305 0 : while (p != u_ptr)
306 : {
307 0 : p--;
308 0 : if (*p != 0)
309 : {
310 0 : limb--; /* make it a ones complement instead */
311 0 : break;
312 : }
313 : }
314 : }
315 :
316 20935468 : return (limb >> (bit_index % GMP_NUMB_BITS)) & 1;
317 : }
|