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
69static inline int
70BN_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
78static void	pollard_pminus1(BIGNUM *); /* print factors for big numbers */
79
80#else
81
82typedef ubig	BIGNUM;
83typedef 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
92static int	BN_dec2bn(BIGNUM **, const char *);
93static int	BN_hex2bn(BIGNUM **, const char *);
94static BN_ULONG BN_div_word(BIGNUM *, BN_ULONG);
95static void	BN_print_fp(FILE *, const BIGNUM *);
96
97#endif
98
99static void	BN_print_dec_fp(FILE *, const BIGNUM *);
100static void	convert_str2bn(BIGNUM **, char *);
101static bool	is_hex_str(char *);
102static void	pr_fact(BIGNUM *);	/* print factors of a value */
103static void	pr_print(BIGNUM *);	/* print a prime */
104static void	usage(void);
105
106static BN_CTX	*ctx;			/* just use a global context */
107static int	hflag;
108
109int
110main(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 */
165static void
166pr_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
225static void
226pr_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
237static void
238usage(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 */
247static void
248pollard_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);
259newbase:
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 */
300static void
301BN_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
314static void
315BN_print_fp(FILE *fp, const BIGNUM *num)
316{
317	fprintf(fp, "%lx", (unsigned long)*num);
318}
319
320static void
321BN_print_dec_fp(FILE *fp, const BIGNUM *num)
322{
323	fprintf(fp, "%lu", (unsigned long)*num);
324}
325
326static int
327BN_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
336static int
337BN_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
346static BN_ULONG
347BN_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 */
362static bool
363is_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.  */
382static void
383convert_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