LLL assertion errors

[robert....@gmail.com]

Hi everyone,

I'm trying to use the new-ish LLL implementation for domain matrices, but it's failing assertions in some of my use cases. Specifically, when the input matrix has really big integers in it.

Here's an example from a recent session:

In [2]: r = sp.sqrt(2).evalf(100)

In [3]: v = [[1, 0, 0, sp.floor(10**100 * r**2)]]

In [4]: v += [[0, 1, 0, sp.floor(10**100 * r**1)]]

In [5]: v += [[0, 0, 1, sp.floor(10**100)]]

In [6]: m = DomainMatrix.from_list(v, ZZ)
...

In [8]: m.lll()

...

assert all([mu_small(i, j) for i in range(m) for j in range(i)])
AssertionError:

Here, m has integers with ~100 digits, and for some reason the implementation does not actually reduce the Gram-Schmidt coefficients mu to be less than 1/2 in absolute value. It does it perfectly if the entries of m have ~5 or ~10 digits.

Am I doing something stupid here, or does the implementation need fixing?

Robert

[Oscar Benjamin]

There is some problem with the implementation in SymPy. If you have
python-flint installed it works because it uses the FLINT
implementation.

>>> print(Matrix(v).to_DM().lll().to_Matrix())
Matrix([[-1, 0, 2, 0],
[27838432316546329948653649573487634649724968645328,
-49211860671581597598021395402360695743160591979141,
13919216158273164974326824786743817324862484322663,
-54612160756961593886340461194352260530027628350507],
[39369488537265278078417116321888556594528473583313,
-69596080791365824871634123933719086624312421613319,
19684744268632639039208558160944278297264236791656,
125969787835055438793681187979620981611206915732087]])

The assertion that fails is here:
https://github.com/sympy/sympy/blob/e65784eb972b9e3a7f62520a8069f5dbb77ca96b/sympy/polys/matrices/lll.py#L85

> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/695ef34e-ea30-495e-af4d-ffe29c5a71a0n%40googlegroups.com.

[Robert Dougherty-Bliss]

Wow, that's great news! I missed the incorporation of flint.

1. Is there up-to-date documentation about flint being used in sympy
anywhere? I see your blog post
https://oscarbenjamin.github.io/blog/czi/post2.html
but it would be nice to have something in the official docs.

2. Is it worth fixing the sympy implementation of LLL? I would be
happy using flint (even making it a required dependency), but I could
investigate fixing the sympy implementation for others.

Robert

> You received this message because you are subscribed to a topic in the Google Groups "sympy" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/topic/sympy/rkoXWBw-oY4/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to sympy+un...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAHVvXxSeUGPkQ2FL5rLA3XoFW78fzD1it3FoJ3c%2B4nOHTAVHQw%40mail.gmail.com.

[Oscar Benjamin]

On Thu, 17 Oct 2024 at 15:59, Robert Dougherty-Bliss
<robert....@gmail.com> wrote:
>
> Wow, that's great news! I missed the incorporation of flint.
>
> 1. Is there up-to-date documentation about flint being used in sympy
> anywhere? I see your blog post
> https://oscarbenjamin.github.io/blog/czi/post2.html
> but it would be nice to have something in the official docs.

There probably should be more in the documentation. It is still a work
in progress though. As of SymPy 1.13 python-flint is used directly for
the ZZ, QQ and GF(p) domain elements and also for the internal
implementation of some things:

- Poly for univariate polynomials over ZZ or QQ
- Algebraic field domains like QQ.algebraic_field(sqrt(2))
- Dense DomainMatrix over ZZ and QQ

Examples:

>>> type(ZZ(2))
flint.types.fmpz.fmpz
>>> type(QQ(2))
flint.types.fmpq.fmpq
>>> type(GF(3)(2))
flint.types.nmod.nmod
>>> type(Poly(x + 1).rep._rep)
flint.types.fmpz_poly.fmpz_poly
>>> type(Matrix([1]).to_DM().to_dense().rep.rep)
flint.types.fmpz_mat.fmpz_mat
>>> type(QQ.algebraic_field(sqrt(2)).one._rep._rep)
flint.types.fmpq_poly.fmpq_poly

Some other things have been added in master since the 1.13 release:

- Univariate polynomials over GF(p)
- Dense DomainMatrix over GF(p)

A lot of work has gone into python-flint itself as well and so the
next release of python-flint (0.7.0) will have multivariate
polynomials over ZZ, QQ and GF(p) and also finite fields like
GF(p**n).

I expect that SymPy 1.14 will support using python-flint's
multivariate polynomials for ZZ, QQ and GF(p) for both Poly and
PolyElement (the implementation of domains like QQ[x,y]). This will be
a more noticeable change for many people I think since functions like
factor, cancel, resultant, groebner etc would get much faster and this
will also make a big difference to DomainMatrix for symbolic matrices.

After that I would say that Arb would be next on the priority list for
features in python-flint that SymPy should make use of.

> 2. Is it worth fixing the sympy implementation of LLL? I would be
> happy using flint (even making it a required dependency), but I could
> investigate fixing the sympy implementation for others.

Yes, I think it is worth fixing. I don't think that SymPy would ever
make python-flint a hard dependency so we still want everything to
work without it.

I can imagine that at some point SymPy will have some features that
are only available when python-flint is installed much like other
optional dependencies (NumPy, matplotlib etc). For core functionality
though we want everything to work exactly the same with or without
python-flint.

In the particular case of LLL there are differences between the SymPy
implementation and the FLINT one. There is apparently the bug that you
have identified here but also the python-flint lll method takes two
parameters delta and eta:

def lll(transform=False, delta=0.99, eta=0.51, rep='zbasis', gram='approx')

The SymPy implementation only takes the delta parameter with a default
value of 3/4. I guess that the FLINT algorithm is a generalisation
with its eta parameter. Ideally these would both do the same thing so
that you get the same results with or without python-flint (this is
the only case I know of where SymPy currently makes use of something
that is not fully equivalent).

The rep parameter to python-flint's lll method can be set to either
zbasis or gram but as far as I can tell rep='gram' just borks the
whole computer consuming all memory so best to stick to zbasis (as
SymPy does when using python-flint for LLL).

Also the implementation in SymPy uses exact rational arithmetic for
the Gram-Schmidt basis which is (asymptotically) slower and is not the
default in python-flint (gram='approx'). I'm not sure if double
precision is sufficient for these things though.

--
Oscar