Skip to content

Commit

Permalink
sagemathgh-35997: Improve getting matrix entries after permutation in…
Browse files Browse the repository at this point in the history
… _palp_PM_max()

<!-- ^^^^^
Please provide a concise, informative and self-explanatory title.
Don't put issue numbers in there, do this in the PR body below.
For example, instead of "Fixes sagemath#1234" use "Introduce new method to
calculate 1+1"
-->
<!-- Describe your changes here in detail --> Modify the _palp_PM_max()
to obtain the matrix entry after permutation without make copy of the
matrix. @mkoeppe
before
```
sage: %timeit -n 10 -r 3 o = lattice_polytope.cross_polytope(3);
o._palp_PM_max();
36.8 ms ± 5.19 ms per loop (mean ± std. dev. of 3 runs, 10 loops each)
```
after
```
sage: %timeit -n 10 -r 3 o = lattice_polytope.cross_polytope(3);
o._palp_PM_max();
4.88 ms ± 2.23 ms per loop (mean ± std. dev. of 3 runs, 10 loops each)

```
<!-- Why is this change required? What problem does it solve? -->
<!-- If this PR resolves an open issue, please link to it here. For
example "Fixes sagemath#12345". -->
<!-- If your change requires a documentation PR, please link it
appropriately. -->

### 📝 Checklist

<!-- Put an `x` in all the boxes that apply. -->
<!-- If your change requires a documentation PR, please link it
appropriately -->
<!-- If you're unsure about any of these, don't hesitate to ask. We're
here to help! -->
<!-- Feel free to remove irrelevant items. -->

- [x] The title is concise, informative, and self-explanatory.
- [x] The description explains in detail what this PR is about.
- [ ] I have linked a relevant issue or discussion.
- [ ] I have created tests covering the changes.
- [ ] I have updated the documentation accordingly.

### ⌛ Dependencies

<!-- List all open PRs that this PR logically depends on
- sagemath#12345: short description why this is a dependency
- sagemath#34567: ...
-->

<!-- If you're unsure about any of these, don't hesitate to ask. We're
here to help! -->

URL: sagemath#35997
Reported by: xuluze
Reviewer(s): Matthias Köppe
  • Loading branch information
Release Manager committed Jul 30, 2023
2 parents 845963c + 2df0d10 commit 79dc351
Show file tree
Hide file tree
Showing 3 changed files with 43 additions and 56 deletions.
6 changes: 3 additions & 3 deletions build/pkgs/configure/checksums.ini
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
tarball=configure-VERSION.tar.gz
sha1=017780512407f8da87a091c485cf8c7162b5adaf
md5=4e3e25464ee850c2593f16364bd6b206
cksum=531215747
sha1=d39abedb0ee7beb3d61b625d4627e10c347320f2
md5=6a0dcd364b964e10b594339094d64a36
cksum=2756767966
2 changes: 1 addition & 1 deletion build/pkgs/configure/package-version.txt
Original file line number Diff line number Diff line change
@@ -1 +1 @@
98520945e069c603400cf2ae902a843f7a1bda13
259402d021bcc64831e1f24b26da34b76820b566
91 changes: 39 additions & 52 deletions src/sage/geometry/lattice_polytope.py
Original file line number Diff line number Diff line change
Expand Up @@ -3188,11 +3188,6 @@ def _palp_PM_max(self, check=False):
sage: all(results) # long time
True
"""
def PGE(S, u, v):
if u == v:
return S.one()
return S((u, v), check=False)

PM = self.vertex_facet_pairing_matrix()
n_v = PM.ncols()
n_f = PM.nrows()
Expand All @@ -3202,53 +3197,45 @@ def PGE(S, u, v):
# and find all the ways of making the first row of PM_max
def index_of_max(iterable):
# returns the index of max of any iterable
m, x = 0, iterable[0]
for k, l in enumerate(iterable):
if l > x:
m, x = k, l
return m
return max(enumerate(iterable), key=lambda x: x[1])[0]

n_s = 1
permutations = {0 : [S_f.one(), S_v.one()]}
permutations = {0: [S_f.one(), S_v.one()]}
for j in range(n_v):
m = index_of_max(
[(PM.with_permuted_columns(permutations[0][1]))[0][i]
for i in range(j, n_v)])
m = index_of_max(PM[0, i] for i in range(j, n_v))
if m > 0:
permutations[0][1] = PGE(S_v, j + 1, m + j + 1) * permutations[0][1]
permutations[0][1] = S_v((j + 1, m + j + 1), check=False) * permutations[0][1]
first_row = list(PM[0])

# Arrange other rows one by one and compare with first row
for k in range(1, n_f):
# Error for k == 1 already!
permutations[n_s] = [S_f.one(), S_v.one()]
m = index_of_max(PM.with_permuted_columns(permutations[n_s][1])[k])
m = index_of_max(PM[k, permutations[n_s][1](j+1) - 1] for j in range(n_v))
if m > 0:
permutations[n_s][1] = PGE(S_v, 1, m+1) * permutations[n_s][1]
d = ((PM.with_permuted_columns(permutations[n_s][1]))[k][0]
permutations[n_s][1] = S_v((1, m + 1), check=False) * permutations[n_s][1]
d = (PM[k, permutations[n_s][1](1) - 1]
- permutations[0][1](first_row)[0])
if d < 0:
# The largest elt of this row is smaller than largest elt
# in 1st row, so nothing to do
continue
# otherwise:
for i in range(1, n_v):
m = index_of_max(
[PM.with_permuted_columns(permutations[n_s][1])[k][j]
for j in range(i, n_v)])
m = index_of_max(PM[k, permutations[n_s][1](j+1) - 1] for j in range(i,n_v))
if m > 0:
permutations[n_s][1] = PGE(S_v, i + 1, m + i + 1) \
permutations[n_s][1] = S_v((i + 1, m + i + 1), check=False) \
* permutations[n_s][1]
if d == 0:
d = (PM.with_permuted_columns(permutations[n_s][1])[k][i]
d = (PM[k, permutations[n_s][1](i+1) - 1]
-permutations[0][1](first_row)[i])
if d < 0:
break
if d < 0:
# This row is smaller than 1st row, so nothing to do
del permutations[n_s]
continue
permutations[n_s][0] = PGE(S_f, 1, k + 1) * permutations[n_s][0]
permutations[n_s][0] = S_f((1, k + 1), check=False) * permutations[n_s][0]
if d == 0:
# This row is the same, so we have a symmetry!
n_s += 1
Expand All @@ -3258,9 +3245,9 @@ def index_of_max(iterable):
first_row = list(PM[k])
permutations = {0: permutations[n_s]}
n_s = 1
permutations = {k:permutations[k] for k in permutations if k < n_s}
permutations = {k: permutations[k] for k in permutations if k < n_s}

b = PM.with_permuted_rows_and_columns(*permutations[0])[0]
b = tuple(PM[permutations[0][0](1) - 1, permutations[0][1](j+1) - 1] for j in range(n_v))
# Work out the restrictions the current permutations
# place on other permutations as a automorphisms
# of the first row
Expand Down Expand Up @@ -3294,34 +3281,35 @@ def index_of_max(iterable):
# between 0 and S(0)
for s in range(l, n_f):
for j in range(1, S[0]):
v = PM.with_permuted_rows_and_columns(
*permutations_bar[n_p])[s]
if v[0] < v[j]:
permutations_bar[n_p][1] = PGE(S_v, 1, j + 1) * permutations_bar[n_p][1]
v0 = PM[permutations_bar[n_p][0](s+1) - 1, permutations_bar[n_p][1](1) - 1]
vj = PM[permutations_bar[n_p][0](s+1) - 1, permutations_bar[n_p][1](j+1) - 1]
if v0 < vj:
permutations_bar[n_p][1] = S_v((1, j + 1), check=False) * permutations_bar[n_p][1]
if ccf == 0:
l_r[0] = PM.with_permuted_rows_and_columns(
*permutations_bar[n_p])[s][0]
permutations_bar[n_p][0] = PGE(S_f, l + 1, s + 1) * permutations_bar[n_p][0]
l_r[0] = PM[permutations_bar[n_p][0](s+1) - 1, permutations_bar[n_p][1](1) - 1]
if s != l:
permutations_bar[n_p][0] = S_f((l + 1, s + 1), check=False) * permutations_bar[n_p][0]
n_p += 1
ccf = 1
permutations_bar[n_p] = copy(permutations[k])
else:
d1 = PM.with_permuted_rows_and_columns(
*permutations_bar[n_p])[s][0]
d1 = PM[permutations_bar[n_p][0](s+1) - 1, permutations_bar[n_p][1](1) - 1]
d = d1 - l_r[0]
if d < 0:
# We move to the next line
continue
elif d==0:
# Maximal values agree, so possible symmetry
permutations_bar[n_p][0] = PGE(S_f, l + 1, s + 1) * permutations_bar[n_p][0]
if s != l:
permutations_bar[n_p][0] = S_f((l + 1, s + 1), check=False) * permutations_bar[n_p][0]
n_p += 1
permutations_bar[n_p] = copy(permutations[k])
else:
# We found a greater maximal value for first entry.
# It becomes our new reference:
l_r[0] = d1
permutations_bar[n_p][0] = PGE(S_f, l + 1, s + 1) * permutations_bar[n_p][0]
if s != l:
permutations_bar[n_p][0] = S_f((l + 1, s + 1), check=False) * permutations_bar[n_p][0]
# Forget previous work done
cf = 0
permutations_bar = {0:copy(permutations_bar[n_p])}
Expand All @@ -3343,18 +3331,16 @@ def index_of_max(iterable):
s -= 1
# Find the largest value in this symmetry block
for j in range(c + 1, h):
v = PM.with_permuted_rows_and_columns(
*permutations_bar[s])[l]
if (v[c] < v[j]):
permutations_bar[s][1] = PGE(S_v, c + 1, j + 1) * permutations_bar[s][1]
vc = PM[(permutations_bar[s][0])(l+1) - 1, (permutations_bar[s][1])(c+1) - 1]
vj = PM[(permutations_bar[s][0])(l+1) - 1, (permutations_bar[s][1])(j+1) - 1]
if (vc < vj):
permutations_bar[s][1] = S_v((c + 1, j + 1), check=False) * permutations_bar[s][1]
if ccf == 0:
# Set reference and carry on to next permutation
l_r[c] = PM.with_permuted_rows_and_columns(
*permutations_bar[s])[l][c]
l_r[c] = PM[(permutations_bar[s][0])(l+1) - 1, (permutations_bar[s][1])(c+1) - 1]
ccf = 1
else:
d1 = PM.with_permuted_rows_and_columns(
*permutations_bar[s])[l][c]
d1 = PM[(permutations_bar[s][0])(l+1) - 1, (permutations_bar[s][1])(c+1) - 1]
d = d1 - l_r[c]
if d < 0:
n_p -= 1
Expand Down Expand Up @@ -3383,7 +3369,7 @@ def index_of_max(iterable):
# the restrictions the last worked out
# row imposes.
c = 0
M = (PM.with_permuted_rows_and_columns(*permutations[0]))[l]
M = tuple(PM[permutations[0][0](l+1) - 1, permutations[0][1](j+1) - 1] for j in range(n_v))
while c < n_v:
s = S[c] + 1
S[c] = c + 1
Expand Down Expand Up @@ -5196,11 +5182,10 @@ def _palp_canonical_order(V, PM_max, permutations):
in 2-d lattice M, (1,3,2,4))
"""
n_v = PM_max.ncols()
n_f = PM_max.nrows()
S_v = SymmetricGroup(n_v)
p_c = S_v.one()
M_max = [max([PM_max[i][j] for i in range(n_f)]) for j in range(n_v)]
S_max = [sum([PM_max[i][j] for i in range(n_f)]) for j in range(n_v)]
M_max = [max(row[j] for row in PM_max.rows()) for j in range(n_v)]
S_max = sum(PM_max)
for i in range(n_v):
k = i
for j in range(i + 1, n_v):
Expand All @@ -5213,13 +5198,15 @@ def _palp_canonical_order(V, PM_max, permutations):
p_c = S_v((1 + i, 1 + k), check=False) * p_c
# Create array of possible NFs.
permutations = [p_c * l[1] for l in permutations.values()]
Vs = [(V.column_matrix().with_permuted_columns(sig).hermite_form(), sig)
Vmatrix = V.column_matrix()
Vs = [(Vmatrix.with_permuted_columns(sig).hermite_form(), sig)
for sig in permutations]
Vmin = min(Vs, key=lambda x:x[0])
vertices = [V.module()(_) for _ in Vmin[0].columns()]
Vmodule = V.module()
vertices = [Vmodule(_) for _ in Vmin[0].columns()]
for v in vertices:
v.set_immutable()
return (PointCollection(vertices, V.module()), Vmin[1])
return (PointCollection(vertices, Vmodule), Vmin[1])


def _palp_convert_permutation(permutation):
Expand Down

0 comments on commit 79dc351

Please sign in to comment.