Opened 9 months ago
Closed 8 months ago
#29866 closed defect (fixed)
closest_vector for IntegerLattice is broken
Reported by: | dimpase | Owned by: | |
---|---|---|---|
Priority: | blocker | Milestone: | sage-9.2 |
Component: | number theory | Keywords: | round |
Cc: | malb, slelievre, chapoton | Merged in: | |
Authors: | Samuel Lelièvre, Dima Pasechnik | Reviewers: | Matthias Koeppe, Samuel Lelièvre |
Report Upstream: | N/A | Work issues: | |
Branch: | 99f0df7 (Commits, GitHub, GitLab) | Commit: | 99f0df74fa4b4b1c836b8b850ceeb516ef012d22 |
Dependencies: | Stopgaps: |
Description (last modified by )
In this example we compute v
, the closest vector to the vector u
that lies in the integer lattice L
. Problem: v
is not in L
.
sage: from sage.modules.free_module_integer import IntegerLattice sage: M = matrix(ZZ, [[20957228, -4966110], [9411844, 19625639]]) sage: L = IntegerLattice(M) sage: u = vector([-423434678248195, -18882583298608161305227077482]) sage: v = L.closest_vector(u) sage: v in L False sage: MGAP = libgap(M) sage: vGAP = libgap(v) sage: libgap.SolutionIntMat(MGAP, vGAP) # independent test that v is not in L fail
This is a simplified version of Taylor Huang's example posted on sage-devel.
This used to work in Sage 8.4, and was broken in 8.5.beta3 by the
change from round(x, 6)
to x.n(digits=6)
in #26648.
The sage-devel thread also reports the membership test to be broken:
from sage.modules.free_module_integer import IntegerLattice coef = Matrix([-44429982080874270968379672793605458, 98931650854481334735580708522902113]) bMat = Matrix([[20957228, -4966110], [9411844, 19625639]]) L = IntegerLattice(bMat) coef*bMat in L # should be true, but returns False bMatGAP = libgap(bMat) v = libgap(coef*bMat)[0] sol = libgap.SolutionIntMat(bMatGAP, v); sol # prints # [ -423434671769860, -18882583298608161305225815339 ] vector(sol.sage())*bMat == vector(coef*bMat) # sanity check - prints True
but this is merely a type "error", in the sense that coef
is not a vector,
but a 1x2 matrix (a different type, and automatic conversion does not happen here).
Not sure whether the latter warrants a fix.
Change History (39)
comment:1 Changed 9 months ago by
- Description modified (diff)
comment:2 Changed 9 months ago by
- Cc slelievre added
- Description modified (diff)
comment:3 Changed 9 months ago by
comment:4 Changed 9 months ago by
Hi, I'm the person who came up with these examples. I can be sure that it is not only about membership testing. I've checked my examples with Smith normal form and indeed some of the closest vectors returned is not in L. So both membership and closest_vector are broken.
comment:5 Changed 9 months ago by
we confirm that closest vector is broken, nothing to do with membership test bug
sage: MGAP=libgap(M) sage: vGAP=libgap(v) sage: libgap.SolutionIntMat(MGAP,vGAP) fail
comment:6 Changed 9 months ago by
- Description modified (diff)
- Priority changed from critical to blocker
- Summary changed from closest_vector for IntegerLattice is broken to closest_vector for IntegerLattice and membership test are broken
comment:7 Changed 9 months ago by
- Description modified (diff)
comment:8 Changed 9 months ago by
- Description modified (diff)
comment:9 Changed 9 months ago by
- Description modified (diff)
comment:10 Changed 9 months ago by
closest vector got broken somewhere between Sage 8.4 (works) and 8.5 (does not). Membership test already broken in 8.2.
comment:11 Changed 9 months ago by
the heavy lifting is done by fplll/fpylll.
Python 3.7.7 (default, May 30 2020, 01:27:43) [GCC 8.3.1 20190518] on linux Type "help", "copyright", "credits" or "license" for more information. >>> from fpylll import CVP, IntegerMatrix >>> A = IntegerMatrix.from_matrix([[20957228, -4966110], [9411844, 19625639]]) >>> t=(-423434678248195, -18882583298608161305227077482) >>> v0 = CVP.closest_vector(A, t) >>> v0 (-423434681577544, -18882583298608161305224294058)
this is the same vector as computed by Sage code in the ticket description, oops. So the bug may be narrowed to fp(y)lll.
Note a comment on https://github.com/fplll/fpylll/issues/124 which says that the basis should be LLL-reduced for this to work correctly.
comment:12 Changed 9 months ago by
Thus, in view of comment:11 (corrected) the closest_vector bug is a Sage interface to fp(y)lll bug.
rather, Sage has a different implementation --- which seems to be unchanged since 2015 (i.e. much before Sage 8.5...)
comment:13 Changed 9 months ago by
- Description modified (diff)
comment:14 Changed 9 months ago by
- Description modified (diff)
- Summary changed from closest_vector for IntegerLattice and membership test are broken to closest_vector for IntegerLattice is broken
The containment is "broken", as the 1x2 matrix coef*bMat
fails to be coersed to a vector.
Namely, with the example in the
sage: vector(coef)*bMat in L True sage: vector(coef*bMat) in L True sage: coef*bMat in L False sage: M=matrix([1,1]) sage: M*bMat in L False
so this is not really a bug, it's just a conversion issue, absence of automatic conversion.
comment:15 Changed 9 months ago by
This seems suspicious.
sage: from sage.modules.free_module_integer import IntegerLattice sage: a = vector(ZZ, [20957228, -4966110]) sage: b = vector(ZZ, [9411844, 19625639]) sage: L = IntegerLattice([a, b]) sage: R = L.voronoi_relevant_vectors() sage: S = [-b, -a, b - a, b, a, a - b] sage: print([r - s for r, s in zip(R, S)]) [(0, 0), (0, 0), (0, -1), (0, 0), (0, 0), (0, 1)] sage: print(R[-1]); print(a - b) (11545384, -24591748) (11545384, -24591749)
comment:16 Changed 9 months ago by
note that on Sage 8.4 (where closest_vector() still worked) this is
... sage: print([r - s for r, s in zip(R, S)]) [(-2133540, 44217388), (11545384, -24591749), (-9411844, -19625639), (2133540, -44217388), (0, 0), (-2133540, 44217388)] sage: print(R[-1]); print(a - b) (9411844, 19625639) (11545384, -24591749)
comment:17 Changed 9 months ago by
If fact, R computed in comment:15 is the same as the one in comment:16 (with Sage 8.4), but in a different order.
comment:18 Changed 9 months ago by
The voronoi cell is slightly off. It should be exactly P
below.
sage: V = L.voronoi_cell() sage: ieqs = [(x^2 + y^2, 2*x, 2*y) for x, y in (-b, -a, a-b, b, a, b-a)] sage: P = Polyhedron(ieqs=ieqs) sage: VP = V.plot(fill='cyan', alpha=0.5) + P.plot(fill='magenta', alpha=0.5) sage: for v in V.vertices(): ....: x, y = v.vector() ....: VP.show(xmin=x-1, xmax=x+1, ymin=y-1, ymax=y+1, figsize=3)
comment:19 Changed 9 months ago by
indeed, on Sage 9.2.beta1
sage: V.vertices() (A vertex at (-483879947381762629327/57254902852288, 315998472693125499957/28627451426144), A vertex at (483879947381762629327/57254902852288, -315998472693125499957/28627451426144), A vertex at (5728193175888805065973/458039243775532, 1390651404779390505975/229019621887766), A vertex at (-708599510920502622251/229019617181844, 1552002451289776620941/114509808590922), A vertex at (708599510920502622251/229019617181844, -1552002451289776620941/114509808590922), A vertex at (-5728193175888805065973/458039243775532, -1390651404779390505975/229019621887766))
but on Sage 8.4
sage: V.vertices() (A vertex at (5728193175888805065973/458039243775532, 1390651404779390505975/229019621887766), A vertex at (3871039688862599879323/458039243775532, -2527988039232444116235/229019621887766), A vertex at (1417199267595526864965/458039243775532, -3104005018306403526499/229019621887766), A vertex at (-5728193175888805065973/458039243775532, -1390651404779390505975/229019621887766), A vertex at (-3871039688862599879323/458039243775532, 2527988039232444116235/229019621887766), A vertex at (-1417199267595526864965/458039243775532, 3104005018306403526499/229019621887766))
quite a bit different
comment:20 Changed 9 months ago by
I'll try bisecting Sage 8.4 and 8.5. Might take a while...
comment:21 Changed 9 months ago by
Could it be #26648 with this change to the function diamond_cut
in src/sage/modules/diamond_cutting.py
?
- hv = [QQ(round(elmt, 6)) for elmt in hv] + hv = [QQ(elmt.n(digits=6)) for elmt in hv]
comment:22 Changed 9 months ago by
good catch - this change is certainly bad. cf.
sage: r=5555.55555555555555555 sage: round(r,6) 5555.555556 sage: r.n(digits=6) 5555.56
why it's always 6 that suffices is another question, though
comment:23 Changed 9 months ago by
yes, this fixes this error:
-
src/sage/modules/diamond_cutting.py
a b def diamond_cut(V, GM, C, verbose=False): 154 154 (A vertex at (2), A vertex at (0)) 155 155 """ 156 156 # coerce to floats 157 from sage.misc.functional import round 157 158 GM = GM.n() 158 159 C = float(C) 159 160 if verbose: … … def diamond_cut(V, GM, C, verbose=False): 223 224 cut_count += 1 224 225 if verbose: 225 226 print("\n%d) Cut using normal vector %s" % (cut_count, hv)) 226 hv = [QQ(elmt.n(digits=6)) for elmt in hv] 227 hv = [QQ(round(elmt,6)) for elmt in hv] 228 # hv = [QQ(elmt.n(digits=6)) for elmt in hv] 227 229 inequalities.append(plane_inequality(hv)) 228 230 229 231 if verbose:
comment:24 Changed 9 months ago by
- Branch set to u/dimpase/modules/corrround
- Cc chapoton added
- Commit set to d734de1617b0ab38a60d7bc2e27f8e77107bede4
- Status changed from new to needs_review
I wonder how many more of these round()
calls were butchered this way, or this was an one-off event.
New commits:
d734de1 | the correct precision of 6 digits after dot
|
comment:25 Changed 9 months ago by
The line
hv = [QQ(round(elmt, 6)) for elmt in hv]
still bothers me a lot. Why 6
there? Is it just waiting for an example with bigger numbers to break it again?
comment:26 Changed 9 months ago by
- Commit changed from d734de1617b0ab38a60d7bc2e27f8e77107bede4 to 1faa6aceeb6799ea31517c9b84b3e4164153e78d
Branch pushed to git repo; I updated commit sha1. New commits:
1faa6ac | get rid of conversion to floats and round()
|
comment:27 Changed 9 months ago by
this should be robust now
comment:28 follow-up: ↓ 29 Changed 9 months ago by
- Description modified (diff)
- Keywords round added
Could it be that the conversion to QQ
was there for inexact vectors,
say vectors over RDF
or RR
?
Should we add the example in the ticket description as a doctest?
comment:29 in reply to: ↑ 28 ; follow-up: ↓ 30 Changed 9 months ago by
Replying to slelievre:
Could it be that the conversion to
RDF
orRR
?
for IntegerLattice this looks irrelevant.
probably it was a dodgy old workaround, from pre-PPL days.
Should we add the example in the ticket description as a doctest?
OK, I can do this.
comment:30 in reply to: ↑ 29 ; follow-up: ↓ 32 Changed 9 months ago by
Replying to dimpase:
Replying to slelievre:
Could it be that the conversion to
RDF
orRR
?for IntegerLattice this looks irrelevant.
Even though we are working with an integer lattice, doesn't this offer to find the closest lattice point to a non-integer vector?
comment:31 Changed 9 months ago by
- Commit changed from 1faa6aceeb6799ea31517c9b84b3e4164153e78d to 685ac9bdb3b9451c97707add3118d63b1006f366
Branch pushed to git repo; I updated commit sha1. New commits:
685ac9b | adding test from trac #29866
|
comment:32 in reply to: ↑ 30 Changed 9 months ago by
Replying to slelievre:
Replying to dimpase:
Replying to slelievre:
Could it be that the conversion to
RDF
orRR
?for IntegerLattice this looks irrelevant.
Even though we are working with an integer lattice, doesn't this offer to find the closest lattice point to a non-integer vector?
the code we're changing only deals with the lattice, not with the vector we are finding a closest lattice point for.
comment:33 Changed 9 months ago by
Right! One last thing, QQ is now imported but not used, can you remove the import?
comment:34 Changed 9 months ago by
- Commit changed from 685ac9bdb3b9451c97707add3118d63b1006f366 to 99f0df74fa4b4b1c836b8b850ceeb516ef012d22
Branch pushed to git repo; I updated commit sha1. New commits:
99f0df7 | remove unused import
|
comment:35 Changed 9 months ago by
done
comment:36 Changed 9 months ago by
please review
comment:37 Changed 9 months ago by
- Reviewers set to Matthias Koeppe
- Status changed from needs_review to positive_review
comment:38 Changed 9 months ago by
- Reviewers changed from Matthias Koeppe to Matthias Koeppe, Samuel Lelièvre
Thanks!
comment:39 Changed 8 months ago by
- Branch changed from u/dimpase/modules/corrround to 99f0df74fa4b4b1c836b8b850ceeb516ef012d22
- Resolution set to fixed
- Status changed from positive_review to closed
this might be just a bug in membership testing
prints
False
. (cf. continutiation of the thread on sage-devel in ticket description)