Merge branch 'main' into load_dataset-added

docs/api.rst

@@ -218,6 +218,8 @@ UGRID1D Topology
 
     Ugrid1d.node_edge_connectivity
     Ugrid1d.node_node_connectivity
+    Ugrid1d.directed_node_node_connectivity
+    Ugrid1d.directed_edge_edge_connectivity
     Ugrid1d.get_connectivity_matrix
 
     Ugrid1d.copy
@@ -300,6 +302,8 @@ UGRID2D Topology
     Ugrid2d.edge_node_connectivity
     Ugrid2d.face_edge_connectivity
     Ugrid2d.face_face_connectivity
+    Ugrid1d.directed_node_node_connectivity
+    Ugrid1d.directed_edge_edge_connectivity
     Ugrid2d.get_connectivity_matrix
 
     Ugrid2d.validate_edge_node_connectivity

docs/changelog.rst

@@ -30,7 +30,8 @@ Added
 - This method is used in :meth:`xugrid.UgridDataArray.from_structured2d` and
   :meth:`xugrid.UgridDataset.from_structured2d` when the optional arguments
   ``x_bounds`` and ``y_bounds`` are provided.
-
+- Added :attr:`xugrid.Ugrid1d.directed_edge_edge_connectivity` and
+  :attr:`xugrid.Ugrid2d.directed_edge_edge_connectivity`.
 
 [0.12.2] 2025-01-31
 -------------------

tests/test_connectivity.py

@@ -492,6 +492,45 @@ def test_directed_node_node_connectivity():
     assert np.array_equal(expected_edges, coo.data)
 
 
+def test_directed_edge_edge_connectivity():
+    r"""
+    Derive connectivity for this network:
+
+                             5
+                          4 --- 6
+                         /
+                        / 3
+       0     1     2   /
+    0 --- 1 --- 2 --- 3 
+                       \
+                        \ 4
+                         \
+                          5 
+    """
+    edge_node_connectivity = np.array(
+        [
+            [0, 1],  # edge 0
+            [1, 2],  # edge 1
+            [2, 3],  # edge 2
+            [3, 4],  # edge 3
+            [3, 5],  # edge 4
+            [4, 6],  # edge 5
+        ]
+    )
+    node_edge_connectivity = connectivity.invert_dense_to_sparse(edge_node_connectivity)
+    csr = connectivity.directed_edge_edge_connectivity(
+        edge_node_connectivity, node_edge_connectivity
+    )
+    assert isinstance(csr, sparse.csr_matrix)
+
+    coo = csr.tocoo()
+    actual = np.column_stack([coo.row, coo.col])
+    expected = np.array([[0, 1], [1, 2], [2, 3], [2, 4], [3, 5]])
+    assert np.array_equal(actual, expected)
+    # Test through which node the connection is formed.
+    assert np.array_equal(coo.data, [1, 2, 3, 3, 4])
+
+
 def test_face_face_connectivity():
     edge_faces = np.array(
         [

tests/test_ugrid1d.py

@@ -95,6 +95,8 @@ def test_ugrid1d_properties():
     assert grid.bounds == (0.0, 0.0, 2.0, 2.0)
     assert isinstance(grid.node_edge_connectivity, sparse.csr_matrix)
     assert isinstance(grid.node_node_connectivity, sparse.csr_matrix)
+    assert isinstance(grid.directed_node_node_connectivity, sparse.csr_matrix)
+    assert isinstance(grid.directed_edge_edge_connectivity, sparse.csr_matrix)
 
     expected_coords = [
         [[0.0, 0.0], [1.0, 1.0]],

tests/test_ugrid2d.py

@@ -151,6 +151,8 @@ def test_ugrid2d_properties():
     assert grid.bounds == (0.0, 0.0, 2.0, 2.0)
     assert isinstance(grid.node_node_connectivity, sparse.csr_matrix)
     assert isinstance(grid.node_edge_connectivity, sparse.csr_matrix)
+    assert isinstance(grid.directed_node_node_connectivity, sparse.csr_matrix)
+    assert isinstance(grid.directed_edge_edge_connectivity, sparse.csr_matrix)
     edge_node_coords = grid.edge_node_coordinates
     face_node_coords = grid.face_node_coordinates
     assert edge_node_coords.shape == (10, 2, 2)

xugrid/ugrid/connectivity.py

@@ -526,6 +526,28 @@ def node_node_connectivity(edge_node_connectivity: IntArray) -> sparse.csr_matri
     return coo_matrix.tocsr()
 
 
+def directed_edge_edge_connectivity(
+    edge_node_connectivity: IntArray,
+    node_edge_connectivity: sparse.csr_matrix,
+) -> sparse.csr_matrix:
+    # - Get the second node of each edge.
+    # - Find which edges are connected to this second node.
+    # - This will include the edge we started with.
+    # - Hence, make room for one more.
+    # - Then filter self -> self away.
+    n_edge = len(edge_node_connectivity)
+    second_node = edge_node_connectivity[:, 1]
+    n_downstream = node_edge_connectivity.getnnz(axis=1)[second_node]
+    upstream_edges = np.repeat(np.arange(n_edge), n_downstream)
+    downstream_edges = node_edge_connectivity[second_node].indices
+    node_index = np.repeat(second_node, n_downstream)
+    valid = downstream_edges != upstream_edges
+    return sparse.csr_matrix(
+        (node_index[valid], (upstream_edges[valid], downstream_edges[valid])),
+        shape=(n_edge, n_edge),
+    )
+
+
 def structured_connectivity(active: IntArray) -> AdjacencyMatrix:
     nrow, ncol = active.shape
     nodes = np.arange(nrow * ncol).reshape(nrow, ncol)

xugrid/ugrid/ugridbase.py

@@ -720,6 +720,25 @@ def directed_node_node_connectivity(self) -> csr_matrix:
         """
         return connectivity.directed_node_node_connectivity(self.edge_node_connectivity)
 
+    @property
+    def directed_edge_edge_connectivity(self) -> csr_matrix:
+        """
+        Directed edge to edge connectivity.
+
+        The connectivity is represented as an adjacency matrix in CSR format,
+        with the row and column indices as a (0-based) edge index. The data of
+        the matrix contains the node index of the common node through which
+        the connection is formed.
+
+        Returns
+        -------
+        connectivity: csr_matrix
+        """
+        return connectivity.directed_edge_edge_connectivity(
+            self.edge_node_connectivity,
+            self.node_edge_connectivity,
+        )
+
     @staticmethod
     def _connectivity_weights(connectivity: csr_matrix, coordinates: FloatArray):
         xy = coordinates