/*
Calculate terms for the following OEIS sequences:
A116893: Numbers k such that gcd(k! + 1, k^k + 1) > 1.
A116892: Values of gcd(k! + 1, k^k + 1), when greater than 1.
A116894: Numbers k such that gcd(k! + 1, k^k + 1) is neither 1 nor 2k + 1.
A134656: Corresponding GCD values in A116894.
This program prints 4 tables, each with columns 'n k',
for 1 <= lower_k <= k < upper_k.
The program can be run multiple times, for increasing ranges of k, and each
output merged.
Compile on Linux as, e.g.: clang -O3 -o a116893 a116893.c -lgmp
Invoke on Linux as, e.g.: ./a116893 1 1000001 1 1
where the first two parameters are, respectively, lower and upper bounds
for k, and the third and fourth parameters are the starting points for n
in, respectively, A116893/A116892 and A116894/A134656.
Created by Nick Hobson (n@nickhobson.com), Feb 19, 2024.
*/

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <gmp.h>
#include <errno.h>

#define COUNTOF(a) (sizeof(a) / sizeof(*(a)))

struct file {
    FILE *ptr;
    char name[12];
};

static void open_f(struct file *f)
{
    errno = 0;
    if (!(f->ptr = fopen(f->name, "w"))) {
        fprintf(stderr, "Error opening %s for write: %s\n", f->name,
                strerror(errno));
        exit(EXIT_FAILURE);
    }
    setvbuf(f->ptr, NULL, _IOLBF, 256);  // Ensure line buffering.
}

static void close_f(struct file *f)
{
    errno = 0;
    if (fclose(f->ptr)) {
        fprintf(stderr, "Error closing %s: %s\n", f->name, strerror(errno));
    }
}

int main(int argc, char *argv[])
{
    // No validation for input parameters, which should be positive integers.
    unsigned long lower_k = (argc > 1) ? strtoul(argv[1], NULL, 10) : 1;
    unsigned long upper_k = (argc > 2) ? strtoul(argv[2], NULL, 10) : 6000;
    unsigned long n[2];
    n[0] = (argc > 3) ? strtoul(argv[3], NULL, 10) : 1;
    n[1] = (argc > 4) ? strtoul(argv[4], NULL, 10) : 1;

    struct file f[] = {{.name = "a116893.txt"},
                       {.name = "a116894.txt"},
                       {.name = "a116892.txt"},
                       {.name = "a134656.txt"}};
    for (size_t i = 0; i < COUNTOF(f); i++) {
        open_f(&f[i]);
    }

    mpz_t fac, pow, gcd;
    mpz_inits(fac, pow, gcd, NULL);
    mpz_fac_ui(fac, lower_k - 1);
    for (unsigned long k = lower_k; k < upper_k; k++) {
        mpz_mul_ui(fac, fac, k);
        mpz_add_ui(fac, fac, 1);
        mpz_ui_pow_ui(pow, k, k);
        mpz_add_ui(pow, pow, 1);
        mpz_gcd(gcd, pow, fac);
        mpz_sub_ui(fac, fac, 1);
        // gcd > 1: write to a116893.txt and a116892.txt.
        if (mpz_cmp_ui(gcd, 1) != 0) {
            fprintf(f[0].ptr, "%lu %lu\n", n[0], k);
            fprintf(f[2].ptr, "%lu ", n[0]++);
            mpz_out_str(f[2].ptr, 10, gcd);
            fputc('\n', f[2].ptr);
            // gcd > 1 and gcd != 2k + 1: write to a116894.txt and a134656.txt.
            if (mpz_cmp_ui(gcd, 2 * k + 1) != 0) {
                fprintf(f[1].ptr, "%lu %lu\n", n[1], k);
                fprintf(f[3].ptr, "%lu ", n[1]++);
                mpz_out_str(f[3].ptr, 10, gcd);
                fputc('\n', f[3].ptr);
            }
        }
    }
    mpz_clears(fac, pow, gcd, NULL);

    for (size_t i = 0; i < COUNTOF(f); i++) {
        close_f(&f[i]);
    }
    return EXIT_SUCCESS;
}