xref: /freebsd/usr.bin/factor/factor.c (revision bdcbfde31e8e9b343f113a1956384bdf30d1ed62)
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   chongo@toad.com,   ...!{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
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
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
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
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
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
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
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
329 BN_print_fp(FILE *fp, const BIGNUM *num)
330 {
331 	fprintf(fp, "%lx", (unsigned long)*num);
332 }
333 
334 static void
335 BN_print_dec_fp(FILE *fp, const BIGNUM *num)
336 {
337 	fprintf(fp, "%lu", (unsigned long)*num);
338 }
339 
340 static int
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
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
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
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
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