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