diff --git a/src/gems/expression/copy.py b/src/gems/expression/copy.py index 4dc984ec..fe87cf04 100644 --- a/src/gems/expression/copy.py +++ b/src/gems/expression/copy.py @@ -18,6 +18,7 @@ AllTimeSumNode, CeilNode, ComparisonNode, + DualNode, ExpressionNode, FloorNode, LiteralNode, @@ -26,6 +27,7 @@ ParameterNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, ScenarioOperatorNode, TimeEvalNode, TimeShiftNode, @@ -95,6 +97,12 @@ def maximum(self, node: MaxNode) -> ExpressionNode: def minimum(self, node: MinNode) -> ExpressionNode: return MinNode([visit(op, self) for op in node.operands]) + def dual(self, node: DualNode) -> ExpressionNode: + return DualNode(node.constraint_id) + + def reduced_cost(self, node: ReducedCostNode) -> ExpressionNode: + return ReducedCostNode(node.variable_id) + def copy_expression(expression: ExpressionNode) -> ExpressionNode: return visit(expression, CopyVisitor()) diff --git a/src/gems/expression/degree.py b/src/gems/expression/degree.py index 99ed7e56..6bd6f4c4 100644 --- a/src/gems/expression/degree.py +++ b/src/gems/expression/degree.py @@ -14,20 +14,24 @@ import gems.expression.scenario_operator from gems.expression.expression import ( + AdditionNode, AllTimeSumNode, + BinaryOperatorNode, CeilNode, + DualNode, FloorNode, MaxNode, MinNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, TimeEvalNode, TimeShiftNode, TimeSumNode, + UnaryOperatorNode, ) from .expression import ( - AdditionNode, ComparisonNode, DivisionNode, ExpressionNode, @@ -112,6 +116,12 @@ def maximum(self, node: MaxNode) -> int | float: def minimum(self, node: MinNode) -> int | float: return 0 if all(visit(op, self) == 0 for op in node.operands) else math.inf + def dual(self, node: DualNode) -> int | float: + return math.inf + + def reduced_cost(self, node: ReducedCostNode) -> int | float: + return math.inf + def compute_degree(expression: ExpressionNode) -> int | float: return visit(expression, ExpressionDegreeVisitor()) @@ -129,3 +139,16 @@ def is_linear(expr: ExpressionNode) -> bool: True if the expression is linear with respect to variables. """ return compute_degree(expr) <= 1 + + +def contains_dual_or_reduced_cost(expr: ExpressionNode) -> bool: + """Return True if expr contains any DualNode or ReducedCostNode.""" + if isinstance(expr, (DualNode, ReducedCostNode)): + return True + if isinstance(expr, (AdditionNode, MaxNode, MinNode)): + return any(contains_dual_or_reduced_cost(o) for o in expr.operands) + if isinstance(expr, BinaryOperatorNode): + return contains_dual_or_reduced_cost(expr.left) or contains_dual_or_reduced_cost(expr.right) + if isinstance(expr, UnaryOperatorNode): + return contains_dual_or_reduced_cost(expr.operand) + return False diff --git a/src/gems/expression/equality.py b/src/gems/expression/equality.py index 625d7d4a..094f5810 100644 --- a/src/gems/expression/equality.py +++ b/src/gems/expression/equality.py @@ -29,11 +29,13 @@ AllTimeSumNode, BinaryOperatorNode, CeilNode, + DualNode, FloorNode, MaxNode, MinNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, ScenarioOperatorNode, TimeEvalNode, TimeShiftNode, @@ -103,6 +105,10 @@ def visit(self, left: ExpressionNode, right: ExpressionNode) -> bool: return self.maximum(left, right) if isinstance(left, MinNode) and isinstance(right, MinNode): return self.minimum(left, right) + if isinstance(left, DualNode) and isinstance(right, DualNode): + return self.dual(left, right) + if isinstance(left, ReducedCostNode) and isinstance(right, ReducedCostNode): + return self.reduced_cost(left, right) raise NotImplementedError(f"Equality not implemented for {left.__class__}") def literal(self, left: LiteralNode, right: LiteralNode) -> bool: @@ -193,6 +199,12 @@ def minimum(self, left: MinNode, right: MinNode) -> bool: self.visit(l, r) for l, r in zip(left.operands, right.operands) ) + def dual(self, left: DualNode, right: DualNode) -> bool: + return left.constraint_id == right.constraint_id + + def reduced_cost(self, left: ReducedCostNode, right: ReducedCostNode) -> bool: + return left.variable_id == right.variable_id + def expressions_equal( left: ExpressionNode, right: ExpressionNode, abs_tol: float = 0, rel_tol: float = 0 diff --git a/src/gems/expression/evaluate.py b/src/gems/expression/evaluate.py index 11d2f8cd..fdc40a34 100644 --- a/src/gems/expression/evaluate.py +++ b/src/gems/expression/evaluate.py @@ -18,11 +18,13 @@ from gems.expression.expression import ( AllTimeSumNode, CeilNode, + DualNode, FloorNode, MaxNode, MinNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, TimeEvalNode, TimeShiftNode, TimeSumNode, @@ -124,6 +126,12 @@ def maximum(self, node: MaxNode) -> float: def minimum(self, node: MinNode) -> float: return min(visit(op, self) for op in node.operands) + def dual(self, node: DualNode) -> float: + raise NotImplementedError("dual() is not statically evaluable.") + + def reduced_cost(self, node: ReducedCostNode) -> float: + raise NotImplementedError("reduced_cost() is not statically evaluable.") + def evaluate(expression: ExpressionNode, value_provider: ValueProvider) -> float: return visit(expression, EvaluationVisitor(value_provider)) diff --git a/src/gems/expression/expression.py b/src/gems/expression/expression.py index c7cb6bea..5dae3244 100644 --- a/src/gems/expression/expression.py +++ b/src/gems/expression/expression.py @@ -329,6 +329,16 @@ def __post_init__(self) -> None: ) +@dataclass(frozen=True, eq=False) +class DualNode(ExpressionNode): + constraint_id: str + + +@dataclass(frozen=True, eq=False) +class ReducedCostNode(ExpressionNode): + variable_id: str + + def sum_expressions(expressions: Sequence[ExpressionNode]) -> ExpressionNode: if len(expressions) == 0: return LiteralNode(0) diff --git a/src/gems/expression/indexing.py b/src/gems/expression/indexing.py index 36ac822c..0460e61c 100644 --- a/src/gems/expression/indexing.py +++ b/src/gems/expression/indexing.py @@ -22,6 +22,7 @@ CeilNode, ComparisonNode, DivisionNode, + DualNode, ExpressionNode, FloorNode, LiteralNode, @@ -32,6 +33,7 @@ ParameterNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, ScenarioOperatorNode, TimeEvalNode, TimeShiftNode, @@ -48,6 +50,9 @@ def get_parameter_structure(self, name: str) -> IndexingStructure: ... @abstractmethod def get_variable_structure(self, name: str) -> IndexingStructure: ... + @abstractmethod + def get_constraint_structure(self, name: str) -> IndexingStructure: ... + @dataclass(frozen=True) class TimeScenarioIndexingVisitor(ExpressionVisitor[IndexingStructure]): @@ -137,6 +142,12 @@ def maximum(self, node: MaxNode) -> IndexingStructure: def minimum(self, node: MinNode) -> IndexingStructure: return self._combine(node.operands) + def dual(self, node: DualNode) -> IndexingStructure: + return self.context.get_constraint_structure(node.constraint_id) + + def reduced_cost(self, node: ReducedCostNode) -> IndexingStructure: + return self.context.get_variable_structure(node.variable_id) + def compute_indexation( expression: ExpressionNode, provider: IndexingStructureProvider diff --git a/src/gems/expression/parsing/parse_expression.py b/src/gems/expression/parsing/parse_expression.py index 47316bf8..57a326a9 100644 --- a/src/gems/expression/parsing/parse_expression.py +++ b/src/gems/expression/parsing/parse_expression.py @@ -9,7 +9,7 @@ # SPDX-License-Identifier: MPL-2.0 # # This file is part of the Antares project. -from dataclasses import dataclass +from dataclasses import dataclass, field from typing import Set from antlr4 import CommonTokenStream, InputStream @@ -20,8 +20,10 @@ from gems.expression.expression import ( Comparator, ComparisonNode, + DualNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, maximum, minimum, ) @@ -33,11 +35,12 @@ @dataclass(frozen=True) class ModelIdentifiers: """ - Allows to distinguish between parameters and variables. + Allows to distinguish between parameters, variables, and constraints. """ variables: Set[str] parameters: Set[str] + constraints: Set[str] = field(default_factory=set) def is_variable(self, identifier: str) -> bool: return identifier in self.variables @@ -175,12 +178,35 @@ def visitAllTimeSum(self, ctx: ExprParser.AllTimeSumContext) -> ExpressionNode: shifted_expr = ctx.expr().accept(self) # type: ignore return shifted_expr.time_sum() + def _visit_dual(self, arg_exprs: list) -> ExpressionNode: + if len(arg_exprs) != 1: + raise ValueError("dual() requires exactly 1 argument.") + cid = arg_exprs[0].getText() # type: ignore + if cid not in self.identifiers.constraints: + raise ValueError(f"'{cid}' is not a constraint of the model.") + return DualNode(cid) + + def _visit_reduced_cost(self, arg_exprs: list) -> ExpressionNode: + if len(arg_exprs) != 1: + raise ValueError("reduced_cost() requires exactly 1 argument.") + vid = arg_exprs[0].getText() # type: ignore + if vid not in self.identifiers.variables: + raise ValueError(f"'{vid}' is not a variable of the model.") + return ReducedCostNode(vid) + # Visit a parse tree produced by ExprParser#function. def visitFunction(self, ctx: ExprParser.FunctionContext) -> ExpressionNode: function_name: str = ctx.IDENTIFIER().getText() # type: ignore arg_list = ctx.argList() # type: ignore + arg_exprs = arg_list.expr() if arg_list is not None else [] # type: ignore + + if function_name == "dual": + return self._visit_dual(arg_exprs) + if function_name == "reduced_cost": + return self._visit_reduced_cost(arg_exprs) + args: list[ExpressionNode] = ( - [expr.accept(self) for expr in arg_list.expr()] # type: ignore + [expr.accept(self) for expr in arg_exprs] # type: ignore if arg_list is not None else [] ) diff --git a/src/gems/expression/print.py b/src/gems/expression/print.py index 54b1f015..a7ca3850 100644 --- a/src/gems/expression/print.py +++ b/src/gems/expression/print.py @@ -16,12 +16,14 @@ from gems.expression.expression import ( AllTimeSumNode, CeilNode, + DualNode, ExpressionNode, FloorNode, MaxNode, MinNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, TimeEvalNode, TimeShiftNode, TimeSumNode, @@ -128,6 +130,12 @@ def maximum(self, node: MaxNode) -> str: def minimum(self, node: MinNode) -> str: return "min(" + ", ".join(visit(op, self) for op in node.operands) + ")" + def dual(self, node: DualNode) -> str: + return f"dual({node.constraint_id})" + + def reduced_cost(self, node: ReducedCostNode) -> str: + return f"reduced_cost({node.variable_id})" + def print_expr(expression: ExpressionNode) -> str: return visit(expression, PrinterVisitor()) diff --git a/src/gems/expression/visitor.py b/src/gems/expression/visitor.py index 3a894648..a2687ed8 100644 --- a/src/gems/expression/visitor.py +++ b/src/gems/expression/visitor.py @@ -24,6 +24,7 @@ CeilNode, ComparisonNode, DivisionNode, + DualNode, ExpressionNode, FloorNode, LiteralNode, @@ -34,6 +35,7 @@ ParameterNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, ScenarioOperatorNode, TimeEvalNode, TimeShiftNode, @@ -110,6 +112,12 @@ def maximum(self, node: MaxNode) -> T: ... @abstractmethod def minimum(self, node: MinNode) -> T: ... + @abstractmethod + def dual(self, node: DualNode) -> T: ... + + @abstractmethod + def reduced_cost(self, node: ReducedCostNode) -> T: ... + def visit(root: ExpressionNode, visitor: ExpressionVisitor[T]) -> T: """ @@ -153,6 +161,10 @@ def visit(root: ExpressionNode, visitor: ExpressionVisitor[T]) -> T: return visitor.maximum(root) elif isinstance(root, MinNode): return visitor.minimum(root) + elif isinstance(root, DualNode): + return visitor.dual(root) + elif isinstance(root, ReducedCostNode): + return visitor.reduced_cost(root) raise ValueError(f"Unknown expression node type {root.__class__}") diff --git a/src/gems/model/model.py b/src/gems/model/model.py index a08e74d2..eb10cc29 100644 --- a/src/gems/model/model.py +++ b/src/gems/model/model.py @@ -35,7 +35,37 @@ def _make_structure_provider( parameters: Dict[str, Parameter], variables: Dict[str, Variable], + constraints: Optional[Dict[str, Constraint]] = None, ) -> IndexingStructureProvider: + # Pre-compute constraint structures using a params/vars-only base provider. + # Constraint expressions cannot contain dual()/reduced_cost(), so the base + # provider's get_constraint_structure is never invoked during this step. + constraint_structures: Dict[str, IndexingStructure] = {} + if constraints: + class _BaseProvider(IndexingStructureProvider): + def get_parameter_structure(self, name: str) -> IndexingStructure: + return parameters[name].structure + + def get_variable_structure(self, name: str) -> IndexingStructure: + return variables[name].structure + + def get_constraint_structure(self, name: str) -> IndexingStructure: + raise NotImplementedError( + f"Constraint structure for '{name}' not available at this stage." + ) + + base = _BaseProvider() + for cname, c in constraints.items(): + try: + constraint_structures[cname] = compute_indexation(c.expression, base) + except ValueError: + # Constraints containing unresolved port fields (sum_connections) + # cannot be indexed before port resolution; fall back to the most + # general structure so callers can still proceed. + constraint_structures[cname] = IndexingStructure( + time=True, scenario=True + ) + class Provider(IndexingStructureProvider): def get_parameter_structure(self, name: str) -> IndexingStructure: return parameters[name].structure @@ -43,6 +73,9 @@ def get_parameter_structure(self, name: str) -> IndexingStructure: def get_variable_structure(self, name: str) -> IndexingStructure: return variables[name].structure + def get_constraint_structure(self, name: str) -> IndexingStructure: + return constraint_structures[name] + return Provider() @@ -94,7 +127,9 @@ def _is_objective_contribution_valid( raise ValueError("Objective contribution must be a linear expression.") data_structure_provider = _make_structure_provider( - model.parameters, model.variables + model.parameters, + model.variables, + {**model.constraints, **model.binding_constraints}, ) objective_structure = compute_indexation( objective_contribution, data_structure_provider diff --git a/src/gems/model/port.py b/src/gems/model/port.py index d6f6db62..718a105f 100644 --- a/src/gems/model/port.py +++ b/src/gems/model/port.py @@ -29,11 +29,13 @@ AllTimeSumNode, BinaryOperatorNode, CeilNode, + DualNode, FloorNode, MaxNode, MinNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, ScenarioOperatorNode, TimeEvalNode, TimeShiftNode, @@ -162,6 +164,12 @@ def minimum(self, node: MinNode) -> None: def port_field_aggregator(self, node: PortFieldAggregatorNode) -> None: raise ValueError("Port definition cannot contain port field aggregation.") + def dual(self, node: DualNode) -> None: + pass # dual() is permitted in port-field definitions + + def reduced_cost(self, node: ReducedCostNode) -> None: + pass # reduced_cost() is permitted in port-field definitions + def _validate_port_field_expression(definition: PortFieldDefinition) -> None: visit(definition.definition, _PortFieldExpressionChecker()) diff --git a/src/gems/model/resolve_library.py b/src/gems/model/resolve_library.py index f079ec41..fc30475e 100644 --- a/src/gems/model/resolve_library.py +++ b/src/gems/model/resolve_library.py @@ -12,6 +12,7 @@ from typing import Dict, List, Optional, Set from gems.expression import ExpressionNode, literal +from gems.expression.degree import contains_dual_or_reduced_cost from gems.expression.indexing_structure import IndexingStructure from gems.expression.parsing.parse_expression import ModelIdentifiers, parse_expression from gems.model import ( @@ -166,20 +167,39 @@ def _convert_port_type(port_type: PortTypeSchema) -> PortType: ) +def _forbid_dual_or_rc(expr: ExpressionNode, context: str) -> None: + if contains_dual_or_reduced_cost(expr): + raise ValueError( + f"Operators dual/reduced_cost are not allowed in {context}." + ) + + def _resolve_model( input_model: ModelSchema, port_types: Dict[str, PortType], library_id: str ) -> Model: identifiers = ModelIdentifiers( variables={v.id for v in input_model.variables}, parameters={p.id for p in input_model.parameters}, + constraints={c.id for c in input_model.binding_constraints} + | {c.id for c in input_model.constraints}, ) + binding_constraints = [ + _to_constraint(c, identifiers) for c in input_model.binding_constraints + ] + constraints = [_to_constraint(c, identifiers) for c in input_model.constraints] + + for c in binding_constraints + constraints: + _forbid_dual_or_rc(c.expression, f"constraint '{c.name}'") + objective_contributions = None if input_model.objective_contributions: objective_contributions = { contrib.id: parse_expression(contrib.expression, identifiers) for contrib in input_model.objective_contributions } + for oid, expr in objective_contributions.items(): + _forbid_dual_or_rc(expr, f"objective contribution '{oid}'") extra_outputs = ( { @@ -198,10 +218,8 @@ def _resolve_model( _resolve_field_definition(d, identifiers) for d in input_model.port_field_definitions ], - binding_constraints=[ - _to_constraint(c, identifiers) for c in input_model.binding_constraints - ], - constraints=[_to_constraint(c, identifiers) for c in input_model.constraints], + binding_constraints=binding_constraints, + constraints=constraints, objective_contributions=objective_contributions, extra_outputs=extra_outputs, ) diff --git a/src/gems/simulation/extra_output.py b/src/gems/simulation/extra_output.py index 176d3917..f3ef08fb 100644 --- a/src/gems/simulation/extra_output.py +++ b/src/gems/simulation/extra_output.py @@ -31,7 +31,7 @@ import numpy as np import xarray as xr -from gems.expression.expression import VariableNode +from gems.expression.expression import DualNode, ReducedCostNode, VariableNode from gems.model.port import PortFieldId from gems.simulation.vectorized_builder import VectorizedBuilderBase from gems.study.system import Component @@ -107,6 +107,12 @@ class VectorizedExtraOutputBuilder(VectorizedBuilderBase[xr.DataArray]): """ var_solution_arrays: Dict[Tuple[str, str], xr.DataArray] + constraint_dual_arrays: Dict[Tuple[str, str], xr.DataArray] = field( + default_factory=dict + ) + var_reduced_cost_arrays: Dict[Tuple[str, str], xr.DataArray] = field( + default_factory=dict + ) def variable(self, node: VariableNode) -> xr.DataArray: key = (self.model_id, node.name) @@ -116,3 +122,21 @@ def variable(self, node: VariableNode) -> xr.DataArray: f"{self.model_id!r}." ) return self.var_solution_arrays[key] + + def dual(self, node: DualNode) -> xr.DataArray: + key = (self.model_id, node.constraint_id) + if key not in self.constraint_dual_arrays: + raise KeyError( + f"Dual of constraint '{node.constraint_id}' not found for model " + f"{self.model_id!r}." + ) + return self.constraint_dual_arrays[key] + + def reduced_cost(self, node: ReducedCostNode) -> xr.DataArray: + key = (self.model_id, node.variable_id) + if key not in self.var_reduced_cost_arrays: + raise KeyError( + f"Reduced cost of variable '{node.variable_id}' not found for model " + f"{self.model_id!r}." + ) + return self.var_reduced_cost_arrays[key] diff --git a/src/gems/simulation/simulation_table.py b/src/gems/simulation/simulation_table.py index 990f500a..0719af56 100644 --- a/src/gems/simulation/simulation_table.py +++ b/src/gems/simulation/simulation_table.py @@ -273,6 +273,9 @@ def _collect_extra_outputs( if lv.name in solution: var_solution_arrays[(mk, vname)] = solution[lv.name] + constraint_dual_arrays = self._collect_constraint_duals(problem) + var_reduced_cost_arrays = self._collect_reduced_costs(problem) + for mk, components in problem.study.model_components.items(): model = problem.study.models[mk] if not model.extra_outputs: @@ -286,6 +289,8 @@ def _collect_extra_outputs( model_id=mk_, param_arrays=problem.param_arrays, var_solution_arrays=var_solution_arrays, + constraint_dual_arrays=constraint_dual_arrays, + var_reduced_cost_arrays=var_reduced_cost_arrays, port_arrays={}, block_length=problem.block_length, ), @@ -296,6 +301,8 @@ def _collect_extra_outputs( model_id=mk, param_arrays=problem.param_arrays, var_solution_arrays=var_solution_arrays, + constraint_dual_arrays=constraint_dual_arrays, + var_reduced_cost_arrays=var_reduced_cost_arrays, port_arrays=port_arrays, block_length=problem.block_length, ) @@ -321,6 +328,59 @@ def _collect_extra_outputs( return dfs + # ------------------------------------------------------------------------- + # Dual / reduced-cost arrays (helpers for _collect_extra_outputs) + # ------------------------------------------------------------------------- + + @staticmethod + def _collect_constraint_duals( + problem: OptimizationProblem, + ) -> Dict[Tuple[str, str], xr.DataArray]: + """Return constraint shadow prices keyed by (model_key, constraint_name).""" + dual_dataset = problem.linopy_model.dual + result: Dict[Tuple[str, str], xr.DataArray] = {} + for mk in problem.study.model_components: + model = problem.study.models[mk] + prefix = mk.replace("-", "_") + all_constraints = {**model.constraints, **model.binding_constraints} + for cname in all_constraints: + safe = cname.replace(" ", "_").replace("-", "_") + dual_val: xr.DataArray = xr.DataArray(0.0) + lb_name = f"{prefix}__{safe}__lb" + ub_name = f"{prefix}__{safe}__ub" + if lb_name in dual_dataset: + dual_val = dual_val + dual_dataset[lb_name] # type: ignore[operator] + if ub_name in dual_dataset: + dual_val = dual_val + dual_dataset[ub_name] # type: ignore[operator] + result[(mk, cname)] = dual_val + return result + + @staticmethod + def _collect_reduced_costs( + problem: OptimizationProblem, + ) -> Dict[Tuple[str, str], xr.DataArray]: + """Return variable reduced costs keyed by (model_key, var_name).""" + solver_model = getattr(problem.linopy_model, "solver_model", None) + if solver_model is None or not hasattr(solver_model, "getSolution"): + return {} + try: + solution = solver_model.getSolution() + col_dual_vals = list(solution.col_dual) + vlabels = problem.linopy_model.matrices.vlabels + rc_series = pd.Series(col_dual_vals, index=vlabels, dtype=float) + rc_series.loc[-1] = float("nan") + + result: Dict[Tuple[str, str], xr.DataArray] = {} + for (mk, vname), lv in problem._linopy_vars.items(): + idx = np.ravel(lv.labels.values) + rc_vals = rc_series.reindex(idx).values.reshape(lv.labels.shape) + result[(mk, vname)] = xr.DataArray( + rc_vals, coords=lv.labels.coords, dims=lv.labels.dims + ) + return result + except Exception: + return {} + # ------------------------------------------------------------------------- # Objective value # ------------------------------------------------------------------------- diff --git a/src/gems/simulation/vectorized_builder.py b/src/gems/simulation/vectorized_builder.py index 0d5d05b5..486a9560 100644 --- a/src/gems/simulation/vectorized_builder.py +++ b/src/gems/simulation/vectorized_builder.py @@ -47,6 +47,7 @@ CeilNode, ComparisonNode, DivisionNode, + DualNode, ExpressionNode, FloorNode, LiteralNode, @@ -57,6 +58,7 @@ ParameterNode, PortFieldAggregatorNode, PortFieldNode, + ReducedCostNode, ScenarioOperatorNode, TimeEvalNode, TimeShiftNode, @@ -382,6 +384,18 @@ def minimum(self, node: MinNode) -> VectorizedExpr: result = xr.where(result <= op, result, op) # type: ignore[no-untyped-call,assignment,operator] return result # type: ignore[return-value] + def dual(self, node: DualNode) -> VectorizedExpr: + raise NotImplementedError( + f"dual() is only available in the extra-output builder, " + f"not in {type(self).__name__}." + ) + + def reduced_cost(self, node: ReducedCostNode) -> VectorizedExpr: + raise NotImplementedError( + f"reduced_cost() is only available in the extra-output builder, " + f"not in {type(self).__name__}." + ) + # ------------------------------------------------------------------ # # Private helpers # # ------------------------------------------------------------------ # @@ -537,6 +551,12 @@ def minimum(self, node: MinNode) -> xr.DataArray: lambda a, b: xr.where(a <= b, a, b), ops # type: ignore[return-value,no-untyped-call] ) + def dual(self, node: DualNode) -> xr.DataArray: + raise NotImplementedError + + def reduced_cost(self, node: ReducedCostNode) -> xr.DataArray: + raise NotImplementedError + def _and_mask( a: Optional[xr.DataArray], b: Optional[xr.DataArray] @@ -693,3 +713,9 @@ def variable(self, node: VariableNode) -> Optional[xr.DataArray]: def port_field(self, node: PortFieldNode) -> Optional[xr.DataArray]: return None + + def dual(self, node: DualNode) -> Optional[xr.DataArray]: + return None + + def reduced_cost(self, node: ReducedCostNode) -> Optional[xr.DataArray]: + return None diff --git a/src/gems/study/system.py b/src/gems/study/system.py index 3fd1bd08..1a808cb4 100644 --- a/src/gems/study/system.py +++ b/src/gems/study/system.py @@ -18,6 +18,8 @@ from dataclasses import dataclass, field, replace from typing import Any, Dict, Iterable, List, Optional +from gems.expression.degree import contains_dual_or_reduced_cost +from gems.expression.expression import ExpressionNode, PortFieldAggregatorNode, PortFieldNode from gems.model import PortField, PortType from gems.model.model import Model from gems.model.port import PortFieldId @@ -53,6 +55,58 @@ class PortRef: port_id: str +def _uses_sum_connections_on(expr: ExpressionNode, port_name: str, field_name: str) -> bool: + """Return True if expr contains sum_connections(port_name.field_name).""" + if ( + isinstance(expr, PortFieldAggregatorNode) + and isinstance(expr.operand, PortFieldNode) + and expr.operand.port_name == port_name + and expr.operand.field_name == field_name + ): + return True + for child in _children(expr): + if _uses_sum_connections_on(child, port_name, field_name): + return True + return False + + +def _check_no_dual_rc_sum_connections( + master_ref: PortRef, slave_ref: PortRef, field_name: str +) -> None: + master_model = master_ref.component.model + slave_model = slave_ref.component.model + master_port_id = PortFieldId(port_name=master_ref.port_id, field_name=field_name) + master_def = master_model.port_fields_definitions.get(master_port_id) + if master_def is None or not contains_dual_or_reduced_cost(master_def.definition): + return + for bc in slave_model.binding_constraints.values(): + if _uses_sum_connections_on(bc.expression, slave_ref.port_id, field_name): + raise ValueError( + f"Port-field definition '{master_port_id}' contains " + f"dual/reduced_cost (non-linear) and cannot be aggregated " + f"via sum_connections in a binding-constraint of model " + f"'{slave_model.id}'." + ) + + +def _children(expr: ExpressionNode) -> list: + from gems.expression.expression import ( + AdditionNode, + BinaryOperatorNode, + MaxNode, + MinNode, + UnaryOperatorNode, + ) + + if isinstance(expr, (AdditionNode, MaxNode, MinNode)): + return list(expr.operands) + if isinstance(expr, BinaryOperatorNode): + return [expr.left, expr.right] + if isinstance(expr, UnaryOperatorNode): + return [expr.operand] + return [] + + @dataclass() class PortsConnection: port1: PortRef @@ -93,9 +147,10 @@ def __validate_ports(self) -> None: f"Port field {field_name} on {port_1.port_name} has 2 definitions." ) - self.master_port[PortField(name=field_name)] = ( - self.port1 if def1 else self.port2 - ) + master_ref = self.port1 if def1 else self.port2 + slave_ref = self.port2 if def1 else self.port1 + self.master_port[PortField(name=field_name)] = master_ref + _check_no_dual_rc_sum_connections(master_ref, slave_ref, field_name) def get_port_type(self) -> PortType: port_1 = self.port1.component.model.ports.get(self.port1.port_id) diff --git a/tests/e2e/functional/studies/10_5/case_description.txt b/tests/e2e/functional/studies/10_5/case_description.txt new file mode 100644 index 00000000..40509d32 --- /dev/null +++ b/tests/e2e/functional/studies/10_5/case_description.txt @@ -0,0 +1,15 @@ +Test 10_5 : one time step + +One area (base_zone) + +One load : 90 MW + +Two generators, with max_power = 100 MW. One generator has got a marginal cost = 10, the other equals to 50. + + +Expected results: + +objective = 900 +base_zone.price = 10 +gas_base_zone.generation_reduced_cost = 0 +oil_base_zone.generation_reduced_cost = 40 \ No newline at end of file diff --git a/tests/e2e/functional/studies/10_5/input/model-libraries/library-eo.yml b/tests/e2e/functional/studies/10_5/input/model-libraries/library-eo.yml new file mode 100644 index 00000000..7ea50986 --- /dev/null +++ b/tests/e2e/functional/studies/10_5/input/model-libraries/library-eo.yml @@ -0,0 +1,95 @@ +# Copyright (c) 2024, RTE (https://www.rte-france.com) +# +# See AUTHORS.txt +# +# This Source Code Form is subject to the terms of the Mozilla Public +# License, v. 2.0. If a copy of the MPL was not distributed with this +# file, You can obtain one at http://mozilla.org/MPL/2.0/. +# +# SPDX-License-Identifier: MPL-2.0 +# +# This file is part of the Antares project. +library: + id: library-eo + description: Library for test cases on extra-outputs + + + port-types: + - id: flow + description: A port which transfers power flow + fields: + - id: flow + + models: + - id: area + parameters: + - id: spillage_cost + - id: ens_cost + variables: + - id: spillage + lower-bound: 0 + upper-bound: 1000000 + variable-type: continuous + - id: unsupplied_energy + lower-bound: 0 + upper-bound: 1000000 + variable-type: continuous + ports: + - id: balance_port + type: flow + binding-constraints: + - id: balance + expression: sum_connections(balance_port.flow) = spillage - unsupplied_energy + objective-contributions: + - id: operational_objective + expression: sum(spillage_cost * spillage + ens_cost * unsupplied_energy) + extra-outputs: + - id: spill_cost_contribution + expression: spillage_cost * spillage + - id: ens_cost_contribution + expression: ens_cost * unsupplied_energy + - id: price + expression: dual(balance) + + - id: load + parameters: + - id: load_value + time-dependent: true + scenario-dependent: false + ports: + - id: balance_port + type: flow + port-field-definitions: + - port: balance_port + field: flow + definition: -load_value + + + - id: simple_gen + parameters: + - id: p_max_cluster # timeseries that takes outages into account + scenario-dependent: false + time-dependent: false + - id: marginal_generation_cost + scenario-dependent: false + time-dependent: false + variables: + - id: generation + lower-bound: 0.0 + upper-bound: p_max_cluster + variable-type: continuous + ports: + - id: balance_port + type: flow + port-field-definitions: + - port: balance_port + field: flow + definition: generation + + objective-contributions: + - id: operational_objective + expression: sum(marginal_generation_cost * generation) + extra-outputs: + - id: generation_reduced_cost + expression: reduced_cost(generation) + diff --git a/tests/e2e/functional/studies/10_5/input/optim-config.yml b/tests/e2e/functional/studies/10_5/input/optim-config.yml new file mode 100644 index 00000000..d6b6556c --- /dev/null +++ b/tests/e2e/functional/studies/10_5/input/optim-config.yml @@ -0,0 +1,11 @@ +time-scope: + first-time-step: 0 + last-time-step: 0 + +solver-options: + name: highs + logs: true + parameters: "THREADS 1" + +scenario-scope: + nb-scenarios: 1 diff --git a/tests/e2e/functional/studies/10_5/input/system.yml b/tests/e2e/functional/studies/10_5/input/system.yml new file mode 100644 index 00000000..db179d4d --- /dev/null +++ b/tests/e2e/functional/studies/10_5/input/system.yml @@ -0,0 +1,82 @@ +# Copyright (c) 2024, RTE (https://www.rte-france.com) +# +# See AUTHORS.txt +# +# This Source Code Form is subject to the terms of the Mozilla Public +# License, v. 2.0. If a copy of the MPL was not distributed with this +# file, You can obtain one at http://mozilla.org/MPL/2.0/. +# +# SPDX-License-Identifier: MPL-2.0 +# +# This file is part of the Antares project. + + # Study case 10_5 : to test extra-output features + # Small test case with one node (1 thermal cluster, 1 one load) + +system: + + id: case_10_5 + components: + - id: base_zone + model: library-eo.area + parameters: + - id: spillage_cost + time-dependent: false + scenario-dependent: false + value: 10000 + - id: ens_cost + time-dependent: false + scenario-dependent: false + value: 20000 + + + - id: load_base_zone + model: library-eo.load + parameters: + - id: load_value + time-dependent: false + scenario-dependent: false + value : 90 + + + - id: gas_base_zone + model: library-eo.simple_gen + parameters: + - id: p_max_cluster # timeseries that takes outages into account + scenario-dependent: false + time-dependent: false + value: 100 + - id: marginal_generation_cost + scenario-dependent: false + time-dependent: false + value: 10 + + - id: oil_base_zone + model: library-eo.simple_gen + parameters: + - id: p_max_cluster # timeseries that takes outages into account + scenario-dependent: false + time-dependent: false + value: 100 + - id: marginal_generation_cost + scenario-dependent: false + time-dependent: false + value: 50 + + + connections: + - component1: base_zone + port1: balance_port + component2: load_base_zone + port2: balance_port + + - component1: base_zone + port1: balance_port + component2: gas_base_zone + port2: balance_port + + - component1: base_zone + port1: balance_port + component2: oil_base_zone + port2: balance_port + diff --git a/tests/e2e/functional/studies/10_5_1/case_description.txt b/tests/e2e/functional/studies/10_5_1/case_description.txt new file mode 100644 index 00000000..dfea30f3 --- /dev/null +++ b/tests/e2e/functional/studies/10_5_1/case_description.txt @@ -0,0 +1,15 @@ +Test 10_5 : one time step + +One area (base_zone) + +One load : 110 MW + +Two generators, with max_power = 100 MW. One generator has got a marginal cost = 10, the other equals to 50. + + +Expected results: + +objective = 1500 +base_zone.price = 50 +gas_base_zone.generation_reduced_cost = -40 +oil_base_zone.generation_reduced_cost = 0 diff --git a/tests/e2e/functional/studies/10_5_1/input/model-libraries/library-eo.yml b/tests/e2e/functional/studies/10_5_1/input/model-libraries/library-eo.yml new file mode 100644 index 00000000..7ea50986 --- /dev/null +++ b/tests/e2e/functional/studies/10_5_1/input/model-libraries/library-eo.yml @@ -0,0 +1,95 @@ +# Copyright (c) 2024, RTE (https://www.rte-france.com) +# +# See AUTHORS.txt +# +# This Source Code Form is subject to the terms of the Mozilla Public +# License, v. 2.0. If a copy of the MPL was not distributed with this +# file, You can obtain one at http://mozilla.org/MPL/2.0/. +# +# SPDX-License-Identifier: MPL-2.0 +# +# This file is part of the Antares project. +library: + id: library-eo + description: Library for test cases on extra-outputs + + + port-types: + - id: flow + description: A port which transfers power flow + fields: + - id: flow + + models: + - id: area + parameters: + - id: spillage_cost + - id: ens_cost + variables: + - id: spillage + lower-bound: 0 + upper-bound: 1000000 + variable-type: continuous + - id: unsupplied_energy + lower-bound: 0 + upper-bound: 1000000 + variable-type: continuous + ports: + - id: balance_port + type: flow + binding-constraints: + - id: balance + expression: sum_connections(balance_port.flow) = spillage - unsupplied_energy + objective-contributions: + - id: operational_objective + expression: sum(spillage_cost * spillage + ens_cost * unsupplied_energy) + extra-outputs: + - id: spill_cost_contribution + expression: spillage_cost * spillage + - id: ens_cost_contribution + expression: ens_cost * unsupplied_energy + - id: price + expression: dual(balance) + + - id: load + parameters: + - id: load_value + time-dependent: true + scenario-dependent: false + ports: + - id: balance_port + type: flow + port-field-definitions: + - port: balance_port + field: flow + definition: -load_value + + + - id: simple_gen + parameters: + - id: p_max_cluster # timeseries that takes outages into account + scenario-dependent: false + time-dependent: false + - id: marginal_generation_cost + scenario-dependent: false + time-dependent: false + variables: + - id: generation + lower-bound: 0.0 + upper-bound: p_max_cluster + variable-type: continuous + ports: + - id: balance_port + type: flow + port-field-definitions: + - port: balance_port + field: flow + definition: generation + + objective-contributions: + - id: operational_objective + expression: sum(marginal_generation_cost * generation) + extra-outputs: + - id: generation_reduced_cost + expression: reduced_cost(generation) + diff --git a/tests/e2e/functional/studies/10_5_1/input/optim-config.yml b/tests/e2e/functional/studies/10_5_1/input/optim-config.yml new file mode 100644 index 00000000..d6b6556c --- /dev/null +++ b/tests/e2e/functional/studies/10_5_1/input/optim-config.yml @@ -0,0 +1,11 @@ +time-scope: + first-time-step: 0 + last-time-step: 0 + +solver-options: + name: highs + logs: true + parameters: "THREADS 1" + +scenario-scope: + nb-scenarios: 1 diff --git a/tests/e2e/functional/studies/10_5_1/input/system.yml b/tests/e2e/functional/studies/10_5_1/input/system.yml new file mode 100644 index 00000000..ff888815 --- /dev/null +++ b/tests/e2e/functional/studies/10_5_1/input/system.yml @@ -0,0 +1,82 @@ +# Copyright (c) 2024, RTE (https://www.rte-france.com) +# +# See AUTHORS.txt +# +# This Source Code Form is subject to the terms of the Mozilla Public +# License, v. 2.0. If a copy of the MPL was not distributed with this +# file, You can obtain one at http://mozilla.org/MPL/2.0/. +# +# SPDX-License-Identifier: MPL-2.0 +# +# This file is part of the Antares project. + + # Study case 10_5 : to test extra-output features + # Small test case with one node (1 thermal cluster, 1 one load) + +system: + + id: case_10_5 + components: + - id: base_zone + model: library-eo.area + parameters: + - id: spillage_cost + time-dependent: false + scenario-dependent: false + value: 10000 + - id: ens_cost + time-dependent: false + scenario-dependent: false + value: 20000 + + + - id: load_base_zone + model: library-eo.load + parameters: + - id: load_value + time-dependent: false + scenario-dependent: false + value : 110 + + + - id: gas_base_zone + model: library-eo.simple_gen + parameters: + - id: p_max_cluster # timeseries that takes outages into account + scenario-dependent: false + time-dependent: false + value: 100 + - id: marginal_generation_cost + scenario-dependent: false + time-dependent: false + value: 10 + + - id: oil_base_zone + model: library-eo.simple_gen + parameters: + - id: p_max_cluster # timeseries that takes outages into account + scenario-dependent: false + time-dependent: false + value: 100 + - id: marginal_generation_cost + scenario-dependent: false + time-dependent: false + value: 50 + + + connections: + - component1: base_zone + port1: balance_port + component2: load_base_zone + port2: balance_port + + - component1: base_zone + port1: balance_port + component2: gas_base_zone + port2: balance_port + + - component1: base_zone + port1: balance_port + component2: oil_base_zone + port2: balance_port + diff --git a/tests/e2e/functional/studies/10_5_2/case_description.txt b/tests/e2e/functional/studies/10_5_2/case_description.txt new file mode 100644 index 00000000..e19a88e7 --- /dev/null +++ b/tests/e2e/functional/studies/10_5_2/case_description.txt @@ -0,0 +1,32 @@ +Test 10_5_2 : three time steps + +One area (base_zone). WIth unsp_cost = 20000 + +One load : +- t= 1 : 80 MW +- t= 2 : 150 MW +- t= 3 : 201 MW + +One generator (gas_base_zone) has got max_power = 100 MW and following marginal cost +- t = 1 : 10 €/MWh, +- t = 2 : 15 €/MWh, +- t = 3 : 40 €/MWh, + +One generator (oil_base_zone) has got max_power = 100 MW and following marginal cost +- t = 1 : 30 €/MWh, +- t = 2 : 10 €/MWh, +- t = 3 : 10 €/MWh, + + +Expected results: + +objective = 27550 +base_zone.price_1 = 10 +base_zone.price_2 = 15 +base_zone.price_3 = 20000 +gas_base_zone.generation_reduced_cost_1 = 0 +gas_base_zone.generation_reduced_cost_2 = 0 +gas_base_zone.generation_reduced_cost_3 = - 19960 +oil_base_zone.generation_reduced_cost_1 = 20 +oil_base_zone.generation_reduced_cost_2 = - 5 +oil_base_zone.generation_reduced_cost_3 = - 19990 \ No newline at end of file diff --git a/tests/e2e/functional/studies/10_5_2/input/data-series/gas_cost_serie.tsv b/tests/e2e/functional/studies/10_5_2/input/data-series/gas_cost_serie.tsv new file mode 100644 index 00000000..0a9c618a --- /dev/null +++ b/tests/e2e/functional/studies/10_5_2/input/data-series/gas_cost_serie.tsv @@ -0,0 +1,3 @@ +10 +15 +40 \ No newline at end of file diff --git a/tests/e2e/functional/studies/10_5_2/input/data-series/load_serie.tsv b/tests/e2e/functional/studies/10_5_2/input/data-series/load_serie.tsv new file mode 100644 index 00000000..35759cd6 --- /dev/null +++ b/tests/e2e/functional/studies/10_5_2/input/data-series/load_serie.tsv @@ -0,0 +1,3 @@ +80 +150 +201 \ No newline at end of file diff --git a/tests/e2e/functional/studies/10_5_2/input/data-series/load_ts.tsv b/tests/e2e/functional/studies/10_5_2/input/data-series/load_ts.tsv new file mode 100644 index 00000000..0179cf77 --- /dev/null +++ b/tests/e2e/functional/studies/10_5_2/input/data-series/load_ts.tsv @@ -0,0 +1,5 @@ +20 +40 +60 +80 +110 \ No newline at end of file diff --git a/tests/e2e/functional/studies/10_5_2/input/data-series/oil_cost_serie.tsv b/tests/e2e/functional/studies/10_5_2/input/data-series/oil_cost_serie.tsv new file mode 100644 index 00000000..bafa962d --- /dev/null +++ b/tests/e2e/functional/studies/10_5_2/input/data-series/oil_cost_serie.tsv @@ -0,0 +1,3 @@ +30 +10 +10 \ No newline at end of file diff --git a/tests/e2e/functional/studies/10_5_2/input/model-libraries/library-eo.yml b/tests/e2e/functional/studies/10_5_2/input/model-libraries/library-eo.yml new file mode 100644 index 00000000..4808c5d3 --- /dev/null +++ b/tests/e2e/functional/studies/10_5_2/input/model-libraries/library-eo.yml @@ -0,0 +1,97 @@ +# Copyright (c) 2024, RTE (https://www.rte-france.com) +# +# See AUTHORS.txt +# +# This Source Code Form is subject to the terms of the Mozilla Public +# License, v. 2.0. If a copy of the MPL was not distributed with this +# file, You can obtain one at http://mozilla.org/MPL/2.0/. +# +# SPDX-License-Identifier: MPL-2.0 +# +# This file is part of the Antares project. +library: + id: library-eo + description: Library for test cases on extra-outputs + + + port-types: + - id: flow + description: A port which transfers power flow + fields: + - id: flow + + models: + - id: area + parameters: + - id: spillage_cost + - id: ens_cost + variables: + - id: spillage + lower-bound: 0 + upper-bound: 1000000 + variable-type: continuous + - id: unsupplied_energy + lower-bound: 0 + upper-bound: 1000000 + variable-type: continuous + ports: + - id: balance_port + type: flow + binding-constraints: + - id: balance + expression: sum_connections(balance_port.flow) = spillage - unsupplied_energy + objective-contributions: + - id: operational_objective + expression: sum(spillage_cost * spillage + ens_cost * unsupplied_energy) + extra-outputs: + - id: spill_cost_contribution + expression: spillage_cost * spillage + - id: ens_cost_contribution + expression: ens_cost * unsupplied_energy + - id: price + expression: dual(balance) + + - id: load + parameters: + - id: load_value + time-dependent: true + scenario-dependent: false + ports: + - id: balance_port + type: flow + port-field-definitions: + - port: balance_port + field: flow + definition: -load_value + + + - id: simple_gen + parameters: + - id: p_max_cluster # timeseries that takes outages into account + scenario-dependent: false + time-dependent: false + - id: marginal_generation_cost + scenario-dependent: false + time-dependent: true + variables: + - id: generation + lower-bound: 0.0 + upper-bound: p_max_cluster + variable-type: continuous + ports: + - id: balance_port + type: flow + port-field-definitions: + - port: balance_port + field: flow + definition: generation + + objective-contributions: + - id: operational_objective + expression: sum(marginal_generation_cost * generation) + extra-outputs: + - id: generation_reduced_cost + expression: reduced_cost(generation) + + + \ No newline at end of file diff --git a/tests/e2e/functional/studies/10_5_2/input/optim-config.yml b/tests/e2e/functional/studies/10_5_2/input/optim-config.yml new file mode 100644 index 00000000..10f6408d --- /dev/null +++ b/tests/e2e/functional/studies/10_5_2/input/optim-config.yml @@ -0,0 +1,11 @@ +time-scope: + first-time-step: 0 + last-time-step: 2 + +solver-options: + name: highs + logs: true + parameters: "THREADS 1" + +scenario-scope: + nb-scenarios: 1 diff --git a/tests/e2e/functional/studies/10_5_2/input/system.yml b/tests/e2e/functional/studies/10_5_2/input/system.yml new file mode 100644 index 00000000..ce65a37f --- /dev/null +++ b/tests/e2e/functional/studies/10_5_2/input/system.yml @@ -0,0 +1,87 @@ +# Copyright (c) 2024, RTE (https://www.rte-france.com) +# +# See AUTHORS.txt +# +# This Source Code Form is subject to the terms of the Mozilla Public +# License, v. 2.0. If a copy of the MPL was not distributed with this +# file, You can obtain one at http://mozilla.org/MPL/2.0/. +# +# SPDX-License-Identifier: MPL-2.0 +# +# This file is part of the Antares project. + +# Study case 10_4 : to test extra-output features +# Small test case with one node (1 thermal cluster, 1 one load) + +system: + + id: case_10_5 + components: + - id: base_zone + model: library-eo.area + parameters: + - id: spillage_cost + time-dependent: false + scenario-dependent: false + value: 1000 + - id: ens_cost + time-dependent: false + scenario-dependent: false + value: 20000 + + + - id: load_base_zone + model: library-eo.load + parameters: + - id: load_value + time-dependent: true + scenario-dependent: false + value: load_serie + + + - id: gas_base_zone + model: library-eo.simple_gen + parameters: + - id: p_max_cluster # timeseries that takes outages into account + scenario-dependent: false + time-dependent: false + value: 100 + - id: marginal_generation_cost + scenario-dependent: false + time-dependent: true + value: gas_cost_serie + + - id: oil_base_zone + model: library-eo.simple_gen + parameters: + - id: p_max_cluster # timeseries that takes outages into account + scenario-dependent: false + time-dependent: false + value: 100 + - id: marginal_generation_cost + scenario-dependent: false + time-dependent: true + value: oil_cost_serie + + + connections: + - component1: base_zone + port1: balance_port + component2: load_base_zone + port2: balance_port + + - component1: base_zone + port1: balance_port + component2: gas_base_zone + port2: balance_port + + - component1: base_zone + port1: balance_port + component2: oil_base_zone + port2: balance_port + + + + + + \ No newline at end of file diff --git a/tests/e2e/functional/test_dual_reduced_cost.py b/tests/e2e/functional/test_dual_reduced_cost.py new file mode 100644 index 00000000..2e47675e --- /dev/null +++ b/tests/e2e/functional/test_dual_reduced_cost.py @@ -0,0 +1,154 @@ +# Copyright (c) 2024, RTE (https://www.rte-france.com) +# +# See AUTHORS.txt +# +# This Source Code Form is subject to the terms of the Mozilla Public +# License, v. 2.0. If a copy of the MPL was not distributed with this +# file, You can obtain one at http://mozilla.org/MPL/2.0/. +# +# SPDX-License-Identifier: MPL-2.0 +# +# This file is part of the Antares project. + +""" +End-to-end tests for the dual() and reduced_cost() extra-output operators. + +Studies 10_5 and 10_5_1: single timestep, two generators, one load. +Study 10_5_2: three timesteps, time-varying costs and loads. + +Negative tests verify that dual/reduced_cost are rejected when used in +contexts where they are not allowed (constraints, objectives). +""" + +from pathlib import Path + +import pytest + +from gems.simulation import TimeBlock, build_problem +from gems.simulation.simulation_table import SimulationTableBuilder +from gems.study.folder import load_study + +STUDIES_DIR = Path(__file__).parent / "studies" + + +@pytest.mark.parametrize( + "study_id, expected", + [ + ( + "10_5", + { + "objective": 900.0, + "base_zone.price": 10.0, + "gas_base_zone.generation_reduced_cost": 0.0, + "oil_base_zone.generation_reduced_cost": 40.0, + }, + ), + ( + "10_5_1", + { + "objective": 1500.0, + "base_zone.price": 50.0, + "gas_base_zone.generation_reduced_cost": -40.0, + "oil_base_zone.generation_reduced_cost": 0.0, + }, + ), + ], +) +def test_dual_reduced_cost_single_timestep( + study_id: str, expected: dict +) -> None: + """Verify nodal price (dual) and reduced costs for single-timestep studies.""" + study = load_study(STUDIES_DIR / study_id) + time_block = TimeBlock(1, [0]) + problem = build_problem(study, time_block, [0]) + problem.solve(solver_name="highs") + + assert problem.termination_condition == "optimal" + assert problem.objective_value == pytest.approx(expected["objective"]) + + st = SimulationTableBuilder().build(problem) + + price = st.component("base_zone").output("price").value( + time_index=0, scenario_index=0 + ) + assert price == pytest.approx(expected["base_zone.price"]) + + gas_rc = st.component("gas_base_zone").output("generation_reduced_cost").value( + time_index=0, scenario_index=0 + ) + assert gas_rc == pytest.approx(expected["gas_base_zone.generation_reduced_cost"]) + + oil_rc = st.component("oil_base_zone").output("generation_reduced_cost").value( + time_index=0, scenario_index=0 + ) + assert oil_rc == pytest.approx(expected["oil_base_zone.generation_reduced_cost"]) + + +def test_dual_reduced_cost_multi_timestep() -> None: + """Verify nodal prices and reduced costs for a 3-timestep study (10_5_2).""" + study = load_study(STUDIES_DIR / "10_5_2") + time_block = TimeBlock(1, [0, 1, 2]) + problem = build_problem(study, time_block, [0]) + problem.solve(solver_name="highs") + + assert problem.termination_condition == "optimal" + assert problem.objective_value == pytest.approx(27550.0) + + st = SimulationTableBuilder().build(problem) + + # Nodal prices at t=0,1,2 + for t, expected_price in enumerate([10.0, 15.0, 20000.0]): + price = st.component("base_zone").output("price").value( + time_index=t, scenario_index=0 + ) + assert price == pytest.approx(expected_price), ( + f"price at t={t}: expected {expected_price}, got {price}" + ) + + # Gas generator reduced costs at t=0,1,2 + for t, expected_rc in enumerate([0.0, 0.0, -19960.0]): + rc = st.component("gas_base_zone").output("generation_reduced_cost").value( + time_index=t, scenario_index=0 + ) + assert rc == pytest.approx(expected_rc, abs=1e-3), ( + f"gas RC at t={t}: expected {expected_rc}, got {rc}" + ) + + # Oil generator reduced costs at t=0,1,2 + for t, expected_rc in enumerate([20.0, -5.0, -19990.0]): + rc = st.component("oil_base_zone").output("generation_reduced_cost").value( + time_index=t, scenario_index=0 + ) + assert rc == pytest.approx(expected_rc, abs=1e-3), ( + f"oil RC at t={t}: expected {expected_rc}, got {rc}" + ) + + +def test_dual_in_constraint_is_rejected() -> None: + """dual() in a constraint expression must be caught by the library resolver.""" + from gems.expression.parsing.parse_expression import ModelIdentifiers, parse_expression + from gems.model.resolve_library import _forbid_dual_or_rc + + ids = ModelIdentifiers( + variables={"x"}, + parameters=set(), + constraints={"balance"}, + ) + expr = parse_expression("dual(balance) + x", ids) + with pytest.raises(ValueError, match="Operators dual/reduced_cost are not allowed"): + _forbid_dual_or_rc(expr, "constraint 'bad'") + + +def test_reduced_cost_in_objective_is_rejected() -> None: + """reduced_cost() in an objective contribution must be caught by the library resolver.""" + from gems.expression.parsing.parse_expression import ModelIdentifiers, parse_expression + from gems.model.resolve_library import _forbid_dual_or_rc + + ids = ModelIdentifiers( + variables={"x"}, + parameters=set(), + constraints=set(), + ) + expr = parse_expression("reduced_cost(x)", ids) + with pytest.raises(ValueError, match="Operators dual/reduced_cost are not allowed"): + _forbid_dual_or_rc(expr, "objective contribution 'obj'") diff --git a/tests/unittests/expressions/visitor/test_indexing.py b/tests/unittests/expressions/visitor/test_indexing.py index 6a40f0d7..28d707d5 100644 --- a/tests/unittests/expressions/visitor/test_indexing.py +++ b/tests/unittests/expressions/visitor/test_indexing.py @@ -33,6 +33,9 @@ def get_parameter_structure(self, name: str) -> IndexingStructure: def get_variable_structure(self, name: str) -> IndexingStructure: return IndexingStructure(True, True) + def get_constraint_structure(self, name: str) -> IndexingStructure: + return IndexingStructure(True, True) + def test_shift() -> None: x = var("x") @@ -102,6 +105,9 @@ def get_parameter_structure(self, name: str) -> IndexingStructure: def get_variable_structure(self, name: str) -> IndexingStructure: return IndexingStructure(True, True) + def get_constraint_structure(self, name: str) -> IndexingStructure: + raise NotImplementedError() + provider = CustomStructureProvider() assert compute_indexation(expr, provider) == IndexingStructure(True, True) diff --git a/tests/unittests/simulation/test_simulation_table_accessor.py b/tests/unittests/simulation/test_simulation_table_accessor.py index 0769f364..ce0d953e 100644 --- a/tests/unittests/simulation/test_simulation_table_accessor.py +++ b/tests/unittests/simulation/test_simulation_table_accessor.py @@ -50,6 +50,12 @@ class FakeStudy: class FakeLinopyModel: solution: dict + @property + def dual(self) -> xr.Dataset: + return xr.Dataset() + + solver_model = None + @dataclass class FakeProblem: diff --git a/tests/unittests/simulation/test_simulation_table_export.py b/tests/unittests/simulation/test_simulation_table_export.py index 83f17631..75c3459d 100644 --- a/tests/unittests/simulation/test_simulation_table_export.py +++ b/tests/unittests/simulation/test_simulation_table_export.py @@ -50,6 +50,12 @@ class FakeStudy: class FakeLinopyModel: solution: dict + @property + def dual(self) -> xr.Dataset: + return xr.Dataset() + + solver_model = None + @dataclass class FakeProblem: diff --git a/tests/unittests/simulation/test_simulation_table_mock.py b/tests/unittests/simulation/test_simulation_table_mock.py index 0da96099..f1886a64 100644 --- a/tests/unittests/simulation/test_simulation_table_mock.py +++ b/tests/unittests/simulation/test_simulation_table_mock.py @@ -48,6 +48,12 @@ class FakeLinopyModel: solution: dict # lv.name -> xr.DataArray + @property + def dual(self) -> xr.Dataset: + return xr.Dataset() + + solver_model = None + @dataclass class FakeProblem: