1 /*
2 Copyright (C) 2010 Fredrik Johansson
3 Copyright (C) 2012 Lina Kulakova
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 <http://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 "fmpz_mod_mat.h"
21 #include "ulong_extras.h"
22
23 void
fmpz_mod_mat_randrank(fmpz_mod_mat_t mat,flint_rand_t state,slong rank)24 fmpz_mod_mat_randrank(fmpz_mod_mat_t mat, flint_rand_t state, slong rank)
25 {
26 slong i;
27 fmpz * diag;
28
29 if (rank < 0 || rank > mat->mat->r || rank > mat->mat->c)
30 {
31 flint_printf("Exception (fmpz_mod_mat_randrank). Impossible rank.\n");
32 flint_abort();
33 }
34
35 diag = _fmpz_vec_init(rank);
36 for (i = 0; i < rank; i++)
37 {
38 fmpz_randtest_mod(&diag[i], state, mat->mod);
39 while (fmpz_is_zero(&diag[i]))
40 {
41 fmpz_randtest_mod(&diag[i], state, mat->mod);
42 }
43 }
44
45 fmpz_mat_randpermdiag(mat->mat, state, diag, rank);
46
47 _fmpz_vec_clear(diag, rank);
48 }
49
50 static void
check_rref(fmpz_mod_mat_t A)51 check_rref(fmpz_mod_mat_t A)
52 {
53 slong i, j, prev_pivot, prev_row_zero;
54
55 prev_row_zero = 0;
56 prev_pivot = -1;
57
58 for (i = 0; i < A->mat->r; i++)
59 {
60 for (j = 0; j < A->mat->c; j++)
61 {
62 /* Found nonzero entry */
63 if (!fmpz_is_zero(A->mat->rows[i] + j))
64 {
65 if (prev_row_zero)
66 {
67 flint_printf("nonzero row after zero row\n");
68 abort();
69 }
70
71 if (j <= prev_pivot)
72 {
73 flint_printf("pivot not strictly to the right of previous\n");
74 abort();
75 }
76
77 prev_pivot = j;
78 break;
79 }
80
81 prev_row_zero = (j + 1 == A->mat->c);
82 }
83 }
84 }
85
86
87 int
main(void)88 main(void)
89 {
90 fmpz_mod_mat_t A;
91 fmpz_t p;
92 slong i, m, n, d, r, rank;
93 slong *perm;
94
95 FLINT_TEST_INIT(state);
96
97 flint_printf("rref....");
98 fflush(stdout);
99
100 /* Maximally sparse matrices of given rank */
101 for (i = 0; i < 10000; i++)
102 {
103 m = n_randint(state, 10);
104 n = n_randint(state, 10);
105 perm = flint_malloc(FLINT_MAX(1, m) * sizeof(slong));
106
107 fmpz_init_set_ui(p, n_randtest_prime(state, 0));
108
109 for (r = 0; r <= FLINT_MIN(m, n); r++)
110 {
111 fmpz_mod_mat_init(A, m, n, p);
112
113 fmpz_mod_mat_randrank(A, state, r);
114
115 rank = fmpz_mod_mat_rref(perm, A);
116
117 if (r < rank)
118 {
119 fmpz_mod_mat_print_pretty(A);
120 flint_printf("FAIL:\n");
121 flint_printf("wrong rank!\n");
122 abort();
123 }
124
125 check_rref(A);
126
127 fmpz_mod_mat_clear(A);
128 }
129
130 fmpz_clear(p);
131 flint_free(perm);
132 }
133
134 /* Dense */
135 for (i = 0; i < 10000; i++)
136 {
137 m = n_randint(state, 5);
138 n = n_randint(state, 4);
139 perm = flint_malloc(FLINT_MAX(1, m) * sizeof(slong));
140
141 fmpz_init_set_ui(p, n_randtest_prime(state, 0));
142
143 for (r = 0; r <= FLINT_MIN(m, n); r++)
144 {
145 d = n_randint(state, 2 * m * n + 1);
146
147 fmpz_mod_mat_init(A, m, n, p);
148 fmpz_mod_mat_randrank(A, state, r);
149
150 fmpz_mat_randops(A->mat, state, d);
151
152 _fmpz_mod_mat_reduce(A);
153
154 rank = fmpz_mod_mat_rref(perm, A);
155
156 if (r < rank)
157 {
158 flint_printf("FAIL:\n");
159 flint_printf("wrong rank!\n");
160 abort();
161 }
162
163 check_rref(A);
164
165 fmpz_mod_mat_clear(A);
166 }
167
168 fmpz_clear(p);
169 flint_free(perm);
170 }
171
172 FLINT_TEST_CLEANUP(state);
173
174 flint_printf("PASS\n");
175 return 0;
176 }
177