1 /*
2 Copyright (C) 2011 Fredrik Johansson
3 Copyright (C) 2012 Sebastian Pancratz
4
5 This file is part of FLINT.
6
7 FLINT is free software: you can redistribute it and/or modify it under
8 the terms of the GNU Lesser General Public License (LGPL) as published
9 by the Free Software Foundation; either version 2.1 of the License, or
10 (at your option) any later version. See <https://www.gnu.org/licenses/>.
11 */
12
13 #include <stdio.h>
14 #include <stdlib.h>
15 #include <gmp.h>
16 #include "flint.h"
17 #include "fmpz.h"
18 #include "fmpz_mat.h"
19 #include "ulong_extras.h"
20
21 int
main(void)22 main(void)
23 {
24 slong m, n, rep, i;
25 FLINT_TEST_INIT(state);
26
27 flint_printf("charpoly_berkowitz....");
28 fflush(stdout);
29
30 for (rep = 0; rep < 1000 * flint_test_multiplier(); rep++)
31 {
32 fmpz_mat_t A, B, C, D;
33 fmpz_poly_t f, g;
34
35 m = n_randint(state, 4);
36 n = m;
37
38 fmpz_mat_init(A, m, n);
39 fmpz_mat_init(B, m, n);
40 fmpz_mat_init(C, m, m);
41 fmpz_mat_init(D, n, n);
42 fmpz_poly_init(f);
43 fmpz_poly_init(g);
44
45 fmpz_mat_randtest(A, state, 10);
46 fmpz_mat_randtest(B, state, 10);
47
48 fmpz_mat_mul(C, A, B);
49 fmpz_mat_mul(D, B, A);
50
51 fmpz_mat_charpoly_berkowitz(f, C);
52 fmpz_mat_charpoly_berkowitz(g, D);
53
54 if (!fmpz_poly_equal(f, g))
55 {
56 flint_printf("FAIL: charpoly(AB) != charpoly(BA).\n");
57 flint_printf("Matrix A:\n"), fmpz_mat_print(A), flint_printf("\n");
58 flint_printf("Matrix B:\n"), fmpz_mat_print(B), flint_printf("\n");
59 flint_printf("cp(AB) = "), fmpz_poly_print_pretty(f, "X"), flint_printf("\n");
60 flint_printf("cp(BA) = "), fmpz_poly_print_pretty(g, "X"), flint_printf("\n");
61 abort();
62 }
63
64 fmpz_mat_clear(A);
65 fmpz_mat_clear(B);
66 fmpz_mat_clear(C);
67 fmpz_mat_clear(D);
68 fmpz_poly_clear(f);
69 fmpz_poly_clear(g);
70 }
71
72 for (rep = 0; rep < 1000 * flint_test_multiplier(); rep++)
73 {
74 fmpz_t c;
75 fmpz_mat_t A, B;
76 fmpz_poly_t f, g;
77
78 m = n_randint(state, 4);
79 n = m;
80
81 fmpz_init(c);
82 fmpz_mat_init(A, m, n);
83 fmpz_mat_init(B, m, n);
84 fmpz_poly_init(f);
85 fmpz_poly_init(g);
86
87 fmpz_mat_randtest(A, state, 10);
88 fmpz_mat_set(B, A);
89
90 for (i = 0; i < 10; i++)
91 {
92 fmpz_randtest(c, state, 5);
93 fmpz_mat_similarity(B, n_randint(state, m), c);
94 }
95
96 fmpz_mat_charpoly_berkowitz(f, A);
97 fmpz_mat_charpoly_berkowitz(g, B);
98
99 if (!fmpz_poly_equal(f, g))
100 {
101 flint_printf("FAIL: charpoly(A) != charpoly(P^{-1}AP).\n");
102 flint_printf("Matrix A:\n"), fmpz_mat_print(A), flint_printf("\n");
103 flint_printf("Matrix P^{-1}AP:\n"), fmpz_mat_print(B), flint_printf("\n");
104 flint_printf("cp(A) = "), fmpz_poly_print_pretty(f, "X"), flint_printf("\n");
105 flint_printf("cp(P^{-1}AP) = "), fmpz_poly_print_pretty(g, "X"), flint_printf("\n");
106 abort();
107 }
108
109 fmpz_clear(c);
110 fmpz_mat_clear(A);
111 fmpz_mat_clear(B);
112 fmpz_poly_clear(f);
113 fmpz_poly_clear(g);
114 }
115
116 FLINT_TEST_CLEANUP(state);
117
118 flint_printf("PASS\n");
119 return 0;
120 }
121