diff --git a/examples/04_Scenarios/scenario_example.py b/examples/04_Scenarios/scenario_example.py index c68d1bbe5..3edb6e7c0 100644 --- a/examples/04_Scenarios/scenario_example.py +++ b/examples/04_Scenarios/scenario_example.py @@ -19,7 +19,7 @@ 'High Demand':[30, 0, 100, 118, 125, 20, 20, 20, 20]}, index=timesteps) power_prices = np.array([0.08, 0.09]) - flow_system = fx.FlowSystem(timesteps=timesteps, scenarios=scenarios) + flow_system = fx.FlowSystem(timesteps=timesteps, scenarios=scenarios, scenario_weights=np.array([0.5, 0.6])) # --- Define Energy Buses --- # These represent nodes, where the used medias are balanced (electricity, heat, and gas) diff --git a/flixopt/components.py b/flixopt/components.py index d471dcf6e..1b745b54c 100644 --- a/flixopt/components.py +++ b/flixopt/components.py @@ -211,7 +211,7 @@ def transform_data(self, flow_system: 'FlowSystem') -> None: self.relative_loss_per_hour = flow_system.create_time_series( f'{self.label_full}|relative_loss_per_hour', self.relative_loss_per_hour ) - if self.initial_charge_state != 'lastValueOfSim': + if not isinstance(self.initial_charge_state, str): self.initial_charge_state = flow_system.create_time_series( f'{self.label_full}|initial_charge_state', self.initial_charge_state, has_time_dim=False ) diff --git a/flixopt/core.py b/flixopt/core.py index e1c3c361d..850e01c04 100644 --- a/flixopt/core.py +++ b/flixopt/core.py @@ -653,7 +653,7 @@ def reset(self) -> None: Reset selections to include all timesteps and scenarios. This is equivalent to clearing all selections. """ - self.clear_selection() + self.set_selection(None, None) def restore_data(self) -> None: """ @@ -755,13 +755,7 @@ def update_stored_data(self, value: xr.DataArray) -> None: return self._stored_data = new_data - self.clear_selection() # Reset selections to full dataset - - def clear_selection(self, timesteps: bool = True, scenarios: bool = True) -> None: - if timesteps: - self._selected_timesteps = None - if scenarios: - self._selected_scenarios = None + self.set_selection(None, None) # Reset selections to full dataset def set_selection(self, timesteps: Optional[pd.DatetimeIndex] = None, scenarios: Optional[pd.Index] = None) -> None: """ @@ -773,15 +767,15 @@ def set_selection(self, timesteps: Optional[pd.DatetimeIndex] = None, scenarios: """ # Only update timesteps if the series has time dimension if self.has_time_dim: - if timesteps is None: - self.clear_selection(timesteps=True, scenarios=False) + if timesteps is None or timesteps.equals(self._stored_data.indexes['time']): + self._selected_timesteps = None else: self._selected_timesteps = timesteps # Only update scenarios if the series has scenario dimension if self.has_scenario_dim: - if scenarios is None: - self.clear_selection(timesteps=False, scenarios=True) + if scenarios is None or scenarios.equals(self._stored_data.indexes['scenario']): + self._selected_scenarios = None else: self._selected_scenarios = scenarios @@ -1077,22 +1071,6 @@ def add_time_series( # Return the TimeSeries object return time_series - def clear_selection(self, timesteps: bool = True, scenarios: bool = True) -> None: - """ - Clear selection for timesteps and/or scenarios. - - Args: - timesteps: Whether to clear timesteps selection - scenarios: Whether to clear scenarios selection - """ - if timesteps: - self._update_selected_timesteps(timesteps=None) - if scenarios: - self._selected_scenarios = None - - for ts in self._time_series.values(): - ts.clear_selection(timesteps=timesteps, scenarios=scenarios) - def set_selection(self, timesteps: Optional[pd.DatetimeIndex] = None, scenarios: Optional[pd.Index] = None) -> None: """ Set active subset for timesteps and scenarios. @@ -1102,35 +1080,30 @@ def set_selection(self, timesteps: Optional[pd.DatetimeIndex] = None, scenarios: scenarios: Scenarios to activate, or None to clear """ if timesteps is None: - self.clear_selection(timesteps=True, scenarios=False) + self._selected_timesteps = None + self._selected_timesteps_extra = None else: - self._update_selected_timesteps(timesteps) + self._selected_timesteps = self._validate_timesteps(timesteps, self._full_timesteps) + self._selected_timesteps_extra = self._create_timesteps_with_extra( + timesteps, self._calculate_hours_of_final_timestep(timesteps, self._full_timesteps) + ) if scenarios is None: - self.clear_selection(timesteps=False, scenarios=True) + self._selected_scenarios = None else: - self._selected_scenarios = self._validate_scenarios(scenarios) + self._selected_scenarios = self._validate_scenarios(scenarios, self._full_scenarios) - # Apply the selection to all TimeSeries objects - self._propagate_selection_to_time_series() + self._selected_hours_per_timestep = self.calculate_hours_per_timestep(self.timesteps_extra, self.scenarios) - def _update_selected_timesteps(self, timesteps: Optional[pd.DatetimeIndex]) -> None: - """ - Updates the timestep and related metrics (timesteps_extra, hours_per_timestep) based on the current selection. - """ - if timesteps is None: - self._selected_timesteps = None - self._selected_timesteps_extra = None - self._selected_hours_per_timestep = None - return + # Apply the selection to all TimeSeries objects + for ts_name, ts in self._time_series.items(): + if ts.has_time_dim: + timesteps = self.timesteps_extra if ts_name in self._has_extra_timestep else self.timesteps + else: + timesteps = None - self._selected_timesteps = self._validate_timesteps(timesteps, self._full_timesteps) - self._selected_timesteps_extra = self._create_timesteps_with_extra( - timesteps, self._calculate_hours_of_final_timestep(timesteps, self._full_timesteps) - ) - self._selected_hours_per_timestep = self.calculate_hours_per_timestep( - self._selected_timesteps_extra, self._selected_scenarios - ) + ts.set_selection(timesteps=timesteps, scenarios=self.scenarios if ts.has_scenario_dim else None) + self._propagate_selection_to_time_series() def as_dataset(self, with_extra_timestep: bool = True, with_constants: bool = True) -> xr.Dataset: """ @@ -1188,7 +1161,7 @@ def _propagate_selection_to_time_series(self) -> None: """Apply the current selection to all TimeSeries objects.""" for ts_name, ts in self._time_series.items(): if ts.has_time_dim: - timesteps = self.timesteps_extra if ts_name in self._has_extra_timestep else self.timesteps + timesteps = self.timesteps_extra if ts_name in self._has_extra_timestep else self.timesteps else: timesteps = None @@ -1482,3 +1455,28 @@ def get_numeric_stats(data: xr.DataArray, decimals: int = 2, padd: int = 10, by_ std = data.std().item() return f'{mean:{format_spec}} (mean), {median:{format_spec}} (median), {min_val:{format_spec}} (min), {max_val:{format_spec}} (max), {std:{format_spec}} (std)' + + +def extract_data( + data: Optional[Union[int, float, xr.DataArray, TimeSeries]], + if_none: Any = None +) -> Any: + """ + Convert data to xr.DataArray. + + Args: + data: The data to convert (scalar, array, or DataArray) + if_none: The value to return if data is None + + Returns: + DataArray with the converted data, or the value specified by if_none + """ + if data is None: + return if_none + if isinstance(data, TimeSeries): + return data.selected_data + if isinstance(data, xr.DataArray): + return data + if isinstance(data, (int, float, np.integer, np.floating)): + return xr.DataArray(data) + raise TypeError(f'Unsupported data type: {type(data).__name__}') \ No newline at end of file diff --git a/flixopt/effects.py b/flixopt/effects.py index 0cf165d66..2da561a36 100644 --- a/flixopt/effects.py +++ b/flixopt/effects.py @@ -12,7 +12,7 @@ import linopy import numpy as np -from .core import NumericDataTS, ScenarioData, TimeSeries, TimeSeriesCollection, TimestepData +from .core import NumericDataTS, ScenarioData, TimeSeries, TimeSeriesCollection, TimestepData, extract_data from .features import ShareAllocationModel from .structure import Element, ElementModel, Interface, Model, SystemModel, register_class_for_io @@ -125,8 +125,8 @@ def __init__(self, model: SystemModel, element: Effect): label_of_element=self.label_of_element, label='invest', label_full=f'{self.label_full}(invest)', - total_max=self.element.maximum_invest, - total_min=self.element.minimum_invest, + total_max=extract_data(self.element.maximum_invest), + total_min=extract_data(self.element.minimum_invest), ) ) @@ -138,14 +138,10 @@ def __init__(self, model: SystemModel, element: Effect): label_of_element=self.label_of_element, label='operation', label_full=f'{self.label_full}(operation)', - total_max=self.element.maximum_operation, - total_min=self.element.minimum_operation, - min_per_hour=self.element.minimum_operation_per_hour.selected_data - if self.element.minimum_operation_per_hour is not None - else None, - max_per_hour=self.element.maximum_operation_per_hour.selected_data - if self.element.maximum_operation_per_hour is not None - else None, + total_max=extract_data(self.element.maximum_operation), + total_min=extract_data(self.element.minimum_operation), + min_per_hour=extract_data(self.element.minimum_operation_per_hour), + max_per_hour=extract_data(self.element.maximum_operation_per_hour), ) ) @@ -155,8 +151,8 @@ def do_modeling(self): self.total = self.add( self._model.add_variables( - lower=self.element.minimum_total if self.element.minimum_total is not None else -np.inf, - upper=self.element.maximum_total if self.element.maximum_total is not None else np.inf, + lower=extract_data(self.element.minimum_total, if_none=-np.inf), + upper=extract_data(self.element.maximum_total, if_none=np.inf), coords=self._model.get_coords(time_dim=False), name=f'{self.label_full}|total', ), diff --git a/flixopt/elements.py b/flixopt/elements.py index 80085cd0c..2ff49567e 100644 --- a/flixopt/elements.py +++ b/flixopt/elements.py @@ -10,7 +10,7 @@ import numpy as np from .config import CONFIG -from .core import PlausibilityError, Scalar, ScenarioData, TimestepData +from .core import PlausibilityError, Scalar, ScenarioData, TimestepData, extract_data from .effects import EffectValuesUserTimestep from .features import InvestmentModel, OnOffModel, PreventSimultaneousUsageModel from .interface import InvestParameters, OnOffParameters @@ -375,8 +375,8 @@ def do_modeling(self): self.total_flow_hours = self.add( self._model.add_variables( - lower=self.element.flow_hours_total_min if self.element.flow_hours_total_min is not None else 0, - upper=self.element.flow_hours_total_max if self.element.flow_hours_total_max is not None else np.inf, + lower=extract_data(self.element.flow_hours_total_min, 0), + upper=extract_data(self.element.flow_hours_total_max, np.inf), coords=self._model.get_coords(time_dim=False), name=f'{self.label_full}|total_flow_hours', ), @@ -456,16 +456,16 @@ def flow_rate_lower_bound_relative(self) -> TimestepData: """Returns the lower bound of the flow_rate relative to its size""" fixed_profile = self.element.fixed_relative_profile if fixed_profile is None: - return self.element.relative_minimum.selected_data - return fixed_profile.selected_data + return extract_data(self.element.relative_minimum) + return extract_data(fixed_profile) @property def flow_rate_upper_bound_relative(self) -> TimestepData: """ Returns the upper bound of the flow_rate relative to its size""" fixed_profile = self.element.fixed_relative_profile if fixed_profile is None: - return self.element.relative_maximum.selected_data - return fixed_profile.selected_data + return extract_data(self.element.relative_maximum) + return extract_data(fixed_profile) @property def flow_rate_lower_bound(self) -> TimestepData: @@ -478,8 +478,8 @@ def flow_rate_lower_bound(self) -> TimestepData: if isinstance(self.element.size, InvestParameters): if self.element.size.optional: return 0 - return self.flow_rate_lower_bound_relative * self.element.size.minimum_size - return self.flow_rate_lower_bound_relative * self.element.size + return self.flow_rate_lower_bound_relative * extract_data(self.element.size.minimum_size) + return self.flow_rate_lower_bound_relative * extract_data(self.element.size) @property def flow_rate_upper_bound(self) -> TimestepData: @@ -488,8 +488,8 @@ def flow_rate_upper_bound(self) -> TimestepData: Further constraining might be done in OnOffModel and InvestmentModel """ if isinstance(self.element.size, InvestParameters): - return self.flow_rate_upper_bound_relative * self.element.size.maximum_size - return self.flow_rate_upper_bound_relative * self.element.size + return self.flow_rate_upper_bound_relative * extract_data(self.element.size.maximum_size) + return self.flow_rate_upper_bound_relative * extract_data(self.element.size) class BusModel(ElementModel): diff --git a/flixopt/features.py b/flixopt/features.py index 31c1dbadb..b8243d794 100644 --- a/flixopt/features.py +++ b/flixopt/features.py @@ -10,7 +10,7 @@ import numpy as np from .config import CONFIG -from .core import Scalar, ScenarioData, TimeSeries, TimestepData +from .core import Scalar, ScenarioData, TimeSeries, TimestepData, extract_data from .interface import InvestParameters, OnOffParameters, Piecewise from .structure import Model, SystemModel @@ -45,8 +45,8 @@ def __init__( def do_modeling(self): self.size = self.add( self._model.add_variables( - lower=0 if self.parameters.optional else self.parameters.minimum_size*1, - upper=self.parameters.maximum_size*1, + lower=0 if self.parameters.optional else extract_data(self.parameters.minimum_size), + upper=extract_data(self.parameters.maximum_size), name=f'{self.label_full}|size', coords=self._model.get_coords(time_dim=False), ), @@ -295,8 +295,8 @@ def do_modeling(self): self.total_on_hours = self.add( self._model.add_variables( - lower=self._on_hours_total_min, - upper=self._on_hours_total_max, + lower=extract_data(self._on_hours_total_min), + upper=extract_data(self._on_hours_total_max), coords=self._model.get_coords(time_dim=False), name=f'{self.label_full}|on_hours_total', ), @@ -440,7 +440,7 @@ def do_modeling(self): # Create count variable for number of switches self.switch_on_nr = self.add( self._model.add_variables( - upper=self._switch_on_max, + upper=extract_data(self._switch_on_max), lower=0, name=f'{self.label_full}|switch_on_nr', ), @@ -534,7 +534,7 @@ def do_modeling(self): self.duration = self.add( self._model.add_variables( lower=0, - upper=self._maximum_duration if self._maximum_duration is not None else mega, + upper=extract_data(self._maximum_duration, mega), coords=self._model.get_coords(), name=f'{self.label_full}|hours', ), @@ -588,7 +588,7 @@ def do_modeling(self): ) # Handle initial condition - if 0 < self.previous_duration < self._minimum_duration.isel(time=0): + if 0 < self.previous_duration < self._minimum_duration.isel(time=0).max(): self.add( self._model.add_constraints( self._state_variable.isel(time=0) == 1, name=f'{self.label_full}|initial_minimum' @@ -613,7 +613,7 @@ def previous_duration(self) -> Scalar: """Computes the previous duration of the state variable""" #TODO: Allow for other/dynamic timestep resolutions return ConsecutiveStateModel.compute_consecutive_hours_in_state( - self._previous_states, self._model.hours_per_step.isel(time=0).item() + self._previous_states, self._model.hours_per_step.isel(time=0).values.flatten()[0] ) @staticmethod @@ -715,8 +715,8 @@ def do_modeling(self): defining_bounds=self._defining_bounds, previous_values=self._previous_values, use_off=self.parameters.use_off, - on_hours_total_min=self.parameters.on_hours_total_min, - on_hours_total_max=self.parameters.on_hours_total_max, + on_hours_total_min=extract_data(self.parameters.on_hours_total_min), + on_hours_total_max=extract_data(self.parameters.on_hours_total_max), effects_per_running_hour=self.parameters.effects_per_running_hour, ) self.add(self.state_model) @@ -965,8 +965,8 @@ def __init__( label_of_element: Optional[str] = None, label: Optional[str] = None, label_full: Optional[str] = None, - total_max: Optional[Scalar] = None, - total_min: Optional[Scalar] = None, + total_max: Optional[ScenarioData] = None, + total_min: Optional[ScenarioData] = None, max_per_hour: Optional[TimestepData] = None, min_per_hour: Optional[TimestepData] = None, ): @@ -1009,8 +1009,8 @@ def do_modeling(self): if self._has_time_dim: self.total_per_timestep = self.add( self._model.add_variables( - lower=-np.inf if (self._min_per_hour is None) else self._min_per_hour * self._model.hours_per_step, - upper=np.inf if (self._max_per_hour is None) else self._max_per_hour * self._model.hours_per_step, + lower=-np.inf if (self._min_per_hour is None) else extract_data(self._min_per_hour) * self._model.hours_per_step, + upper=np.inf if (self._max_per_hour is None) else extract_data(self._max_per_hour) * self._model.hours_per_step, coords=self._model.get_coords(time_dim=True, scenario_dim=self._has_scenario_dim), name=f'{self.label_full}|total_per_timestep', ), diff --git a/tests/test_scenarios.py b/tests/test_scenarios.py new file mode 100644 index 000000000..5b9105a68 --- /dev/null +++ b/tests/test_scenarios.py @@ -0,0 +1,333 @@ +import matplotlib.pyplot as plt +import numpy as np +import pandas as pd +from linopy.testing import assert_linequal +import xarray as xr +import pytest +import flixopt as fx + +from flixopt.commons import Effect, FullCalculation, InvestParameters, Sink, Source, Storage, TimeSeriesData, solvers +from flixopt.elements import Bus, Flow +from flixopt.flow_system import FlowSystem + +from .conftest import create_linopy_model, create_calculation_and_solve + + +@pytest.fixture +def test_system(): + """Create a basic test system with scenarios.""" + # Create a two-day time index with hourly resolution + timesteps = pd.date_range( + "2023-01-01", periods=48, freq="h", name="time" + ) + + # Create two scenarios + scenarios = pd.Index(["Scenario A", "Scenario B"], name="scenario") + + # Create scenario weights as TimeSeriesData + # Using TimeSeriesData to avoid conversion issues + scenario_weights = TimeSeriesData(np.array([0.7, 0.3])) + + # Create a flow system with scenarios + flow_system = FlowSystem( + timesteps=timesteps, + scenarios=scenarios, + scenario_weights=scenario_weights # Use TimeSeriesData for weights + ) + + # Create demand profiles that differ between scenarios + # Scenario A: Higher demand in first day, lower in second day + # Scenario B: Lower demand in first day, higher in second day + demand_profile_a = np.concatenate([ + np.sin(np.linspace(0, 2*np.pi, 24)) * 5 + 10, # Day 1, max ~15 + np.sin(np.linspace(0, 2*np.pi, 24)) * 2 + 5 # Day 2, max ~7 + ]) + + demand_profile_b = np.concatenate([ + np.sin(np.linspace(0, 2*np.pi, 24)) * 2 + 5, # Day 1, max ~7 + np.sin(np.linspace(0, 2*np.pi, 24)) * 5 + 10 # Day 2, max ~15 + ]) + + # Stack the profiles into a 2D array (time, scenario) + demand_profiles = np.column_stack([demand_profile_a, demand_profile_b]) + + # Create the necessary model elements + # Create buses + electricity_bus = Bus("Electricity") + + # Create a demand sink with scenario-dependent profiles + demand = Flow( + label="Demand", + bus=electricity_bus.label_full, + fixed_relative_profile=demand_profiles + ) + demand_sink = Sink("Demand", sink=demand) + + # Create a power source with investment option + power_gen = Flow( + label="Generation", + bus=electricity_bus.label_full, + size=InvestParameters( + minimum_size=0, + maximum_size=20, + specific_effects={"Costs": 100} # €/kW + ), + effects_per_flow_hour={"Costs": 20} # €/MWh + ) + generator = Source("Generator", source=power_gen) + + # Create a storage for electricity + storage_charge = Flow( + label="Charge", + bus=electricity_bus.label_full, + size=10 + ) + storage_discharge = Flow( + label="Discharge", + bus=electricity_bus.label_full, + size=10 + ) + storage = Storage( + label="Battery", + charging=storage_charge, + discharging=storage_discharge, + capacity_in_flow_hours=InvestParameters( + minimum_size=0, + maximum_size=50, + specific_effects={"Costs": 50} # €/kWh + ), + eta_charge=0.95, + eta_discharge=0.95, + initial_charge_state="lastValueOfSim" + ) + + # Create effects and objective + cost_effect = Effect( + label="Costs", + unit="€", + description="Total costs", + is_standard=True, + is_objective=True + ) + + # Add all elements to the flow system + flow_system.add_elements( + electricity_bus, + generator, + demand_sink, + storage, + cost_effect + ) + + # Return the created system and its components + return { + "flow_system": flow_system, + "timesteps": timesteps, + "scenarios": scenarios, + "electricity_bus": electricity_bus, + "demand": demand, + "demand_sink": demand_sink, + "generator": generator, + "power_gen": power_gen, + "storage": storage, + "storage_charge": storage_charge, + "storage_discharge": storage_discharge, + "cost_effect": cost_effect + } + +@pytest.fixture +def flow_system_complex_scenarios() -> fx.FlowSystem: + """ + Helper method to create a base model with configurable parameters + """ + thermal_load = np.array([30, 0, 90, 110, 110, 20, 20, 20, 20]) + electrical_load = np.array([40, 40, 40, 40, 40, 40, 40, 40, 40]) + flow_system = fx.FlowSystem(pd.date_range('2020-01-01', periods=9, freq='h', name='time'), + pd.Index(['A', 'B', 'C'], name='scenario')) + # Define the components and flow_system + flow_system.add_elements( + fx.Effect('costs', '€', 'Kosten', is_standard=True, is_objective=True), + fx.Effect('CO2', 'kg', 'CO2_e-Emissionen', specific_share_to_other_effects_operation={'costs': 0.2}), + fx.Effect('PE', 'kWh_PE', 'Primärenergie', maximum_total=3.5e3), + fx.Bus('Strom'), + fx.Bus('Fernwärme'), + fx.Bus('Gas'), + fx.Sink('Wärmelast', sink=fx.Flow('Q_th_Last', 'Fernwärme', size=1, fixed_relative_profile=thermal_load)), + fx.Source( + 'Gastarif', source=fx.Flow('Q_Gas', 'Gas', size=1000, effects_per_flow_hour={'costs': 0.04, 'CO2': 0.3}) + ), + fx.Sink('Einspeisung', sink=fx.Flow('P_el', 'Strom', effects_per_flow_hour=-1 * electrical_load)), + ) + + boiler = fx.linear_converters.Boiler( + 'Kessel', + eta=0.5, + on_off_parameters=fx.OnOffParameters(effects_per_running_hour={'costs': 0, 'CO2': 1000}), + Q_th=fx.Flow( + 'Q_th', + bus='Fernwärme', + load_factor_max=1.0, + load_factor_min=0.1, + relative_minimum=5 / 50, + relative_maximum=1, + previous_flow_rate=50, + size=fx.InvestParameters( + fix_effects=1000, fixed_size=50, optional=False, specific_effects={'costs': 10, 'PE': 2} + ), + on_off_parameters=fx.OnOffParameters( + on_hours_total_min=0, + on_hours_total_max=1000, + consecutive_on_hours_max=10, + consecutive_on_hours_min=1, + consecutive_off_hours_max=10, + effects_per_switch_on=0.01, + switch_on_total_max=1000, + ), + flow_hours_total_max=1e6, + ), + Q_fu=fx.Flow('Q_fu', bus='Gas', size=200, relative_minimum=0, relative_maximum=1), + ) + + invest_speicher = fx.InvestParameters( + fix_effects=0, + piecewise_effects=fx.PiecewiseEffects( + piecewise_origin=fx.Piecewise([fx.Piece(5, 25), fx.Piece(25, 100)]), + piecewise_shares={ + 'costs': fx.Piecewise([fx.Piece(50, 250), fx.Piece(250, 800)]), + 'PE': fx.Piecewise([fx.Piece(5, 25), fx.Piece(25, 100)]), + }, + ), + optional=False, + specific_effects={'costs': 0.01, 'CO2': 0.01}, + minimum_size=0, + maximum_size=1000, + ) + speicher = fx.Storage( + 'Speicher', + charging=fx.Flow('Q_th_load', bus='Fernwärme', size=1e4), + discharging=fx.Flow('Q_th_unload', bus='Fernwärme', size=1e4), + capacity_in_flow_hours=invest_speicher, + initial_charge_state=0, + maximal_final_charge_state=10, + eta_charge=0.9, + eta_discharge=1, + relative_loss_per_hour=0.08, + prevent_simultaneous_charge_and_discharge=True, + ) + + flow_system.add_elements(boiler, speicher) + + return flow_system + + +@pytest.fixture +def flow_system_piecewise_conversion_scenarios(flow_system_complex_scenarios) -> fx.FlowSystem: + """ + Use segments/Piecewise with numeric data + """ + flow_system = flow_system_complex_scenarios + + flow_system.add_elements( + fx.LinearConverter( + 'KWK', + inputs=[fx.Flow('Q_fu', bus='Gas')], + outputs=[ + fx.Flow('P_el', bus='Strom', size=60, relative_maximum=55, previous_flow_rate=10), + fx.Flow('Q_th', bus='Fernwärme'), + ], + piecewise_conversion=fx.PiecewiseConversion( + { + 'P_el': fx.Piecewise( + [ + fx.Piece(np.linspace(5, 6, len(flow_system.time_series_collection.timesteps)), 30), + fx.Piece(40, np.linspace(60, 70, len(flow_system.time_series_collection.timesteps))), + ] + ), + 'Q_th': fx.Piecewise([fx.Piece(6, 35), fx.Piece(45, 100)]), + 'Q_fu': fx.Piecewise([fx.Piece(12, 70), fx.Piece(90, 200)]), + } + ), + on_off_parameters=fx.OnOffParameters(effects_per_switch_on=0.01), + ) + ) + + return flow_system + + +def test_scenario_weights(flow_system_piecewise_conversion_scenarios): + """Test that scenario weights are correctly used in the model.""" + scenarios = flow_system_piecewise_conversion_scenarios.time_series_collection.scenarios + weights = np.linspace(0.5, 1, len(scenarios)) / np.sum(np.linspace(0.5, 1, len(scenarios))) + flow_system_piecewise_conversion_scenarios.scenario_weights = weights + model = create_linopy_model(flow_system_piecewise_conversion_scenarios) + np.testing.assert_allclose(model.scenario_weights.values, weights) + assert_linequal(model.objective.expression, + (model.variables['costs|total'] * weights).sum() + model.variables['Penalty|total']) + assert np.isclose(model.scenario_weights.sum().item(), 1.0) + +def test_scenario_dimensions_in_variables(flow_system_piecewise_conversion_scenarios): + """Test that all time variables are correctly broadcasted to scenario dimensions.""" + model = create_linopy_model(flow_system_piecewise_conversion_scenarios) + for var in model.variables: + assert model.variables[var].dims in [('time', 'scenario'), ('scenario',), ()] + +def test_full_scenario_optimization(flow_system_piecewise_conversion_scenarios): + """Test a full optimization with scenarios and verify results.""" + scenarios = flow_system_piecewise_conversion_scenarios.time_series_collection.scenarios + weights = np.linspace(0.5, 1, len(scenarios)) / np.sum(np.linspace(0.5, 1, len(scenarios))) + flow_system_piecewise_conversion_scenarios.scenario_weights = weights + calc = create_calculation_and_solve(flow_system_piecewise_conversion_scenarios, + solver=fx.solvers.GurobiSolver(mip_gap=0.01, time_limit_seconds=60), + name='test_full_scenario') + calc.results.to_file() + + res = fx.results.CalculationResults.from_file('results', 'test_full_scenario') + fx.FlowSystem.from_dataset(res.flow_system) + calc = create_calculation_and_solve( + flow_system_piecewise_conversion_scenarios, + solver=fx.solvers.GurobiSolver(mip_gap=0.01, time_limit_seconds=60), + name='test_full_scenario', + ) + +@pytest.mark.slow +def test_io_persistance(flow_system_piecewise_conversion_scenarios): + """Test a full optimization with scenarios and verify results.""" + scenarios = flow_system_piecewise_conversion_scenarios.time_series_collection.scenarios + weights = np.linspace(0.5, 1, len(scenarios)) / np.sum(np.linspace(0.5, 1, len(scenarios))) + flow_system_piecewise_conversion_scenarios.scenario_weights = weights + calc = create_calculation_and_solve(flow_system_piecewise_conversion_scenarios, + solver=fx.solvers.HighsSolver(mip_gap=0.001, time_limit_seconds=60), + name='test_full_scenario') + calc.results.to_file() + + res = fx.results.CalculationResults.from_file('results', 'test_full_scenario') + flow_system_2 = fx.FlowSystem.from_dataset(res.flow_system) + calc_2 = create_calculation_and_solve( + flow_system_2, + solver=fx.solvers.HighsSolver(mip_gap=0.001, time_limit_seconds=60), + name='test_full_scenario_2', + ) + + np.testing.assert_allclose(calc.results.objective, calc_2.results.objective, rtol=0.001) + + +def test_scenarios_selection(flow_system_piecewise_conversion_scenarios): + flow_system = flow_system_piecewise_conversion_scenarios + scenarios = flow_system_piecewise_conversion_scenarios.time_series_collection.scenarios + weights = np.linspace(0.5, 1, len(scenarios)) / np.sum(np.linspace(0.5, 1, len(scenarios))) + flow_system_piecewise_conversion_scenarios.scenario_weights = weights + calc = fx.FullCalculation(flow_system=flow_system_piecewise_conversion_scenarios, + selected_scenarios=flow_system.time_series_collection.scenarios[0:2], + name='test_full_scenario') + calc.do_modeling() + calc.solve(fx.solvers.GurobiSolver(mip_gap=0.01, time_limit_seconds=60)) + + calc.results.to_file() + flow_system_2 = fx.FlowSystem.from_dataset(calc.results.flow_system) + + assert calc.results.solution.indexes['scenario'].equals(flow_system.time_series_collection.scenarios[0:2]) + + assert flow_system_2.time_series_collection.scenarios.equals(flow_system.time_series_collection.scenarios[0:2]) + + np.testing.assert_allclose(flow_system_2.scenario_weights.selected_data.values, weights[0:2]) + diff --git a/tests/test_timeseries.py b/tests/test_timeseries.py index e2432e784..237935e59 100644 --- a/tests/test_timeseries.py +++ b/tests/test_timeseries.py @@ -91,7 +91,7 @@ def test_selection_methods(self, sample_timeseries, sample_timesteps): assert sample_timeseries.selected_data.equals(sample_timeseries.stored_data.sel(time=subset_index)) # Clear selection - sample_timeseries.clear_selection() + sample_timeseries.set_selection() assert sample_timeseries._selected_timesteps is None assert sample_timeseries.selected_data.equals(sample_timeseries.stored_data) @@ -408,7 +408,7 @@ def test_scenario_selection(self, sample_scenario_timeseries, sample_scenario_in ) # Clear selection - sample_scenario_timeseries.clear_selection(timesteps=False, scenarios=True) + sample_scenario_timeseries.set_selection() assert sample_scenario_timeseries._selected_scenarios is None def test_all_equal_with_scenarios(self, sample_timesteps, sample_scenario_index): @@ -561,7 +561,7 @@ def test_selection_propagation(self, sample_allocator, sample_timesteps): assert len(ts3._selected_timesteps) == len(subset_timesteps) + 1 # Clear selection - sample_allocator.clear_selection() + sample_allocator.set_selection() # Check selection cleared assert ts1._selected_timesteps is None @@ -673,7 +673,7 @@ def test_selection_propagation_with_scenarios( assert ts1.selected_data.shape == (2, 2) # 2 timesteps, 2 scenarios # Clear selections - sample_scenario_allocator.clear_selection() + sample_scenario_allocator.set_selection() assert ts1._selected_timesteps is None assert ts1.active_timesteps.equals(sample_scenario_allocator.timesteps) assert ts1._selected_scenarios is None