Ian Grant

/*
 * cnr.c
 *
 * (c) Ian Douglas Grant 2010.  All rights reserved.
 *
 * A command line interface tool for calculation of and displaying an iteration
 * of combinations of n choose r elements.
 *
 */



/* preprocessor statements */

#include <stdlib.h>     /* standard C includes first */
#include <stdio.h>
#include <ctype.h>
#include <string.h>

#include <gmp.h>

/*
 * Multiple precision data types and operations.
 * ---------------------------------------------
 *
 * This program requires arbitrary precision. This is achieved through the use
 * of the GNU Multiple Precision library.
 *
 * gmp.h :  needs to be included after stdio.h for certain prototypes required
 * for the multiple precision arithmetic calls mpz_* etc.  this pretty much
 * limits this code to unix-like OSs but I am not sure, maybe cygwin and
 * knowing what you are doing on win32 might help.
 *
 * http://gmplib.org/
 *
 * As far as I know this library is clearly the best option around. NetBSD
 * pkgsrc should have devel/gmp and FreeBSD ports collection has
 * math/libgmp4.
 */




/*
 * function prototypes
 *
 */

void        usage     (void);
void        factorial (const unsigned long, mpz_t);




/*
 * function: main() descr:    program start point & main functionality.
 */

int
main(int argc, char *argv[])
{
    mpz_t       term     [4];   /* local declarations */
    mpz_t       cnr;
    mpz_t       i;
    int            *comb = 0;
    int     cmp = 0;
    unsigned long   n = 0, r = 0;
    unsigned long   m = 0, max = 0, j = 0, x = 0;
    char        ch = ' ';
    char        str       [80];
    char           *p;
    int     loop = 1;

    if (argc != 3) {    /* input arg sanity checking &
                 * initialisations */
        usage();
        return 1;
    }
    for (x = 0; x < strlen(argv[1]); x++) {
        if (isdigit(argv[1][x]) == 0) {
            usage();
            return 1;
        }
    }
    for (x = 0; x < strlen(argv[2]); x++) {
        if (isdigit(argv[2][x]) == 0) {
            usage();
            return 1;
        }
    }
    if ((n = atoi(argv[1])) < 0 || n > 99) {
        usage();
        return 1;
    }
    if ((r = atoi(argv[2])) < 0 || r > n) {
        usage();
        return 1;
    }
    if ((comb = malloc(r * sizeof(int))) == 0) {
        printf("malloc() failure\n");
        return 2;
    }
    mpz_init(term[0]);
    mpz_init(term[1]);
    mpz_init(term[2]);
    mpz_init(term[3]);
    mpz_init(cnr);
    mpz_init(i);

    factorial(n, term[0]);  /* main functionality begins */
    factorial(n - r, term[1]);
    factorial(r, term[2]);
    mpz_mul(term[3], term[1], term[2]);
    mpz_fdiv_q(cnr, term[0], term[3]);  /* C(n,r) = n!/(n-r)!r!  */
    printf("C(%d,%d) = ", n, r);
    mpz_out_str(NULL, 10, cnr);
    while (loop) {      /* user interaction loop */
        printf("\nIterate combinations ? > ");
        p = str;
        while ((*p++ = getchar()) != '\n') {
            if ((str - p) > 50)
                break;
        }
        *p = '\0';
        ch = str[0];
        if (ch == 'n' || ch == 'N') {
            printf("OK done.\n");   /* handle choice made */
            return 0;
        } else if (ch == 'y' || ch == 'Y')
            loop = 0;
    }
    for (x = 0; x < r; x++) /* C(n, r) iterator code block */
        comb[x] = x + 1;
    for (x = 0; x < r; x++)
        printf("%d ", comb[x]);
    printf("\n");
    for (mpz_set_ui(i, 1); cmp = mpz_cmp(i, cnr); mpz_add_ui(i, i, 1)) {
        m = r - 1, max = n;
        while (comb[m] == max)
            m--, max--; /* this algorithm was written */
        comb[m]++;  /* in C from the pseudocode notation */
        for (j = m + 1; j < r; j++) /* within Discrete
                         * Mathematics 5th ed. */
            comb[j] = comb[j - 1] + 1;  /* by Richard
                             * Johnsonbaugh */
        for (x = 0; x < r; x++)
            printf("%d ", comb[x]);
        printf("\n");
    }
    free(comb);
    return 0;
}




/*
 * function: usage() descr:  prints a response to incorrect input arguments
 */
void
usage(void)
{
    printf("usage: cnr [n] [r]\n");
    printf("       Where n is an integer\n");
    printf("       larger than -1 and smaller\n");
    printf("       than 100;  where r is an\n");
    printf("       integer larger than 0 and smaller\n");
    printf("       than or equal to n.\n");
}




/*
 * function: factorial() descr: n(n-1)! note: make sure you have called
 * mpz_init() on the  argument mpz_t rop. factorial() follows the form of gmp
 * functions and is return void. You may or may not wish to keep this snippet
 * in such a form. returning a mpz_t would not be a problem.
 */

void
factorial(const unsigned long n, mpz_t rop)
{
    mpz_t       answer;
    unsigned long   x = 0;

    mpz_init(answer);

    if (!n) {       /* A very important fact as */
        mpz_set_ui(rop, 1); /* C(0, 0) = 1 */
        return;
    }
    x = n - 1;
    mpz_set_ui(answer, n);
    while (x) {     /* Iteration avoids function call */
        mpz_mul_ui(answer, answer, x);  /* overhead */
        x--;
    }
    mpz_set(rop, answer);
    return;
}