usr.bin/factor/factor.c

   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 }