feat(clearing): support directional transfer capacities in complex clearing

assume/common/utils.py

@@ -550,6 +550,130 @@ def create_incidence_matrix(lines, buses, zones_id=None):
     return incidence_matrix
 
 
+def aggregate_line_capacities(
+    lines: pd.DataFrame,
+    incidence_matrix: pd.DataFrame,
+    zones_id: str = None,
+    node_mapping: dict = None,
+) -> pd.DataFrame:
+    """
+    Compute forward and reverse capacities for each line (or aggregated edge) used
+    by the transport-model clearing (`complex_clearing`).
+
+    The function returns a DataFrame indexed by the columns of `incidence_matrix`
+    with columns `cap_forward` and `cap_reverse` (absolute MW). If a column in
+    `incidence_matrix` matches a line index in `lines`, the capacities are taken
+    from that physical line. Otherwise the function attempts a zone-pair style
+    aggregation: physical lines are mapped to node pairs using `node_mapping`
+    (or by treating buses as nodes), and capacities are summed per zone-pair.
+
+    Directional columns in `lines` take precedence:
+      - `s_nom_forward` used for forward (bus0 -> bus1)
+      - `s_nom_reverse` used for reverse (bus1 -> bus0)
+    If missing, fallback to `s_nom * s_max_pu` for that direction.
+
+    Args:
+        lines: DataFrame of lines (indexed by line id).
+        incidence_matrix: Incidence matrix whose columns identify edges/lines.
+        zones_id: Optional zones identifier (unused here, kept for API compatibility).
+        node_mapping: Optional mapping from bus id -> node/zone id.
+
+    Returns:
+        pd.DataFrame: indexed by `incidence_matrix.columns` with columns
+            ['cap_forward', 'cap_reverse'].
+    """
+
+    # prepare defaults for each physical line
+    per_line_caps = {}
+    for line_idx, line in lines.iterrows():
+        s_max_pu = (
+            lines.at[line_idx, "s_max_pu"]
+            if "s_max_pu" in lines.columns
+            and not pd.isna(lines.at[line_idx, "s_max_pu"])
+            else 1.0
+        )
+        default_capacity = lines.at[line_idx, "s_nom"] * s_max_pu
+
+        if "s_nom_forward" in lines.columns and not pd.isna(
+            lines.at[line_idx, "s_nom_forward"]
+        ):
+            cap_f = lines.at[line_idx, "s_nom_forward"]
+        else:
+            cap_f = default_capacity
+
+        if "s_nom_reverse" in lines.columns and not pd.isna(
+            lines.at[line_idx, "s_nom_reverse"]
+        ):
+            cap_r = lines.at[line_idx, "s_nom_reverse"]
+        else:
+            cap_r = default_capacity
+
+        per_line_caps[line_idx] = {
+            "cap_forward": float(cap_f),
+            "cap_reverse": float(cap_r),
+        }
+
+    # If all incidence columns directly match physical lines, return per-line caps
+    cols = list(incidence_matrix.columns)
+    if all(col in per_line_caps for col in cols):
+        df = pd.DataFrame.from_dict(per_line_caps, orient="index")
+        # Ensure ordering matches incidence_matrix.columns
+        return df.reindex(cols)
+
+    # Otherwise, attempt to aggregate by node-pair keys (zone-pair aggregation)
+    # Build mapping from physical line -> node pair key
+    if node_mapping is None:
+        # identity mapping: bus id -> bus id
+        node_mapping = {}
+        for _, row in lines.iterrows():
+            node_mapping[row["bus0"]] = row["bus0"]
+            node_mapping[row["bus1"]] = row["bus1"]
+
+    agg_caps = {col: {"cap_forward": 0.0, "cap_reverse": 0.0} for col in cols}
+
+    for line_idx, line in lines.iterrows():
+        bus0 = line["bus0"]
+        bus1 = line["bus1"]
+        node0 = node_mapping.get(bus0, bus0)
+        node1 = node_mapping.get(bus1, bus1)
+
+        # Determine forward/reverse capacities for this physical line
+        caps = per_line_caps[line_idx]
+
+        # Try matching a column that corresponds to the node0->node1 direction
+        key_f = f"{node0}_{node1}"
+        key_r = f"{node1}_{node0}"
+
+        if key_f in agg_caps:
+            agg_caps[key_f]["cap_forward"] += caps["cap_forward"]
+            agg_caps[key_f]["cap_reverse"] += caps["cap_reverse"]
+        elif key_r in agg_caps:
+            # If the aggregated column uses reversed ordering, still add capacities
+            agg_caps[key_r]["cap_forward"] += caps["cap_forward"]
+            agg_caps[key_r]["cap_reverse"] += caps["cap_reverse"]
+        else:
+            # final fallback: if no matching aggregated key, try to add to any column
+            # that contains either node name (best-effort)
+            matched = False
+            for col in cols:
+                if str(node0) in str(col) and str(node1) in str(col):
+                    agg_caps[col]["cap_forward"] += caps["cap_forward"]
+                    agg_caps[col]["cap_reverse"] += caps["cap_reverse"]
+                    matched = True
+                    break
+            if not matched:
+                # give up and skip mapping this physical line
+                logger.debug(
+                    f"aggregate_line_capacities: could not map line {line_idx} to incidence column"
+                )
+
+    df = pd.DataFrame.from_dict(agg_caps, orient="index")
+    # ensure numeric types
+    df["cap_forward"] = df["cap_forward"].astype(float)
+    df["cap_reverse"] = df["cap_reverse"].astype(float)
+    return df
+
+
 def normalize_availability(powerplants_df, availability_df):
     # Create a copy of the availability dataframe to avoid modifying the original
     normalized_df = availability_df.copy()

assume/markets/clearing_algorithms/complex_clearing.py

@@ -12,7 +12,11 @@
 from pyomo.opt import OptSolver, SolverFactory, TerminationCondition
 
 from assume.common.market_objects import MarketConfig, MarketProduct, Orderbook
-from assume.common.utils import create_incidence_matrix, get_supported_solver_pyomo
+from assume.common.utils import (
+    aggregate_line_capacities,
+    create_incidence_matrix,
+    get_supported_solver_pyomo,
+)
 from assume.markets.base_market import MarketRole
 
 # Set the log level to WARNING
@@ -32,6 +36,7 @@ def market_clearing_opt_constraints(
     with_linked_bids: bool,
     incidence_matrix: pd.DataFrame,
     lines: pd.DataFrame,
+    directional_capacities: pd.DataFrame = None,
 ):
     """
     Adds the constraints to the model.
@@ -184,17 +189,27 @@ def energy_balance_rule(model, node, t):
         model.transmission_constr = pyo.ConstraintList()
         for t in model.T:
             for line in model.lines:
-                # s_max_pu might also be time variant. but for now we assume it is static
-                s_max_pu = (
-                    lines.at[line, "s_max_pu"]
-                    if "s_max_pu" in lines.columns
-                    and not pd.isna(lines.at[line, "s_max_pu"])
-                    else 1.0
-                )
-                capacity = lines.at[line, "s_nom"] * s_max_pu
-                # Limit the flow on each line
-                model.transmission_constr.add(model.flows[t, line] <= capacity)
-                model.transmission_constr.add(model.flows[t, line] >= -capacity)
+                # If precomputed directional capacities are provided, use them
+                if (
+                    directional_capacities is not None
+                    and line in directional_capacities.index
+                ):
+                    cap_forward = directional_capacities.at[line, "cap_forward"]
+                    cap_reverse = directional_capacities.at[line, "cap_reverse"]
+                    model.transmission_constr.add(model.flows[t, line] <= cap_forward)
+                    model.transmission_constr.add(model.flows[t, line] >= -cap_reverse)
+                else:
+                    # s_max_pu might also be time variant. but for now we assume it is static
+                    s_max_pu = (
+                        lines.at[line, "s_max_pu"]
+                        if "s_max_pu" in lines.columns
+                        and not pd.isna(lines.at[line, "s_max_pu"])
+                        else 1.0
+                    )
+                    capacity = lines.at[line, "s_nom"] * s_max_pu
+                    # Limit the flow on each line (symmetric fallback)
+                    model.transmission_constr.add(model.flows[t, line] <= capacity)
+                    model.transmission_constr.add(model.flows[t, line] >= -capacity)
 
 
 def market_clearing_opt_objective(model: pyo.ConcreteModel, orders: Orderbook):
@@ -220,6 +235,7 @@ def market_clearing_opt(
     with_linked_bids: bool,
     incidence_matrix: pd.DataFrame = None,
     lines: pd.DataFrame = None,
+    directional_capacities: pd.DataFrame = None,
     solver: OptSolver = None,
     solver_options: dict = {},
     func_constraints=market_clearing_opt_constraints,
@@ -267,7 +283,14 @@ def market_clearing_opt(
     model = pyo.ConcreteModel()
 
     func_constraints(
-        model, orders, market_products, mode, with_linked_bids, incidence_matrix, lines
+        model,
+        orders,
+        market_products,
+        mode,
+        with_linked_bids,
+        incidence_matrix,
+        lines,
+        directional_capacities,
     )
 
     func_objective(model, orders)
@@ -382,7 +405,27 @@ def __init__(self, marketconfig: MarketConfig):
                 # Nodal Case
                 self.incidence_matrix = create_incidence_matrix(self.lines, buses)
                 self.nodes = buses.index.values
-
+            # Pre-compute directional capacities for use in the clearing constraints
+            try:
+                self.directional_capacities = aggregate_line_capacities(
+                    self.lines,
+                    self.incidence_matrix,
+                    zones_id=self.zones_id,
+                    node_mapping=self.node_to_zone,
+                )
+            except Exception:
+                self.directional_capacities = None
+
+            # Informational log if input contains directional columns
+            if self.lines is not None:
+                has_directional = (
+                    "s_nom_forward" in self.lines.columns
+                    or "s_nom_reverse" in self.lines.columns
+                )
+                if has_directional:
+                    logger.info(
+                        "Directional NTC columns detected in lines data. Asymmetric transfer limits will be applied."
+                    )
         self.log_flows = self.marketconfig.param_dict.get("log_flows", False)
         self.pricing_mechanism = self.marketconfig.param_dict.get(
             "pricing_mechanism", "pay_as_clear"
@@ -518,6 +561,7 @@ def clear(
                 with_linked_bids=with_linked_bids,
                 incidence_matrix=self.incidence_matrix,
                 lines=self.lines,
+                directional_capacities=getattr(self, "directional_capacities", None),
                 solver=self.solver,
                 solver_options=self.solver_options,
             )

docs/source/example_simulations.rst

@@ -50,6 +50,9 @@ Overview of Example Simulations
    * - small_with_zonal_clearing‡
      - example_01d
      - Implements zonal market clearing.
+   * - small_with_directional_zonal_clearing‡
+     - example_01j
+     - Demonstrates complex clearing with directional transfer capacities between zones.
    * - market_study_eom
      - example_01f
      - Showcases comparison of single market to multi market. Case 1 in [3]_

docs/source/market_mechanism.rst

@@ -79,7 +79,13 @@ is given by: :math:`\mathbf{a}_{c, p} \: u_c \leq u_{p} \quad \forall \: c, p \i
 
 with the incidence matrix :math:`\mathbf{a}_{c, p}` defining the links between bids as 1, if c is linked as child to p, 0 else.
 
-Flows in the network are limited by the Net Transfer Capacity ('s_nom' * 's_max_pu') of each line l: :math:`\quad -NTC_{l} \leq F_{l, t} \leq NTC_{l} \quad \forall \: l \in \mathcal{L}, t \in \mathcal{T}`,
+Flows in the network are limited by the Net Transfer Capacity of each line :math:`l`. By default ASSUME derives a symmetric limit from ``s_nom * s_max_pu`` and applies it in both directions:
+
+.. math::
+
+  -NTC_{l} \leq F_{l, t} \leq NTC_{l} \quad \forall \: l \in \mathcal{L}, t \in \mathcal{T}
+
+If ``lines.csv`` additionally provides ``s_nom_forward`` and/or ``s_nom_reverse``, the complex clearing uses these directional limits instead. This allows transport constraints to differ by flow direction, which is useful for zonal representations with asymmetric commercial transfer capacities.
 
 Because with this algorithm, paradoxically accepted bids (PABs) can occur, the objective is solved in an iterative manner:
 

docs/source/release_notes.rst

@@ -15,6 +15,7 @@ Upcoming Release
 
 **New Features:**
   - **Generic Forecasting Interface**: This interface enables to specify different forecast algorithms for preprocess, initialization and update during runtime. They can be specified in the config.yaml or unit csv files. For more information about currently implemented algorithms and how to specify them please read the documentation on Unit forecasts.
+  - **Directional transfer capacities in complex clearing**: ``complex_clearing`` can now use asymmetric line limits from ``s_nom_forward`` and ``s_nom_reverse`` instead of assuming the same transfer capacity in both directions.
 
 **Improvements:**
   - **In complex clearing, the solver instance is now created once during initialization of the clearing role and reused for each market clearing**. This improves performance for e.g. year-long simulations.

examples/inputs/example_01j/buses.csv

@@ -0,0 +1,4 @@
+name,v_nom,zone_id,x,y
+north_1,380.0,north_1,10.0,54.0
+north_2,380.0,north_2,9.5,53.5
+south,380.0,south,11.6,48.1

examples/inputs/example_01j/buses.csv.license

@@ -0,0 +1,3 @@
+SPDX-FileCopyrightText: ASSUME Developers
+
+SPDX-License-Identifier: AGPL-3.0-or-later