1 /*
2     Copyright (C) 2011 Fredrik Johansson
3 
4     This file is part of FLINT.
5 
6     FLINT is free software: you can redistribute it and/or modify it under
7     the terms of the GNU Lesser General Public License (LGPL) as published
8     by the Free Software Foundation; either version 2.1 of the License, or
9     (at your option) any later version.  See <http://www.gnu.org/licenses/>.
10 */
11 
12 #include <stdio.h>
13 #include <stdlib.h>
14 #include <gmp.h>
15 #include <mpfr.h>
16 #include "flint.h"
17 #include "arith.h"
18 #include "profiler.h"
19 #include "fmpz.h"
20 #include "fmpz_mat.h"
21 #include "fmpq_poly.h"
22 
23 
main()24 int main()
25 {
26     fmpz * num1;
27     fmpz * num2;
28     fmpz * num3;
29     fmpz * den1;
30     fmpz * den2;
31     fmpz * den3;
32     slong i, n, N;
33 
34     FLINT_TEST_INIT(state);
35 
36     flint_printf("bernoulli_number_vec....");
37     fflush(stdout);
38 
39     N = 2000;
40 
41     num1 = _fmpz_vec_init(N);
42     num2 = _fmpz_vec_init(N);
43     num3 = _fmpz_vec_init(N);
44     den1 = _fmpz_vec_init(N);
45     den2 = _fmpz_vec_init(N);
46     den3 = _fmpz_vec_init(N);
47 
48     for (n = 0; n < N; n += (n<100) ? 1 : n/3)
49     {
50         _arith_bernoulli_number_vec_recursive(num1, den1, n);
51         _arith_bernoulli_number_vec_multi_mod(num2, den2, n);
52         _arith_bernoulli_number_vec_zeta(num3, den3, n);
53 
54         for (i = 0; i < n; i++)
55         {
56             if (!fmpz_equal(num1 + i, num2 + i) ||
57                 !fmpz_equal(num1 + i, num3 + i))
58             {
59                 flint_printf("FAIL: n = %wd, numerator of B_%wd\n", n, i);
60                 flint_printf("recursive: "); fmpz_print(num1 + i); flint_printf("\n");
61                 flint_printf("multi_mod: "); fmpz_print(num2 + i); flint_printf("\n");
62                 flint_printf("zeta:      "); fmpz_print(num3 + i); flint_printf("\n");
63                 abort();
64             }
65 
66             if (!fmpz_equal(den1 + i, den2 + i) ||
67                 !fmpz_equal(den1 + i, den3 + i))
68             {
69                 flint_printf("FAIL: n = %wd, denominator of B_%wd\n", n, i);
70                 flint_printf("recursive: "); fmpz_print(den1 + i); flint_printf("\n");
71                 flint_printf("multi_mod: "); fmpz_print(den2 + i); flint_printf("\n");
72                 flint_printf("zeta:      "); fmpz_print(den3 + i); flint_printf("\n");
73                 abort();
74             }
75         }
76     }
77 
78     _fmpz_vec_clear(num1, N);
79     _fmpz_vec_clear(num2, N);
80     _fmpz_vec_clear(num3, N);
81     _fmpz_vec_clear(den1, N);
82     _fmpz_vec_clear(den2, N);
83     _fmpz_vec_clear(den3, N);
84 
85     FLINT_TEST_CLEANUP(state);
86     flint_printf("PASS\n");
87     return 0;
88 }
89