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 <stdlib.h>
13 #include "flint.h"
14 #include "fmpz.h"
15 #include "fmpz_poly.h"
16 #include "fmpz_poly_mat.h"
17
18 slong
fmpz_poly_mat_nullspace(fmpz_poly_mat_t res,const fmpz_poly_mat_t mat)19 fmpz_poly_mat_nullspace(fmpz_poly_mat_t res, const fmpz_poly_mat_t mat)
20 {
21 slong i, j, k, n, rank, nullity;
22 slong * pivots;
23 slong * nonpivots;
24 fmpz_poly_mat_t tmp;
25 fmpz_poly_t den;
26
27 n = mat->c;
28
29 fmpz_poly_init(den);
30 fmpz_poly_mat_init_set(tmp, mat);
31 rank = fmpz_poly_mat_rref(tmp, den, tmp);
32 nullity = n - rank;
33
34 fmpz_poly_mat_zero(res);
35
36 if (rank == 0)
37 {
38 for (i = 0; i < nullity; i++)
39 fmpz_poly_set_ui(res->rows[i] + i, UWORD(1));
40 }
41 else if (nullity)
42 {
43 pivots = flint_malloc(rank * sizeof(slong));
44 nonpivots = flint_malloc(nullity * sizeof(slong));
45
46 for (i = j = k = 0; i < rank; i++)
47 {
48 while (fmpz_poly_is_zero(tmp->rows[i] + j))
49 {
50 nonpivots[k] = j;
51 k++;
52 j++;
53 }
54 pivots[i] = j;
55 j++;
56 }
57 while (k < nullity)
58 {
59 nonpivots[k] = j;
60 k++;
61 j++;
62 }
63
64 fmpz_poly_set(den, tmp->rows[0] + pivots[0]);
65
66 for (i = 0; i < nullity; i++)
67 {
68 for (j = 0; j < rank; j++)
69 fmpz_poly_set(res->rows[pivots[j]] + i,
70 tmp->rows[j] + nonpivots[i]);
71 fmpz_poly_neg(res->rows[nonpivots[i]] + i, den);
72 }
73
74 flint_free(pivots);
75 flint_free(nonpivots);
76 }
77
78 fmpz_poly_clear(den);
79 fmpz_poly_mat_clear(tmp);
80 return nullity;
81 }
82