Signed-off-by: Mattias Andrée <maand...@kth.se>
---
LICENSE | 1 +
Makefile | 4 +
README | 1 +
factor.1 | 62 ++++++
factor.c | 667 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
5 files changed, 735 insertions(+)
create mode 100644 factor.1
create mode 100644 factor.c
diff --git a/LICENSE b/LICENSE
index 1ff215b..bd75908 100644
--- a/LICENSE
+++ b/LICENSE
@@ -61,3 +61,4 @@ Authors/contributors include:
© 2015 Wolfgang Corcoran-Mathe <first.lord.of.t...@gmail.com>
© 2016 Mattias Andrée <maand...@kth.se>
© 2016 Eivind Uggedal <eiv...@uggedal.com>
+© 2016 Joel Mickelin <joel.micke...@mykolab.ch>
diff --git a/Makefile b/Makefile
index 5fd0b50..50f30cd 100644
--- a/Makefile
+++ b/Makefile
@@ -104,6 +104,7 @@ BIN =\
env\
expand\
expr\
+ factor\
false\
find\
flock\
@@ -189,6 +190,9 @@ $(BIN): $(LIB) $(@:=.o)
$(OBJ): $(HDR) config.mk
+sbase: LDFLAGS += -lgmp -lpthread
+factor: LDFLAGS += -lgmp -lpthread
+
.o:
$(CC) $(LDFLAGS) -o $@ $< $(LIB)
diff --git a/README b/README
index aa18ede..c5834c9 100644
--- a/README
+++ b/README
@@ -34,6 +34,7 @@ The following tools are implemented:
=*|o env .
#*|o expand .
#*|o expr .
+ * factor .
=*|o false .
= find .
=* x flock .
diff --git a/factor.1 b/factor.1
new file mode 100644
index 0000000..7386d39
--- /dev/null
+++ b/factor.1
@@ -0,0 +1,62 @@
+.Dd 2016-02-25
+.Dt FACTOR 1
+.Os sbase
+.Sh NAME
+.Nm factor
+.Nd prime factor numbers
+.Sh SYNOPSIS
+.Nm
+.Op Fl c Ar N
+.Op Fl p Ar N
+.Op Fl q
+.Op Ar number ...
+.Sh DESCRIPTION
+.Nm
+prints the prime factors of
+.Ar number .
+If no
+.Ar number
+has been specified, the numbers to factor are
+line by line from stdin.
+.Sh OPTIONS
+.Bl -tag -width Ds
+.It Fl c Ar N
+Select the certainty for prime factors. The probably
+that a composite number will be listed as a prime
+factor will be less than 4^(-\fIN\fP).
+.It Fl p Ar N
+Use at most
+.Ar N
+threads to factor each
+.Ar number .
+This is seldom useful.
+Note, the input numbers are not factored concurrently,
+only there factors.
+.It Fl q
+Mark factors with \(aq?\(aq if it is not certain that
+they are prime.
+.El
+.Sh EXIT STATUS
+.Bl -tag -width Ds
+.It 0
+All numbers were successfully factorized.
+.It > 1
+An error occurred.
+.El
+.Sh NOTES
+The factors are not output in ascending order,
+but in order of discovery.
+.Pp
+The prime test is two-stage. The first stage
+is unaffected by the
+.Fl c
+flag. In this stage it is possible that the
+number is designated as a prime with complete
+certainty. In second stage, which is affected
+by the
+.Fl c
+flag, the number cannot be designated as a
+prime with complete certainty.
+.Sh SEE ALSO
+.Xr bc 1
+
diff --git a/factor.c b/factor.c
new file mode 100644
index 0000000..350dfc1
--- /dev/null
+++ b/factor.c
@@ -0,0 +1,667 @@
+/* See LICENSE file for copyright and license details. */
+#include <stdlib.h>
+#include <stdint.h>
+#include <stdio.h>
+#include <ctype.h>
+#include <errno.h>
+#include <semaphore.h>
+#include <pthread.h>
+
+#include "util.h"
+
+/*
+ * Optimisations that have been excluded.
+ *
+ * For composite N, find R, B ∈ ℕ : B↑R = N. If found, continue by
+ * multiplying root_order by R and factor B instead of N. This is
+ * done by calculating N↑−P for all primes P : P ≤ log₃ N, if the
+ * result is an integer, replace N with N↑−P, multiply R with
+ * P and try P again. The final N is B. This is useful in some
+ * cases, but not overall.
+ */
+
+static int certainty = 5;
+static int parallel = 1;
+static const char *questioned = "";
+
+#if !defined(USE_GMP) && !defined(USE_TOMMATH)
+# define USE_GMP /* GMP is significantly fast than libtommath */
+#endif
+
+#if defined(USE_GMP)
+# include <gmp.h>
+# define lowest_bit(x) mpz_scan1(x, 0)
+# define shift_right(r, x, steps) mpz_tdiv_q_2exp(r, x, steps)
+# define prime_test(x) mpz_probab_prime_p(x, certainty)
+# define to_string(x) mpz_get_str(strbuf, 10, x)
+# define div_mod(q, r, n, d) mpz_tdiv_qr(q, r, n, d)
+# define bit_count(x) mpz_sizeinbase(x, 2)
+# define gcd(r, a, b) mpz_gcd(r, a, b)
+# define zabs(r, x) mpz_abs(r, x)
+# define zmul(r, a, b) mpz_mul(r, a, b)
+# define zmul_i(r, a, b) mpz_mul_ui(r, a, b)
+# define zadd(r, a, b) mpz_add(r, a, b)
+# define zadd_i(r, a, b) mpz_add_ui(r, a, b)
+# define zsub(r, a, b) mpz_sub(r, a, b)
+# define zmod(r, a, b) mpz_mod(r, a, b)
+# define zsqrt(r, x) mpz_sqrt(r, x)
+# define zsqrt_test(r, x) mpz_root(r, x, 2)
+# define zinit(x) mpz_init(x)
+# define zfree(x) mpz_clear(x)
+# define zset(r, x) mpz_set(r, x)
+# define zset_i(r, x) mpz_set_ui(r, x)
+# define zcmp(a, b) mpz_cmp(a, b)
+# define zcmp_i(a, b) mpz_cmp_ui(a, b)
+# define zparse(r, s) mpz_init_set_str(r, s, 10);
+typedef mpz_t bigint_t;
+#elif defined(USE_TOMMATH)
+# include <tommath.h>
+# define lowest_bit(x) mp_cnt_lsb(x)
+# define shift_right(r, x, steps) mp_div_2d(x, steps, r, 0)
+static int prime_test(mp_int *x) { int ret; mp_prime_is_prime(x,
certainty, &ret); return ret;}
+# define to_string(x) (mp_todecimal(x, strbuf), strbuf)
+# define div_mod(q, r, n, d) mp_div(n, d, q, r)
+# define bit_count(x) mp_count_bits(x)
+# define gcd(r, a, b) mp_gcd(a, b, r)
+# define zabs(r, x) mp_abs(x, r)
+# define zmul(r, a, b) mp_mul(a, b, r)
+# define zmul_i(r, a, b) (zset_i(ctx->tmp, b), zmul(r, a,
ctx->tmp))
+# define zadd(r, a, b) mp_add(a, b, r)
+# define zadd_i(r, a, b) (zset_i(ctx->tmp, b), zadd(r, a,
ctx->tmp))
+# define zsub(r, a, b) mp_sub(a, b, r)
+# define zmod(r, a, b) mp_mod(a, b, r)
+# define zsqrt(r, x) mp_sqrt(x, r)
+static int zsqrt_test(mp_int *r, mp_int *x) { int ret; mp_is_square(x, &ret);
zsqrt(r, x); return ret; }
+# define zinit(x) mp_init(x)
+# define zfree(x) mp_clear(x)
+/*# define zset(r, x) mp_copy(x, r) /\* TODO Is this really
correct? */
+static void zset(mp_int *r, mp_int *x) { mp_zero(r); zadd(r, r, x); }
+# define zset_i(r, x) mp_set_int(r, x)
+# define zcmp(a, b) mp_cmp(a, b)
+# define zcmp_i(a, b) (zset_i(ctx->tmp, b), zcmp(a, ctx->tmp))
+# define zparse(r, s) (zinit(r), mp_read_radix(r, s, 10))
+typedef mp_int bigint_t[1];
+#endif
+
+#define elementsof(x) (sizeof(x) / sizeof(*x))
+#define is_factorised(x) (!zcmp_i(x, 1))
+
+enum { NO = 0, MAYBE, YES };
+
+enum { POLLARDS_RHO_INITIAL_SEED = 0 };
+enum { POLLARDS_RHO_SEED_INCREASEMENT = 500 };
+
+struct context {
+ bigint_t div_n, div_q, div_r, div_d;
+ bigint_t *div_stack;
+ size_t div_stack_size;
+ bigint_t factor, d, x, y, conj_a, conj_b;
+#ifdef USE_TOMMATH
+ bigint_t tmp;
+#endif
+};
+
+struct thread_data {
+ bigint_t integer;
+ size_t root_order;
+};
+
+#define LIST_CONSTANTS X(3) X(5) X(7) X(11) X(13) X(17)
+static bigint_t constants[18];
+
+#define _5(x) x, x, x, x, x
+#define _25(x) _5(_5(x))
+#define _50(x) _25(x), _25(x)
+static const long pollards_rho_seeds[] = {
+ [0] = _50(2),
+ [4] = _50(11),
+ [8] = _50(101),
+ [16] = _50(503),
+ [54] = _25(4993),
+ [64] = _5(6029),
+ [71] = _5(6500),
+ [72] = _5(7001),
+ [74] = _5(7559),
+ [78] = _5(8017),
+ [82] = _5(7559),
+ [83] = _5(7000),
+ [84] = _5(6500),
+ [85] = _5(6029),
+ [86] = _5(10711),
+ [89] = _5(17041),
+ [92] = _5(34511)
+};
+
+static char *strbuf;
+static void (*output_primes)(bigint_t, int, size_t);
+static void (*subfactor)(struct context *, bigint_t, size_t, bigint_t,
bigint_t);
+static sem_t semaphore;
+static pthread_mutex_t print_mutex;
+static pthread_mutex_t counter_mutex;
+static pthread_cond_t counter_condition;
+static int running = 0;
+#ifdef DEBUG
+static bigint_t result;
+static bigint_t expected;
+#endif
+
+static void
+context_init(struct context *ctx, bigint_t integer)
+{
+ size_t n;
+
+ if (!integer) {
+ ctx->div_stack_size = 0;
+ ctx->div_stack = 0;
+ } else {
+ n = ctx->div_stack_size = bit_count(integer);
+ ctx->div_stack = emalloc(n * sizeof(bigint_t));
+ while (n--)
+ zinit(ctx->div_stack[n]);
+ }
+
+ zinit(ctx->div_n);
+ zinit(ctx->div_q);
+ zinit(ctx->div_r);
+ zinit(ctx->div_d);
+
+ zinit(ctx->factor);
+ zinit(ctx->d);
+ zinit(ctx->x);
+ zinit(ctx->y);
+ zinit(ctx->conj_a);
+ zinit(ctx->conj_b);
+
+#ifdef USE_TOMMATH
+ zinit(ctx->tmp);
+#endif
+}
+
+static void
+context_reinit(struct context *ctx, bigint_t integer)
+{
+ size_t i, n = bit_count(integer);
+
+ if (n > ctx->div_stack_size) {
+ i = ctx->div_stack_size;
+ ctx->div_stack_size = n;
+ ctx->div_stack = erealloc(ctx->div_stack, n * sizeof(bigint_t));
+ while (i < n)
+ zinit(ctx->div_stack[i++]);
+ }
+}
+
+static void
+context_free(struct context *ctx)
+{
+ size_t n;
+
+ for (n = ctx->div_stack_size; n--;)
+ zfree(ctx->div_stack[n]);
+ free(ctx->div_stack);
+
+ zfree(ctx->div_n);
+ zfree(ctx->div_q);
+ zfree(ctx->div_r);
+ zfree(ctx->div_d);
+
+ zfree(ctx->factor);
+ zfree(ctx->d);
+ zfree(ctx->x);
+ zfree(ctx->y);
+ zfree(ctx->conj_a);
+ zfree(ctx->conj_b);
+
+#ifdef USE_TOMMATH
+ zfree(ctx->tmp);
+#endif
+}
+
+static void
+parallel_init(void)
+{
+ if (sem_init(&semaphore, 0, parallel))
+ eprintf("sem_init:");
+ if ((errno = pthread_mutex_init(&print_mutex, 0)))
+ eprintf("pthread_mutex_init:");
+ if ((errno = pthread_mutex_init(&counter_mutex, 0)))
+ eprintf("pthread_mutex_init:");
+ if ((errno = pthread_cond_init(&counter_condition, 0)))
+ eprintf("pthread_cond_init:");
+}
+
+static void
+parallel_term(void)
+{
+ if (sem_destroy(&semaphore))
+ eprintf("sem_destroy:");
+ if ((errno = pthread_mutex_destroy(&print_mutex)))
+ eprintf("pthread_mutex_destroy:");
+ if ((errno = pthread_mutex_destroy(&counter_mutex)))
+ eprintf("pthread_mutex_destroy:");
+ if ((errno = pthread_cond_destroy(&counter_condition)))
+ eprintf("pthread_cond_destroy:");
+}
+
+static void
+output_primes_nonparallel(bigint_t factor, int is_prime, size_t power)
+{
+ const char *fstr = to_string(factor);
+ const char *qstr = is_prime == MAYBE ? questioned : "";
+ while (power--) {
+ printf(" %s%s", fstr, qstr);
+#ifdef DEBUG
+ zmul(result, result, factor);
+#endif
+ }
+}
+
+static void
+output_primes_parallel(bigint_t factor, int is_prime, size_t power)
+{
+ const char *fstr, *qstr = is_prime == MAYBE ? questioned : "";
+ if (pthread_mutex_lock(&print_mutex))
+ eprintf("pthread_mutex_lock:");
+ fstr = to_string(factor);
+ while (power--) {
+ printf(" %s%s", fstr, qstr);
+#ifdef DEBUG
+ zmul(result, result, factor);
+#endif
+ }
+ if (pthread_mutex_unlock(&print_mutex))
+ eprintf("pthread_mutex_unlock:");
+}
+
+static ssize_t
+iterated_division(struct context *ctx, bigint_t remainder, bigint_t numerator,
+ bigint_t denominator, size_t root_order, int is_prime)
+{
+ /*
+ * Just like n↑m by squaring, excepted this is iterated division.
+ */
+
+ const char *dstr = root_order ? to_string(denominator) : 0;
+ const char *nstr = is_prime == MAYBE ? questioned : "";
+ size_t partial_times = 1, times = 0, out, i;
+ bigint_t *n = &ctx->div_n, *q = &ctx->div_q, *r = &ctx->div_r, *d =
&ctx->div_d;
+ bigint_t *div_stack = ctx->div_stack;
+
+ zset(*n, numerator);
+ zset(*d, denominator);
+ zset(*div_stack++, denominator);
+
+ for (;;) {
+ zmul(*d, *d, *d);
+ if (zcmp(*d, *n) <= 0) {
+ zset(*div_stack++, *d);
+ partial_times <<= 1;
+ } else {
+ out = root_order * partial_times;
+ for (; partial_times; out >>= 1, partial_times >>= 1) {
+ div_mod(*q, *r, *n, *--div_stack);
+ if (!zcmp_i(*r, 0)) {
+ for (i = 0; i < out; i++)
+ printf(" %s%s", dstr, nstr);
+ times |= partial_times;
+ zset(*n, *q);
+ }
+ }
+ zset(remainder, *n);
+ return times;
+ }
+ }
+}
+
+static void
+subfactor_proper(struct context *ctx, bigint_t factor, bigint_t c, bigint_t
next, size_t root_order)
+{
+ size_t bits, cd;
+ bigint_t *d = &ctx->d, *x = &ctx->x, *y = &ctx->y;
+ bigint_t *conj_a = &ctx->conj_a, *conj_b = &ctx->conj_b;
+ int is_prime;
+ long seed;
+
+beginning:
+ bits = bit_count(factor);
+
+ if (bits < elementsof(pollards_rho_seeds))
+ seed = pollards_rho_seeds[bits];
+ else
+ seed = pollards_rho_seeds[elementsof(pollards_rho_seeds) - 1];
+
+ zadd_i(*x, c, seed);
+ zset(*y, *x);
+
+ for (;;) {
+ /*
+ * There may exist a number b = (A = ⌊√n⌋ + 1)² − n such that B
= √b is an integer
+ * in which case n = (A − B)(A + B) [n is(!) odd composite].
If not, the only the
+ * trivial iteration of A (A += 1) seems to be the one not
consuming your entire
+ * CPU and it is also guaranteed to find the factors, but it is
just so slow.
+ */
+ zsqrt(*conj_a, factor);
+ zadd_i(*conj_a, *conj_a, 1);
+ zmul(*conj_b, *conj_a, *conj_a);
+ zsub(*conj_b, *conj_b, factor);
+
+ if (zsqrt_test(*conj_b, *conj_b)) {
+ zadd(next, *conj_a, *conj_b);
+ zsub(factor, *conj_a, *conj_b);
+ subfactor(ctx, next, root_order, 0, 0);
+ subfactor(ctx, factor, root_order, c, next);
+ break;
+ }
+
+ /* Pollard's rho algorithm with Floyd's cycle-finding algorithm
and seed.
+ * A special-purpose integer factorisation algorithm used for
factoring
+ * integers with small factors. */
+ do {
+ zmul(*x, *x, *x), zadd_i(*x, *x, seed);
+ zmul(*y, *y, *y), zadd_i(*y, *y, seed);
+ zmul(*y, *y, *y), zadd_i(*y, *y, seed);
+ zmod(*x, *x, factor);
+ zmod(*y, *y, factor);
+
+ zsub(*d, *x, *y);
+ zabs(*d, *d);
+ gcd(*d, *d, factor);
+ } while (!zcmp_i(*d, 1));
+
+ if (!zcmp(factor, *d)) {
+ if ((is_prime = prime_test(factor))) {
+ output_primes(factor, is_prime, root_order);
+ break;
+ } else {
+ zadd_i(c, c, POLLARDS_RHO_SEED_INCREASEMENT);
+ goto beginning;
+ }
+ }
+
+ if ((is_prime = prime_test(factor))) {
+ iterated_division(ctx, factor, factor, *d, root_order,
is_prime);
+ if (is_factorised(factor))
+ break;
+ } else {
+ cd = iterated_division(ctx, factor, factor, *d, 0, 0);
+ zset(next, *d);
+ subfactor(ctx, next, root_order * cd, 0, 0);
+ if (is_factorised(factor))
+ break;
+ }
+
+ if ((is_prime = prime_test(factor))) {
+ output_primes(factor, is_prime, root_order);
+ break;
+ }
+ }
+}
+
+static void
+subfactor_nonparallel(struct context *ctx, bigint_t integer, size_t root_order,
+ bigint_t reuse_seed, bigint_t reuse_next)
+{
+ if (reuse_seed) {
+ zset_i(reuse_seed, POLLARDS_RHO_INITIAL_SEED);
+ subfactor_proper(ctx, integer, reuse_seed, reuse_next,
root_order);
+ } else {
+ bigint_t seed, next;
+ zinit(seed);
+ zset_i(seed, POLLARDS_RHO_INITIAL_SEED);
+ zinit(next);
+ subfactor_proper(ctx, integer, seed, next, root_order);
+ zfree(seed);
+ zfree(next);
+ }
+}
+
+static void *
+subfactor_parallel_threaded(void *data_)
+{
+ struct thread_data *data = data_;
+ struct context ctx;
+
+ if (sem_wait(&semaphore))
+ eprintf("sem_wait:");
+
+ output_primes(data->integer, 2, 1);
+
+ context_init(&ctx, data->integer);
+ subfactor_nonparallel(&ctx, data->integer, data->root_order, 0, 0);
+ context_free(&ctx);
+ zfree(data->integer);
+ free(data);
+
+ if (sem_post(&semaphore))
+ eprintf("sem_post:");
+
+ if ((errno = pthread_mutex_lock(&counter_mutex)))
+ eprintf("pthread_mutex_lock:");
+ if (--running == 0) {
+ if ((errno = pthread_cond_signal(&counter_condition)))
+ eprintf("pthread_cond_signal:");
+ }
+ if ((errno = pthread_mutex_unlock(&counter_mutex)))
+ eprintf("pthread_mutex_unlock:");
+
+ return NULL;
+}
+
+static void
+subfactor_parallel(struct context *ctx, bigint_t integer, size_t root_order,
+ bigint_t reuse_seed, bigint_t reuse_next)
+{
+ if (reuse_seed) {
+ subfactor_nonparallel(ctx, integer, root_order, reuse_seed,
reuse_next);
+ } else {
+ struct thread_data *data = emalloc(sizeof(*data));
+ pthread_t thread;
+ zinit(data->integer);
+ zset(data->integer, integer);
+ data->root_order = root_order;
+
+ if ((errno = pthread_mutex_lock(&counter_mutex)))
+ eprintf("pthread_mutex_lock:");
+ running++;
+
+ if ((errno = pthread_mutex_unlock(&counter_mutex)))
+ eprintf("pthread_mutex_unlock:");
+
+ if ((errno = pthread_create(&thread, 0,
subfactor_parallel_threaded, data)))
+ eprintf("pthread_create:");
+ if ((errno = pthread_detach(thread)))
+ eprintf("pthread_detach:");
+ }
+}
+
+static int
+factor(struct context *ctx, char *integer_str, bigint_t reusable_seed,
bigint_t reusable_next)
+{
+ bigint_t integer;
+ size_t i, power;
+ int is_prime;
+
+ if (!*integer_str)
+ goto invalid;
+ for (i = 0; integer_str[i]; i++)
+ if (!isdigit(integer_str[i]))
+ goto invalid;
+
+ zparse(integer, integer_str);
+#ifdef DEBUG
+ zset_i(result, 1);
+ zset(expected, integer);
+#endif
+
+ strbuf = integer_str;
+
+ while (*integer_str == '0' && *integer_str != 0) integer_str++;
+ printf("%s:", integer_str);
+
+ /* Behave like GNU factor: print empty set for 0 and 1, and pretend 0
is positive. */
+ if (zcmp_i(integer, 1) <= 0)
+ goto done;
+
+ /* Remove factors of two. */
+ power = lowest_bit(integer);
+ if (power > 0) {
+ shift_right(integer, integer, power);
+ while (power--) {
+ printf(" 2");
+#ifdef DEBUG
+ zmul_i(result, result, 2);
+#endif
+ }
+ if (is_factorised(integer))
+ goto done;
+ }
+
+ context_reinit(ctx, integer);
+
+ /* Remove factors of other tiny primes. */
+#ifdef DEBUG
+# define print_prime(factor) printf(" "#factor), zmul_i(result, result,
factor);
+#else
+# define print_prime(factor) printf(" "#factor);
+#endif
+#define X(factor)\
+ power = iterated_division(ctx, integer, integer, constants[factor], 0,
YES);\
+ if (power > 0) {\
+ while (power--)\
+ print_prime(factor);\
+ if (is_factorised(integer))\
+ goto done;\
+ }
+ LIST_CONSTANTS;
+#undef X
+
+ if ((is_prime = prime_test(integer))) {
+ output_primes(integer, is_prime, 1);
+ goto done;
+ }
+
+ if (parallel < 2) {
+ subfactor(ctx, integer, 1, reusable_seed, reusable_next);
+ } else {
+ if (sem_trywait(&semaphore))
+ eprintf("sem_trywait:");
+
+ running++;
+ subfactor(ctx, integer, 1, reusable_seed, reusable_next);
+
+ if (sem_post(&semaphore))
+ eprintf("sem_post:");
+
+ if ((errno = pthread_mutex_lock(&counter_mutex)))
+ eprintf("pthread_mutex_lock:");
+ if (--running != 0) {
+ if ((errno = pthread_cond_wait(&counter_condition,
&counter_mutex)))
+ eprintf("pthread_cond_wait:");
+ }
+ if ((errno = pthread_mutex_unlock(&counter_mutex)))
+ eprintf("pthread_mutex_unlock:");
+ }
+
+#ifdef DEBUG
+ if (zcmp(result, expected))
+ fprintf(stderr, "\033[1;31mIncorrect factorisation of
%s\033[m\n", to_string(expected));
+#endif
+
+done:
+ zfree(integer);
+ printf("\n");
+ return 0;
+invalid:
+ weprintf("%s is not a valid non-negative integer\n", integer_str);
+ return 1;
+}
+
+static void
+usage(void)
+{
+ eprintf("usage: %s [-c N] [-p N] [-q] [number ...]\n", argv0);
+}
+
+int
+main(int argc, char *argv[])
+{
+ long temp;
+ int ret = 0;
+ struct context ctx;
+ bigint_t reusable_seed, reusable_next;
+
+ ARGBEGIN {
+ case 'c':
+ temp = strtol(EARGF(usage()), NULL, 10);
+ if (temp < 1)
+ eprintf("value of -c must be a positive integer\n");
+ certainty = temp > INT_MAX ? INT_MAX : (int)temp;
+ break;
+ case 'p':
+ temp = atoi(EARGF(usage()));
+ if (temp < 1)
+ eprintf("value of -p must be a positive integer\n");
+ parallel = temp > INT_MAX ? INT_MAX : (int)temp;
+ parallel = temp > SEM_VALUE_MAX ? SEM_VALUE_MAX : (int)temp;
+ break;
+ case 'q':
+ questioned = "?";
+ break;
+ default:
+ usage();
+ } ARGEND;
+
+#define X(x) zinit(constants[x]), zset_i(constants[x], x);
+ LIST_CONSTANTS;
+#undef X
+
+ if (parallel > 1) {
+ output_primes = output_primes_parallel;
+ subfactor = subfactor_parallel;
+ parallel_init();
+ } else {
+ output_primes = output_primes_nonparallel;
+ subfactor = subfactor_nonparallel;
+ }
+
+ context_init(&ctx, 0);
+ zinit(reusable_seed);
+ zinit(reusable_next);
+#ifdef DEBUG
+ zinit(result);
+ zinit(expected);
+#endif
+
+ if (*argv) {
+ while (*argv)
+ ret |= factor(&ctx, *argv++, reusable_seed,
reusable_next);
+ } else {
+ ssize_t n;
+ size_t size = 0;
+ char *line = 0;
+ while ((n = getline(&line, &size, stdin)) >= 0) {
+ if (n > 0 && line[n - 1] == '\n')
+ n--;
+ line[n] = 0;
+ ret |= *line ? factor(&ctx, line, reusable_seed,
reusable_next) : 0;
+ }
+ free(line);
+ }
+
+ context_free(&ctx);
+ zfree(reusable_seed);
+ zfree(reusable_next);
+#ifdef DEBUG
+ zfree(result);
+ zfree(expected);
+#endif
+
+ if (parallel > 1)
+ parallel_term();
+
+#define X(x) zfree(constants[x]);
+ LIST_CONSTANTS;
+#undef X
+
+ return fshut(stdout, "<stdout>") || ret;
+}
--
2.7.1
