Add HarmonicDifferentials

benchmark/geometry/harmonic_differentials.py

@@ -0,0 +1,43 @@
+# ********************************************************************
+#  This file is part of sage-flatsurf.
+#
+#        Copyright (C) 2022 Julian Rüth
+#
+#  sage-flatsurf is free software: you can redistribute it and/or modify
+#  it under the terms of the GNU General Public License as published by
+#  the Free Software Foundation, either version 2 of the License, or
+#  (at your option) any later version.
+#
+#  sage-flatsurf is distributed in the hope that it will be useful,
+#  but WITHOUT ANY WARRANTY; without even the implied warranty of
+#  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+#  GNU General Public License for more details.
+#
+#  You should have received a copy of the GNU General Public License
+#  along with sage-flatsurf. If not, see <https://www.gnu.org/licenses/>.
+# ********************************************************************
+import flatsurf
+
+
+def time_harmonic_differential(surface):
+    if surface == "TORUS":
+        surface = flatsurf.translation_surfaces.torus((1, 0), (0, 1))
+    elif surface == "3413":
+        E = flatsurf.EquiangularPolygons(3, 4, 13)
+        P = E.an_element()
+        surface = flatsurf.similarity_surfaces.billiard(P, rational=True).minimal_cover(cover_type="translation")
+    else:
+        raise NotImplementedError
+
+    surface = surface.delaunay_triangulation()
+    surface.set_immutable()
+    Ω = flatsurf.HarmonicDifferentials(surface)
+    a = flatsurf.SimplicialHomology(surface).gens()[0]
+    H = flatsurf.SimplicialCohomology(surface)
+    Ω(H({a: 1}), check=False)
+
+
+time_harmonic_differential.params = ([
+    "TORUS",
+    "3413",
+])

doc/environment.yml

@@ -1,4 +1,4 @@
-# This file lists the minimial dependencies needed to build the documentation
+# This file lists the minimal dependencies needed to build the documentation
 # for sage-flatsurf. Create a conda environment with these dependencies
 # preinstalled with:
 # conda env create -n sage-flatsurf-docbuild -f ../environment.yml

flatsurf/__init__.py

@@ -3,6 +3,10 @@
 """
 from flatsurf.version import version as __version__
 
+from flatsurf.geometry.homology import SimplicialHomology
+from flatsurf.geometry.cohomology import SimplicialCohomology
+from flatsurf.geometry.harmonic_differentials import HarmonicDifferentials
+
 from flatsurf.geometry.polygon import (
     Polygon,
     polygons,
@@ -34,5 +38,5 @@
     ConeSurface,
     RationalConeSurface,
     HalfTranslationSurface,
-    TranslationSurface,
+    TranslationSurface
 )

flatsurf/features.py

@@ -4,7 +4,7 @@
 # ####################################################################
 #  This file is part of sage-flatsurf.
 #
-#        Copyright (C) 2021 Julian Rüth
+#        Copyright (C) 2021-2023 Julian Rüth
 #
 #  sage-flatsurf is free software: you can redistribute it and/or modify
 #  it under the terms of the GNU General Public License as published by
@@ -28,3 +28,7 @@
 pyflatsurf_feature = PythonModule(
     "pyflatsurf", url="https://github.com/flatsurf/flatsurf/#install-with-conda"
 )
+
+pyeantic_feature = PythonModule(
+    "pyeantic", url="https://github.com/flatsurf/e-antic/#install-with-conda"
+)

flatsurf/geometry/categories/cone_surfaces.py

@@ -382,6 +382,9 @@ def angles(self, numerical=False, return_adjacent_edges=False):
 
                         return angles
 
+                    def singularities(self):
+                        return [self.point(edges[0][0], self.polygon(edges[0][0]).vertex(edges[0][1])) for (angle, edges) in self.angles(return_adjacent_edges=True) if angle > 1]
+
                 class Connected(SurfaceCategoryWithAxiom):
                     r"""
                     The category of oriented connected cone surfaces without boundary.

flatsurf/geometry/categories/dilation_surfaces.py

@@ -338,6 +338,7 @@ def apply_matrix(self, m, in_place=True, mapping=False):
                 if not us.is_mutable():
                     raise ValueError("in-place changes only work for mutable surfaces")
                 for label in self.labels():
+                    # TODO: replace_polygon() does not exist anymore.
                     us.replace_polygon(label, m * self.polygon(label))
                 if m.det() < self.base_ring().zero():
                     # Polygons were all reversed orientation. Need to redo gluings.

flatsurf/geometry/categories/euclidean_polygons.py

@@ -1777,11 +1777,12 @@ def circumscribing_circle(self):
                     circle = circle_from_three_points(
                         self.vertex(0), self.vertex(1), self.vertex(2), self.base_ring()
                     )
-                    for i in range(3, len(self.vertices())):
-                        if not circle.point_position(self.vertex(i)) == 0:
-                            raise ValueError(
-                                "Vertex " + str(i) + " is not on the circle."
-                            )
+                    # TODO: Temporarily disabled because it fails over inexact rings.
+                    # for i in range(3, len(self.vertices())):
+                    #     if not circle.point_position(self.vertex(i)) == 0:
+                    #         raise ValueError(
+                    #             "Vertex " + str(i) + " is not on the circle."
+                    #         )
                     return circle
 
                 def subdivide(self):
@@ -1870,7 +1871,7 @@ def subdivide_edges(self, parts=2):
                     steps = [e / parts for e in self.edges()]
                     from flatsurf import Polygon
 
-                    return Polygon(edges=[e for e in steps for p in range(parts)])
+                    return Polygon(edges=[e for e in steps for p in range(parts)]).translate(self.vertex(0))
 
                 def j_invariant(self):
                     r"""

flatsurf/geometry/categories/polygonal_surfaces.py

@@ -1075,6 +1075,10 @@ def is_with_boundary(self):
                     False
 
                 """
+                if "WithoutBoundary" in self.category().axioms():
+                    return False
+                if "WithBoundary" in self.category().axioms():
+                    return True
                 for label in self.labels():
                     for edge in range(len(self.polygon(label).vertices())):
                         cross = self.opposite_edge(label, edge)

flatsurf/geometry/categories/similarity_surfaces.py

@@ -2073,7 +2073,7 @@ def subdivide(self):
 
                 Subdivision of this surface yields a surface with three triangles::
 
-                    sage: T = S.subdivide()
+                    sage: T = S.subdivide().codomain()
                     sage: T.labels()
                     (('Δ', 0), ('Δ', 1), ('Δ', 2))
 
@@ -2097,7 +2097,7 @@ def subdivide(self):
                     sage: S.glue(("Δ", 0), ("□", 2))
                     sage: S.glue(("□", 1), ("□", 3))
 
-                    sage: T = S.subdivide()
+                    sage: T = S.subdivide().codomain()
 
                     sage: T.labels()
                     (('Δ', 0), ('□', 2), ('Δ', 1), ('Δ', 2), ('□', 3), ('□', 1), ('□', 0))
@@ -2152,7 +2152,11 @@ def subdivide(self):
                         if opposite is not None:
                             surface.glue(((label, p), 0), (opposite, 0))
 
-                return surface
+                surface.set_immutable()
+
+                from flatsurf.geometry.deformation import SubdivideDeformation
+                return SubdivideDeformation(self, surface)
+
 
             def subdivide_edges(self, parts=2):
                 r"""
@@ -2188,7 +2192,7 @@ def subdivide_edges(self, parts=2):
                     sage: S.glue(("Δ", 0), ("□", 2))
                     sage: S.glue(("□", 1), ("□", 3))
 
-                    sage: T = S.subdivide_edges()
+                    sage: T = S.subdivide_edges().codomain()
                     sage: list(sorted(T.gluings()))
                     [(('Δ', 0), ('□', 5)),
                      (('Δ', 1), ('□', 4)),
@@ -2229,7 +2233,11 @@ def subdivide_edges(self, parts=2):
                                     ),
                                 )
 
-                return surface
+                surface.set_immutable()
+
+                from flatsurf.geometry.deformation import SubdivideEdgesDeformation
+
+                return SubdivideEdgesDeformation(self, surface, parts=parts)
 
     class Rational(SurfaceCategoryWithAxiom):
         r"""

flatsurf/geometry/categories/translation_surfaces.py

@@ -714,3 +714,26 @@ def canonicalize(self, in_place=None):
                     s.relabel(labels, in_place=True)
                     s.set_immutable()
                     return s
+
+                def _pyflatsurf(self):
+                    r"""
+                    Return a surface backed by libflatsurf that is isomorphic to this
+                    surface and an isomorphism to that surface.
+
+                    EXAMPLES::
+
+                        sage: from flatsurf import MutableOrientedSimilaritySurface, Polygon
+                        sage: S = MutableOrientedSimilaritySurface(QQ)
+                        sage: S.add_polygon(Polygon(vertices=[(0, 0), (1, 0), (1, 1)]), label=0)
+                        0
+                        sage: S.add_polygon(Polygon(vertices=[(0, 0), (1, 1), (0, 1)]), label=1)
+                        1
+                        sage: S.glue((0, 0), (1, 1))
+                        sage: S.glue((0, 1), (1, 2))
+                        sage: S.glue((0, 2), (1, 0))
+                        sage: S._to_pyflatsurf()  # optional: pyflatsurf
+
+                    """
+                    from flatsurf.geometry.pyflatsurf.surface import Surface_pyflatsurf
+
+                    return Surface_pyflatsurf._from_flatsurf(self)

flatsurf/geometry/cohomology.py

@@ -0,0 +1,166 @@
+r"""
+TODO: Document this module.
+"""
+######################################################################
+#  This file is part of sage-flatsurf.
+#
+#        Copyright (C) 2022 Julian Rüth
+#
+#  sage-flatsurf is free software: you can redistribute it and/or modify
+#  it under the terms of the GNU General Public License as published by
+#  the Free Software Foundation, either version 2 of the License, or
+#  (at your option) any later version.
+#
+#  sage-flatsurf is distributed in the hope that it will be useful,
+#  but WITHOUT ANY WARRANTY; without even the implied warranty of
+#  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+#  GNU General Public License for more details.
+#
+#  You should have received a copy of the GNU General Public License
+#  along with sage-flatsurf. If not, see <https://www.gnu.org/licenses/>.
+######################################################################
+
+from sage.structure.parent import Parent
+from sage.structure.element import Element
+from sage.categories.all import SetsWithPartialMaps
+from sage.structure.unique_representation import UniqueRepresentation
+from sage.misc.cachefunc import cached_method
+
+
+class SimplicialCohomologyClass(Element):
+    def __init__(self, parent, values):
+        super().__init__(parent)
+
+        self._values = values
+
+    def _repr_(self):
+        return repr(self._values)
+
+    def __call__(self, homology):
+        r"""
+        Evaluate this class at an element of homology.
+
+        EXAMPLES::
+
+            sage: from flatsurf import translation_surfaces, SimplicialCohomology
+            sage: T = translation_surfaces.torus((1, 0), (0, 1)).delaunay_triangulation()
+            sage: T.set_immutable()
+            sage: H = SimplicialCohomology(T)
+
+            sage: γ = H.homology().gens()[0]
+            sage: f = H({γ: 1.337})
+            sage: f(γ)
+            1.33700000000000
+
+        """
+        return sum(
+            self._values.get(gen, 0) * homology.coefficient(gen)
+            for gen in self.parent().homology().gens()
+        )
+
+    def _add_(self, other):
+        r"""
+        Return the pointwise sum of two homology classes.
+
+        EXAMPLES::
+
+            sage: from flatsurf import translation_surfaces, SimplicialCohomology
+            sage: T = translation_surfaces.torus((1, 0), (0, 1)).delaunay_triangulation()
+            sage: T.set_immutable()
+            sage: H = SimplicialCohomology(T)
+
+            sage: γ = H.homology().gens()[0]
+            sage: f = H({γ: 1.337})
+            sage: (f + f)(γ) == 2*f(γ)
+            True
+
+        """
+        values = {}
+        for gen in self.parent().homology().gens():
+            value = self(gen) + other(gen)
+            if value:
+                values[gen] = value
+
+        return self.parent()(values)
+
+
+class SimplicialCohomology(UniqueRepresentation, Parent):
+    r"""
+    Absolute simplicial cohomology of the ``surface`` with ``coefficients``.
+
+    EXAMPLES::
+
+        sage: from flatsurf import translation_surfaces, SimplicialCohomology
+        sage: T = translation_surfaces.torus((1, 0), (0, 1)).delaunay_triangulation()
+        sage: T.set_immutable()
+
+    Currently, surfaces must be Delaunay triangulated to compute their homology::
+
+        sage: SimplicialCohomology(T)
+        H¹(Translation Surface in H_1(0) built from 2 isosceles triangles; Real Field with 53 bits of precision)
+
+    """
+    Element = SimplicialCohomologyClass
+
+    @staticmethod
+    def __classcall__(cls, surface, coefficients=None, homology=None, category=None):
+        r"""
+        Normalize parameters used to construct cohomology.
+
+        TESTS:
+
+        Cohomology is unique and cached::
+
+            sage: from flatsurf import translation_surfaces, SimplicialCohomology
+            sage: T = translation_surfaces.torus((1, 0), (0, 1)).delaunay_triangulation()
+            sage: T.set_immutable()
+            sage: SimplicialCohomology(T) is SimplicialCohomology(T)
+            True
+
+        """
+        from sage.all import RR
+        from flatsurf.geometry.homology import SimplicialHomology
+        return super().__classcall__(cls, surface, coefficients or RR, homology or SimplicialHomology(surface), category or SetsWithPartialMaps())
+
+    def __init__(self, surface, coefficients, homology, category):
+        # TODO: Not checking this anymore. Do we need it?
+        # if surface != surface.delaunay_triangulation():
+        #     # TODO: This is a silly limitation in here.
+        #     raise NotImplementedError("Surface must be Delaunay triangulated")
+
+        Parent.__init__(self, category=category)
+
+        self._surface = surface
+        self._homology = homology
+        self._coefficients = coefficients
+
+    def _repr_(self):
+        return f"H¹({self._surface}; {self._coefficients})"
+
+    def _element_constructor_(self, x):
+        if not x:
+            x = {}
+
+        if isinstance(x, dict):
+            assert all(k in self.homology().gens() for k in x.keys())
+            return self.element_class(self, x)
+
+        raise NotImplementedError
+
+    @cached_method
+    def homology(self):
+        r"""
+        Return the homology of the underlying space (with integer
+        coefficients.)
+
+        EXAMPLES::
+
+        sage: from flatsurf import translation_surfaces, SimplicialCohomology
+        sage: T = translation_surfaces.torus((1, 0), (0, 1)).delaunay_triangulation()
+        sage: T.set_immutable()
+        sage: H = SimplicialCohomology(T)
+        sage: H.homology()
+        H₁(Translation Surface in H_1(0) built from 2 isosceles triangles; Integer Ring)
+
+        """
+        return self._homology