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 #ifdef __FBSDID
46 __FBSDID("$FreeBSD$");
47 #endif
48 #endif /* not lint */
49
50 /*
51 * factor - factor a number into primes
52 *
53 * By: Landon Curt Noll [email protected], ...!{sun,tolsoft}!hoptoad!chongo
54 *
55 * chongo <for a good prime call: 391581 * 2^216193 - 1> /\oo/\
56 *
57 * usage:
58 * factor [-h] [number] ...
59 *
60 * The form of the output is:
61 *
62 * number: factor1 factor1 factor2 factor3 factor3 factor3 ...
63 *
64 * where factor1 <= factor2 <= factor3 <= ...
65 *
66 * If no args are given, the list of numbers are read from stdin.
67 */
68
69 #include <ctype.h>
70 #include <err.h>
71 #include <errno.h>
72 #include <inttypes.h>
73 #include <limits.h>
74 #include <stdio.h>
75 #include <stdlib.h>
76 #include <unistd.h>
77
78 #include "primes.h"
79
80 #ifdef HAVE_OPENSSL
81
82 #include <openssl/bn.h>
83
84 #define PRIME_CHECKS 5
85
86 static void pollard_pminus1(BIGNUM *); /* print factors for big numbers */
87
88 #else
89
90 typedef ubig BIGNUM;
91 typedef u_long BN_ULONG;
92
93 #define BN_CTX int
94 #define BN_CTX_new() NULL
95 #define BN_new() ((BIGNUM *)calloc(sizeof(BIGNUM), 1))
96 #define BN_is_zero(v) (*(v) == 0)
97 #define BN_is_one(v) (*(v) == 1)
98 #define BN_mod_word(a, b) (*(a) % (b))
99
100 static int BN_dec2bn(BIGNUM **a, const char *str);
101 static int BN_hex2bn(BIGNUM **a, const char *str);
102 static BN_ULONG BN_div_word(BIGNUM *, BN_ULONG);
103 static void BN_print_fp(FILE *, const BIGNUM *);
104
105 #endif
106
107 static void BN_print_dec_fp(FILE *, const BIGNUM *);
108
109 static void pr_fact(BIGNUM *); /* print factors of a value */
110 static void pr_print(BIGNUM *); /* print a prime */
111 static void usage(void);
112
113 static BN_CTX *ctx; /* just use a global context */
114 static int hflag;
115
116 int
main(int argc,char * argv[])117 main(int argc, char *argv[])
118 {
119 BIGNUM *val;
120 int ch;
121 char *p, buf[LINE_MAX]; /* > max number of digits. */
122
123 ctx = BN_CTX_new();
124 val = BN_new();
125 if (val == NULL)
126 errx(1, "can't initialise bignum");
127
128 while ((ch = getopt(argc, argv, "h")) != -1)
129 switch (ch) {
130 case 'h':
131 hflag++;
132 break;
133 case '?':
134 default:
135 usage();
136 }
137 argc -= optind;
138 argv += optind;
139
140 /* No args supplied, read numbers from stdin. */
141 if (argc == 0)
142 for (;;) {
143 if (fgets(buf, sizeof(buf), stdin) == NULL) {
144 if (ferror(stdin))
145 err(1, "stdin");
146 exit (0);
147 }
148 for (p = buf; isblank(*p); ++p);
149 if (*p == '\n' || *p == '\0')
150 continue;
151 if (*p == '-')
152 errx(1, "negative numbers aren't permitted.");
153 if (BN_dec2bn(&val, buf) == 0 &&
154 BN_hex2bn(&val, buf) == 0)
155 errx(1, "%s: illegal numeric format.", buf);
156 pr_fact(val);
157 }
158 /* Factor the arguments. */
159 else
160 for (; *argv != NULL; ++argv) {
161 if (argv[0][0] == '-')
162 errx(1, "negative numbers aren't permitted.");
163 if (BN_dec2bn(&val, argv[0]) == 0 &&
164 BN_hex2bn(&val, argv[0]) == 0)
165 errx(1, "%s: illegal numeric format.", argv[0]);
166 pr_fact(val);
167 }
168 exit(0);
169 }
170
171 /*
172 * pr_fact - print the factors of a number
173 *
174 * Print the factors of the number, from the lowest to the highest.
175 * A factor will be printed multiple times if it divides the value
176 * multiple times.
177 *
178 * Factors are printed with leading tabs.
179 */
180 static void
pr_fact(BIGNUM * val)181 pr_fact(BIGNUM *val)
182 {
183 const ubig *fact; /* The factor found. */
184
185 /* Firewall - catch 0 and 1. */
186 if (BN_is_zero(val)) /* Historical practice; 0 just exits. */
187 exit(0);
188 if (BN_is_one(val)) {
189 printf("1: 1\n");
190 return;
191 }
192
193 /* Factor value. */
194
195 if (hflag) {
196 fputs("0x", stdout);
197 BN_print_fp(stdout, val);
198 } else
199 BN_print_dec_fp(stdout, val);
200 putchar(':');
201 for (fact = &prime[0]; !BN_is_one(val); ++fact) {
202 /* Look for the smallest factor. */
203 do {
204 if (BN_mod_word(val, (BN_ULONG)*fact) == 0)
205 break;
206 } while (++fact <= pr_limit);
207
208 /* Watch for primes larger than the table. */
209 if (fact > pr_limit) {
210 #ifdef HAVE_OPENSSL
211 BIGNUM *bnfact;
212
213 bnfact = BN_new();
214 BN_set_word(bnfact, *(fact - 1));
215 if (!BN_sqr(bnfact, bnfact, ctx))
216 errx(1, "error in BN_sqr()");
217 if (BN_cmp(bnfact, val) > 0 ||
218 BN_is_prime_ex(val, PRIME_CHECKS, NULL, NULL) == 1)
219 pr_print(val);
220 else
221 pollard_pminus1(val);
222 #else
223 pr_print(val);
224 #endif
225 break;
226 }
227
228 /* Divide factor out until none are left. */
229 do {
230 printf(hflag ? " 0x%" PRIx64 "" : " %" PRIu64 "", *fact);
231 BN_div_word(val, (BN_ULONG)*fact);
232 } while (BN_mod_word(val, (BN_ULONG)*fact) == 0);
233
234 /* Let the user know we're doing something. */
235 fflush(stdout);
236 }
237 putchar('\n');
238 }
239
240 static void
pr_print(BIGNUM * val)241 pr_print(BIGNUM *val)
242 {
243 if (hflag) {
244 fputs(" 0x", stdout);
245 BN_print_fp(stdout, val);
246 } else {
247 putchar(' ');
248 BN_print_dec_fp(stdout, val);
249 }
250 }
251
252 static void
usage(void)253 usage(void)
254 {
255 fprintf(stderr, "usage: factor [-h] [value ...]\n");
256 exit(1);
257 }
258
259 #ifdef HAVE_OPENSSL
260
261 /* pollard p-1, algorithm from Jim Gillogly, May 2000 */
262 static void
pollard_pminus1(BIGNUM * val)263 pollard_pminus1(BIGNUM *val)
264 {
265 BIGNUM *base, *rbase, *num, *i, *x;
266
267 base = BN_new();
268 rbase = BN_new();
269 num = BN_new();
270 i = BN_new();
271 x = BN_new();
272
273 BN_set_word(rbase, 1);
274 newbase:
275 if (!BN_add_word(rbase, 1))
276 errx(1, "error in BN_add_word()");
277 BN_set_word(i, 2);
278 BN_copy(base, rbase);
279
280 for (;;) {
281 BN_mod_exp(base, base, i, val, ctx);
282 if (BN_is_one(base))
283 goto newbase;
284
285 BN_copy(x, base);
286 BN_sub_word(x, 1);
287 if (!BN_gcd(x, x, val, ctx))
288 errx(1, "error in BN_gcd()");
289
290 if (!BN_is_one(x)) {
291 if (BN_is_prime_ex(x, PRIME_CHECKS, NULL, NULL) == 1)
292 pr_print(x);
293 else
294 pollard_pminus1(x);
295 fflush(stdout);
296
297 BN_div(num, NULL, val, x, ctx);
298 if (BN_is_one(num))
299 return;
300 if (BN_is_prime_ex(num, PRIME_CHECKS, NULL,
301 NULL) == 1) {
302 pr_print(num);
303 fflush(stdout);
304 return;
305 }
306 BN_copy(val, num);
307 }
308 if (!BN_add_word(i, 1))
309 errx(1, "error in BN_add_word()");
310 }
311 }
312
313 /*
314 * Sigh.. No _decimal_ output to file functions in BN.
315 */
316 static void
BN_print_dec_fp(FILE * fp,const BIGNUM * num)317 BN_print_dec_fp(FILE *fp, const BIGNUM *num)
318 {
319 char *buf;
320
321 buf = BN_bn2dec(num);
322 if (buf == NULL)
323 return; /* XXX do anything here? */
324 fprintf(fp, "%s", buf);
325 free(buf);
326 }
327
328 #else
329
330 static void
BN_print_fp(FILE * fp,const BIGNUM * num)331 BN_print_fp(FILE *fp, const BIGNUM *num)
332 {
333 fprintf(fp, "%lx", (unsigned long)*num);
334 }
335
336 static void
BN_print_dec_fp(FILE * fp,const BIGNUM * num)337 BN_print_dec_fp(FILE *fp, const BIGNUM *num)
338 {
339 fprintf(fp, "%lu", (unsigned long)*num);
340 }
341
342 static int
BN_dec2bn(BIGNUM ** a,const char * str)343 BN_dec2bn(BIGNUM **a, const char *str)
344 {
345 char *p;
346
347 errno = 0;
348 **a = strtoul(str, &p, 10);
349 return (errno == 0 && (*p == '\n' || *p == '\0'));
350 }
351
352 static int
BN_hex2bn(BIGNUM ** a,const char * str)353 BN_hex2bn(BIGNUM **a, const char *str)
354 {
355 char *p;
356
357 errno = 0;
358 **a = strtoul(str, &p, 16);
359 return (errno == 0 && (*p == '\n' || *p == '\0'));
360 }
361
362 static BN_ULONG
BN_div_word(BIGNUM * a,BN_ULONG b)363 BN_div_word(BIGNUM *a, BN_ULONG b)
364 {
365 BN_ULONG mod;
366
367 mod = *a % b;
368 *a /= b;
369 return mod;
370 }
371
372 #endif
373