# polynomial kernel

Is there any particular reason why computing the kernel of a matrix over a polynomial ring might fail? Consider the following example:

sage: PR.<x1, x2, x3, y1, y2, y3> = QQ[]
sage: def P(t):
...       u = 1 - t
...       x = x1*t*t + x2*t*u*2 + x3*u*u
...       y = y1*t*t + y2*t*u*2 + y3*u*u
...       return (x, y)
sage: def C(t):
...       x, y = P(t)
...       return vector(PR, [x^2, y^2, x*y, x, y, 1])
sage: M = Matrix([C(t) for t in [0, 1, -1, 2, 1/2]])
sage: M.right_kernel()
Traceback (most recent call last):
...
ArithmeticError: Ideal Ideal (x3^2, x1^2) of Multivariate Polynomial Ring
in x1, x2, x3, y1, y2, y3 over Rational Field not principal


It seems to me as if being able to compute the reduced echelon form (using the frac method) should be enough to actually compute a kernel as well, without any additional requirements. I will demonstrate this in an answer below. So why does computing the echelon form work, but computing the kernel not?

edit retag close merge delete

Sort by » oldest newest most voted

As stated in the question, one can obtain a pretty good result using the echelon form:

def mykernel(m):
r = m.base_ring()
m = m.echelon_form('frac')
res = []
for i in range(m.nrows(), m.ncols()):
v = m.column(i)
vr = v.base_ring()
v = v.list()
v = v + [0]*(m.ncols()-m.nrows())
v[i] = -1
v = vector(vr, v)
v = v*v.denominator()
v = v.change_ring(r)
v = v/gcd(flatten([p.coefficients() for p in v]))
res.append(v)
return res


The code above makes assumptions which will be valid in most but not in all cases. In particular, it assumes that the left square of the echelon form will be an identity matrix. This is not neccessarily true, so a proper implementation would have to check for the actual location of the ones. But the general concept remains valid: every column in the echelon form which is not a unit vector corresponds to an element of the kernel.

In the example case, this yields the following vector:

[(-y1^2 + 4*y1*y2 - 4*y2^2 - 2*y1*y3 + 4*y2*y3 - y3^2, -x1^2 + 4*x1*x2 -
4*x2^2 - 2*x1*x3 + 4*x2*x3 - x3^2, 2*x1*y1 - 4*x2*y1 + 2*x3*y1 - 4*x1*y2
+ 8*x2*y2 - 4*x3*y2 + 2*x1*y3 - 4*x2*y3 + 2*x3*y3, 2*x3*y1^2 -
4*x2*y1*y2 - 4*x3*y1*y2 + 4*x1*y2^2 + 4*x3*y2^2 - 2*x1*y1*y3 +
8*x2*y1*y3 - 2*x3*y1*y3 - 4*x1*y2*y3 - 4*x2*y2*y3 + 2*x1*y3^2, 4*x2^2*y1
- 2*x1*x3*y1 - 4*x2*x3*y1 + 2*x3^2*y1 - 4*x1*x2*y2 + 8*x1*x3*y2 -
4*x2*x3*y2 + 2*x1^2*y3 - 4*x1*x2*y3 + 4*x2^2*y3 - 2*x1*x3*y3, -x3^2*y1^2
+ 4*x2*x3*y1*y2 - 4*x1*x3*y2^2 - 4*x2^2*y1*y3 + 2*x1*x3*y1*y3 +
4*x1*x2*y2*y3 - x1^2*y3^2)]

more