Raised as a by-product of #31
In approximant basis computations which work with increasing approximation orders (currently: mbasis incrementally with order += 1, or pmbasis with larger steps), once the processed approximation order is large enough so that a bunch of rows of a kernel basis have already been found in the approximant basis, these rows are not affected by subsequent computations, i.e. these only affect the remaining non-kernel rows. This also means that the degrees in the basis will grow unevenly: the overall sum of row degrees continues to grow at the same rate but this growth is only supported by the non-kernel rows.
Here is an example with extreme parameter (single column) which makes this quite visible. An m x 1 vector F of degree d typically has a kernel basis of degree close to d/m, which will appear in an approximant basis of F at order about d + d/m:
sage: xring.<x> = GF(9001)[]
sage: F = matrix.random(xring, 6, 1, degree=999)
sage: P = F.minimal_approximant_basis(1200)
sage: P[1:6,:].is_minimal_kernel_basis(F)
True
sage: P.degree_matrix()
[201 200 200 200 200 200]
[200 200 199 199 199 199]
[200 200 200 199 199 199]
[200 200 200 200 199 199]
[200 200 200 200 200 199]
[199 199 199 199 199 199]
This means that any approximant basis at larger order, say 2d for example, will simply be the same approximant basis with its single non-kernel row multiplied by a large power of the variable.
sage: P = F.minimal_approximant_basis(2000)
sage: P.degree_matrix()
[1001 1000 1000 1000 1000 1000]
[ 200 200 199 199 199 199]
[ 200 200 200 199 199 199]
[ 200 200 200 200 199 199]
[ 200 200 200 200 200 199]
[ 199 199 199 199 199 199]
Obtaining the latter, with just multiplication by x, should be particularly cheap. However, in the current code, that is far from being the case. The first reason is that there is no detection that we are only working on a single row (mbasis continues to compute nullspaces of constant m x 1 vectors and to update the bases via multiplication, which ignores the fact that these nullspaces are very close to an identity matrix). A second reason is that the overall degree of the basis increases quite fast despite its average degree remaining small, making this situation not well suited to an nmod_mat_poly_t representation (list of m x m matrix coefficient for each degree) as used in mbasis, and also not well suited to generic implementations of evaluation-interpolation multiplications which will only depend on this overall degree, despite the very special product that should be done here.
See #31 for an example where doubling the approximation order multiplies the running time by 36 (instead of a factor between 2 and 4).
This issue typically does not arise on input whose degree is close to the approximation order (but might still happen if it has a low degree kernel basis, which is not expected generically). However, as soon as the order is somewhat larger than this input matrix degree, there may be substiantal gains to be obtained, by trying to detect that a kernel basis (or a good part of it) has been found: then it would remain to solve an approximation problem of smaller matrix dimension (focusing on non-kernel rows) and eventually using a basis multiplication to combine its result into a global approximant basis. This could be based on comparing the input order with the expected order that provides a kernel basis, but we should make sure that adding this will not slow down some frequent cases.
Raised as a by-product of #31
In approximant basis computations which work with increasing approximation orders (currently: mbasis incrementally with
order += 1, or pmbasis with larger steps), once the processed approximation order is large enough so that a bunch of rows of a kernel basis have already been found in the approximant basis, these rows are not affected by subsequent computations, i.e. these only affect the remaining non-kernel rows. This also means that the degrees in the basis will grow unevenly: the overall sum of row degrees continues to grow at the same rate but this growth is only supported by the non-kernel rows.Here is an example with extreme parameter (single column) which makes this quite visible. An
m x 1vectorFof degreedtypically has a kernel basis of degree close tod/m, which will appear in an approximant basis ofFat order aboutd + d/m:This means that any approximant basis at larger order, say
2dfor example, will simply be the same approximant basis with its single non-kernel row multiplied by a large power of the variable.Obtaining the latter, with just multiplication by
x, should be particularly cheap. However, in the current code, that is far from being the case. The first reason is that there is no detection that we are only working on a single row (mbasis continues to compute nullspaces of constantm x 1vectors and to update the bases via multiplication, which ignores the fact that these nullspaces are very close to an identity matrix). A second reason is that the overall degree of the basis increases quite fast despite its average degree remaining small, making this situation not well suited to annmod_mat_poly_trepresentation (list ofm x mmatrix coefficient for each degree) as used in mbasis, and also not well suited to generic implementations of evaluation-interpolation multiplications which will only depend on this overall degree, despite the very special product that should be done here.See #31 for an example where doubling the approximation order multiplies the running time by 36 (instead of a factor between 2 and 4).
This issue typically does not arise on input whose degree is close to the approximation order (but might still happen if it has a low degree kernel basis, which is not expected generically). However, as soon as the order is somewhat larger than this input matrix degree, there may be substiantal gains to be obtained, by trying to detect that a kernel basis (or a good part of it) has been found: then it would remain to solve an approximation problem of smaller matrix dimension (focusing on non-kernel rows) and eventually using a basis multiplication to combine its result into a global approximant basis. This could be based on comparing the input order with the expected order that provides a kernel basis, but we should make sure that adding this will not slow down some frequent cases.