reduced_basis for number field multiples wrong

[haraldschilly]

This was reported by the "report a problem" link:

The reduced_basis for number fields applies LLL to basis of the maximal order in order to get a reduced basis. The problem is that after applying LLL, it multiples the transformation matrix by the vector [1,a,a²,...], where 'a' is generator of the field - instead of using the vector it actually used for the maximal order.

How to reproduce:

x = QQ['x'].0
k = NumberField(x^6 + 2218926655879913714112*x^4 - \
32507675650290949030789018433536*x^3 + \
4923635504174417014460581055002374467948544*x^2 - \
36066074010564497464129951249279114076897746988630016*x + \
264187244046129768986806800244258952598300346857154900812365824,'a')
new_basis = k.reduced_basis()
print new_basis[0].minpoly()

This prints a crazy big polynomial. Looking a bit at what is returned it is clear that reduced_basis is multiplying the transformation matrix by the wrong vector.

To solve you just need to multiply by the vector of generators of the maximal order:

mp = pari( k.Minkowski_embedding(k.maximal_order().gens(), prec=120) )
ml = sage.libs.pari.gen.gen.qflll(mp)
T = matrix([[ZZ(y) for y in list(x)] for x in list(ml)]).transpose()
r = Matrix(O.gens())* T
for rr in r[0]:
print rr.minpoly()

This works fine.

Component: number fields

Stopgaps: todo

Author: John Cremona

Branch/Commit: 8f1f51d

Reviewer: Jeroen Demeyer

Issue created by migration from https://trac.sagemath.org/ticket/10017

[sagetrac-fwclarke]

comment:1

Looks like the definition of O was omitted. The following slightly simplified code works ok for me:

sage: O = k.ring_of_integers()
sage: mp = pari( k.Minkowski_embedding(O.gens(), prec=120) )
sage: ml = sage.libs.pari.gen.gen.qflll(mp)
sage: T = ml.sage()
sage: r = vector(O.gens())*T
sage: for rr in r:
....:     rr.minpoly()

[sagetrac-jakobkroeker]

Stopgaps: todo

[JohnCremona]

comment:7

While the original diagnosis is not wrong, there is a more basic problem as this smaller example shows:

sage: K.<a> = NumberField(x^3-10)
sage: ZK = K.ring_of_integers()
sage: ZK.basis()
[1/3*a^2 + 1/3*a + 1/3, a, a^2]
sage: ZK([1,2,3])
3*a^2 + 2*a + 1

Although the ring of integers is an order which knows its own Z-basis, when you construct an element of the order from an integer vector it ignores that and uses the power_basis of the field and not its integral basis.

I think the problem is in these lines of _element_constructor() (order.py lines 1142-1143):

       if not is_Element(x) or x.parent() is not self._K:
            x = self._K(x)

There appears to be a related problem in the constructor of OrderElementAbsolute in number_field_element:

       K = order.number_field()
        NumberFieldElement_absolute.__init__(self, K, f)
        self._number_field = K

where the parameter f is not documented and there are no relevant tests, but this seems to be saying "take the data in f and construct a field element from it".

Perhaps the simplest solution is to go to those lines in order.py and insert code to handle the case where the input data is an integer vector of the right length.

[JohnCremona]

comment:8

After doing the above "simplest solution", which did fix the x^3-10 example above, I saw that it's not the fix needed for the original which is to take the vectors output by the LLL reduction and push them into the ring of integers rather than the field. Now I get

sage: [c.minpoly() for c in new_basis]

[x - 1,
 x^6 + 3*x^5 + 258*x^4 + 511*x^3 + 3564*x^2 + 3309*x + 2347,
 x^6 - 24*x^5 + 6126*x^4 - 312664*x^3 + 5407566*x^2 - 33643572*x + 95921443,
 x^6 + 18*x^5 + 3366*x^4 + 82008*x^3 + 886962*x^2 - 840726*x + 5521647,
 x^6 + 27*x^5 + 9579*x^4 + 623358*x^3 + 5060091*x^2 - 139224285*x + 880944177,
 x^6 - 72*x^5 + 65286*x^4 - 10762768*x^3 + 473072922*x^2 - 2502686322*x + 54227921641]

which is a lot better.

Patch up when I have added doctests...

[JohnCremona]

comment:9

The doctest based on the original post works fine for me when run manually, but when run as a doctest something weird happens -- it takes ages and then gives an enormous python crash. I have no idea what this is but hope that a human can help before this breaks all the patchbots out there...

---

New commits:

6bee3fc

#10017: fix construction of order elements from list, and number field reduced basis

[JohnCremona]

Author: John Cremona

[JohnCremona]

Branch: u/cremona/10017

[JohnCremona]

Commit: 6bee3fc

[sagetrac-git]

Changed commit from 6bee3fc to 5d51b32

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

0000be4	Merge branch 'develop' into 10017
5d51b32	#10017 added 2 lines to new doctest

[jdemeyer]

comment:11

Quick comments (more might come):

1. Replace See (:trac:`10017`):: by Check that :trac:`10017` is fixed:: or something along those lines.

2. Replace

sage: assert R.discriminant()==k.discriminant()

by

sage: R.discriminant() == k.discriminant()
True

[sagetrac-git]

Changed commit from 5d51b32 to 8b3eb02

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

8b3eb02

#10017 minor changes to docstring

[jdemeyer]

comment:13

In this block of code:

            M = self.Minkowski_embedding(self.integral_basis(), prec=prec)
            T = pari(M).qflll().sage()
            ZK = self.ring_of_integers()
            self.__reduced_basis = [ ZK(v.list()) for v in T.columns() ]

you are assuming that self.integral_basis() corresponds to the basis chosen in ZK(v.list()). After all, that's the bug on this ticket, right? To me, this looks way too implicit. I would prefer to see the code written in such a way that you are actually using self.integral_basis() in this computation, instead of the ZK() call. I.e. something like

ZK = self.integral_basis()
reduced_basis = [sum(x * g for x, g in zip(v.list(), ZK)) for v in T.columns()]

[JohnCremona]

comment:14

Yes, you are right, and perhaps being explicit like this is safest. I still think (do you agree?) that being able to do ZK(list) is useful and must give the linear combination of ZK's basis with coefficients from the list.

A related different point is that after computing the reduced basis, it is stored as self.__reduced_basis, but not actually used anywhere. Perhaps we should implement the ability to change the integral basis itself; but that might be too dangerous since various functions using pari directly will assume that the integral basis does not change.

[videlec]

comment:15

(fixing exponentiation in ticket description)

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

39fc227

#10017: small stylistic improvement

[jdemeyer]

comment:23

The patchbot reports a failure on 32-bit:

sage -t --long src/sage/rings/number_field/number_field.py
**********************************************************************
File "src/sage/rings/number_field/number_field.py", line 5599, in sage.rings.number_field.number_field.NumberField_generic.reduced_basis
Failed example:
    [c.minpoly() for c in new_basis]
Expected:
    [x - 1,
     x^6 + 3*x^5 + 258*x^4 + 511*x^3 + 3564*x^2 + 3309*x + 2347,
     x^6 - 24*x^5 + 6126*x^4 - 312664*x^3 + 5407566*x^2 - 33643572*x + 95921443,
     x^6 + 18*x^5 + 3366*x^4 + 82008*x^3 + 886962*x^2 - 840726*x + 5521647,
     x^6 + 27*x^5 + 9579*x^4 + 623358*x^3 + 5060091*x^2 - 139224285*x + 880944177,
     x^6 - 72*x^5 + 65286*x^4 - 10762768*x^3 + 473072922*x^2 - 2502686322*x + 54227921641]
Got:
    [x - 1,
     x^6 - 45*x^5 + 6993*x^4 - 519048*x^3 + 8139204*x^2 + 64936080*x + 394386192,
     x^6 + 3*x^5 + 258*x^4 + 511*x^3 + 3564*x^2 + 3309*x + 2347,
     x^6 - 12*x^5 + 459*x^4 - 13068*x^3 + 171423*x^2 - 1008720*x + 2218941,
     x^6 - 72*x^5 + 21618*x^4 - 2725154*x^3 + 96844941*x^2 - 367602102*x + 5299367596,
     x^6 - 36*x^5 + 4104*x^4 - 225746*x^3 + 4831632*x^2 - 39603546*x + 130361299]
**********************************************************************

[jdemeyer]

Reviewer: Jeroen Demeyer

[JohnCremona]

comment:24

Would it be OK to tag the current output as 64-bit and this alternative as 32-bit? Of course this "reduced basis" is in no way canonical.

[jdemeyer]

comment:25

let me check if the prec argument makes a difference.

[jdemeyer]

comment:26

On 64-bit, the output stabilizes for prec >= 117; the result it

[x - 1,
 x^2 - x + 1,
 x^6 + 3*x^5 - 102*x^4 - 103*x^3 + 10572*x^2 - 59919*x + 127657,
 x^6 - 3*x^5 - 102*x^4 + 315*x^3 + 10254*x^2 - 80955*x + 198147,
 x^3 - 171*x + 848,
 x^6 + 171*x^4 + 1696*x^3 + 29241*x^2 + 145008*x + 719104]

[JohnCremona]

comment:27

That's a lot better -- can you check if it's the same on 32-bit using the same prec (say 120)?

[jdemeyer]

comment:28

Yes, on 32-bit I get the same result for prec >= 68.

So the best solution is increasing the prec in the doctest. While you're at it, could you please replace k = NumberField(...,'a') by k.<a> = NumberField(...)?

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

22746f8

#10017: adjust doctest for 32-bit

[sagetrac-git]

Changed commit from 39fc227 to 22746f8

[JohnCremona]

comment:30

Apologies for the delay -- I pushed to u/cremona/20017 by mistake 4 days ago. It's now where it should be.

[jdemeyer]

comment:31

OK, looks good but let's wait for the patchbot.

[jdemeyer]

comment:32

Sorry to say this, but the branch doesn't merge cleanly...

[JohnCremona]

comment:33

Replying to @jdemeyer:

Sorry to say this, but the branch doesn't merge cleanly...

Not your fault! I am just building the latest beta and will merge it with my branch shortly.

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

f56a7e8	Merge branch 'develop' into 10017
8f1f51d	#10017: fix deprecation warning for Minkowski_embedding

[sagetrac-git]

Changed commit from 22746f8 to 8f1f51d

[JohnCremona]

comment:35

OK, I merged, fixed the small conflict and a deprecation warning. Let's see if the patchbots like this one better.

[vbraun]

Changed branch from u/cremona/10017 to 8f1f51d