sagemath / sage

Main repository of SageMath
https://www.sagemath.org
Other
1.44k stars 481 forks source link

Polyhedron.affine_hull: more output options #27366

Closed dkrenn closed 3 years ago

dkrenn commented 5 years ago

At the moment when calling .affine_hull, either the polyhedron or the affine map can be returned. If both needed, parts need to be recomputed, so we extend the parameters to allow returning both at the same time.

Moreover, we also allow to additionally return the section map, i.e., the right inverse of the projection map. This is a preparation for #27365 and #31659.

Depends on #30551

CC: @jplab @videlec @kliem

Component: geometry

Keywords: polytope

Author: Daniel Krenn, Matthias Koeppe, Jonathan Kliem

Branch/Commit: eee1aad

Reviewer: Matthias Koeppe, Jonathan Kliem

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

mkoeppe commented 3 years ago
comment:36

Thanks a lot! Great idea to add the test method

kliem commented 3 years ago
comment:37

I started adding doctests and I realized that they are extremely long and nobody is actually going to check them...

mkoeppe commented 3 years ago

Changed branch from u/gh-kliem/affine-hull-more to u/mkoeppe/affine-hull-more

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from 69ceebc to 9629620

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

9629620Polyhedron_base.affine_hull_projection: Fix up use of echelong form
mkoeppe commented 3 years ago

Reviewer: Matthias Koeppe, ...

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from 9629620 to 8c24c5c

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

8c24c5cPolyhedron_base.affine_hull_projection: Remove last occurrence of 'parametric_form'
kliem commented 3 years ago
comment:43

There is still a doctest that fails terribly. (After fixing an incorrect doctest that I will push in a minute).

sage: p = data['polyhedron'].linear_transformation(data['section_map'][0].matrix().transpose()) + data['section_map'][1]                                                            
sage: p = Polyhedron([(0,0)], base_ring=RDF)                                                                                                                                        
sage: data = p.affine_hull_projection(return_all_data=True, orthogonal=True, extend=True)                                                                                           
sage: p1 = data['polyhedron'].linear_transformation(data['section_map'][0].matrix().transpose()) + data['section_map'][1]                                                           
sage: p1 == p                                                                                                                                                                       
True
sage: p1 = data['polyhedron'].linear_transformation(data['section_map'][0].matrix().transpose()) + data['section_map'][1]                
sage: p1 == p                                                         
True

Somehow linear transformation applied to a 0-dimensional polyhedron isn't stable. The problem really seems to be connected to underlying matrices:

sage: p = Polyhedron([(0,0)], base_ring=RDF)                                                                                                                                        
sage: data = p.affine_hull_projection(return_all_data=True, orthogonal=True, extend=True)                                                                                           
sage: M = data['section_map'][0].matrix().transpose()                                                                                                                               
sage: V = data['polyhedron'].vertices_matrix()                                                                                                                                      
sage: M*V                                                                                                                                                                           
[3.445383e-316]
[          0.0]
sage: M*V                                                                                                                                                                           
[3.445383e-316]
[          0.0]
sage: M*V                                                                                                                                                                           
[3.17525635e-316]
[            0.0]
sage: p1 = data['polyhedron'].linear_transformation(data['section_map'][0].matrix().transpose()) + data['section_map'][1]                                                           
sage: M*V                                                                                                                                                                           
[0.0]
[0.0]
sage: M*V                                                                                                                                                                           
[1.0]
[0.0]
sage: M*V                                                                                                                                                                           
[0.0]
[1.0]
sage: M*V                                                                                                                                                                           
[0.0]
[0.0]
sage: M*V                                                                                                                                                                           
[0.0]
[0.0]

Apparently, the content of the matrix is never initialized. You might say that multiplication of those matrices isn't well-defined. However, if you think of this as linear maps, there is only one linear map from RLF^0 to RLF^2.

kliem commented 3 years ago
comment:44

Ok, seems to work now.

Fixing the problem with real matrices was actually quite easy. The matrix was not initialized as I thought. Should I open a seperate ticket for this?


New commits:

9f5560ainitialize empty matrix after trivial multiplication
f9faa02minimal extension only avoid AA if the base ring is not already AA
kliem commented 3 years ago

Changed branch from u/mkoeppe/affine-hull-more to u/gh-kliem/affine-hull-more

kliem commented 3 years ago

Changed commit from 8c24c5c to f9faa02

mkoeppe commented 3 years ago
comment:45

Thanks for fixing this! I think it's fine to keep it on this ticket

mkoeppe commented 3 years ago
comment:46

Green patchbot.

But we may want to consider making the result a dataclass (https://docs.python.org/3/library/dataclasses.html) instead of a dictionary

mkoeppe commented 3 years ago
comment:47

(depends on #30551 - Drop Python 3.6 support - which is planned for 9.4 anyway)

mkoeppe commented 3 years ago
comment:48
from dataclasses import dataclass
from typing import Any
@dataclass
class AffineHullProjectionData:
   polyhedron: Any
   projection_linear_map: Any
   projection_translation: Any
   section_linear_map: Any
   section_translation: Any

Note, unpacking the pairs -- this would be future-proof in case we ever decide to add proper affine linear maps to Sage

mkoeppe commented 3 years ago

Changed branch from u/gh-kliem/affine-hull-more to u/mkoeppe/affine-hull-more

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from f9faa02 to 57fd3e1

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

57fd3e1Fixup doctest formatting
kliem commented 3 years ago
comment:51

This ticket should depend now on #30551, right?

Apparently sage.geometry.polyhedron.base is a startup module, which is terrible. I chased it down and it seems one needs to do a series of lazy/more careful imports.

Adding some lazy imports above, I can avoid the import of `sage.geometry.polyhedron.base'.

mkoeppe commented 3 years ago

Dependencies: #30551

mkoeppe commented 3 years ago
comment:53

Yes, I've added the dependency. But I don't think we need to merge in the branch of that.

mkoeppe commented 3 years ago
comment:54

I have opened #31705 for the lazy import improvements

kliem commented 3 years ago
comment:55
                L_section = linear_transformation(matrix(len(pivots), self.ambient_dim(),
                                                         [E[i] for i in range(len(pivots))]).transpose(),
                                                  side='right')

This is not always an inverse to A:

sage: set_random_seed(0)                                                                                                              
sage: M = random_matrix(ZZ, 5, 5, distribution='uniform')                                                                             
sage: while True: 
....:     M = random_matrix(ZZ, 5, 5, distribution='uniform')   
....:     if M.rank() != 5: 
....:         break 
....:                                                                                                                                 
sage: P = Polyhedron(M)                                                                                                               
sage: P._test_affine_hull_projection()                                                                                                
---------------------------------------------------------------------------
AssertionError                            Traceback (most recent call last)
<ipython-input-23-f3c92d590b2b> in <module>
----> 1 P._test_affine_hull_projection()

~/Applications/sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/base.py in _test_affine_hull_projection(self, tester, verbose, **options)
  10487             else:
  10488                 self_extend = self
> 10489             tester.assertEqual(data.polyhedron.linear_transformation(M)
  10490                                + data.section_translation,
  10491                                self_extend)

/usr/lib/python3.8/unittest/case.py in assertEqual(self, first, second, msg)
    910         """
    911         assertion_func = self._getAssertEqualityFunc(first, second)
--> 912         assertion_func(first, second, msg=msg)
    913 
    914     def assertNotEqual(self, first, second, msg=None):

/usr/lib/python3.8/unittest/case.py in _baseAssertEqual(self, first, second, msg)
    903             standardMsg = '%s != %s' % _common_shorten_repr(first, second)
    904             msg = self._formatMessage(msg, standardMsg)
--> 905             raise self.failureException(msg)
    906 
    907     def assertEqual(self, first, second, msg=None):

AssertionError: A 4-dimensional polyhedron in QQ^5 defined as the convex hull of 5 vertices != A 4-dimensional polyhedron in ZZ^5 defined as the convex hull of 5 vertices
kliem commented 3 years ago

Changed branch from u/mkoeppe/affine-hull-more to u/gh-kliem/affine-hull-more

kliem commented 3 years ago

Changed commit from 57fd3e1 to b17aa8c

kliem commented 3 years ago

New commits:

b17aa8cadd another doctest
7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

Changed commit from b17aa8c to eee1aad

7ed8c4ca-6d56-4ae9-953a-41e42b4ed313 commented 3 years ago

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

eee1aadfix section map
kliem commented 3 years ago

Changed reviewer from Matthias Koeppe, ... to Matthias Koeppe, Jonathan Kliem

kliem commented 3 years ago
comment:59

I'm happy with it now.

mkoeppe commented 3 years ago
comment:60

Thanks a lot for finding this elegant fix.

vbraun commented 3 years ago

Changed branch from u/gh-kliem/affine-hull-more to eee1aad