I276 enabled new layering functionality

CHANGELOG.md

@@ -5,6 +5,7 @@ This project adheres to [Semantic Versioning](https://semver.org/).
 ## Unreleased
 
 ### Added
+- [#305](https://github.com/equinor/flownet/pull/305) Functionality of generating seperate flownets per layer is enabled, additional nodes are split between the layers according to the volume of concave hull of the layers.
 - [#298](https://github.com/equinor/flownet/pull/298) Connections between (well)nodes that go through non-reservoir are now removed. Angle threshold export to user.
 - [#296](https://github.com/equinor/flownet/pull/296) Adapted perforation strategy to allow for layering + bug fixes in the 'multiple' and 'multiple_based_on_workovers' perforation strategies.
 - [#284](https://github.com/equinor/flownet/pull/284) Added the option to specify cumulative phase rates as observations (WOPT, WWPT, WGPT, WGIT, WWIT)

src/flownet/ahm/_run_ahm.py

@@ -416,14 +416,20 @@ def run_flownet_history_matching(
         field_data.faults if config.model_parameters.fault_mult else None
     )
 
-    concave_hull_bounding_boxes: Optional[np.ndarray] = None
+    concave_hull_list: Optional[List[np.ndarray]] = None
     if config.flownet.data_source.concave_hull:
-        concave_hull_bounding_boxes = field_data.grid_cell_bounding_boxes
+        concave_hull_list = []
+        for layer_id in df_well_connections["LAYER_ID"].unique():
+            concave_hull_list.append(
+                field_data.grid_cell_bounding_boxes(layer_id=layer_id)
+            )
 
     df_entity_connections: pd.DataFrame = create_connections(
-        df_well_connections[["WELL_NAME", "X", "Y", "Z"]].drop_duplicates(keep="first"),
+        df_well_connections[["WELL_NAME", "X", "Y", "Z", "LAYER_ID"]].drop_duplicates(
+            keep="first"
+        ),
         config,
-        concave_hull_bounding_boxes=concave_hull_bounding_boxes,
+        concave_hull_list=concave_hull_list,
     )
 
     network = NetworkModel(

src/flownet/config_parser/_config_parser.py

@@ -1621,6 +1621,12 @@ def parse_config(
                 f"Valid way to define the layering is as nested list with "
                 f"layer intervals. e.g. ((1, 5),(6, 9),(10, 15))."
             )
+    if layers and not config.flownet.data_source.concave_hull:
+        raise ValueError(
+            "Concave_hull should be True when layers are defined. "
+            "The concave hulls of the layers are used to split "
+            "the number of additional nodes between the layers."
+        )
 
     req_relp_parameters: List[str] = []
     if (

src/flownet/data/from_flow.py

@@ -1,7 +1,7 @@
 import warnings
 import operator
 from pathlib import Path
-from typing import Union, List, Optional, Tuple
+from typing import Union, List, Tuple
 
 import numpy as np
 import pandas as pd
@@ -317,7 +317,7 @@ def _faults(self) -> pd.DataFrame:
 
         return df_faults.drop(["I", "J", "K"], axis=1)
 
-    def _grid_cell_bounding_boxes(self, layer_id: Optional[int] = None) -> np.ndarray:
+    def grid_cell_bounding_boxes(self, layer_id: int) -> np.ndarray:
         """
         Function to get the bounding box (x, y and z min + max) for all grid cells
 
@@ -328,7 +328,7 @@ def _grid_cell_bounding_boxes(self, layer_id: Optional[int] = None) -> np.ndarra
             A (active grid cells x 6) numpy array with columns [ xmin, xmax, ymin, ymax, zmin, zmax ]
             filtered on layer_id if not None.
         """
-        if layer_id is not None:
+        if self._layers:
             (k_min, k_max) = tuple(map(operator.sub, self._layers[layer_id], (1, 1)))
         else:
             (k_min, k_max) = (0, self._grid.nz)
@@ -362,11 +362,6 @@ def get_unique_regions(self, name: str) -> np.ndarray:
         """array with unique 'name' regions"""
         return np.unique(self._init[name][0])
 
-    @property
-    def grid_cell_bounding_boxes(self) -> np.ndarray:
-        """Boundingboxes for all gridcells"""
-        return self._grid_cell_bounding_boxes()
-
     @property
     def faults(self) -> pd.DataFrame:
         """dataframe with all fault data"""

src/flownet/network_model/_generate_connections.py

@@ -94,9 +94,51 @@ def _remove_recorded_connections(
     return conn_matrix
 
 
+def _split_additional_flow_nodes(
+    total_additional_nodes: int, concave_hull_list: List[np.ndarray]
+) -> List[int]:
+    """
+    This function splits the additional_flow_nodes defined in the config over the layers.
+    The division is based on the sum of the volume of the boundingboxes in the layer.
+
+    Args:
+        total_additional_nodes: The total number of additional nodes to add to the model
+        concave_hull_list: List of boundingboxes per layer, i.e., numpy array with x, y, z min/max
+            boundingboxes for each grid block
+
+    Returns:
+        List of additional nodes per layer (len of list, same as number of layers)
+
+    """
+    volumes = []
+    for xyz in concave_hull_list:
+        if len(xyz.shape) == 1:
+            num_bounding_boxes = 1
+            xyz = xyz[None, :]
+        else:
+            num_bounding_boxes = xyz.shape[0]
+        volume = sum(
+            [
+                (
+                    (xyz[i, 1] - xyz[i, 0])
+                    * (xyz[i, 3] - xyz[i, 2])
+                    * (xyz[i, 5] - xyz[i, 4])
+                )
+                for i in range(0, num_bounding_boxes)
+            ]
+        )
+        volumes.append(volume)
+
+    fractions = [v / sum(volumes) for v in volumes]
+    sep_add_flownodes = [round(frac * total_additional_nodes) for frac in fractions]
+
+    return sep_add_flownodes
+
+
 def _generate_connections(
     df_coordinates: pd.DataFrame,
     configuration: Any,
+    additional_flow_nodes: int,
     concave_hull_bounding_boxes: Optional[np.ndarray] = None,
 ) -> Tuple[List[Coordinate], List[Coordinate]]:
     """
@@ -106,7 +148,8 @@ def _generate_connections(
 
     Args:
         df_coordinates: coordinates on original DataFrame format
-        configuration: FlowNet configuration yaml
+        configuration: Flownet configuration yaml,
+        additional_flow_nodes: Number of additional flow nodes to generate
         concave_hull_bounding_boxes: Numpy array with x, y, z min/max boundingboxes for each grid block
 
     Returns:
@@ -128,7 +171,7 @@ def _generate_connections(
         tuple(elem)
         for elem in mitchell_best_candidate_modified_3d(
             well_perforations,
-            num_added_flow_nodes=configuration.flownet.additional_flow_nodes,
+            num_added_flow_nodes=additional_flow_nodes,
             num_candidates=configuration.flownet.additional_node_candidates,
             hull_factor=configuration.flownet.hull_factor,
             concave_hull_bounding_boxes=concave_hull_bounding_boxes,
@@ -319,7 +362,7 @@ def _create_entity_connection_matrix(
     aquifer_ends: List[Coordinate],
     max_distance_fraction: float,
     max_distance: float,
-    concave_hull_bounding_boxes: Optional[np.ndarray] = None,
+    concave_hull_list: Optional[List[np.ndarray]] = None,
     n_non_reservoir_evaluation: Optional[int] = 10,
 ) -> pd.DataFrame:
     """
@@ -333,7 +376,8 @@ def _create_entity_connection_matrix(
         aquifer_ends: List of coordinates of all aquifer ends
         max_distance_fraction: Fraction of longest connection distance to be removed
         max_distance: Maximum distance between nodes, removed otherwise
-        concave_hull_bounding_boxes: Numpy array with x, y, z min/max boundingboxes for each grid block
+        concave_hull_list: List of boundingboxes per layer, i.e., numpy array with x, y, z min/max
+            boundingboxes for each grid block
         n_non_reservoir_evaluation: Number of equally spaced points along a connection to check fornon-reservoir.
 
     Returns:
@@ -357,7 +401,7 @@ def _create_entity_connection_matrix(
         str_start_entity = __get_entity_str(df_coordinates, start)
         str_end_entity = __get_entity_str(df_coordinates, end)
 
-        if concave_hull_bounding_boxes is not None:
+        if concave_hull_list is not None:
             tube_coordinates = linspace(
                 start=start,
                 stop=end,
@@ -367,7 +411,12 @@ def _create_entity_connection_matrix(
                 axis=1,
             ).T
 
-            if not all(check_in_hull(concave_hull_bounding_boxes, tube_coordinates)):
+            if not any(
+                [
+                    all(check_in_hull(concave_hull, tube_coordinates))
+                    for concave_hull in concave_hull_list
+                ]
+            ):
                 continue
 
         df_out = df_out.append(
@@ -522,7 +571,7 @@ def _generate_aquifer_connections(
 def create_connections(
     df_coordinates: pd.DataFrame,
     configuration: Any,
-    concave_hull_bounding_boxes: Optional[np.ndarray] = None,
+    concave_hull_list: Optional[List[np.ndarray]] = None,
 ) -> pd.DataFrame:
     """
     Creates additional flow nodes to increase complexity of field simulation structure so that history-matching can
@@ -534,17 +583,45 @@ def create_connections(
     Args:
         df_coordinates: Original structure of entity and X, Y, Z coords
         configuration: FlowNet configuration yaml as dictionary
-        concave_hull_bounding_boxes: Numpy array with x, y, z min/max boundingboxes for each grid block
+        concave_hull_list: List of boundingboxes per layer, i.e., numpy array with x, y, z min/max
+            boundingboxes for each grid block
 
     Returns:
         Desired restructuring of start-end coordinates into separate columns, as per Flow needs.
 
     """
-    starts, ends = _generate_connections(
-        df_coordinates=df_coordinates,
-        configuration=configuration,
-        concave_hull_bounding_boxes=concave_hull_bounding_boxes,
-    )
+    if df_coordinates["LAYER_ID"].nunique() > 1 and concave_hull_list is not None:
+        additional_flow_nodes_list = _split_additional_flow_nodes(
+            total_additional_nodes=configuration.flownet.additional_flow_nodes,
+            concave_hull_list=concave_hull_list,
+        )
+    else:
+        additional_flow_nodes_list = [configuration.flownet.additional_flow_nodes]
+
+    starts: List[Coordinate] = []
+    ends: List[Coordinate] = []
+
+    if concave_hull_list is not None:
+        for i, layer_id in enumerate(df_coordinates["LAYER_ID"].unique()):
+            starts_append, ends_append = _generate_connections(
+                df_coordinates=df_coordinates[df_coordinates["LAYER_ID"] == layer_id],
+                configuration=configuration,
+                additional_flow_nodes=additional_flow_nodes_list[i],
+                concave_hull_bounding_boxes=concave_hull_list[i],
+            )
+            starts.extend(starts_append)
+            ends.extend(ends_append)
+    else:
+        for i, layer_id in enumerate(df_coordinates["LAYER_ID"].unique()):
+            starts_append, ends_append = _generate_connections(
+                df_coordinates=df_coordinates[df_coordinates["LAYER_ID"] == layer_id],
+                configuration=configuration,
+                additional_flow_nodes=additional_flow_nodes_list[i],
+                concave_hull_bounding_boxes=concave_hull_list,
+            )
+            starts.extend(starts_append)
+            ends.extend(ends_append)
+
     aquifer_starts: List[Coordinate] = []
     aquifer_ends: List[Coordinate] = []
 
@@ -566,6 +643,6 @@ def create_connections(
         aquifer_ends,
         configuration.flownet.max_distance_fraction,
         configuration.flownet.max_distance,
-        concave_hull_bounding_boxes=concave_hull_bounding_boxes,
+        concave_hull_list=concave_hull_list,
         n_non_reservoir_evaluation=configuration.flownet.n_non_reservoir_evaluation,
     ).reset_index()

tests/test_data_flow.py

@@ -56,15 +56,17 @@ def test_grid_cell_bounding_boxes() -> None:
         "multiple_based_on_workovers",
     )
 
-    # Test no argument and entire field being equal
+    # Test one layer for the whole field and no layers equal
     flowdata._layers = ((1, flowdata.grid.nz),)
-    assert_almost_equal(
-        flowdata._grid_cell_bounding_boxes(), flowdata._grid_cell_bounding_boxes(0)
-    )
+    field_one_layer = flowdata.grid_cell_bounding_boxes(0)
+
+    flowdata._layers = ()
+    field_no_layer = flowdata.grid_cell_bounding_boxes(0)
+    assert_almost_equal(field_one_layer, field_no_layer)
 
     # Test zero'th layer id
     flowdata._layers = ((1, 2), (3, 4))
-    result = flowdata._grid_cell_bounding_boxes(0)
+    result = flowdata.grid_cell_bounding_boxes(0)
     active_cells = EclRegion(flowdata.grid, True)
     active_cells.select_kslice(
         *tuple(map(operator.sub, flowdata._layers[0], (1, 1))), intersect=True
@@ -75,7 +77,7 @@ def test_grid_cell_bounding_boxes() -> None:
 
     # Test last layer id
     flowdata._layers = ((1, 2), (3, 4))
-    result = flowdata._grid_cell_bounding_boxes(1)
+    result = flowdata.grid_cell_bounding_boxes(1)
     active_cells = EclRegion(flowdata.grid, True)
     active_cells.select_kslice(
         *tuple(map(operator.sub, flowdata._layers[-1], (1, 1))), intersect=True

tests/test_generate_connections.py

@@ -8,6 +8,7 @@
     _generate_connections,
     _create_entity_connection_matrix,
     _is_angle_too_large,
+    _split_additional_flow_nodes,
 )
 
 DATA = {
@@ -104,8 +105,11 @@ def test_generate_connections() -> None:
     config.flownet.random_seed = 1
     config.flownet.angle_threshold = None
 
+    # pylint: disable=no-member
     starts, ends = _generate_connections(
-        df_coordinates=DF_COORDINATES, configuration=config
+        df_coordinates=DF_COORDINATES,
+        configuration=config,
+        additional_flow_nodes=config.flownet.additional_flow_nodes,
     )
 
     assert len(starts) == len(ends)
@@ -115,6 +119,7 @@ def test_generate_connections() -> None:
     starts, ends = _generate_connections(
         df_coordinates=DF_COORDINATES,
         configuration=config,
+        additional_flow_nodes=config.flownet.additional_flow_nodes,
         concave_hull_bounding_boxes=np.array([0, 2, 0, 2, 0, 2]).reshape(-1, 6),
     )
 
@@ -126,7 +131,9 @@ def test_generate_connections() -> None:
     config.flownet.angle_threshold = 150
 
     starts, ends = _generate_connections(
-        df_coordinates=DF_COORDINATES, configuration=config
+        df_coordinates=DF_COORDINATES,
+        configuration=config,
+        additional_flow_nodes=config.flownet.additional_flow_nodes,
     )
 
     assert len(starts) == len(ends)
@@ -136,12 +143,32 @@ def test_generate_connections() -> None:
     # Test removal of all connections
     config.flownet.angle_threshold = 1
     starts, ends = _generate_connections(
-        df_coordinates=DF_COORDINATES, configuration=config
+        df_coordinates=DF_COORDINATES,
+        configuration=config,
+        additional_flow_nodes=config.flownet.additional_flow_nodes,
     )
 
     assert len(starts) == len(ends) == 0
 
 
+def test_split_additional_flow_nodes() -> None:
+
+    total_additional_nodes = 100
+    concave_hull_list = [
+        np.array([[0, 2, 0, 2, 0, 2], [2, 4, 0, 2, 0, 2]]),
+        np.array([0, 2, 0, 2, 2, 4]),
+    ]
+
+    sep_add_flownodes = _split_additional_flow_nodes(
+        total_additional_nodes=total_additional_nodes,
+        concave_hull_list=concave_hull_list,
+    )
+
+    assert sum(sep_add_flownodes) == total_additional_nodes
+    assert sep_add_flownodes[0] == 67
+    assert sep_add_flownodes[1] == 33
+
+
 def test_create_entity_connection_matrix() -> None:
 
     df = _create_entity_connection_matrix(