Bareiss: Division-free Linear Algebra for computing Nullspace

Python file for an implementation of Bareiss’ algorithm

This module offers an implementation of Bareiss algorithm (https://www.ams.org/journals/mcom/1968-22-103/S0025-5718-1968-0226829-0/). Such algorithm computes (fraction-free) an echelon form of a matrix whose coefficients are in an integral domain.

AUTHORS:

  • Antonio Jimenez-Pastor (2016-10-01): initial version

class ajpastor.misc.bareiss.BareissAlgorithm(parent, matrix, is_zero=<function BareissAlgorithm.<lambda> at 0x7f193190c0d0>, relations=[])

Bases: ajpastor.misc.linear_solver.LinearSystemSolver

This class represents the application of the Bareiss’ algorithm over a matrix with polynomial coefficients.

Bareiss’ algorithm is a division-free algorithm to compute an echelon form of a matrix over an integral domain. The main idea from Bareiss algorithm is to perform some type of Gauss-Jordan elimination on the matrix but keeping track that we never get denominators and ensuring that, in the end, the main diagonal has the extra property that \(d_{i+1,i+1}\ |\ d_{i,i}\).

This implementation only works with polynomial coefficients, however we allow to provide a method that solves the membership problem for the ideal \(I = \{p(X) \in R[X]\ |\ p(X) = 0\}\). If such method is not provided, the algorithm work over the polynomials as no relation between the variables exist.

This algorithm allows to solve homogeneous system (no inhomogeneous term is allowed).

INPUT:

  • parent: the ring where the coefficients will be treated. It has to be a polynomial ring or its field of fractions.

  • matrix: matrix for the homogeneous system.

  • method: method for the membership problem for the ideal \(I\).

  • relations: list of relation known for the variables of parent. This, together with method is use to check the membership to the ideal \(I\).

EXAMPLES:

sage: from ajpastor.misc.bareiss import *
sage: R.<a,b> = PolynomialRing(QQ)
sage: M = Matrix([[a^2, b], [-b, 1]])
sage: BA = BareissAlgorithm(R, M)
sage: BA.parent()
Multivariate Polynomial Ring in a, b over Rational Field
sage: BA.echelon_form()
[1 0]
[0 1]
sage: BA.rank()
2
sage: BA = BareissAlgorithm(FractionField(R), M)
sage: BA.parent()
Multivariate Polynomial Ring in a, b over Rational Field

We can see how the relations involved on the coefficients are important:

sage: BA = BareissAlgorithm(R, M, lambda p : p.reduce([a^2+b^2]) == 0)
sage: BA.echelon_form()
[-b 1]
[ 0 0]
sage: BA.rank()
1

This class also works with non-square matrices:

sage: M = Matrix([[a, b, a],[-b,a,-b]])
sage: BA = BareissAlgorithm(R, M)
sage: BA.echelon_form()
[-1  0 -1]
[ 0  1  0]
sage: BA.rank()
2
sage: BA.syzygy()
[ 1]
[ 0]
[-1]
sage: # Example where M is not square and we have relations
sage: BA = BareissAlgorithm(R, M, relations=[a^2+b^2])
sage: BA.echelon_form()
[a b a]
[0 0 0]
sage: BA.rank()
1
sage: BA.syzygy()
[-b -a]
[ a  0]
[ 0  a]
sage: # Another example with M being non-square
sage: M = Matrix([[a,b],[-b,a],[a^3 - 3*b^4, 1+b^5]])
sage: BA = BareissAlgorithm(R, M)
sage: BA.echelon_form()
[1 0]
[0 1]
[0 0]
sage: BA.rank()
2

This algorithm can find the relations on the fly thanks to the input is_zero:

sage: M = Matrix([[a, b, a],[-b,a,-b]])
sage: BA = BareissAlgorithm(R, M, lambda p : p.reduce([a^2+b^2]) == 0)
sage: BA.relations()
[0]
sage: H = BA.echelon_form()
sage: BA.relations()
[a^2 + b^2]

The algorithm can improve the given relations on the fly:

sage: M = Matrix([[a^2, b^2+2*a^2,b],[-b^2, a^2, 3],[1, 1, 1]])
sage: BA = BareissAlgorithm(R, M, relations=[a^4 + 2*a^2*b^2 + b^4])
sage: BA.echelon_form()
[1 0 0]
[0 1 0]
[0 0 1]
sage: BA.relations()
[a^4 + 2*a^2*b^2 + b^4]
sage: BA = BareissAlgorithm(R, M, lambda p : p.reduce([a^2+b^2]) == 0, [a^4 + 2*a^2*b^2 + b^4])
sage: BA.echelon_form()
[1 1 0]
[0 0 1]
[0 0 0]
sage: BA.relations()
[a^2 + b^2]