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
| 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: head/games/factor/factor.c 199815 2009-11-26 00:38:13Z fanf $");
| 46__FBSDID("$FreeBSD: head/games/factor/factor.c 216598 2010-12-20 19:07:56Z uqs $");
|
47#endif
48#endif /* not lint */
49
50/*
51 * factor - factor a number into primes
52 *
53 * By: Landon Curt Noll chongo@toad.com, ...!{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 <limits.h>
73#include <stdio.h>
74#include <stdlib.h>
75#include <unistd.h>
76
77#include "primes.h"
78
79#ifdef HAVE_OPENSSL
80
81#include <openssl/bn.h>
82
83#define PRIME_CHECKS 5
84
85static void pollard_pminus1(BIGNUM *); /* print factors for big numbers */
86
87#else
88
89typedef ubig BIGNUM;
90typedef u_long BN_ULONG;
91
92#define BN_CTX int
93#define BN_CTX_new() NULL
94#define BN_new() ((BIGNUM *)calloc(sizeof(BIGNUM), 1))
95#define BN_is_zero(v) (*(v) == 0)
96#define BN_is_one(v) (*(v) == 1)
97#define BN_mod_word(a, b) (*(a) % (b))
98
99static int BN_dec2bn(BIGNUM **a, const char *str);
100static int BN_hex2bn(BIGNUM **a, const char *str);
101static BN_ULONG BN_div_word(BIGNUM *, BN_ULONG);
102static void BN_print_fp(FILE *, const BIGNUM *);
103
104#endif
105
106static void BN_print_dec_fp(FILE *, const BIGNUM *);
107
108static void pr_fact(BIGNUM *); /* print factors of a value */
109static void pr_print(BIGNUM *); /* print a prime */
110static void usage(void);
111
112static BN_CTX *ctx; /* just use a global context */
113static int hflag;
114
115int
116main(int argc, char *argv[])
117{
118 BIGNUM *val;
119 int ch;
120 char *p, buf[LINE_MAX]; /* > max number of digits. */
121
122 ctx = BN_CTX_new();
123 val = BN_new();
124 if (val == NULL)
125 errx(1, "can't initialise bignum");
126
127 while ((ch = getopt(argc, argv, "h")) != -1)
128 switch (ch) {
129 case 'h':
130 hflag++;
131 break;
132 case '?':
133 default:
134 usage();
135 }
136 argc -= optind;
137 argv += optind;
138
139 /* No args supplied, read numbers from stdin. */
140 if (argc == 0)
141 for (;;) {
142 if (fgets(buf, sizeof(buf), stdin) == NULL) {
143 if (ferror(stdin))
144 err(1, "stdin");
145 exit (0);
146 }
147 for (p = buf; isblank(*p); ++p);
148 if (*p == '\n' || *p == '\0')
149 continue;
150 if (*p == '-')
151 errx(1, "negative numbers aren't permitted.");
152 if (BN_dec2bn(&val, buf) == 0 &&
153 BN_hex2bn(&val, buf) == 0)
154 errx(1, "%s: illegal numeric format.", buf);
155 pr_fact(val);
156 }
157 /* Factor the arguments. */
158 else
159 for (; *argv != NULL; ++argv) {
160 if (argv[0][0] == '-')
161 errx(1, "negative numbers aren't permitted.");
162 if (BN_dec2bn(&val, argv[0]) == 0 &&
163 BN_hex2bn(&val, argv[0]) == 0)
164 errx(1, "%s: illegal numeric format.", argv[0]);
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 */
179static void
180pr_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));
| 47#endif
48#endif /* not lint */
49
50/*
51 * factor - factor a number into primes
52 *
53 * By: Landon Curt Noll chongo@toad.com, ...!{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 <limits.h>
73#include <stdio.h>
74#include <stdlib.h>
75#include <unistd.h>
76
77#include "primes.h"
78
79#ifdef HAVE_OPENSSL
80
81#include <openssl/bn.h>
82
83#define PRIME_CHECKS 5
84
85static void pollard_pminus1(BIGNUM *); /* print factors for big numbers */
86
87#else
88
89typedef ubig BIGNUM;
90typedef u_long BN_ULONG;
91
92#define BN_CTX int
93#define BN_CTX_new() NULL
94#define BN_new() ((BIGNUM *)calloc(sizeof(BIGNUM), 1))
95#define BN_is_zero(v) (*(v) == 0)
96#define BN_is_one(v) (*(v) == 1)
97#define BN_mod_word(a, b) (*(a) % (b))
98
99static int BN_dec2bn(BIGNUM **a, const char *str);
100static int BN_hex2bn(BIGNUM **a, const char *str);
101static BN_ULONG BN_div_word(BIGNUM *, BN_ULONG);
102static void BN_print_fp(FILE *, const BIGNUM *);
103
104#endif
105
106static void BN_print_dec_fp(FILE *, const BIGNUM *);
107
108static void pr_fact(BIGNUM *); /* print factors of a value */
109static void pr_print(BIGNUM *); /* print a prime */
110static void usage(void);
111
112static BN_CTX *ctx; /* just use a global context */
113static int hflag;
114
115int
116main(int argc, char *argv[])
117{
118 BIGNUM *val;
119 int ch;
120 char *p, buf[LINE_MAX]; /* > max number of digits. */
121
122 ctx = BN_CTX_new();
123 val = BN_new();
124 if (val == NULL)
125 errx(1, "can't initialise bignum");
126
127 while ((ch = getopt(argc, argv, "h")) != -1)
128 switch (ch) {
129 case 'h':
130 hflag++;
131 break;
132 case '?':
133 default:
134 usage();
135 }
136 argc -= optind;
137 argv += optind;
138
139 /* No args supplied, read numbers from stdin. */
140 if (argc == 0)
141 for (;;) {
142 if (fgets(buf, sizeof(buf), stdin) == NULL) {
143 if (ferror(stdin))
144 err(1, "stdin");
145 exit (0);
146 }
147 for (p = buf; isblank(*p); ++p);
148 if (*p == '\n' || *p == '\0')
149 continue;
150 if (*p == '-')
151 errx(1, "negative numbers aren't permitted.");
152 if (BN_dec2bn(&val, buf) == 0 &&
153 BN_hex2bn(&val, buf) == 0)
154 errx(1, "%s: illegal numeric format.", buf);
155 pr_fact(val);
156 }
157 /* Factor the arguments. */
158 else
159 for (; *argv != NULL; ++argv) {
160 if (argv[0][0] == '-')
161 errx(1, "negative numbers aren't permitted.");
162 if (BN_dec2bn(&val, argv[0]) == 0 &&
163 BN_hex2bn(&val, argv[0]) == 0)
164 errx(1, "%s: illegal numeric format.", argv[0]);
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 */
179static void
180pr_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 BN_sqr(bnfact, bnfact, ctx);
| 214 if (!BN_sqr(bnfact, bnfact, ctx))
215 errx(1, "error in BN_sqr()");
|
215 if (BN_cmp(bnfact, val) > 0 ||
216 BN_is_prime(val, PRIME_CHECKS,
217 NULL, 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%lx" : " %lu", *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
239static void
240pr_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
251static void
252usage(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 */
261static void
262pollard_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);
273newbase:
| 216 if (BN_cmp(bnfact, val) > 0 ||
217 BN_is_prime(val, PRIME_CHECKS,
218 NULL, 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%lx" : " %lu", *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
240static void
241pr_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
252static void
253usage(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 */
262static void
263pollard_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);
274newbase:
|
274 BN_add_word(rbase, 1);
| 275 if (!BN_add_word(rbase, 1))
276 errx(1, "error in BN_add_word()");
|
275 BN_set_word(i, 2);
276 BN_copy(base, rbase);
277
278 for (;;) {
279 BN_mod_exp(base, base, i, val, ctx);
280 if (BN_is_one(base))
281 goto newbase;
282
283 BN_copy(x, base);
284 BN_sub_word(x, 1);
| 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);
|
285 BN_gcd(x, x, val, ctx);
| 287 if (!BN_gcd(x, x, val, ctx))
288 errx(1, "error in BN_gcd()");
|
286
287 if (!BN_is_one(x)) {
288 if (BN_is_prime(x, PRIME_CHECKS, NULL, NULL,
289 NULL) == 1)
290 pr_print(x);
291 else
292 pollard_pminus1(x);
293 fflush(stdout);
294
295 BN_div(num, NULL, val, x, ctx);
296 if (BN_is_one(num))
297 return;
298 if (BN_is_prime(num, PRIME_CHECKS, NULL, NULL,
299 NULL) == 1) {
300 pr_print(num);
301 fflush(stdout);
302 return;
303 }
304 BN_copy(val, num);
305 }
| 289
290 if (!BN_is_one(x)) {
291 if (BN_is_prime(x, PRIME_CHECKS, NULL, NULL,
292 NULL) == 1)
293 pr_print(x);
294 else
295 pollard_pminus1(x);
296 fflush(stdout);
297
298 BN_div(num, NULL, val, x, ctx);
299 if (BN_is_one(num))
300 return;
301 if (BN_is_prime(num, PRIME_CHECKS, NULL, NULL,
302 NULL) == 1) {
303 pr_print(num);
304 fflush(stdout);
305 return;
306 }
307 BN_copy(val, num);
308 }
|
306 BN_add_word(i, 1);
| 309 if (!BN_add_word(i, 1))
310 errx(1, "error in BN_add_word()");
|
307 }
308}
309
310/*
311 * Sigh.. No _decimal_ output to file functions in BN.
312 */
313static void
314BN_print_dec_fp(FILE *fp, const BIGNUM *num)
315{
316 char *buf;
317
318 buf = BN_bn2dec(num);
319 if (buf == NULL)
320 return; /* XXX do anything here? */
321 fprintf(fp, buf);
322 free(buf);
323}
324
325#else
326
327static void
328BN_print_fp(FILE *fp, const BIGNUM *num)
329{
330 fprintf(fp, "%lx", (unsigned long)*num);
331}
332
333static void
334BN_print_dec_fp(FILE *fp, const BIGNUM *num)
335{
336 fprintf(fp, "%lu", (unsigned long)*num);
337}
338
339static int
340BN_dec2bn(BIGNUM **a, const char *str)
341{
342 char *p;
343
344 errno = 0;
345 **a = strtoul(str, &p, 10);
346 return (errno == 0 && (*p == '\n' || *p == '\0'));
347}
348
349static int
350BN_hex2bn(BIGNUM **a, const char *str)
351{
352 char *p;
353
354 errno = 0;
355 **a = strtoul(str, &p, 16);
356 return (errno == 0 && (*p == '\n' || *p == '\0'));
357}
358
359static BN_ULONG
360BN_div_word(BIGNUM *a, BN_ULONG b)
361{
362 BN_ULONG mod;
363
364 mod = *a % b;
365 *a /= b;
366 return mod;
367}
368
369#endif
| 311 }
312}
313
314/*
315 * Sigh.. No _decimal_ output to file functions in BN.
316 */
317static void
318BN_print_dec_fp(FILE *fp, const BIGNUM *num)
319{
320 char *buf;
321
322 buf = BN_bn2dec(num);
323 if (buf == NULL)
324 return; /* XXX do anything here? */
325 fprintf(fp, buf);
326 free(buf);
327}
328
329#else
330
331static void
332BN_print_fp(FILE *fp, const BIGNUM *num)
333{
334 fprintf(fp, "%lx", (unsigned long)*num);
335}
336
337static void
338BN_print_dec_fp(FILE *fp, const BIGNUM *num)
339{
340 fprintf(fp, "%lu", (unsigned long)*num);
341}
342
343static int
344BN_dec2bn(BIGNUM **a, const char *str)
345{
346 char *p;
347
348 errno = 0;
349 **a = strtoul(str, &p, 10);
350 return (errno == 0 && (*p == '\n' || *p == '\0'));
351}
352
353static int
354BN_hex2bn(BIGNUM **a, const char *str)
355{
356 char *p;
357
358 errno = 0;
359 **a = strtoul(str, &p, 16);
360 return (errno == 0 && (*p == '\n' || *p == '\0'));
361}
362
363static BN_ULONG
364BN_div_word(BIGNUM *a, BN_ULONG b)
365{
366 BN_ULONG mod;
367
368 mod = *a % b;
369 *a /= b;
370 return mod;
371}
372
373#endif
|