Characteristic and minimal polynomials

fmpz_mat_charpoly computes the characteristic polynomial $\chi_A(x) = \det(xI - A)$ of a square integer matrix, and fmpz_mat_minpoly its minimal polynomial, the monic polynomial $\mu_A$ of least degree with $\mu_A(A) = 0$. Both have integer coefficients, and $\mu_A$ divides $\chi_A$. For all but tiny matrices both are computed modulo word-size primes and reconstructed by the Chinese remainder theorem; the results are always proven.

Characteristic polynomial#

fmpz_mat_charpoly chooses the method from the dimension $n$ and the bit size $b$ of the largest entry:

fmpz_mat_charpoly_berkowitz is a wrapper of the generic _gr_mat_charpoly_berkowitz, which uses $O(n^4)$ additions and multiplications and no divisions (Berkowitz's algorithm); its intermediate results are coefficients of characteristic polynomials of submatrices, so their size stays controlled.

fmpz_mat_charpoly_modular computes the bound $C$ of fmpz_mat_charpoly_bound on the absolute values of the coefficients, chooses $\lceil (\log_2 C + 1)/63 \rceil$ primes above $2^{63}$, reduces $A$ modulo all of them at once with fmpz_mat_multi_mod_ui, computes the characteristic polynomials modulo the primes in parallel with nmod_mat_charpoly ($O(n^\omega)$ or $O(n^3)$ operations each), and reconstructs the coefficients in the symmetric range with _fmpz_vec_multi_CRT_ui. Since the number of primes is fixed by the a priori bound, there is no early termination. With $\log_2 C \approx n(b + \log_2 n)$, about $n(b + \log_2 n)/63$ characteristic polynomials modulo primes are needed, each costing $O(n^3)$ word operations, so the total is $O(n^4 (b + \log n))$ word operations, plus the multimodular reduction and reconstruction.

The bound: the coefficient of $x^{n-k}$ is $\pm$ the sum of the $\binom{n}{k}$ principal $k \times k$ minors of $A$, each of which is at most the product of the $k$ largest row norms (and of the $k$ largest column norms) by Hadamard's inequality.

Minimal polynomial#

fmpz_mat_minpoly and fmpz_mat_minpoly_modular compute minimal polynomials modulo primes $p$ above $2^{63}$ with nmod_mat_minpoly_with_gens, which also records which standard basis vectors were used as generators of Krylov subspaces. For each prime:

• if the degree modulo $p$ is smaller than the degree seen before, $p$ is discarded: the minimal polynomial modulo $p$ divides the reduction of $\mu_A$, so a smaller degree identifies a bad prime;
• if it is larger, all previous primes were bad and are discarded;
• otherwise the coefficients are combined by the Chinese remainder theorem.

When the reconstruction is unchanged by a new prime, the candidate $f$ is tested: $f(A) e_i = 0$ is checked exactly for every recorded generator $e_i$, by $\deg f$ matrix-vector products. The Krylov spaces of these generators span $\mathbb{Q}^n$ (they span $\mathbb{F}_p^n$), so the test proves $\mu_A \mid f$; and $\deg f$, the degree modulo a prime, is at most $\deg \mu_A$. Hence $f = \mu_A$. A bound on the coefficients of $\mu_A$ [Dum2007], from a bound $\rho$ on the spectral radius obtained with ovals of Cassini, caps the number of primes.

Related functions#

• nmod_mat_charpoly and nmod_mat_minpoly work modulo a word-size modulus.
• fmpq_mat_charpoly and fmpq_mat_minpoly scale a rational matrix to an integer matrix by the common denominator $q$ of its entries, call these functions, and rescale the coefficient of $x^{n-k}$ by $q^{-k}$.
• fmpz_poly_factor factors the result, for example to find the rational eigenvalues.

Example#

A $10 \times 10$ matrix with a repeated eigenvalue structure: the characteristic polynomial has degree 10 while the minimal polynomial has degree 4. The Berkowitz and modular algorithms agree:

#include <flint/flint.h>
#include <flint/fmpz.h>
#include <flint/fmpz_poly.h>
#include <flint/fmpz_mat.h>

int main(void)
{
    fmpz_mat_t A;
    fmpz_poly_t c, c2, m;
    slong i, n = 10;

    fmpz_mat_init(A, n, n);
    fmpz_poly_init(c);
    fmpz_poly_init(c2);
    fmpz_poly_init(m);

    /* block diagonal: eigenvalue 2 (blocks of sizes 2, 2, 1),
       eigenvalue -1 (blocks of sizes 2, 1, 1, 1) */
    for (i = 0; i < 5; i++)
        fmpz_set_si(fmpz_mat_entry(A, i, i), 2);
    for (i = 5; i < 10; i++)
        fmpz_set_si(fmpz_mat_entry(A, i, i), -1);
    fmpz_one(fmpz_mat_entry(A, 0, 1));
    fmpz_one(fmpz_mat_entry(A, 2, 3));
    fmpz_one(fmpz_mat_entry(A, 5, 6));

    /* conjugate by a unimodular matrix to hide the structure */
    {
        fmpz_mat_t U, Ui, T;
        fmpz_t den;
        flint_rand_t state;

        flint_rand_init(state);
        fmpz_init(den);
        fmpz_mat_init(U, n, n);
        fmpz_mat_init(Ui, n, n);
        fmpz_mat_init(T, n, n);
        fmpz_mat_one(U);
        fmpz_mat_randops(U, state, 30);
        fmpz_mat_inv(Ui, den, U);   /* den = +-1 */
        fmpz_mat_mul(T, U, A);
        fmpz_mat_mul(A, T, Ui);
        fmpz_mat_scalar_divexact_fmpz(A, A, den);
        fmpz_mat_clear(U);
        fmpz_mat_clear(Ui);
        fmpz_mat_clear(T);
        fmpz_clear(den);
        flint_rand_clear(state);
    }

    fmpz_mat_charpoly(c, A);
    flint_printf("charpoly: %{fmpz_poly}\n", c);
    fmpz_mat_charpoly_berkowitz(c2, A);
    flint_printf("berkowitz agrees: %d\n", fmpz_poly_equal(c, c2));
    fmpz_mat_charpoly_modular(c2, A);
    flint_printf("modular agrees: %d\n", fmpz_poly_equal(c, c2));

    fmpz_mat_minpoly(m, A);
    flint_printf("minpoly: %{fmpz_poly}\n", m);

    fmpz_mat_clear(A);
    fmpz_poly_clear(c);
    fmpz_poly_clear(c2);
    fmpz_poly_clear(m);
    return 0;
}

charpoly: x^10 - 5 * x^9 + 30 * x^7 - 15 * x^6 - 81 * x^5 + 30 * x^4 + 120 * x^3 - 80 * x - 32
berkowitz agrees: 1
modular agrees: 1
minpoly: x^4 - 2 * x^3 - 3 * x^2 + 4 * x + 4

The characteristic polynomial is $(x - 2)^5 (x + 1)^5$ and the minimal polynomial $(x - 2)^2 (x + 1)^2$, the exponents being the sizes of the largest Jordan blocks.