diff --git a/CHANGELOG.md b/CHANGELOG.md index f59b76ec6..d16443eeb 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -30,12 +30,14 @@ Please remove all irrelevant sections before releasing. Until here --> ## [Unreleased] - ????-??-?? -Multi-Period and stochastic modeling is coming to flixopt in this release. +This Release brings Multi-year-investments and stochastic modeling to flixopt. +Further, IO methods were improved and resampling and selection of parts of the FlowSystem is now possible. +Several internal improvements were made to the codebase. -In this release, we introduce the following new features: -#### Multi-period-support + +#### Multi-year-investments A flixopt model might be modeled with a "year" dimension. -This enables to model transformation pathways over multiple years. +This enables to model transformation pathways over multiple years with several investment decisions #### Stochastic modeling A flixopt model can be modeled with a scenario dimension. @@ -67,17 +69,17 @@ The weighted sum of the total objective effect of each scenario is used as the o ### Added -* FlowSystem Restoring: The used FlowSystem will now get restired from the results (lazily). ALll Parameters can be safely acessed anytime after the solve. -* FLowResults added as a new class to store the results of Flows. They can now be accessed directly. +* FlowSystem Restoring: The used FlowSystem is now accessible directly form the results without manual restoring (lazily). All Parameters can be safely accessed anytime after the solve. +* FlowResults added as a new class to store the results of Flows. They can now be accessed directly. * Added precomputed DataArrays for `size`s, `flow_rate`s and `flow_hour`s. * Added `effects_per_component()`-Dataset to Results that stores the direct (and indirect) effects of each component. This greatly improves the evaluation of the impact of individual Components, even with many and complex effects. -* Improved filter methods for Results -* Balanced storage - Storage charging and discharging sizes can now be forced to be equal in when optimizing their size. +* Improved filter methods in `results.py` +* Balanced storage - Storage charging and discharging sizes can now be forced to be equal when optimizing their size by the `balanced` parameter. * Added Example for 2-stage Investment decisions leveraging the resampling of a FlowSystem -* New Storage Parameter: `relative_minimum_final_charge_state` and `relative_maximum_final_charge_state` parameter for final state control +* New Storage Parameter: `relative_minimum_final_charge_state` and `relative_maximum_final_charge_state` parameter for final state control. Default to last value of `relative_minimum_charge_state` and `relative_maximum_charge_state`, which will prevent change of behaviour for most users. ### Changed -* **BREAKING**: `relative_minimum_charge_state` and `relative_maximum_charge_state` don't have an extra timestep anymore. The final charge state can now be constrained by parameters `relative_minimum_final_charge_state` and `relative_maximum_final_charge_state` instead +* **BREAKING**: `relative_minimum_charge_state` and `relative_maximum_charge_state` don't have an extra timestep anymore. * **BREAKING**: Renamed class `SystemModel` to `FlowSystemModel` * **BREAKING**: Renamed class `Model` to `Submodel` * **BREAKING**: Renamed `mode` parameter in plotting methods to `style` @@ -87,7 +89,7 @@ The weighted sum of the total objective effect of each scenario is used as the o * Enhanced FlowSystem interface with improved `__repr__()` and `__str__()` methods * Improved Model Structure - Views and organisation is now divided into: * Model: The main Model (linopy.Model) that is used to create and store the variables and constraints for the flow_system. - * Submodel: The base class for all submodels. Each is a subset of the Model, for simpler acess and clearer code. + * Submodel: The base class for all submodels. Each is a subset of the Model, for simpler access and clearer code. ### Deprecated * The `agg_group` and `agg_weight` parameters of `TimeSeriesData` are deprecated and will be removed in a future version. Use `aggregation_group` and `aggregation_weight` instead. @@ -103,7 +105,7 @@ The weighted sum of the total objective effect of each scenario is used as the o * Better type consistency across all framework components ### Known issues -* IO for single Interfaces/Elemenets to Datasets might not work properly if the Interface/Element is not part of a fully transformed and connected FlowSystem. This arrises from Numeric Data not being stored as xr.DataArray by the user. To avoid this, always use the `to_dataset()` on Elements inside a FlowSystem thats connected and transformed. +* IO for single Interfaces/Elements to Datasets might not work properly if the Interface/Element is not part of a fully transformed and connected FlowSystem. This arises from Numeric Data not being stored as xr.DataArray by the user. To avoid this, always use the `to_dataset()` on Elements inside a FlowSystem that's connected and transformed. ### *Development* * **BREAKING**: Calculation.do_modeling() now returns the Calculation object instead of its linopy.Model diff --git a/flixopt/components.py b/flixopt/components.py index 34eac2cfe..4b6f3699e 100644 --- a/flixopt/components.py +++ b/flixopt/components.py @@ -204,13 +204,14 @@ def _plausibility_checks(self) -> None: f'(in {self.label_full}) and variable size is uncommon. Please check if this is intended!' ) - def transform_data(self, flow_system: FlowSystem): - super().transform_data(flow_system) + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: + super().transform_data(flow_system, name_prefix) if self.conversion_factors: self.conversion_factors = self._transform_conversion_factors(flow_system) if self.piecewise_conversion: self.piecewise_conversion.has_time_dim = True - self.piecewise_conversion.transform_data(flow_system, f'{self.label_full}|PiecewiseConversion') + prefix = '|'.join(filter(None, [name_prefix, self.label_full])) + self.piecewise_conversion.transform_data(flow_system, f'{prefix}|PiecewiseConversion') def _transform_conversion_factors(self, flow_system: FlowSystem) -> list[dict[str, xr.DataArray]]: """Converts all conversion factors to internal datatypes""" @@ -261,10 +262,12 @@ class Storage(Component): initial_charge_state: Storage charge state at the beginning of the time horizon. Can be numeric value or 'lastValueOfSim', which is recommended for if the initial start state is not known. Default is 0. - minimal_final_charge_state: Minimum absolute charge state required at the end - of the time horizon. Useful for ensuring energy security or meeting contracts. - maximal_final_charge_state: Maximum absolute charge state allowed at the end - of the time horizon. Useful for preventing overcharge or managing inventory. + minimal_final_charge_state: Minimum absolute charge state required at the end of the time horizon. + maximal_final_charge_state: Maximum absolute charge state allowed at the end of the time horizon. + relative_minimum_final_charge_state: Minimum relative charge state required at the end of the time horizon. + Defaults to the last value of the relative_minimum_charge_state. + relative_maximum_final_charge_state: Maximum relative charge state allowed at the end of the time horizon. + Defaults to the last value of the relative_maximum_charge_state. eta_charge: Charging efficiency factor (0-1 range). Accounts for conversion losses during charging. Default is 1 (perfect efficiency). eta_discharge: Discharging efficiency factor (0-1 range). Accounts for @@ -420,46 +423,47 @@ def create_model(self, model: FlowSystemModel) -> StorageModel: self.submodel = StorageModel(model, self) return self.submodel - def transform_data(self, flow_system: FlowSystem) -> None: - super().transform_data(flow_system) + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: + super().transform_data(flow_system, name_prefix) + base = '|'.join(filter(None, [name_prefix, self.label_full])) self.relative_minimum_charge_state = flow_system.fit_to_model_coords( - f'{self.label_full}|relative_minimum_charge_state', + f'{base}|relative_minimum_charge_state', self.relative_minimum_charge_state, ) self.relative_maximum_charge_state = flow_system.fit_to_model_coords( - f'{self.label_full}|relative_maximum_charge_state', + f'{base}|relative_maximum_charge_state', self.relative_maximum_charge_state, ) - self.eta_charge = flow_system.fit_to_model_coords(f'{self.label_full}|eta_charge', self.eta_charge) - self.eta_discharge = flow_system.fit_to_model_coords(f'{self.label_full}|eta_discharge', self.eta_discharge) + self.eta_charge = flow_system.fit_to_model_coords(f'{base}|eta_charge', self.eta_charge) + self.eta_discharge = flow_system.fit_to_model_coords(f'{base}|eta_discharge', self.eta_discharge) self.relative_loss_per_hour = flow_system.fit_to_model_coords( - f'{self.label_full}|relative_loss_per_hour', self.relative_loss_per_hour + f'{base}|relative_loss_per_hour', self.relative_loss_per_hour ) if not isinstance(self.initial_charge_state, str): self.initial_charge_state = flow_system.fit_to_model_coords( - f'{self.label_full}|initial_charge_state', self.initial_charge_state, dims=['year', 'scenario'] + f'{base}|initial_charge_state', self.initial_charge_state, dims=['year', 'scenario'] ) self.minimal_final_charge_state = flow_system.fit_to_model_coords( - f'{self.label_full}|minimal_final_charge_state', self.minimal_final_charge_state, dims=['year', 'scenario'] + f'{base}|minimal_final_charge_state', self.minimal_final_charge_state, dims=['year', 'scenario'] ) self.maximal_final_charge_state = flow_system.fit_to_model_coords( - f'{self.label_full}|maximal_final_charge_state', self.maximal_final_charge_state, dims=['year', 'scenario'] + f'{base}|maximal_final_charge_state', self.maximal_final_charge_state, dims=['year', 'scenario'] ) self.relative_minimum_final_charge_state = flow_system.fit_to_model_coords( - f'{self.label_full}|relative_minimum_final_charge_state', + f'{base}|relative_minimum_final_charge_state', self.relative_minimum_final_charge_state, dims=['year', 'scenario'], ) self.relative_maximum_final_charge_state = flow_system.fit_to_model_coords( - f'{self.label_full}|relative_maximum_final_charge_state', + f'{base}|relative_maximum_final_charge_state', self.relative_maximum_final_charge_state, dims=['year', 'scenario'], ) if isinstance(self.capacity_in_flow_hours, InvestParameters): - self.capacity_in_flow_hours.transform_data(flow_system, f'{self.label_full}|InvestParameters') + self.capacity_in_flow_hours.transform_data(flow_system, f'{base}|InvestParameters') else: self.capacity_in_flow_hours = flow_system.fit_to_model_coords( - f'{self.label_full}|capacity_in_flow_hours', self.capacity_in_flow_hours, dims=['year', 'scenario'] + f'{base}|capacity_in_flow_hours', self.capacity_in_flow_hours, dims=['year', 'scenario'] ) def _plausibility_checks(self) -> None: @@ -691,14 +695,11 @@ def create_model(self, model) -> TransmissionModel: self.submodel = TransmissionModel(model, self) return self.submodel - def transform_data(self, flow_system: FlowSystem) -> None: - super().transform_data(flow_system) - self.relative_losses = flow_system.fit_to_model_coords( - f'{self.label_full}|relative_losses', self.relative_losses - ) - self.absolute_losses = flow_system.fit_to_model_coords( - f'{self.label_full}|absolute_losses', self.absolute_losses - ) + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: + super().transform_data(flow_system, name_prefix) + base = '|'.join(filter(None, [name_prefix, self.label_full])) + self.relative_losses = flow_system.fit_to_model_coords(f'{base}|relative_losses', self.relative_losses) + self.absolute_losses = flow_system.fit_to_model_coords(f'{base}|absolute_losses', self.absolute_losses) class TransmissionModel(ComponentModel): diff --git a/flixopt/effects.py b/flixopt/effects.py index a0f67b873..f465dcc3b 100644 --- a/flixopt/effects.py +++ b/flixopt/effects.py @@ -16,7 +16,7 @@ import numpy as np import xarray as xr -from .core import Scalar, TemporalData, TemporalDataUser +from .core import NonTemporalDataUser, Scalar, TemporalData, TemporalDataUser from .features import ShareAllocationModel from .structure import Element, ElementModel, FlowSystemModel, Submodel, register_class_for_io @@ -174,41 +174,42 @@ def __init__( self.minimum_total = minimum_total self.maximum_total = maximum_total - def transform_data(self, flow_system: FlowSystem): + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: + base = '|'.join(filter(None, [name_prefix, self.label_full])) self.minimum_operation_per_hour = flow_system.fit_to_model_coords( - f'{self.label_full}|minimum_operation_per_hour', self.minimum_operation_per_hour + f'{base}|minimum_operation_per_hour', self.minimum_operation_per_hour ) self.maximum_operation_per_hour = flow_system.fit_to_model_coords( - f'{self.label_full}|maximum_operation_per_hour', self.maximum_operation_per_hour + f'{base}|maximum_operation_per_hour', self.maximum_operation_per_hour ) self.specific_share_to_other_effects_operation = flow_system.fit_effects_to_model_coords( - f'{self.label_full}|operation->', self.specific_share_to_other_effects_operation, 'operation' + f'{base}|operation->', self.specific_share_to_other_effects_operation, 'operation' ) self.minimum_operation = flow_system.fit_to_model_coords( - f'{self.label_full}|minimum_operation', self.minimum_operation, dims=['year', 'scenario'] + f'{base}|minimum_operation', self.minimum_operation, dims=['year', 'scenario'] ) self.maximum_operation = flow_system.fit_to_model_coords( - f'{self.label_full}|maximum_operation', self.maximum_operation, dims=['year', 'scenario'] + f'{base}|maximum_operation', self.maximum_operation, dims=['year', 'scenario'] ) self.minimum_invest = flow_system.fit_to_model_coords( - f'{self.label_full}|minimum_invest', self.minimum_invest, dims=['year', 'scenario'] + f'{base}|minimum_invest', self.minimum_invest, dims=['year', 'scenario'] ) self.maximum_invest = flow_system.fit_to_model_coords( - f'{self.label_full}|maximum_invest', self.maximum_invest, dims=['year', 'scenario'] + f'{base}|maximum_invest', self.maximum_invest, dims=['year', 'scenario'] ) self.minimum_total = flow_system.fit_to_model_coords( - f'{self.label_full}|minimum_total', + f'{base}|minimum_total', self.minimum_total, dims=['year', 'scenario'], ) self.maximum_total = flow_system.fit_to_model_coords( - f'{self.label_full}|maximum_total', self.maximum_total, dims=['year', 'scenario'] + f'{base}|maximum_total', self.maximum_total, dims=['year', 'scenario'] ) self.specific_share_to_other_effects_invest = flow_system.fit_effects_to_model_coords( - f'{self.label_full}|invest->', + f'{base}|invest->', self.specific_share_to_other_effects_invest, 'invest', dims=['year', 'scenario'], @@ -275,7 +276,7 @@ def _do_modeling(self): TemporalEffectsUser = TemporalDataUser | dict[str, TemporalDataUser] # User-specified Shares to Effects """ This datatype is used to define a temporal share to an effect by a certain attribute. """ -NonTemporalEffectsUser = Scalar | dict[str, Scalar] # User-specified Shares to Effects +NonTemporalEffectsUser = NonTemporalDataUser | dict[str, NonTemporalDataUser] # User-specified Shares to Effects """ This datatype is used to define a scalar share to an effect by a certain attribute. """ TemporalEffects = dict[str, TemporalData] # User-specified Shares to Effects diff --git a/flixopt/elements.py b/flixopt/elements.py index b592e90ba..e67f6bd8f 100644 --- a/flixopt/elements.py +++ b/flixopt/elements.py @@ -97,12 +97,13 @@ def create_model(self, model: FlowSystemModel) -> ComponentModel: self.submodel = ComponentModel(model, self) return self.submodel - def transform_data(self, flow_system: FlowSystem) -> None: + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: if self.on_off_parameters is not None: - self.on_off_parameters.transform_data(flow_system, self.label_full) + prefix = '|'.join(filter(None, [name_prefix, self.label_full])) + self.on_off_parameters.transform_data(flow_system, prefix) for flow in self.inputs + self.outputs: - flow.transform_data(flow_system) + flow.transform_data(flow_system, name_prefix) def _check_unique_flow_labels(self): all_flow_labels = [flow.label for flow in self.inputs + self.outputs] @@ -189,9 +190,10 @@ def create_model(self, model: FlowSystemModel) -> BusModel: self.submodel = BusModel(model, self) return self.submodel - def transform_data(self, flow_system: FlowSystem): + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: + base = '|'.join(filter(None, [name_prefix, self.label_full])) self.excess_penalty_per_flow_hour = flow_system.fit_to_model_coords( - f'{self.label_full}|excess_penalty_per_flow_hour', self.excess_penalty_per_flow_hour + f'{base}|excess_penalty_per_flow_hour', self.excess_penalty_per_flow_hour ) def _plausibility_checks(self) -> None: @@ -417,38 +419,35 @@ def create_model(self, model: FlowSystemModel) -> FlowModel: self.submodel = FlowModel(model, self) return self.submodel - def transform_data(self, flow_system: FlowSystem): - self.relative_minimum = flow_system.fit_to_model_coords( - f'{self.label_full}|relative_minimum', self.relative_minimum - ) - self.relative_maximum = flow_system.fit_to_model_coords( - f'{self.label_full}|relative_maximum', self.relative_maximum - ) + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: + base = '|'.join(filter(None, [name_prefix, self.label_full])) + self.relative_minimum = flow_system.fit_to_model_coords(f'{base}|relative_minimum', self.relative_minimum) + self.relative_maximum = flow_system.fit_to_model_coords(f'{base}|relative_maximum', self.relative_maximum) self.fixed_relative_profile = flow_system.fit_to_model_coords( - f'{self.label_full}|fixed_relative_profile', self.fixed_relative_profile + f'{base}|fixed_relative_profile', self.fixed_relative_profile ) self.effects_per_flow_hour = flow_system.fit_effects_to_model_coords( - self.label_full, self.effects_per_flow_hour, 'per_flow_hour' + base, self.effects_per_flow_hour, 'per_flow_hour' ) self.flow_hours_total_max = flow_system.fit_to_model_coords( - f'{self.label_full}|flow_hours_total_max', self.flow_hours_total_max, dims=['year', 'scenario'] + f'{base}|flow_hours_total_max', self.flow_hours_total_max, dims=['year', 'scenario'] ) self.flow_hours_total_min = flow_system.fit_to_model_coords( - f'{self.label_full}|flow_hours_total_min', self.flow_hours_total_min, dims=['year', 'scenario'] + f'{base}|flow_hours_total_min', self.flow_hours_total_min, dims=['year', 'scenario'] ) self.load_factor_max = flow_system.fit_to_model_coords( - f'{self.label_full}|load_factor_max', self.load_factor_max, dims=['year', 'scenario'] + f'{base}|load_factor_max', self.load_factor_max, dims=['year', 'scenario'] ) self.load_factor_min = flow_system.fit_to_model_coords( - f'{self.label_full}|load_factor_min', self.load_factor_min, dims=['year', 'scenario'] + f'{base}|load_factor_min', self.load_factor_min, dims=['year', 'scenario'] ) if self.on_off_parameters is not None: - self.on_off_parameters.transform_data(flow_system, self.label_full) + self.on_off_parameters.transform_data(flow_system, base) if isinstance(self.size, InvestParameters): - self.size.transform_data(flow_system, self.label_full) + self.size.transform_data(flow_system, base) else: - self.size = flow_system.fit_to_model_coords(f'{self.label_full}|size', self.size, dims=['year', 'scenario']) + self.size = flow_system.fit_to_model_coords(f'{base}|size', self.size, dims=['year', 'scenario']) def _plausibility_checks(self) -> None: # TODO: Incorporate into Variable? (Lower_bound can not be greater than upper bound @@ -687,7 +686,7 @@ def absolute_flow_rate_bounds(self) -> tuple[TemporalData, TemporalData]: if not self.with_investment: # Basic case without investment and without OnOff lb = lb_relative * self.element.size - elif not self.element.size.optional: + elif isinstance(self.element.size, InvestParameters) and not self.element.size.optional: # With non-optional Investment lb = lb_relative * self.element.size.minimum_or_fixed_size diff --git a/flixopt/features.py b/flixopt/features.py index 734e0c3ef..51a4832b7 100644 --- a/flixopt/features.py +++ b/flixopt/features.py @@ -11,7 +11,6 @@ import linopy import numpy as np -from .interface import InvestParameters, OnOffParameters, Piecewise from .modeling import BoundingPatterns, ModelingPrimitives, ModelingUtilities from .structure import FlowSystemModel, Submodel @@ -23,7 +22,19 @@ class InvestmentModel(Submodel): - """Investment model using factory patterns but keeping old interface""" + """ + This feature model is used to model the investment of a variable. + It applies the corresponding bounds to the variable and the on/off state of the variable. + + Args: + model: The optimization model instance + label_of_element: The label of the parent (Element). Used to construct the full label of the model. + parameters: The parameters of the feature model. + label_of_model: The label of the model. This is needed to construct the full label of the model. + + """ + + parameters: InvestParameters def __init__( self, @@ -32,17 +43,6 @@ def __init__( parameters: InvestParameters, label_of_model: str | None = None, ): - """ - This feature model is used to model the investment of a variable. - It applies the corresponding bounds to the variable and the on/off state of the variable. - - Args: - model: The optimization model instance - label_of_element: The label of the parent (Element). Used to construct the full label of the model. - parameters: The parameters of the feature model. - label_of_model: The label of the model. This is needed to construct the full label of the model. - - """ self.piecewise_effects: PiecewiseEffectsModel | None = None self.parameters = parameters super().__init__(model, label_of_element=label_of_element, label_of_model=label_of_model) @@ -210,6 +210,8 @@ def _do_modeling(self): short_name='consecutive_on_hours', minimum_duration=self.parameters.consecutive_on_hours_min, maximum_duration=self.parameters.consecutive_on_hours_max, + duration_per_step=self.hours_per_step, + duration_dim='time', previous_duration=self._get_previous_on_duration(), ) @@ -221,6 +223,8 @@ def _do_modeling(self): short_name='consecutive_off_hours', minimum_duration=self.parameters.consecutive_off_hours_min, maximum_duration=self.parameters.consecutive_off_hours_max, + duration_per_step=self.hours_per_step, + duration_dim='time', previous_duration=self._get_previous_off_duration(), ) # TODO: @@ -251,9 +255,9 @@ def _add_effects(self): # Properties access variables from Submodel's tracking system @property - def total_on_hours(self) -> linopy.Variable | None: + def on_hours_total(self) -> linopy.Variable: """Total on hours variable""" - return self['total_on_hours'] + return self['on_hours_total'] @property def off(self) -> linopy.Variable | None: diff --git a/flixopt/flow_system.py b/flixopt/flow_system.py index 10651b5ed..5a91dc36c 100644 --- a/flixopt/flow_system.py +++ b/flixopt/flow_system.py @@ -50,15 +50,16 @@ class FlowSystem(Interface): This is the main container class that users work with to build and manage their System. Args: - timesteps: The timesteps of the model. - years: The years of the model. - scenarios: The scenarios of the model. - hours_of_last_timestep: The duration of the last time step. Uses the last time interval if not specified - hours_of_previous_timesteps: The duration of previous timesteps. - If None, the first time increment of time_series is used. - This is needed to calculate previous durations (for example consecutive_on_hours). - If you use an array, take care that its long enough to cover all previous values! - weights: The weights of each year and scenario. If None, all have the same weight (normalized to 1). Its recommended to scale the weights to sum up to 1. + timesteps: The timesteps of the model. + years: The years of the model. + scenarios: The scenarios of the model. + hours_of_last_timestep: The duration of the last time step. Uses the last time interval if not specified + hours_of_previous_timesteps: The duration of previous timesteps. + If None, the first time increment of time_series is used. + This is needed to calculate previous durations (for example consecutive_on_hours). + If you use an array, take care that its long enough to cover all previous values! + years_of_last_year: The duration of the last year. Uses the last year interval if not specified + weights: The weights of each year and scenario. If None, all scenarios have the same weight, while the years have the weight of their represented year (all normalized to 1). Its recommended to scale the weights to sum up to 1. Notes: - Creates an empty registry for components and buses, an empty EffectCollection, and a placeholder for a SystemModel. @@ -72,16 +73,21 @@ def __init__( scenarios: pd.Index | None = None, hours_of_last_timestep: float | None = None, hours_of_previous_timesteps: int | float | np.ndarray | None = None, + years_of_last_year: int | None = None, weights: NonTemporalDataUser | None = None, ): self.timesteps = self._validate_timesteps(timesteps) - self.timesteps_extra = self._create_timesteps_with_extra(timesteps, hours_of_last_timestep) self.hours_of_previous_timesteps = self._calculate_hours_of_previous_timesteps( timesteps, hours_of_previous_timesteps ) - self.years = None if years is None else self._validate_years(years) + self.years_of_last_year = years_of_last_year + if years is None: + self.years, self.years_per_year = None, None + else: + self.years = self._validate_years(years) + self.years_per_year = self.calculate_years_per_year(self.years, years_of_last_year) self.scenarios = None if scenarios is None else self._validate_scenarios(scenarios) @@ -175,6 +181,17 @@ def calculate_hours_per_timestep(timesteps_extra: pd.DatetimeIndex) -> xr.DataAr hours_per_step, coords={'time': timesteps_extra[:-1]}, dims='time', name='hours_per_timestep' ) + @staticmethod + def calculate_years_per_year(years: pd.Index, years_of_last_year: int | None = None) -> xr.DataArray: + """Calculate duration of each timestep as a 1D DataArray.""" + years_per_year = np.diff(years) + return xr.DataArray( + np.append(years_per_year, years_of_last_year or years_per_year[-1]), + coords={'year': years}, + dims='year', + name='years_per_year', + ) + @staticmethod def _calculate_hours_of_previous_timesteps( timesteps: pd.DatetimeIndex, hours_of_previous_timesteps: float | np.ndarray | None @@ -263,6 +280,7 @@ def from_dataset(cls, ds: xr.Dataset) -> FlowSystem: timesteps=ds.indexes['time'], years=ds.indexes.get('year'), scenarios=ds.indexes.get('scenario'), + years_of_last_year=reference_structure.get('years_of_last_year'), weights=cls._resolve_dataarray_reference(reference_structure['weights'], arrays_dict) if 'weights' in reference_structure else None, @@ -347,7 +365,6 @@ def fit_to_model_coords( self, name: str, data: TemporalDataUser | NonTemporalDataUser | None, - has_time_dim: bool = True, dims: Collection[FlowSystemDimensions] | None = None, ) -> TemporalData | NonTemporalData | None: """ @@ -356,7 +373,6 @@ def fit_to_model_coords( Args: name: Name of the data data: Data to fit to model coordinates - has_time_dim: Wether to use the time dimension or not dims: Collection of dimension names to use for fitting. If None, all dimensions are used. Returns: @@ -365,18 +381,9 @@ def fit_to_model_coords( if data is None: return None - if dims is None: - coords = self.coords + coords = self.coords - if not has_time_dim: - warnings.warn( - 'has_time_dim is deprecated. Please pass dims to fit_to_model_coords instead.', - DeprecationWarning, - stacklevel=2, - ) - coords.pop('time') - else: - coords = self.coords + if dims is not None: coords = {k: coords[k] for k in dims if k in coords} # Rest of your method stays the same, just pass coords @@ -399,7 +406,6 @@ def fit_effects_to_model_coords( label_prefix: str | None, effect_values: TemporalEffectsUser | NonTemporalEffectsUser | None, label_suffix: str | None = None, - has_time_dim: bool = True, dims: Collection[FlowSystemDimensions] | None = None, ) -> TemporalEffects | NonTemporalEffects | None: """ @@ -414,7 +420,6 @@ def fit_effects_to_model_coords( effect: self.fit_to_model_coords( '|'.join(filter(None, [label_prefix, effect, label_suffix])), value, - has_time_dim=has_time_dim, dims=dims, ) for effect, value in effect_values_dict.items() @@ -772,6 +777,8 @@ def sel( if not self.connected_and_transformed: self.connect_and_transform() + ds = self.to_dataset() + # Build indexers dict from non-None parameters indexers = {} if time is not None: @@ -784,7 +791,7 @@ def sel( if not indexers: return self.copy() # Return a copy when no selection - selected_dataset = self.to_dataset().sel(**indexers) + selected_dataset = ds.sel(**indexers) return self.__class__.from_dataset(selected_dataset) def isel( @@ -807,6 +814,8 @@ def isel( if not self.connected_and_transformed: self.connect_and_transform() + ds = self.to_dataset() + # Build indexers dict from non-None parameters indexers = {} if time is not None: @@ -819,7 +828,7 @@ def isel( if not indexers: return self.copy() # Return a copy when no selection - selected_dataset = self.to_dataset().isel(**indexers) + selected_dataset = ds.isel(**indexers) return self.__class__.from_dataset(selected_dataset) def resample( diff --git a/flixopt/interface.py b/flixopt/interface.py index 416c6ad3a..fea55b17e 100644 --- a/flixopt/interface.py +++ b/flixopt/interface.py @@ -73,7 +73,7 @@ def __init__(self, start: TemporalDataUser, end: TemporalDataUser): self.end = end self.has_time_dim = False - def transform_data(self, flow_system: FlowSystem, name_prefix: str): + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: dims = None if self.has_time_dim else ['year', 'scenario'] self.start = flow_system.fit_to_model_coords(f'{name_prefix}|start', self.start, dims=dims) self.end = flow_system.fit_to_model_coords(f'{name_prefix}|end', self.end, dims=dims) @@ -219,7 +219,7 @@ def __getitem__(self, index) -> Piece: def __iter__(self) -> Iterator[Piece]: return iter(self.pieces) # Enables iteration like for piece in piecewise: ... - def transform_data(self, flow_system: FlowSystem, name_prefix: str): + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: for i, piece in enumerate(self.pieces): piece.transform_data(flow_system, f'{name_prefix}|Piece{i}') @@ -441,7 +441,7 @@ def items(self): """ return self.piecewises.items() - def transform_data(self, flow_system: FlowSystem, name_prefix: str): + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: for name, piecewise in self.piecewises.items(): piecewise.transform_data(flow_system, f'{name_prefix}|{name}') @@ -653,7 +653,7 @@ def has_time_dim(self, value): for piecewise in self.piecewise_shares.values(): piecewise.has_time_dim = value - def transform_data(self, flow_system: FlowSystem, name_prefix: str): + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: self.piecewise_origin.transform_data(flow_system, f'{name_prefix}|PiecewiseEffects|origin') for effect, piecewise in self.piecewise_shares.items(): piecewise.transform_data(flow_system, f'{name_prefix}|PiecewiseEffects|{effect}') @@ -861,7 +861,6 @@ def __init__( specific_effects: NonTemporalEffectsUser | None = None, # costs per Flow-Unit/Storage-Size/... piecewise_effects: PiecewiseEffects | None = None, divest_effects: NonTemporalEffectsUser | None = None, - investment_scenarios: Literal['individual'] | list[int | str] | None = None, ): self.fix_effects: NonTemporalEffectsUser = fix_effects or {} self.divest_effects: NonTemporalEffectsUser = divest_effects or {} @@ -871,10 +870,8 @@ def __init__( self.piecewise_effects = piecewise_effects self.minimum_size = minimum_size if minimum_size is not None else CONFIG.modeling.EPSILON self.maximum_size = maximum_size if maximum_size is not None else CONFIG.modeling.BIG # default maximum - self.investment_scenarios = investment_scenarios - def transform_data(self, flow_system: FlowSystem, name_prefix: str): - self._plausibility_checks(flow_system) + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: self.fix_effects = flow_system.fit_effects_to_model_coords( label_prefix=name_prefix, effect_values=self.fix_effects, @@ -908,21 +905,6 @@ def transform_data(self, flow_system: FlowSystem, name_prefix: str): f'{name_prefix}|fixed_size', self.fixed_size, dims=['year', 'scenario'] ) - def _plausibility_checks(self, flow_system): - if isinstance(self.investment_scenarios, list): - if not set(self.investment_scenarios).issubset(flow_system.scenarios): - raise ValueError( - f'Some scenarios in investment_scenarios are not present in the time_series_collection: ' - f'{set(self.investment_scenarios) - set(flow_system.scenarios)}' - ) - if self.investment_scenarios is not None: - if not self.optional: - if self.minimum_size is not None or self.fixed_size is not None: - logger.warning( - 'When using investment_scenarios, minimum_size and fixed_size should only ne used if optional is True.' - 'Otherwise the investment cannot be 0 incertain scenarios while being non-zero in others.' - ) - @property def minimum_or_fixed_size(self) -> NonTemporalData: return self.fixed_size if self.fixed_size is not None else self.minimum_size @@ -1138,7 +1120,7 @@ def __init__( self.switch_on_total_max: Scalar = switch_on_total_max self.force_switch_on: bool = force_switch_on - def transform_data(self, flow_system: FlowSystem, name_prefix: str): + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: self.effects_per_switch_on = flow_system.fit_effects_to_model_coords( name_prefix, self.effects_per_switch_on, 'per_switch_on' ) diff --git a/flixopt/modeling.py b/flixopt/modeling.py index b29ba7475..4e368d7ae 100644 --- a/flixopt/modeling.py +++ b/flixopt/modeling.py @@ -226,6 +226,8 @@ def consecutive_duration_tracking( short_name: str = None, minimum_duration: TemporalData | None = None, maximum_duration: TemporalData | None = None, + duration_dim: str = 'time', + duration_per_step: int | float | TemporalData = None, previous_duration: TemporalData = 0, ) -> tuple[linopy.Variable, tuple[linopy.Constraint, linopy.Constraint, linopy.Constraint]]: """ @@ -233,9 +235,9 @@ def consecutive_duration_tracking( Mathematical formulation: duration[t] ≤ state[t] * M ∀t - duration[t+1] ≤ duration[t] + hours_per_step[t] ∀t - duration[t+1] ≥ duration[t] + hours_per_step[t] + (state[t+1] - 1) * M ∀t - duration[0] = (hours_per_step[0] + previous_duration) * state[0] + duration[t+1] ≤ duration[t] + duration_per_step[t] ∀t + duration[t+1] ≥ duration[t] + duration_per_step[t] + (state[t+1] - 1) * M ∀t + duration[0] = (duration_per_step[0] + previous_duration) * state[0] If minimum_duration provided: duration[t] ≥ (state[t-1] - state[t]) * minimum_duration[t-1] ∀t > 0 @@ -254,14 +256,13 @@ def consecutive_duration_tracking( if not isinstance(model, Submodel): raise ValueError('ModelingPrimitives.sum_up_variable() can only be used with a Submodel') - hours_per_step = model.hours_per_step - mega = hours_per_step.sum('time') + previous_duration # Big-M value + mega = duration_per_step.sum(duration_dim) + previous_duration # Big-M value # Duration variable duration = model.add_variables( lower=0, upper=maximum_duration if maximum_duration is not None else mega, - coords=model.get_coords(), + coords=state_variable.coords, name=name, short_name=short_name, ) @@ -271,25 +272,26 @@ def consecutive_duration_tracking( # Upper bound: duration[t] ≤ state[t] * M constraints['ub'] = model.add_constraints(duration <= state_variable * mega, name=f'{duration.name}|ub') - # Forward constraint: duration[t+1] ≤ duration[t] + hours_per_step[t] + # Forward constraint: duration[t+1] ≤ duration[t] + duration_per_step[t] constraints['forward'] = model.add_constraints( - duration.isel(time=slice(1, None)) - <= duration.isel(time=slice(None, -1)) + hours_per_step.isel(time=slice(None, -1)), + duration.isel({duration_dim: slice(1, None)}) + <= duration.isel({duration_dim: slice(None, -1)}) + duration_per_step.isel({duration_dim: slice(None, -1)}), name=f'{duration.name}|forward', ) - # Backward constraint: duration[t+1] ≥ duration[t] + hours_per_step[t] + (state[t+1] - 1) * M + # Backward constraint: duration[t+1] ≥ duration[t] + duration_per_step[t] + (state[t+1] - 1) * M constraints['backward'] = model.add_constraints( - duration.isel(time=slice(1, None)) - >= duration.isel(time=slice(None, -1)) - + hours_per_step.isel(time=slice(None, -1)) - + (state_variable.isel(time=slice(1, None)) - 1) * mega, + duration.isel({duration_dim: slice(1, None)}) + >= duration.isel({duration_dim: slice(None, -1)}) + + duration_per_step.isel({duration_dim: slice(None, -1)}) + + (state_variable.isel({duration_dim: slice(1, None)}) - 1) * mega, name=f'{duration.name}|backward', ) - # Initial condition: duration[0] = (hours_per_step[0] + previous_duration) * state[0] + # Initial condition: duration[0] = (duration_per_step[0] + previous_duration) * state[0] constraints['initial'] = model.add_constraints( - duration.isel(time=0) == (hours_per_step.isel(time=0) + previous_duration) * state_variable.isel(time=0), + duration.isel({duration_dim: 0}) + == (duration_per_step.isel({duration_dim: 0}) + previous_duration) * state_variable.isel({duration_dim: 0}), name=f'{duration.name}|initial', ) @@ -297,15 +299,18 @@ def consecutive_duration_tracking( if minimum_duration is not None: constraints['lb'] = model.add_constraints( duration - >= (state_variable.isel(time=slice(None, -1)) - state_variable.isel(time=slice(1, None))) - * minimum_duration.isel(time=slice(None, -1)), + >= ( + state_variable.isel({duration_dim: slice(None, -1)}) + - state_variable.isel({duration_dim: slice(1, None)}) + ) + * minimum_duration.isel({duration_dim: slice(None, -1)}), name=f'{duration.name}|lb', ) # Handle initial condition for minimum duration - if previous_duration > 0 and previous_duration < minimum_duration.isel(time=0).max(): + if previous_duration > 0 and previous_duration < minimum_duration.isel({duration_dim: 0}).max(): constraints['initial_lb'] = model.add_constraints( - state_variable.isel(time=0) == 1, name=f'{duration.name}|initial_lb' + state_variable.isel({duration_dim: 0}) == 1, name=f'{duration.name}|initial_lb' ) variables = {'duration': duration} @@ -585,3 +590,146 @@ def state_transition_bounds( mutex = model.add_constraints(switch_on + switch_off <= 1, name=f'{name}|mutex') return transition, initial, mutex + + @staticmethod + def continuous_transition_bounds( + model: Submodel, + continuous_variable: linopy.Variable, + switch_on: linopy.Variable, + switch_off: linopy.Variable, + name: str, + max_change: float | xr.DataArray, + previous_value: float | xr.DataArray = 0.0, + coord: str = 'time', + ) -> tuple[linopy.Constraint, linopy.Constraint, linopy.Constraint, linopy.Constraint]: + """ + Constrains a continuous variable to only change when switch variables are active. + + Mathematical formulation: + -max_change * (switch_on[t] + switch_off[t]) <= continuous[t] - continuous[t-1] <= max_change * (switch_on[t] + switch_off[t]) ∀t > 0 + -max_change * (switch_on[0] + switch_off[0]) <= continuous[0] - previous_value <= max_change * (switch_on[0] + switch_off[0]) + switch_on[t], switch_off[t] ∈ {0, 1} + + This ensures the continuous variable can only change when switch_on or switch_off is 1. + When both switches are 0, the variable must stay exactly constant. + + Args: + model: The submodel to add constraints to + continuous_variable: The continuous variable to constrain + switch_on: Binary variable indicating when changes are allowed (typically transitions to active state) + switch_off: Binary variable indicating when changes are allowed (typically transitions to inactive state) + name: Base name for the constraints + max_change: Maximum possible change in the continuous variable (Big-M value) + previous_value: Initial value of the continuous variable before first period + coord: Coordinate name for time dimension + + Returns: + Tuple of constraints: (transition_upper, transition_lower, initial_upper, initial_lower) + """ + if not isinstance(model, Submodel): + raise ValueError('ModelingPrimitives.continuous_transition_bounds() can only be used with a Submodel') + + # Transition constraints for t > 0: continuous variable can only change when switches are active + transition_upper = model.add_constraints( + continuous_variable.isel({coord: slice(1, None)}) - continuous_variable.isel({coord: slice(None, -1)}) + <= max_change * (switch_on.isel({coord: slice(1, None)}) + switch_off.isel({coord: slice(1, None)})), + name=f'{name}|transition_ub', + ) + + transition_lower = model.add_constraints( + continuous_variable.isel({coord: slice(None, -1)}) - continuous_variable.isel({coord: slice(1, None)}) + <= max_change * (switch_on.isel({coord: slice(1, None)}) + switch_off.isel({coord: slice(1, None)})), + name=f'{name}|transition_lb', + ) + + # Initial constraints for t = 0 + initial_upper = model.add_constraints( + continuous_variable.isel({coord: 0}) - previous_value + <= max_change * (switch_on.isel({coord: 0}) + switch_off.isel({coord: 0})), + name=f'{name}|initial_ub', + ) + + initial_lower = model.add_constraints( + -continuous_variable.isel({coord: 0}) + previous_value + <= max_change * (switch_on.isel({coord: 0}) + switch_off.isel({coord: 0})), + name=f'{name}|initial_lb', + ) + + return transition_upper, transition_lower, initial_upper, initial_lower + + @staticmethod + def link_changes_to_level_with_binaries( + model: Submodel, + level_variable: linopy.Variable, + increase_variable: linopy.Variable, + decrease_variable: linopy.Variable, + increase_binary: linopy.Variable, + decrease_binary: linopy.Variable, + name: str, + max_change: float | xr.DataArray, + initial_level: float | xr.DataArray = 0.0, + coord: str = 'year', + ) -> tuple[linopy.Constraint, linopy.Constraint, linopy.Constraint, linopy.Constraint, linopy.Constraint]: + """ + Link changes to level evolution with binary control and mutual exclusivity. + + Creates the complete constraint system for ALL time periods: + 1. level[0] = initial_level + increase[0] - decrease[0] + 2. level[t] = level[t-1] + increase[t] - decrease[t] ∀t > 0 + 3. increase[t] <= max_change * increase_binary[t] ∀t + 4. decrease[t] <= max_change * decrease_binary[t] ∀t + 5. increase_binary[t] + decrease_binary[t] <= 1 ∀t + + Args: + model: The submodel to add constraints to + increase_variable: Incremental additions for ALL periods (>= 0) + decrease_variable: Incremental reductions for ALL periods (>= 0) + increase_binary: Binary indicators for increases for ALL periods + decrease_binary: Binary indicators for decreases for ALL periods + level_variable: Level variable for ALL periods + name: Base name for constraints + max_change: Maximum change per period + initial_level: Starting level before first period + coord: Time coordinate name + + Returns: + Tuple of (initial_constraint, transition_constraints, increase_bounds, decrease_bounds, mutual_exclusion) + """ + if not isinstance(model, Submodel): + raise ValueError('BoundingPatterns.link_changes_to_level_with_binaries() can only be used with a Submodel') + + # 1. Initial period: level[0] - initial_level = increase[0] - decrease[0] + initial_constraint = model.add_constraints( + level_variable.isel({coord: 0}) - initial_level + == increase_variable.isel({coord: 0}) - decrease_variable.isel({coord: 0}), + name=f'{name}|initial_level', + ) + + # 2. Transition periods: level[t] = level[t-1] + increase[t] - decrease[t] for t > 0 + transition_constraints = model.add_constraints( + level_variable.isel({coord: slice(1, None)}) + == level_variable.isel({coord: slice(None, -1)}) + + increase_variable.isel({coord: slice(1, None)}) + - decrease_variable.isel({coord: slice(1, None)}), + name=f'{name}|transitions', + ) + + # 3. Increase bounds: increase[t] <= max_change * increase_binary[t] for all t + increase_bounds = model.add_constraints( + increase_variable <= increase_binary * max_change, + name=f'{name}|increase_bounds', + ) + + # 4. Decrease bounds: decrease[t] <= max_change * decrease_binary[t] for all t + decrease_bounds = model.add_constraints( + decrease_variable <= decrease_binary * max_change, + name=f'{name}|decrease_bounds', + ) + + # 5. Mutual exclusivity: increase_binary[t] + decrease_binary[t] <= 1 for all t + mutual_exclusion = model.add_constraints( + increase_binary + decrease_binary <= 1, + name=f'{name}|mutual_exclusion', + ) + + return initial_constraint, transition_constraints, increase_bounds, decrease_bounds, mutual_exclusion diff --git a/flixopt/structure.py b/flixopt/structure.py index b728a711c..01529aba0 100644 --- a/flixopt/structure.py +++ b/flixopt/structure.py @@ -171,7 +171,11 @@ def get_coords( def weights(self) -> int | xr.DataArray: """Returns the scenario weights of the FlowSystem. If None, return weights that are normalized to 1 (one)""" if self.flow_system.weights is None: - weights = self.flow_system.fit_to_model_coords('weights', 1, dims=['year', 'scenario']) + weights = self.flow_system.fit_to_model_coords( + 'weights', + 1 if self.flow_system.years is None else self.flow_system.years_per_year, + dims=['year', 'scenario'], + ) return weights / weights.sum() @@ -218,11 +222,12 @@ class Interface: transform_data(flow_system): Transform data to match FlowSystem dimensions """ - def transform_data(self, flow_system: FlowSystem): + def transform_data(self, flow_system: FlowSystem, name_prefix: str = '') -> None: """Transform the data of the interface to match the FlowSystem's dimensions. Args: flow_system: The FlowSystem containing timing and dimensional information + name_prefix: The prefix to use for the names of the variables. Defaults to '', which results in no prefix. Raises: NotImplementedError: Must be implemented by subclasses diff --git a/tests/todos.txt b/tests/todos.txt deleted file mode 100644 index d4628c259..000000000 --- a/tests/todos.txt +++ /dev/null @@ -1,5 +0,0 @@ -# testing of - # abschnittsweise linear testen - # Komponenten mit offenen Flows - # Binärvariablen ohne max-Wert-Vorgabe des Flows (Binärungenauigkeitsproblem) - # Medien-zulässigkeit