1 /*
2 Copyright (C) 2010 Fredrik Johansson
3 Copyright (C) 2020 William Hart
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_vec.h"
19 #include "fmpz_mat.h"
20 #include "ulong_extras.h"
21
22
23 int
main(void)24 main(void)
25 {
26 fmpz_mat_t A, X, B, AX;
27 fmpz_t den;
28 slong i, m, n, k;
29 int success;
30
31 FLINT_TEST_INIT(state);
32
33 flint_printf("can_solve_fflu....");
34 fflush(stdout);
35
36 /* test random systems (likely not soluble) */
37 for (i = 0; i < 1000 * flint_test_multiplier(); i++)
38 {
39 m = n_randint(state, 20);
40 n = n_randint(state, 20);
41 k = n_randint(state, 20);
42
43 fmpz_mat_init(A, m, k);
44 fmpz_mat_init(B, m, n);
45 fmpz_mat_init(X, k, n);
46 fmpz_mat_init(AX, m, n);
47 fmpz_init(den);
48
49 fmpz_mat_randrank(A, state, n_randint(state, FLINT_MIN(m, k) + 1), 1+n_randint(state, 2)*n_randint(state, 100));
50 fmpz_mat_randtest(B, state, 1+n_randint(state, 2)*n_randint(state, 100));
51
52 /* Dense */
53 if (n_randint(state, 2))
54 fmpz_mat_randops(A, state, 1+n_randint(state, 1 + m*m));
55
56 success = fmpz_mat_can_solve_fflu(X, den, A, B);
57
58 if (success)
59 {
60 fmpz_mat_mul(AX, A, X);
61 fmpz_mat_scalar_divexact_fmpz(AX, AX, den);
62 }
63
64 if (success && !fmpz_mat_equal(AX, B))
65 {
66 flint_printf("FAIL:\n");
67 flint_printf("AX != B!\n");
68 flint_printf("A:\n"), fmpz_mat_print_pretty(A), flint_printf("\n");
69 flint_printf("B:\n"), fmpz_mat_print_pretty(B), flint_printf("\n");
70 flint_printf("X:\n"), fmpz_mat_print_pretty(X), flint_printf("\n");
71 flint_printf("den(X) = "), fmpz_print(den), flint_printf("\n");
72 flint_printf("AX:\n"), fmpz_mat_print_pretty(AX), flint_printf("\n");
73 abort();
74 }
75
76 fmpz_mat_clear(A);
77 fmpz_mat_clear(B);
78 fmpz_mat_clear(X);
79 fmpz_mat_clear(AX);
80 fmpz_clear(den);
81 }
82
83 /* test random soluble systems */
84 for (i = 0; i < 1000 * flint_test_multiplier(); i++)
85 {
86 m = n_randint(state, 20);
87 n = n_randint(state, 20);
88 k = n_randint(state, 20);
89
90 fmpz_mat_init(A, m, k);
91 fmpz_mat_init(B, m, n);
92 fmpz_mat_init(X, k, n);
93 fmpz_mat_init(AX, m, n);
94 fmpz_init(den);
95
96 fmpz_mat_randrank(A, state, n_randint(state, FLINT_MIN(m, k) + 1), 1+n_randint(state, 2)*n_randint(state, 100));
97 fmpz_mat_randtest(X, state, 1+n_randint(state, 2)*n_randint(state, 100));
98
99 /* Dense */
100 if (n_randint(state, 2))
101 fmpz_mat_randops(A, state, 1+n_randint(state, 1 + m*m));
102
103 fmpz_mat_mul(B, A, X);
104
105 fmpz_randtest_not_zero(den, state, 20);
106 fmpz_mat_scalar_mul_fmpz(A, A, den);
107
108 success = fmpz_mat_can_solve_fflu(X, den, A, B);
109
110 if (success)
111 {
112 fmpz_mat_mul(AX, A, X);
113 fmpz_mat_scalar_divexact_fmpz(AX, AX, den);
114 }
115
116 if (!success || !fmpz_mat_equal(AX, B))
117 {
118 flint_printf("FAIL:\n");
119 flint_printf("AX != B!\n");
120 flint_printf("A:\n"), fmpz_mat_print_pretty(A), flint_printf("\n");
121 flint_printf("B:\n"), fmpz_mat_print_pretty(B), flint_printf("\n");
122 flint_printf("X:\n"), fmpz_mat_print_pretty(X), flint_printf("\n");
123 flint_printf("den(X) = "), fmpz_print(den), flint_printf("\n");
124 flint_printf("AX:\n"), fmpz_mat_print_pretty(AX), flint_printf("\n");
125 abort();
126 }
127
128 fmpz_mat_clear(A);
129 fmpz_mat_clear(B);
130 fmpz_mat_clear(X);
131 fmpz_mat_clear(AX);
132 fmpz_clear(den);
133 }
134
135
136 FLINT_TEST_CLEANUP(state);
137
138 flint_printf("PASS\n");
139 return 0;
140 }
141