Cod sursa(job #3236631)

Utilizator mgntMarius B mgnt Data 29 iunie 2024 19:44:44
Problema Dirichlet Scor 100
Compilator c-64 Status done
Runda Arhiva de probleme Marime 3.53 kb


#include <stdio.h>

#include <stdlib.h>

#include <inttypes.h>

#include <stddef.h>



void dirichlet_program__extended_euclid (uint64_t  a, uint64_t  b, uint64_t  * s2, uint64_t  * t2, uint64_t  * d2, uint64_t  * hint) {


    /*      T   Given a, b ∈ ℤ, 0 ≤ a, b, 0 < a + b, it returns d, s, t ∈ ℤ, s.t. d = a·s + b·t.
    **      T
    **      T   The correctness, the time complexity, the overflow issues, all on Wikipedia.
    **      T
    **      T   https://en.wikipedia.org/wiki/Extended_Euclidean_algorithm
    */


    uint64_t  old_r = a;    uint64_t  r = b;

    uint64_t  old_s = 1;    uint64_t  s = 0;

    uint64_t  old_t = 0;    uint64_t  t = 1;

                            uint64_t  q = 0;

                            uint64_t  x = 0;

                            uint64_t  f = 1;

    while (r != 0) {

        q = old_r / r;

        f = 1 ^ f;

        x = old_r;       old_r = r;       r = x - q * r;

        x = old_s;       old_s = s;       s = x - q * s;

        x = old_t;       old_t = t;       t = x - q * t;


    }



    *s2 = old_s;

    *t2 = old_t;

    *d2 = old_r;

    *hint = f;

}


uint64_t dirichlet_program__modular_inverse (uint64_t  a, uint64_t  m) {

    /*      T   Returns t, such that  0 ≤ t < m, and  t·a ≡ 1 (mod m), on condition that  a, m, two positive integers, relatively prime.  */


    uint64_t  const b = m;

    uint64_t  s = 0;        uint64_t  t = 0;        uint64_t  d = 0;        uint64_t  hint = 0;

    dirichlet_program__extended_euclid(a, b, &s, &t, &d, &hint);

    return (hint ? s : (b + s));

}


uint64_t  dirichlet_program__factorial (uint64_t  n, uint64_t  m, uint64_t  n0, uint64_t  fact_n0_mod_m) {

    /*      T   Returns  n! (mod m), given n, m, n₀, and n₀! mod m, 1 ≤ n₀ ≤ n, etc.  */

    uint64_t  product = fact_n0_mod_m;

    uint64_t  factor  = n0 + 1;

    uint64_t  last_factor = n;


    while (factor <= last_factor) {

        product = (product * factor) % m;

        factor += 1;

    }

    return product;

}


int main (void) {


    /*      T   The number of the counted sequences of size n  is the  n -th catalan number Cₙ = (1 / (n + 1))·(²ⁿ ₙ) = (1 / (n + 1))·(((2·n)!)/((n!)²));
    **      T
    **      T   m := 9999991 is a prime,  n ≤ 1.000.000,  and  (2·n)² < m;   modular arithmetic with uint64 will do just fine!
    */

    int  exit_status = 0;

    uint64_t  const m = 9999991;

    FILE  * file_f = NULL;

    uint64_t  n = 0;

    uint64_t  fact_n_mod_m = 0;

    uint64_t  fact_2n_mod_m = 0;

    uint64_t  inv_fact_n_mod_m = 0;

    uint64_t  inv_n_plus_1_mod_m = 0;

    uint64_t  total_count_mod_m = 0;

    FILE  * file_g = NULL;


    /*  read the input */

    file_f = fopen("dirichlet.in", "r");

    (void) fscanf(file_f, "%" SCNu64, &n);

    (void) fclose(file_f);

    file_f = NULL;


    /* solve */

    fact_n_mod_m = dirichlet_program__factorial(n, m, 1, 1);

    fact_2n_mod_m = dirichlet_program__factorial(n + n, m, n, fact_n_mod_m);

    inv_fact_n_mod_m = dirichlet_program__modular_inverse(fact_n_mod_m, m);

    inv_n_plus_1_mod_m = dirichlet_program__modular_inverse(n + 1, m);


    total_count_mod_m = (fact_2n_mod_m * inv_fact_n_mod_m) % m;

    total_count_mod_m = (total_count_mod_m * inv_fact_n_mod_m) % m;

    total_count_mod_m = (total_count_mod_m * inv_n_plus_1_mod_m) % m;


    /*  write output */

    file_g = fopen("dirichlet.out", "w");

    (void) fprintf(file_g, "%" PRIu64 "\n", total_count_mod_m);

    (void) fclose(file_g);

    file_g = NULL;


    /*  pass the baton  */

    return exit_status;

}