LLL and sympy

[T-Rex]

Greetings from a sympy newbie.  I have used other symboic algebra packages like pari-gp, and I'll teaching myself python and sympy now.  

According to   https://oscarbenjamin.github.io/blog/czi/post2.html    sympy now has LLL, but when I consulted the sympy documentation I could not find actual commands that produces an LLL-reduced basis.  I'm a newbie so mostly I overlooked something obvious; I'd be most appreciative for your help and references.  And if this is not the correct avenue for posting such questions, apologies and please let me know where I should post my questions.

I know that I can do this in e.g. sage, but for now I'd like to stick to standard python/sympy.

Thanks for your help and for making available this amazing package!

[Oscar Benjamin]

On Sun, 5 Nov 2023 at 15:04, T-Rex <call...@gmail.com> wrote:
>
> According to https://oscarbenjamin.github.io/blog/czi/post2.html sympy now has LLL, but when I consulted the sympy documentation I could not find actual commands that produces an LLL-reduced basis.

SymPy has a new matrix type called DomainMatrix for which the
LLL-algorithm is implemented. Ideally the LLL-algorithm would be
exposed through the ordinary Matrix type but it is not. The
DomainMatrix.lll method is documented here:
https://docs.sympy.org/latest/modules/polys/domainmatrix.html#sympy.polys.matrices.domainmatrix.DomainMatrix.lll

Starting with an ordinary Matrix M the command to get an LLL-reduced
basis in SymPy 1.12 is:

In [1]: from sympy.polys.matrices import DomainMatrix

In [2]: M = randMatrix(5) # make a random 5x5 Matrix

In [3]: M
Out[3]:
⎡74 3 71 3 27⎤
⎢ ⎥
⎢86 74 14 89 26⎥
⎢ ⎥
⎢58 2 14 54 17⎥
⎢ ⎥
⎢35 53 17 58 31⎥
⎢ ⎥
⎣72 83 33 61 85⎦

In [4]: dM = DomainMatrix.from_Matrix(M) # Convert to DomainMatrix

In [5]: dM.lll().to_Matrix() # Compute LLL-basis and convert back to Matrix
Out[5]:
⎡ 7 -19 17 23 22⎤
⎢ ⎥
⎢51 21 -3 31 -5⎥
⎢ ⎥
⎢-23 51 3 4 14⎥
⎢ ⎥
⎢53 -2 -4 -24 18⎥
⎢ ⎥
⎣-14 -24 -58 -4 13⎦

In SymPy 1.13 there will be a to_DM() method to make this a bit easier
so it is just:

M.to_DM().lll().to_Matrix()

Ideally we should add an .lll() method directly to Matrix to do this
without users needing to perform the conversion explicitly.

Also SymPy 1.13 will be able to use python-flint for this. If you have
python-flint installed (pip install python-flint) and then set the
environment variable SYMPY_GROUND_TYPES=flint then the
DomainMatrix.lll method will automatically use Flint which is a lot
faster for this operation.

--
Oscar

[T-Rex]

Thanks @oscar!  I didn't think I would get a reply from the blog's owner :-)

Re setting of the environment variable SYMPY_GROUND_TYPES=flint:   (1) Would that be in conflict with other parts of sympy/numpy/scipy?  (2) After it is set to flint, is there a way to turn it back to the default?  I can see scenarios where I need to bench mark various options, and it would be really nice to be able to change the flag as I go along in my jupyter notebook, as oppsed to restarting the program from scratch.

Taboo question:  When might Sympy 1.13 be coming out?  I am know that such questions are taboo, but my project has various stages/steps, and depending on the availability of 1.1.3 I might want to work on the LLL-parts first while waiting for 1.1.3 to become available.

MANY THANKS for your help and advice.

[Aaron Meurer]

On Sun, Nov 5, 2023 at 8:22 PM T-Rex <call...@gmail.com> wrote:
>
> Thanks @oscar! I didn't think I would get a reply from the blog's owner :-)
>
> Re setting of the environment variable SYMPY_GROUND_TYPES=flint: (1) Would that be in conflict with other parts of sympy/numpy/scipy? (2) After it is set to flint, is there a way to turn it back to the default? I can see scenarios where I need to bench mark various options, and it would be really nice to be able to change the flag as I go along in my jupyter notebook, as oppsed to restarting the program from scratch.

NumPy, SciPy, and other libraries do not use flint, only SymPy. Other
parts of SymPy do use flint, but it's all completely transparent. The
only difference is the output.

Right now, there's no way to change the option after SymPy is
imported, which is why it is set as an environment variable. So if you
want to test different things, you'll have to restart the notebook.

>
> Taboo question: When might Sympy 1.13 be coming out? I am know that such questions are taboo, but my project has various stages/steps, and depending on the availability of 1.1.3 I might want to work on the LLL-parts first while waiting for 1.1.3 to become available.

There's no reason to not ask this question. And
anyways, Oscar mentioned to me that he wants to make a release soon so
(hopefully) it should be out before the end of the year, if not month.


Aaron Meurer

> --
> 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/54713a7c-886f-46dd-965e-0e288e61d901n%40googlegroups.com.

[Oscar Benjamin]

On Mon, 6 Nov 2023 at 03:22, T-Rex <call...@gmail.com> wrote:
>
> Thanks @oscar! I didn't think I would get a reply from the blog's owner :-)
>
> Re setting of the environment variable SYMPY_GROUND_TYPES=flint: (1) Would that be in conflict with other parts of sympy/numpy/scipy? (2) After it is set to flint, is there a way to turn it back to the default? I can see scenarios where I need to bench mark various options, and it would be really nice to be able to change the flag as I go along in my jupyter notebook, as oppsed to restarting the program from scratch.

It is not possible to do this. It would not be easy to make it work
because using flint for the ground types changes many things starting
with the representation of the lowest level domain types like ZZ, QQ
and GF(p).

There are ways to test the speed of the different representations
though by using explicit conversions between the low-level SDM, DDM
and DFM matrix types:

In [23]: M = randMatrix(50) # 50x50 matrix of 2-digit integers

In [24]: dM = M.to_DM() # DomainMatrix

In [25]: dMd = dM.rep.to_ddm() # internal dense domain matrix rep

In [26]: dMf = dM.rep.to_dfm() # internal dense flint matrix rep

In [27]: %time ok = dMd.lll() # time for SymPy implementation
CPU times: user 1.01 s, sys: 6.95 ms, total: 1.02 s
Wall time: 1.13 s

In [28]: %time ok = dMf.lll() # time when using flint
CPU times: user 18.4 ms, sys: 3.93 ms, total: 22.3 ms
Wall time: 184 ms

The two are not exactly equal because the Flint algorithm is different:

In [34]: dMd.lll() == dMf.lll().to_ddm()
Out[34]: False

There is discussion about the differences in the original pull request
that added LLL:

https://github.com/sympy/sympy/pull/24572

The SymPy implementation uses "classic" LLL with all exact arithmetic
and only has the parameter delta. The Flint version I assume uses
partially approximate arithmetic and has two parameters delta and eta.
I am not sure how exactly to make the output of the Flint version
precisely equal to what the SymPy implementation would return.

Note that you can also use python-flint directly to test what it returns:

In [35]: import flint

In [36]: M = flint.fmpz_mat([[1, 2], [3, 4]])

In [38]: M.lll()
Out[38]:
[1, 0]
[0, 2]

You can see the python-flint code for this method here although it is
calling into Flint for the actual algorithms:

https://github.com/flintlib/python-flint/blob/33e5485bda16339eb334c816639e22e57c4658e2/src/flint/types/fmpz_mat.pyx#L613-L654

The SymPy implementation of the algorithm is here:
https://github.com/sympy/sympy/blob/master/sympy/polys/matrices/lll.py

> Taboo question: When might Sympy 1.13 be coming out?

Soon, but I am not sure exactly when. There are some issues to be
resolved such as whether or not to use python-flint automatically
rather than requiring the SYMPY_GROUND_TYPES=flint to be set
explicitly.

> I am know that such questions are taboo, but my project has various stages/steps, and depending on the availability of 1.1.3 I might want to work on the LLL-parts first while waiting for 1.1.3 to become available.

It is very easy to install the current development version of SymPy
(future 1.13) right now if you want to test it out:

pip install git+https://github.com/sympy/sympy.git

I don't anticipate any changes to the LLL code between now and the
1.13 release except that I might try to add the lll method to Matrix
so that it is not necessary to convert to DomainMatrix explicitly.

Anyone else is also welcome to send a PR to that effect. It would also
be good to add the interface for Hermite and Smith normal forms. Those
are implemented for DomainMatrix but not provided as methods and the
code is not there to call through to flint when available.

--
Oscar

[Oscar Benjamin]

On Mon, 6 Nov 2023 at 11:46, Oscar Benjamin <oscar.j....@gmail.com> wrote:
>
> I don't anticipate any changes to the LLL code between now and the
> 1.13 release except that I might try to add the lll method to Matrix
> so that it is not necessary to convert to DomainMatrix explicitly.

I just did this:
https://github.com/sympy/sympy/pull/25880

As of SymPy 1.13 there will be Matrix methods lll and lll_transform:

In [5]: r = sqrt(2).evalf(10)

In [6]: r
Out[6]: 1.414213562

In [7]: M = Matrix.hstack(eye(3), Matrix([round(10**10*r**2),
round(10**10*r), 10**10]))

In [8]: M
Out[8]:
⎡1 0 0 20000000000⎤
⎢ ⎥
⎢0 1 0 14142135620⎥
⎢ ⎥
⎣0 0 1 10000000000⎦

In [9]: M.lll()
Out[9]:
⎡ -1 0 2 0 ⎤
⎢ ⎥
⎢18928 -33461 9465 19180 ⎥
⎢ ⎥
⎣-48945 86523 -24472 249260⎦

Top row gives coefficients for the minpoly of sqrt(2):

In [11]: print(minpoly(sqrt(2)))
_x**2 - 2

In [12]: print(Poly(M.lll()[0,:-1], x))
Poly(-x**2 + 2, x, domain='ZZ')

The lll_transform method also gives the transformation matrix that
sends the original basis to the reduced basis:

In [25]: B, T = M.lll_transform()

In [26]: T*M == B
Out[26]: True

--
Oscar