1 /*
2 * Copyright (c) 1989, 1993
3 * The Regents of the University of California. All rights reserved.
4 *
5 * This code is derived from software contributed to Berkeley by
6 * Landon Curt Noll.
7 *
8 * Redistribution and use in source and binary forms, with or without
9 * modification, are permitted provided that the following conditions
10 * are met:
11 * 1. Redistributions of source code must retain the above copyright
12 * notice, this list of conditions and the following disclaimer.
13 * 2. Redistributions in binary form must reproduce the above copyright
14 * notice, this list of conditions and the following disclaimer in the
15 * documentation and/or other materials provided with the distribution.
16 * 3. Neither the name of the University nor the names of its contributors
17 * may be used to endorse or promote products derived from this software
18 * without specific prior written permission.
19 *
20 * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
21 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
22 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
23 * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
24 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
25 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
26 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
27 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
28 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
29 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
30 * SUCH DAMAGE.
31 */
32
33 #ifndef lint
34 #include <sys/cdefs.h>
35 #ifdef __COPYRIGHT
36 __COPYRIGHT("@(#) Copyright (c) 1989, 1993\
37 The Regents of the University of California. All rights reserved.");
38 #endif
39 #ifdef __SCCSID
40 __SCCSID("@(#)factor.c 8.4 (Berkeley) 5/4/95");
41 #endif
42 #ifdef __RCSID
43 __RCSID("$NetBSD: factor.c,v 1.19 2009/08/12 05:54:31 dholland Exp $");
44 #endif
45 #endif /* not lint */
46
47 /*
48 * factor - factor a number into primes
49 *
50 * By: Landon Curt Noll [email protected], ...!{sun,tolsoft}!hoptoad!chongo
51 *
52 * chongo <for a good prime call: 391581 * 2^216193 - 1> /\oo/\
53 *
54 * usage:
55 * factor [-h] [number] ...
56 *
57 * The form of the output is:
58 *
59 * number: factor1 factor1 factor2 factor3 factor3 factor3 ...
60 *
61 * where factor1 <= factor2 <= factor3 <= ...
62 *
63 * If no args are given, the list of numbers are read from stdin.
64 */
65
66 #include <ctype.h>
67 #include <err.h>
68 #include <errno.h>
69 #include <inttypes.h>
70 #include <limits.h>
71 #include <stdbool.h>
72 #include <stdio.h>
73 #include <stdlib.h>
74 #include <unistd.h>
75
76 #include "primes.h"
77
78 #ifdef HAVE_OPENSSL
79
80 #include <openssl/bn.h>
81
82 #if OPENSSL_VERSION_NUMBER < 0x30000000L
83 static inline int
BN_check_prime(BIGNUM * p,BN_CTX * ctx,BN_GENCB * cb)84 BN_check_prime(BIGNUM *p, BN_CTX *ctx, BN_GENCB *cb)
85 {
86 const int nchecks = 5;
87
88 return BN_is_prime_ex(p, nchecks, ctx, cb);
89 }
90 #endif
91
92 static void pollard_pminus1(BIGNUM *); /* print factors for big numbers */
93
94 #else
95
96 typedef ubig BIGNUM;
97 typedef u_long BN_ULONG;
98
99 #define BN_CTX int
100 #define BN_CTX_new() NULL
101 #define BN_new() ((BIGNUM *)calloc(sizeof(BIGNUM), 1))
102 #define BN_is_zero(v) (*(v) == 0)
103 #define BN_is_one(v) (*(v) == 1)
104 #define BN_mod_word(a, b) (*(a) % (b))
105
106 static int BN_dec2bn(BIGNUM **, const char *);
107 static int BN_hex2bn(BIGNUM **, const char *);
108 static BN_ULONG BN_div_word(BIGNUM *, BN_ULONG);
109 static void BN_print_fp(FILE *, const BIGNUM *);
110
111 #endif
112
113 static void BN_print_dec_fp(FILE *, const BIGNUM *);
114 static void convert_str2bn(BIGNUM **, char *);
115 static bool is_hex_str(char *);
116 static void pr_fact(BIGNUM *); /* print factors of a value */
117 static void pr_print(BIGNUM *); /* print a prime */
118 static void usage(void);
119
120 static BN_CTX *ctx; /* just use a global context */
121 static int hflag;
122
123 int
main(int argc,char * argv[])124 main(int argc, char *argv[])
125 {
126 BIGNUM *val;
127 int ch;
128 char *p, buf[LINE_MAX]; /* > max number of digits. */
129
130 ctx = BN_CTX_new();
131 val = BN_new();
132 if (val == NULL)
133 errx(1, "can't initialise bignum");
134
135 while ((ch = getopt(argc, argv, "h")) != -1)
136 switch (ch) {
137 case 'h':
138 hflag++;
139 break;
140 case '?':
141 default:
142 usage();
143 }
144 argc -= optind;
145 argv += optind;
146
147 /* No args supplied, read numbers from stdin. */
148 if (argc == 0)
149 for (;;) {
150 if (fgets(buf, sizeof(buf), stdin) == NULL) {
151 if (ferror(stdin))
152 err(1, "stdin");
153 exit (0);
154 }
155 for (p = buf; isblank(*p); ++p);
156 if (*p == '\n' || *p == '\0')
157 continue;
158 convert_str2bn(&val, p);
159 pr_fact(val);
160 }
161 /* Factor the arguments. */
162 else
163 for (p = *argv; p != NULL; p = *++argv) {
164 convert_str2bn(&val, p);
165 pr_fact(val);
166 }
167 exit(0);
168 }
169
170 /*
171 * pr_fact - print the factors of a number
172 *
173 * Print the factors of the number, from the lowest to the highest.
174 * A factor will be printed multiple times if it divides the value
175 * multiple times.
176 *
177 * Factors are printed with leading tabs.
178 */
179 static void
pr_fact(BIGNUM * val)180 pr_fact(BIGNUM *val)
181 {
182 const ubig *fact; /* The factor found. */
183
184 /* Firewall - catch 0 and 1. */
185 if (BN_is_zero(val)) /* Historical practice; 0 just exits. */
186 exit(0);
187 if (BN_is_one(val)) {
188 printf("1: 1\n");
189 return;
190 }
191
192 /* Factor value. */
193
194 if (hflag) {
195 fputs("0x", stdout);
196 BN_print_fp(stdout, val);
197 } else
198 BN_print_dec_fp(stdout, val);
199 putchar(':');
200 for (fact = &prime[0]; !BN_is_one(val); ++fact) {
201 /* Look for the smallest factor. */
202 do {
203 if (BN_mod_word(val, (BN_ULONG)*fact) == 0)
204 break;
205 } while (++fact <= pr_limit);
206
207 /* Watch for primes larger than the table. */
208 if (fact > pr_limit) {
209 #ifdef HAVE_OPENSSL
210 BIGNUM *bnfact;
211
212 bnfact = BN_new();
213 BN_set_word(bnfact, *(fact - 1));
214 if (!BN_sqr(bnfact, bnfact, ctx))
215 errx(1, "error in BN_sqr()");
216 if (BN_cmp(bnfact, val) > 0 ||
217 BN_check_prime(val, NULL, NULL) == 1)
218 pr_print(val);
219 else
220 pollard_pminus1(val);
221 #else
222 pr_print(val);
223 #endif
224 break;
225 }
226
227 /* Divide factor out until none are left. */
228 do {
229 printf(hflag ? " 0x%" PRIx64 "" : " %" PRIu64 "", *fact);
230 BN_div_word(val, (BN_ULONG)*fact);
231 } while (BN_mod_word(val, (BN_ULONG)*fact) == 0);
232
233 /* Let the user know we're doing something. */
234 fflush(stdout);
235 }
236 putchar('\n');
237 }
238
239 static void
pr_print(BIGNUM * val)240 pr_print(BIGNUM *val)
241 {
242 if (hflag) {
243 fputs(" 0x", stdout);
244 BN_print_fp(stdout, val);
245 } else {
246 putchar(' ');
247 BN_print_dec_fp(stdout, val);
248 }
249 }
250
251 static void
usage(void)252 usage(void)
253 {
254 fprintf(stderr, "usage: factor [-h] [value ...]\n");
255 exit(1);
256 }
257
258 #ifdef HAVE_OPENSSL
259
260 /* pollard p-1, algorithm from Jim Gillogly, May 2000 */
261 static void
pollard_pminus1(BIGNUM * val)262 pollard_pminus1(BIGNUM *val)
263 {
264 BIGNUM *base, *rbase, *num, *i, *x;
265
266 base = BN_new();
267 rbase = BN_new();
268 num = BN_new();
269 i = BN_new();
270 x = BN_new();
271
272 BN_set_word(rbase, 1);
273 newbase:
274 if (!BN_add_word(rbase, 1))
275 errx(1, "error in BN_add_word()");
276 BN_set_word(i, 2);
277 BN_copy(base, rbase);
278
279 for (;;) {
280 BN_mod_exp(base, base, i, val, ctx);
281 if (BN_is_one(base))
282 goto newbase;
283
284 BN_copy(x, base);
285 BN_sub_word(x, 1);
286 if (!BN_gcd(x, x, val, ctx))
287 errx(1, "error in BN_gcd()");
288
289 if (!BN_is_one(x)) {
290 if (BN_check_prime(x, NULL, NULL) == 1)
291 pr_print(x);
292 else
293 pollard_pminus1(x);
294 fflush(stdout);
295
296 BN_div(num, NULL, val, x, ctx);
297 if (BN_is_one(num))
298 return;
299 if (BN_check_prime(num, NULL, NULL) == 1) {
300 pr_print(num);
301 fflush(stdout);
302 return;
303 }
304 BN_copy(val, num);
305 }
306 if (!BN_add_word(i, 1))
307 errx(1, "error in BN_add_word()");
308 }
309 }
310
311 /*
312 * Sigh.. No _decimal_ output to file functions in BN.
313 */
314 static void
BN_print_dec_fp(FILE * fp,const BIGNUM * num)315 BN_print_dec_fp(FILE *fp, const BIGNUM *num)
316 {
317 char *buf;
318
319 buf = BN_bn2dec(num);
320 if (buf == NULL)
321 return; /* XXX do anything here? */
322 fprintf(fp, "%s", buf);
323 free(buf);
324 }
325
326 #else
327
328 static void
BN_print_fp(FILE * fp,const BIGNUM * num)329 BN_print_fp(FILE *fp, const BIGNUM *num)
330 {
331 fprintf(fp, "%lx", (unsigned long)*num);
332 }
333
334 static void
BN_print_dec_fp(FILE * fp,const BIGNUM * num)335 BN_print_dec_fp(FILE *fp, const BIGNUM *num)
336 {
337 fprintf(fp, "%lu", (unsigned long)*num);
338 }
339
340 static int
BN_dec2bn(BIGNUM ** a,const char * str)341 BN_dec2bn(BIGNUM **a, const char *str)
342 {
343 char *p;
344
345 errno = 0;
346 **a = strtoul(str, &p, 10);
347 return (errno == 0 ? 1 : 0); /* OpenSSL returns 0 on error! */
348 }
349
350 static int
BN_hex2bn(BIGNUM ** a,const char * str)351 BN_hex2bn(BIGNUM **a, const char *str)
352 {
353 char *p;
354
355 errno = 0;
356 **a = strtoul(str, &p, 16);
357 return (errno == 0 ? 1 : 0); /* OpenSSL returns 0 on error! */
358 }
359
360 static BN_ULONG
BN_div_word(BIGNUM * a,BN_ULONG b)361 BN_div_word(BIGNUM *a, BN_ULONG b)
362 {
363 BN_ULONG mod;
364
365 mod = *a % b;
366 *a /= b;
367 return mod;
368 }
369
370 #endif
371
372 /*
373 * Scan the string from left-to-right to see if the longest substring
374 * is a valid hexadecimal number.
375 */
376 static bool
is_hex_str(char * str)377 is_hex_str(char *str)
378 {
379 char c, *p;
380 bool saw_hex = false;
381
382 for (p = str; *p; p++) {
383 if (isdigit(*p))
384 continue;
385 c = tolower(*p);
386 if (c >= 'a' && c <= 'f') {
387 saw_hex = true;
388 continue;
389 }
390 break; /* Not a hexadecimal digit. */
391 }
392 return saw_hex;
393 }
394
395 /* Convert string pointed to by *str to a bignum. */
396 static void
convert_str2bn(BIGNUM ** val,char * p)397 convert_str2bn(BIGNUM **val, char *p)
398 {
399 int n = 0;
400
401 if (*p == '+') p++;
402 if (*p == '-')
403 errx(1, "negative numbers aren't permitted.");
404 if (*p == '0') {
405 p++;
406 if (*p == 'x' || *p == 'X')
407 n = BN_hex2bn(val, ++p);
408 } else {
409 n = is_hex_str(p) ? BN_hex2bn(val, p) : BN_dec2bn(val, p);
410 }
411 if (n == 0)
412 errx(1, "%s: illegal numeric format.", p);
413 }
414