Examples of half-edge surfaces

This page contains examples of mathematical applications that were implemented/visualized using the half-edge module. It simultaneously extends the half-edge users guide as for some functions and modules it only makes sense to investigate them with given examples.

Computing Laplace invariants of a surface

In this example we compute the Laplace invariants of a given quad-surface, its Christoffel Dual, its face-circle-center dual (fcc_dual, build from centers of incircles of faces) and its Doliwa dual (build from intersection points of diagonals of a face). To ensure the existence of the duals our example surface is an S-CMC.

Note

An S-CMC surface is an S-Isothermic surface of constant mean curvature. It’s quad faces possess incircles and its vertices possess spheres. Four spheres that belong to vertices of the same face intersect the corresponding incircle orthogonally. In this example it is dealt with type 1 S-CMC surfaces that are characterized by the property that two vertex spheres that belong to a common edge touch at one point. S-CMC surfaces are discrete Koenigs nets.

An interior edge of a quad-surface is adjacent to two quadrilaterals. For each quadrilateral the edge opposite to the given one (w.r.t. the quadrilateral) intersects the given edge in a point. Denoting these intersections by \(I_1\) and \(I_2\), the Laplace invariant of an edge \(ij\) can be computed as

\[cr(v_1, I_1, v_2, I_2)\]

where the edge is incident to the vertices \(v_1\) and \(v_2\). A simple way of computing this value is using the coefficients of the planes induces by the adjacent quadrilaterals, see laplace_invariant().

Denoting by \(L_{ij}\) the Laplace invariants of edges of a discrete Koeniegs net, it holds

\[(L_{i0} * L_{i2}) / (L_{i1} * L_{i3}) = 1\]

for the edges in a vertex star around the same vertex (enumerated successively).

Similarly for the Doliwa dual of a Koenigs net (build from the intersection points of the diagonals of the quadrilaterals) this invariant holds for the edges around a face.

Furthermore for an edge and its dual edge of the Doliwa dual, the Laplace invariants coincide.

>>> import os
>>> import ddg
>>> import numpy as np

>>> # Function reads given s-cmc surfaces consisting of a primal surface and its dual build from the centers of face
>>> # circles (fcc dual).
>>> # The primal surface also stores christoffel dual coordinates.
>>> # The test builds the doliwa dual (dual build from intersection points of diagonals).
>>> # The following is tested:
>>> # - product of laplace invariants around a vertex star == 1 of the primal surface (due to being a Koenigs net)
>>> # - product of laplace invariants around a vertex star == 1 of the christoffel dual (due to being a Koenigs net)
>>> # - product of laplace invariants around a face == 1 of the doliwa dual (due to being Doliwa dual of a Koenigs net)
>>> # - same laplace invariant of an edge and its Doliwa dual edge

>>> ###################################
>>> # INITIALIZE
>>> ###################################

>>> # read hds from json
>>> dirname = "testing/tests/ddg/datastructures/halfedge/data"
>>> primal = ddg.halfedge.io.read_json_file(f"{dirname}/primal_14x14.json")
>>> fcc_dual = ddg.halfedge.io.read_json_file(f"{dirname}/dual_14x14.json")
>>>
>>> # assign cell attributes from index
>>> def index_attr_to_value(cell_class, attribute, value_list):
...     for cell in cell_class:
...         attribute[cell] = (
...             value_list[attribute[cell]] if attribute[cell] is not -1 else None
...         )
...
>>> index_attr_to_value(
...     primal.faces, primal.faces.dual_vertex, list([v for v in fcc_dual.verts])
... )
>>> index_attr_to_value(
...     fcc_dual.faces, fcc_dual.faces.dual_vertex, list([v for v in primal.verts])
... )
>>> index_attr_to_value(
...     primal.edges, primal.edges.dual_edge, list([v for v in fcc_dual.edges])
... )
>>> index_attr_to_value(
...     fcc_dual.edges, fcc_dual.edges.dual_edge, list([v for v in primal.edges])
... )

>>> ###################################
>>> # BUILD DOLIWA DUAL OF PRIMAL
>>> ###################################

>>> # SET DIAGONALS INTERSECTION PTS AS FACE ATTRIBUTES
>>> def func_faces(face):
...     Q = np.array(
...         [getattr(v, "co") for v in ddg.halfedge.get.get_vertices(face)]
...     )
...     return ddg.math.euclidean.intersect_diags(*Q)
...
>>> ddg.halfedge.set.set_attr_by_function(
...     primal,
...     func_faces,
...     cell_type="faces",
...     attr_name=f"diagonals_intersection_co",
... )
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>

>>> # BUILD DOLIWA DUAL
>>> # in the sense that the doliwa coordinates are stored in the vertices of the fcc dual surface
>>> # the fcc surface only is needed for the combinatorics
>>> doliwa_co_attr = fcc_dual.verts.add_attribute(f"doliwa_dual_co")
>>> for f in primal.faces:
...     doliwa_co_attr[f.dual_vertex] = getattr(f, f"diagonals_intersection_co")
...

>>> ###############################
>>> # COMPUTE LAPLACE INVARIANTS FOR PRIMAL, CHRISTOFFEL DUAL, FCC DUAL AND DOLIWA DUAL
>>> ###############################
>>> for surf, co_attr in [
...     (primal, "co"),
...     (primal, "christoffel_dual_co"),
...     (fcc_dual, "co"),
...     (fcc_dual, "doliwa_dual_co"),
... ]:
...     # bicolor edges in horizontal and vertical family with True/False
...     func_line_family = (
...         lambda edge: False if edge.radius_and_label_of_head[0] == "R" else True
...     )
...     line_family_attr = ddg.halfedge.set.set_attr_by_function(
...         surf, func_line_family, cell_type="edges", attr_name=f"line_family"
...     )
...     # set laplace invariants as edge attributes
...     func_edges = lambda edge: ddg.halfedge.math.laplace_invariant(
...         edge, co_attr=co_attr, line_family_attr=line_family_attr
...     )
...     ddg.halfedge.set.set_attr_by_function(
...         surf, func_edges, cell_type="edges", attr_name=f"laplace_inv_{co_attr}"
...     )
...
...     # set laplace invariants as vertex attributes (==1 for koenigs)
...     func_verts = lambda vert: ddg.halfedge.math.laplace_invariant_cross(
...         vert, getattr(surf.edges, f"laplace_inv_{co_attr}")
...     )
...     ddg.halfedge.set.set_attr_by_function(
...         surf,
...         func_verts,
...         cell_type="verts",
...         attr_name=f"laplace_inv_cross_{co_attr}",
...     )
...     # set laplace invariants as face attributes (!= 1 for s-cmc)
...     func_faces = lambda face: ddg.halfedge.math.laplace_invariant_quad(
...         face, getattr(surf.edges, f"laplace_inv_{co_attr}")
...     )
...     ddg.halfedge.set.set_attr_by_function(
...         surf,
...         func_faces,
...         cell_type="faces",
...         attr_name=f"laplace_inv_quad_{co_attr}",
...     )
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>
<ddg.datastructures.halfedge._node.NodeAttribute object at 0x...>

>>> ###################################
>>> # ASSERT LAPLACE INVARIANTS
>>> #####################################
>>> for v in ddg.halfedge.get.interior_vertices(primal):
...     assert np.isclose(getattr(v, f"laplace_inv_cross_co"), 1)
...     assert np.isclose(getattr(v, f"laplace_inv_cross_christoffel_dual_co"), 1)
...

>>> for f in fcc_dual.faces:
...     if getattr(f, f"laplace_inv_quad_doliwa_dual_co"):
...         assert np.isclose(getattr(f, f"laplace_inv_quad_doliwa_dual_co"), 1)
...

>>> for e in ddg.halfedge.get.single_edges(primal):
...     if e.dual_edge:
...         if getattr(e, f"laplace_inv_co") and getattr(
...             e.dual_edge, f"laplace_inv_doliwa_dual_co"
...         ):
...             assert np.isclose(
...                 getattr(e, f"laplace_inv_co"),
...                 getattr(e.dual_edge, f"laplace_inv_doliwa_dual_co"),
...             )
...

Note

The function laplace_invariant() takes an keyword argument called line_family_attr that stores an attribute on edges. This may be set to True/False for horizontal/vertical edges. If it is given and set to False simply the value will be returned otherwise 1 / value. In the example below the edges are split into families by the attribute edge.radius_and_label_of_head (given from previous implementations) that will be overwritten in the first step of the code below.