muse2/simulation/optimisation/
constraints.rs

1//! Code for adding constraints to the dispatch optimisation problem.
2use super::VariableMap;
3use crate::asset::{AssetIterator, AssetRef};
4use crate::commodity::{CommodityID, CommodityType};
5use crate::model::Model;
6use crate::region::RegionID;
7use crate::time_slice::{TimeSliceID, TimeSliceSelection};
8use highs::RowProblem as Problem;
9
10/// Corresponding variables for a constraint along with the row offset in the solution
11pub struct KeysWithOffset<T> {
12    offset: usize,
13    keys: Vec<T>,
14}
15
16impl<T> KeysWithOffset<T> {
17    /// Zip the keys with the corresponding dual values in the solution, accounting for the offset
18    pub fn zip_duals<'a>(&'a self, duals: &'a [f64]) -> impl Iterator<Item = (&'a T, f64)> {
19        assert!(
20            self.offset + self.keys.len() <= duals.len(),
21            "Bad constraint keys: dual rows out of range"
22        );
23
24        self.keys.iter().zip(duals[self.offset..].iter().copied())
25    }
26}
27
28/// Indicates the commodity ID and time slice selection covered by each commodity balance constraint
29pub type CommodityBalanceKeys = KeysWithOffset<(CommodityID, RegionID, TimeSliceSelection)>;
30
31/// Indicates the asset ID and time slice covered by each activity constraint
32pub type ActivityKeys = KeysWithOffset<(AssetRef, TimeSliceID)>;
33
34/// The keys for different constraints
35pub struct ConstraintKeys {
36    /// Keys for commodity balance constraints
37    pub commodity_balance_keys: CommodityBalanceKeys,
38    /// Keys for activity constraints
39    pub activity_keys: ActivityKeys,
40}
41
42/// Add asset-level constraints
43///
44/// Note: the ordering of constraints is important, as the dual values of the constraints must later
45/// be retrieved to calculate commodity prices.
46///
47/// # Arguments
48///
49/// * `problem` - The optimisation problem
50/// * `variables` - The variables in the problem
51/// * `model` - The model
52/// * `assets` - The asset pool
53/// * `year` - Current milestone year
54///
55/// # Returns
56///
57/// Keys for the different constraints.
58pub fn add_asset_constraints<'a, I>(
59    problem: &mut Problem,
60    variables: &VariableMap,
61    model: &'a Model,
62    assets: I,
63    year: u32,
64) -> ConstraintKeys
65where
66    I: Iterator<Item = &'a AssetRef> + Clone + 'a,
67{
68    let commodity_balance_keys =
69        add_commodity_balance_constraints(problem, variables, model, assets.clone(), year);
70
71    let activity_keys = add_activity_constraints(problem, variables);
72
73    // Return constraint keys
74    ConstraintKeys {
75        commodity_balance_keys,
76        activity_keys,
77    }
78}
79
80/// Add asset-level input-output commodity balances.
81///
82/// These constraints fix the supply-demand balance for the whole system.
83///
84/// See description in [the dispatch optimisation documentation][1].
85///
86/// [1]: https://energysystemsmodellinglab.github.io/MUSE_2.0/dispatch_optimisation.html#commodity-balance-constraints
87fn add_commodity_balance_constraints<'a, I>(
88    problem: &mut Problem,
89    variables: &VariableMap,
90    model: &'a Model,
91    assets: I,
92    year: u32,
93) -> CommodityBalanceKeys
94where
95    I: Iterator<Item = &'a AssetRef> + Clone + 'a,
96{
97    // Row offset in problem. This line **must** come before we add more constraints.
98    let offset = problem.num_rows();
99
100    let mut keys = Vec::new();
101    let mut terms = Vec::new();
102    for (commodity_id, commodity) in model.commodities.iter() {
103        if !matches!(
104            commodity.kind,
105            CommodityType::SupplyEqualsDemand | CommodityType::ServiceDemand
106        ) {
107            continue;
108        }
109
110        for region_id in model.iter_regions() {
111            for ts_selection in model
112                .time_slice_info
113                .iter_selections_at_level(commodity.time_slice_level)
114            {
115                for (asset, flow) in assets
116                    .clone()
117                    .filter_region(region_id)
118                    .flows_for_commodity(commodity_id)
119                {
120                    // If the commodity has a time slice level of season/annual, the constraint will
121                    // cover multiple time slices
122                    for (time_slice, _) in ts_selection.iter(&model.time_slice_info) {
123                        let var = variables.get(asset, time_slice);
124                        terms.push((var, flow.coeff.value()));
125                    }
126                }
127
128                // Add constraint. For SED commodities, the RHS is zero and for SVD commodities it
129                // is the exogenous demand supplied by the user.
130                let rhs = if commodity.kind == CommodityType::ServiceDemand {
131                    commodity
132                        .demand
133                        .get(&(region_id.clone(), year, ts_selection.clone()))
134                        .unwrap()
135                        .value()
136                } else {
137                    0.0
138                };
139                problem.add_row(rhs..=rhs, terms.drain(..));
140                keys.push((
141                    commodity_id.clone(),
142                    region_id.clone(),
143                    ts_selection.clone(),
144                ))
145            }
146        }
147    }
148
149    CommodityBalanceKeys { offset, keys }
150}
151
152/// Add constraints on the activity of different assets.
153///
154/// This ensures that assets do not exceed their specified capacity and availability for each time
155/// slice.
156fn add_activity_constraints(problem: &mut Problem, variables: &VariableMap) -> ActivityKeys {
157    // Row offset in problem. This line **must** come before we add more constraints.
158    let offset = problem.num_rows();
159
160    let mut keys = Vec::new();
161    for (asset, time_slice, var) in variables.iter() {
162        let limits = asset.get_activity_limits(time_slice);
163        let limits = limits.start().value()..=limits.end().value();
164
165        problem.add_row(limits, [(var, 1.0)]);
166        keys.push((asset.clone(), time_slice.clone()))
167    }
168
169    ActivityKeys { offset, keys }
170}
muse2/simulation/optimisation/constraints.rs

muse2/simulation/optimisation/
constraints.rs
