/*
   Copyright (C) 2022 Stephen Hewitt, cambridgeclarion.org

   This file is part of a programme to research fast Fermat factorisation algorithms.

   This is free software: you can redistribute it and/or modify
   it under the terms of version 3 of the GNU General Public License as
   published by the Free Software Foundation.

   It is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   GNU General Public License for more details.

   You should have received a copy of the GNU General Public Licence
   along with this programme.  If not, see <https://www.gnu.org/licenses/>.

   A copy of the licence is also available at
   https://www.cambridgeclarion.org/download/GPLv3.txt
   It has the following cryptographic hash
   SHA256 3972dc9744f6499f0f9b2dbf76696f2ae7ad8af9b23dde66d6af86c9dfb36986

*/

/*
    This programme is for research on Fermat factorisation algorithms.

    This module for generating numbers of the form pq suitable for testing a Fermat factorisation algorithm.
    It categorises each pq it generates, in terms of which reduction ratio series A/B
    the number will produce in each of the first few prime bases.
    It prints the category designation and the number pq in hex.

    It generates p and q as (pseudo)random primes and relatively close together, as follows.

    It generates a number of length (for example) 500 bits, with the most significant bit set and all the other bits random.
    Then it finds the next highest prime, which it takes as q.
    It generates another number of (for example) 272 bits, again with the most significant bit set and all the other bits random.

    It then finds the next highest prime after (q + this number).  It accepts this prime as p.
    It outputs pq and discards p and q.

    It uses the GMP library for probabilistic identification of primes.

    For more information see the articles on www.cambridgeclarion.org
    "Empirical explorations of faster Fermat factorisation" by Stephen Hewitt
*/

#include <cstdio>
#include <cstring>
#include <vector>
#include <algorithm>

#include "gmp.h"

const char* lsd2series[] = {
#include "lsd2series.gh"
};

// calculate the category value which is a sort of pun in hex eg 0xABB3

int characterise(mpz_t pq, unsigned requested_primes_n, unsigned long* category)
{
    unsigned long lsb3 = mpz_fdiv_ui(pq, 8); // first base 2 which is special case
    unsigned long designator;

    if (!lsb3 & 1)
        return 89; // illegal binary ending: N is even
    else if ((lsb3 & 3) == 3)
        designator = 3;
    else
        designator = lsb3;

    unsigned long A_designator = 0xa;
    unsigned long B_designator = 0xb;

    for (unsigned prime_index = 1; prime_index < requested_primes_n; prime_index++)
    {
        unsigned long base = strlen(lsd2series[prime_index]); // the total number of possible lsd values, which is indeed the base
        unsigned lsd = mpz_fdiv_ui(pq, base);

        A_designator <<= 4;
        B_designator <<= 4;
        designator |= (lsd2series[prime_index][lsd] == 'A') ? A_designator: B_designator;
    }
    *category = designator;
    return 0;
}

// needs to be shared between main and functions that need a random number for now
static gmp_randstate_t state;

void print_random_pq(unsigned p_bits_n, unsigned pq_diff_bits_n)
{
    mpz_t pp;
    mpz_init(pp);
    mpz_t qq;
    mpz_init(qq);
    mpz_t temp;
    mpz_init(temp);

    mpz_urandomb(qq, state, p_bits_n);
    mpz_setbit(qq, p_bits_n - 1); // ensure leading bit is set.
    mpz_nextprime(qq, qq);

    mpz_urandomb(temp, state, pq_diff_bits_n);
    mpz_setbit(temp, pq_diff_bits_n - 1); // ensure leading bit is set.
    mpz_add(pp, temp, qq);
    mpz_nextprime(pp, pp);

    mpz_mul(temp, pp, qq);

    unsigned long category;

    characterise(temp, 4, &category);
    printf("0x%lX, ", category);
    putchar('"');
    mpz_out_str(0, 16, temp);
    putchar('"');
    putchar(',');
    putchar('\n');
}

int main(int argc, char* argv[])
{
    mpz_t pp;
    mpz_init(pp);
    mpz_t qq;
    mpz_init(qq);
    mpz_t temp;
    mpz_init(temp);

    const char* random_seed_string = "49478398549667362729830029384765444327361293891766278172656988327532392307475";

    mp_limb_t tmp;
    if (mpz_set_str(temp, random_seed_string, 10))
        return 29;
    gmp_randinit_mt(state);
    gmp_randseed(state, temp);

    for (unsigned i = 0; i < 200; i++)
        print_random_pq(500, 272);

    return 0;
}