bernoulli/bernoulli.cpp at main · JL2014/bernoulli · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
// g++ -o bernoulli bernoulli.cpp -lflint -L/usr/local/lib
#include <iostream>
#include <string>
#include <flint/fmpz.h>
#include <flint/fmpq.h>
#include <flint/arith.h>
#include <flint/fmpz_factor.h>

using namespace std;

int main(int argc, char* argv[])
{
    if (argc < 2) {
        cerr << "Usage: " << argv[0] << " n" << endl;
        return 1;
    }

    unsigned long n = stoul(argv[1]);
    fmpq_t Bn;
    fmpq_init(Bn);

    arith_bernoulli_number(Bn, n);

    cout << "B(" << n << ") = ";
    fmpz_print(fmpq_numref(Bn));
    cout << "/";
    fmpz_print(fmpq_denref(Bn));
    cout << endl;

    fmpz_t denom;
    fmpz_init(denom);
    fmpz_set(denom, fmpq_denref(Bn));

    fmpz_factor_t fac;
    fmpz_factor_init(fac);
    fmpz_factor(fac, denom);

    slong num_factors = fac->num;
    cout << endl << "Facteurs premiers du dénominateur (" << num_factors << "):" << endl;

    for (slong i = 0; i < num_factors; i++) {
        fmpz_print(fac->p + i);
        cout << ",";
    }
    cout << endl;

    /*for (slong i = 0; i < num_factors; i++) {
        fmpz_t p, pm1, two_n, rem;
        fmpz_init(p);
        fmpz_init(pm1);
        fmpz_init(two_n);
        fmpz_init(rem);

        fmpz_set(p, fac->p + i);
        fmpz_sub_ui(pm1, p, 1); // pm1 = p - 1
        fmpz_set_ui(two_n, 2 * n);

        fmpz_mod(rem, two_n, pm1);
        if (fmpz_is_zero(rem)) {
            fmpz_print(p);
            cout << "  (car " << (2 * n) << " divisible par (p-1)=";
            fmpz_print(pm1);
            cout << ")" << endl;
        }

        fmpz_clear(p);
        fmpz_clear(pm1);
        fmpz_clear(two_n);
        fmpz_clear(rem);
    }*/

    fmpz_clear(denom);
    fmpz_factor_clear(fac);
    fmpq_clear(Bn);

    return 0;
}