diff --git a/src/graph/investment.rs b/src/graph/investment.rs index 06e4fc0d7..5d6ade42d 100644 --- a/src/graph/investment.rs +++ b/src/graph/investment.rs @@ -2,7 +2,7 @@ use super::{CommoditiesGraph, GraphEdge, GraphNode}; use crate::commodity::{CommodityMap, CommodityType}; use crate::region::RegionID; -use crate::simulation::investment::InvestmentSet; +use crate::simulation::market::MarketSet; use highs::{Col, HighsModelStatus, RowProblem, Sense}; use indexmap::IndexMap; use log::warn; @@ -13,7 +13,7 @@ use petgraph::visit::EdgeRef; use petgraph::{Directed, Direction}; use std::collections::HashMap; -type InvestmentGraph = Graph; +type InvestmentGraph = Graph; /// Analyse the commodity graphs for a given year to determine the order in which investment /// decisions should be made. @@ -24,10 +24,10 @@ type InvestmentGraph = Graph; /// all regions are combined into a single `InvestmentGraph`. TODO: at present there can be no /// edges between regions; in future we will want to implement trade as edges between regions, /// but this will have no impact on the following steps. -/// 2. Condense strongly connected components (cycles) into `InvestmentSet::Cycle` nodes. +/// 2. Condense strongly connected components (cycles) into `MarketSet::Cycle` nodes. /// 3. Perform a topological sort on the condensed graph. /// 4. Compute layers for investment based on the topological order, grouping independent sets into -/// `InvestmentSet::Layer`s. +/// `MarketSet::Layer`s. /// /// Arguments: /// * `graphs` - Commodity graphs for each region and year, outputted from `build_commodity_graphs_for_model` @@ -35,13 +35,13 @@ type InvestmentGraph = Graph; /// * `year` - The year to solve the investment order for /// /// # Returns -/// A Vec of `InvestmentSet`s in the order they should be solved, with cycles grouped into -/// `InvestmentSet::Cycle`s and independent sets grouped into `InvestmentSet::Layer`s. +/// A Vec of `MarketSet`s in the order they should be solved, with cycles grouped into +/// `MarketSet::Cycle`s and independent sets grouped into `MarketSet::Layer`s. fn solve_investment_order_for_year( graphs: &IndexMap<(RegionID, u32), CommoditiesGraph>, commodities: &CommodityMap, year: u32, -) -> Vec { +) -> Vec { // Initialise InvestmentGraph for this year from the set of original `CommodityGraph`s let mut investment_graph = init_investment_graph_for_year(graphs, year, commodities); @@ -60,7 +60,7 @@ fn solve_investment_order_for_year( /// Initialise an `InvestmentGraph` for the given year from a set of `CommodityGraph`s /// /// Commodity graphs for each region are first filtered to only include SVD/SED commodities. Each -/// commodity node is then added to a global investment graph as an `InvestmentSet::Single`, with +/// commodity node is then added to a global investment graph as an `MarketSet::Single`, with /// edges preserved from the original commodity graphs. fn init_investment_graph_for_year( graphs: &IndexMap<(RegionID, u32), CommoditiesGraph>, @@ -96,7 +96,7 @@ fn init_investment_graph_for_year( }; ( ni, - combined.add_node(InvestmentSet::Single((cid.clone(), region_id.clone()))), + combined.add_node(MarketSet::Single((cid.clone(), region_id.clone()))), ) }) .collect(); @@ -114,7 +114,7 @@ fn init_investment_graph_for_year( combined } -/// Compresses cycles into `InvestmentSet::Cycle` nodes +/// Compresses cycles into `MarketSet::Cycle` nodes fn compress_cycles(graph: &InvestmentGraph) -> InvestmentGraph { // Detect strongly connected components let mut condensed_graph = condensation(graph.clone(), true); @@ -124,12 +124,12 @@ fn compress_cycles(graph: &InvestmentGraph) -> InvestmentGraph { // Map to a new InvestmentGraph condensed_graph.map( - // Map nodes to InvestmentSet + // Map nodes to MarketSet // If only one member, keep as-is; if multiple members, create Cycle |_, node_weight| match node_weight.len() { 0 => unreachable!("Condensed graph node must have at least one member"), 1 => node_weight[0].clone(), - _ => InvestmentSet::Cycle( + _ => MarketSet::Cycle( node_weight .iter() .flat_map(|s| s.iter_markets()) @@ -145,7 +145,7 @@ fn compress_cycles(graph: &InvestmentGraph) -> InvestmentGraph { /// Order the members of each strongly connected component using a mixed-integer linear program. /// /// `condensed_graph` contains the SCCs detected in the original investment graph, stored as -/// `Vec` node weights. Single-element components are already acyclic, but components +/// `Vec` node weights. Single-element components are already acyclic, but components /// with multiple members require an internal ordering so that the investment algorithm can treat /// them as near-acyclic chains, minimising potential disruption. /// @@ -229,22 +229,22 @@ fn compress_cycles(graph: &InvestmentGraph) -> InvestmentGraph { /// * As with any SCC, at least one pairwise violation is guaranteed. In this ordering, the only /// pairwise violation is between B and C, as C is solved before B, but B may consume C. /// -/// The resulting order replaces the original `InvestmentSet::Cycle` entry inside the condensed +/// The resulting order replaces the original `MarketSet::Cycle` entry inside the condensed /// graph, providing a deterministic processing sequence for downstream logic. #[allow(clippy::too_many_lines)] fn order_sccs( - condensed_graph: &mut Graph, GraphEdge>, + condensed_graph: &mut Graph, GraphEdge>, original_graph: &InvestmentGraph, ) { const EXTERNAL_BIAS: f64 = 0.1; - // Map each investment set back to the node index in the original graph so we can inspect edges. - let node_lookup: HashMap = original_graph + // Map each market set back to the node index in the original graph so we can inspect edges. + let node_lookup: HashMap = original_graph .node_indices() .map(|idx| (original_graph.node_weight(idx).unwrap().clone(), idx)) .collect(); - // Work through each SCC; groups with just one investment set don't need to be ordered. + // Work through each SCC; groups with just one market set don't need to be ordered. for group in condensed_graph.node_indices() { let scc = condensed_graph.node_weight_mut(group).unwrap(); let n = scc.len(); @@ -252,7 +252,7 @@ fn order_sccs( continue; } - // Capture current order and resolve each investment set back to its original graph index. + // Capture current order and resolve each market set back to its original graph index. let original_order = scc.clone(); let original_indices = original_order .iter() @@ -387,7 +387,7 @@ fn order_sccs( } } -/// Compute layers of investment sets from the topological order +/// Compute layers of market sets from the topological order /// /// This function works by computing the rank of each node in the graph based on the longest path /// from any root node to that node. Any nodes with the same rank are independent and can be solved @@ -395,11 +395,11 @@ fn order_sccs( /// to lowest rank (root nodes). /// /// This function computes the ranks of each node, groups nodes by rank, and then produces a final -/// ordered Vec of `InvestmentSet`s which gives the order in which to solve the investment decisions. +/// ordered Vec of `MarketSet`s which gives the order in which to solve the investment decisions. /// -/// Investment sets with the same rank (i.e., can be solved in parallel) are grouped into -/// `InvestmentSet::Layer`. Investment sets that are alone in their rank remain as-is (i.e. either -/// `Single` or `Cycle`). `Layer`s can contain a mix of `Single` and `Cycle` investment sets. +/// Market sets with the same rank (i.e., can be solved in parallel) are grouped into +/// `MarketSet::Layer`. Market sets that are alone in their rank remain as-is (i.e. either +/// `Single` or `Cycle`). `Layer`s can contain a mix of `Single` and `Cycle` market sets. /// /// For example, given the following graph: /// @@ -411,11 +411,11 @@ fn order_sccs( /// D E F /// ``` /// -/// Rank 0: A -> `InvestmentSet::Single` -/// Rank 1: B, C -> `InvestmentSet::Layer` -/// Rank 2: D, E, F -> `InvestmentSet::Layer` +/// Rank 0: A -> `MarketSet::Single` +/// Rank 1: B, C -> `MarketSet::Layer` +/// Rank 2: D, E, F -> `MarketSet::Layer` /// -/// These are returned as a `Vec` from highest rank to lowest (i.e. the D, E, F layer +/// These are returned as a `Vec` from highest rank to lowest (i.e. the D, E, F layer /// first, then the B, C layer, then the singleton A). /// /// Arguments: @@ -425,9 +425,9 @@ fn order_sccs( /// * `order` - The topological order of the graph nodes. Computed using `petgraph::algo::toposort`. /// /// Returns: -/// A Vec of `InvestmentSet`s in the order they should be solved, with independent sets grouped into -/// `InvestmentSet::Layer`s. -fn compute_layers(graph: &InvestmentGraph, order: &[NodeIndex]) -> Vec { +/// A Vec of `MarketSet`s in the order they should be solved, with independent sets grouped into +/// `MarketSet::Layer`s. +fn compute_layers(graph: &InvestmentGraph, order: &[NodeIndex]) -> Vec { // Initialize all ranks to 0 let mut ranks: HashMap<_, usize> = graph.node_indices().map(|n| (n, 0)).collect(); @@ -445,27 +445,27 @@ fn compute_layers(graph: &InvestmentGraph, order: &[NodeIndex]) -> Vec> = vec![Vec::new(); max_rank + 1]; + let mut groups: Vec> = vec![Vec::new(); max_rank + 1]; for node_idx in order { let rank = ranks[node_idx]; let w = graph.node_weight(*node_idx).unwrap().clone(); groups[rank].push(w); } - // Produce final ordered Vec: ranks descending (leaf-first), - // compressing equal-rank nodes into an InvestmentSet::Layer. + // Produce final ordered Vec: ranks descending (leaf-first), + // compressing equal-rank nodes into an MarketSet::Layer. let mut result = Vec::new(); for mut items in groups.into_iter().rev() { if items.is_empty() { unreachable!("Should be no gaps in the ranking") } - // If only one InvestmentSet in the group, we do not need to compress into a layer, so just + // If only one MarketSet in the group, we do not need to compress into a layer, so just // push the single item (this item may be a `Single` or `Cycle`). if items.len() == 1 { result.push(items.remove(0)); // Otherwise, create a layer. The items within the layer may be a mix of `Single` or `Cycle`. } else { - result.push(InvestmentSet::Layer(items)); + result.push(MarketSet::Layer(items)); } } @@ -481,14 +481,14 @@ fn compute_layers(graph: &InvestmentGraph, order: &[NodeIndex]) -> Vec, commodities: &CommodityMap, years: &[u32], -) -> HashMap> { +) -> HashMap> { let mut investment_orders = HashMap::new(); for year in years { let order = solve_investment_order_for_year(commodity_graphs, commodities, *year); @@ -508,7 +508,7 @@ mod tests { #[test] fn order_sccs_simple_cycle() { - let markets = ["A", "B", "C"].map(|id| InvestmentSet::Single((id.into(), "GBR".into()))); + let markets = ["A", "B", "C"].map(|id| MarketSet::Single((id.into(), "GBR".into()))); // Create graph with cycle edges plus an extra dependency B ← D (see doc comment) let mut original = InvestmentGraph::new(); @@ -524,7 +524,7 @@ mod tests { ); } // External market receiving exports from C; encourages C to appear early. - let external = original.add_node(InvestmentSet::Single(("X".into(), "GBR".into()))); + let external = original.add_node(MarketSet::Single(("X".into(), "GBR".into()))); original.add_edge( node_indices[2], external, @@ -532,7 +532,7 @@ mod tests { ); // Single SCC containing all markets. - let mut condensed: Graph, GraphEdge> = Graph::new(); + let mut condensed: Graph, GraphEdge> = Graph::new(); let component = condensed.add_node(markets.to_vec()); order_sccs(&mut condensed, &original); @@ -540,7 +540,7 @@ mod tests { // Expected order corresponds to the example in the doc comment. // Note that C should be first, as it has an outgoing edge to the external market. let expected = ["C", "A", "B"] - .map(|id| InvestmentSet::Single((id.into(), "GBR".into()))) + .map(|id| MarketSet::Single((id.into(), "GBR".into()))) .to_vec(); assert_eq!(condensed.node_weight(component).unwrap(), &expected); @@ -569,11 +569,11 @@ mod tests { let result = solve_investment_order_for_year(&graphs, &commodities, 2020); // Expected order: C, B, A (leaf nodes first) - // No cycles or layers, so all investment sets should be `Single` + // No cycles or layers, so all market sets should be `Single` assert_eq!(result.len(), 3); - assert_eq!(result[0], InvestmentSet::Single(("C".into(), "GBR".into()))); - assert_eq!(result[1], InvestmentSet::Single(("B".into(), "GBR".into()))); - assert_eq!(result[2], InvestmentSet::Single(("A".into(), "GBR".into()))); + assert_eq!(result[0], MarketSet::Single(("C".into(), "GBR".into()))); + assert_eq!(result[1], MarketSet::Single(("B".into(), "GBR".into()))); + assert_eq!(result[2], MarketSet::Single(("A".into(), "GBR".into()))); } #[rstest] @@ -596,11 +596,11 @@ mod tests { let graphs = IndexMap::from([(("GBR".into(), 2020), graph)]); let result = solve_investment_order_for_year(&graphs, &commodities, 2020); - // Should be a single `Cycle` investment set containing both commodities + // Should be a single `Cycle` market set containing both commodities assert_eq!(result.len(), 1); assert_eq!( result[0], - InvestmentSet::Cycle(vec![("A".into(), "GBR".into()), ("B".into(), "GBR".into())]) + MarketSet::Cycle(vec![("A".into(), "GBR".into()), ("B".into(), "GBR".into())]) ); } @@ -637,15 +637,15 @@ mod tests { // Expected order: D, Layer(B, C), A assert_eq!(result.len(), 3); - assert_eq!(result[0], InvestmentSet::Single(("D".into(), "GBR".into()))); + assert_eq!(result[0], MarketSet::Single(("D".into(), "GBR".into()))); assert_eq!( result[1], - InvestmentSet::Layer(vec![ - InvestmentSet::Single(("B".into(), "GBR".into())), - InvestmentSet::Single(("C".into(), "GBR".into())) + MarketSet::Layer(vec![ + MarketSet::Single(("B".into(), "GBR".into())), + MarketSet::Single(("C".into(), "GBR".into())) ]) ); - assert_eq!(result[2], InvestmentSet::Single(("A".into(), "GBR".into()))); + assert_eq!(result[2], MarketSet::Single(("A".into(), "GBR".into()))); } #[rstest] @@ -678,23 +678,23 @@ mod tests { assert_eq!(result.len(), 3); assert_eq!( result[0], - InvestmentSet::Layer(vec![ - InvestmentSet::Single(("C".into(), "GBR".into())), - InvestmentSet::Single(("C".into(), "FRA".into())) + MarketSet::Layer(vec![ + MarketSet::Single(("C".into(), "GBR".into())), + MarketSet::Single(("C".into(), "FRA".into())) ]) ); assert_eq!( result[1], - InvestmentSet::Layer(vec![ - InvestmentSet::Single(("B".into(), "GBR".into())), - InvestmentSet::Single(("B".into(), "FRA".into())) + MarketSet::Layer(vec![ + MarketSet::Single(("B".into(), "GBR".into())), + MarketSet::Single(("B".into(), "FRA".into())) ]) ); assert_eq!( result[2], - InvestmentSet::Layer(vec![ - InvestmentSet::Single(("A".into(), "GBR".into())), - InvestmentSet::Single(("A".into(), "FRA".into())) + MarketSet::Layer(vec![ + MarketSet::Single(("A".into(), "GBR".into())), + MarketSet::Single(("A".into(), "FRA".into())) ]) ); } diff --git a/src/model.rs b/src/model.rs index 4dddaaa13..43d1ef0ce 100644 --- a/src/model.rs +++ b/src/model.rs @@ -4,7 +4,7 @@ use crate::asset::UserAsset; use crate::commodity::{CommodityID, CommodityMap}; use crate::process::ProcessMap; use crate::region::{Region, RegionID, RegionMap}; -use crate::simulation::investment::InvestmentSet; +use crate::simulation::market::MarketSet; use crate::time_slice::TimeSliceInfo; use std::collections::HashMap; use std::path::PathBuf; @@ -33,7 +33,7 @@ pub struct Model { /// User-defined assets pub user_assets: Vec, /// Commodity ordering for each milestone year - pub investment_order: HashMap>, + pub investment_order: HashMap>, } impl Model { diff --git a/src/simulation.rs b/src/simulation.rs index dbaaa4e9c..53f195db4 100644 --- a/src/simulation.rs +++ b/src/simulation.rs @@ -14,6 +14,7 @@ pub mod optimisation; use optimisation::{DispatchRun, FlowMap}; pub mod investment; use investment::perform_agent_investment; +pub mod market; pub mod prices; pub use prices::PriceMap; diff --git a/src/simulation/investment.rs b/src/simulation/investment.rs index d121f408e..5aaf76a13 100644 --- a/src/simulation/investment.rs +++ b/src/simulation/investment.rs @@ -1,7 +1,7 @@ //! Code for performing agent investment. use super::optimisation::{DispatchRun, FlowMap}; use crate::agent::{Agent, AgentID}; -use crate::asset::{Asset, AssetCapacity, AssetIterator, AssetRef, AssetState}; +use crate::asset::{Asset, AssetCapacity, AssetRef, AssetState}; use crate::commodity::{Commodity, CommodityID, CommodityMap}; use crate::model::Model; use crate::output::DataWriter; @@ -9,12 +9,11 @@ use crate::region::RegionID; use crate::simulation::prices::Prices; use crate::time_slice::{TimeSliceID, TimeSliceInfo, TimeSliceLevel}; use crate::units::{ActivityPerCapacity, Capacity, Dimensionless, Flow, FlowPerCapacity}; -use anyhow::{Context, Result, ensure}; +use anyhow::{Result, ensure}; use indexmap::IndexMap; -use itertools::{Itertools, chain}; +use itertools::Itertools; use log::{debug, warn}; use std::collections::{HashMap, HashSet}; -use std::fmt::Display; use strum::IntoEnumIterator; pub mod appraisal; @@ -25,136 +24,10 @@ use appraisal::{ }; /// A map of demand across time slices for a specific market -type DemandMap = IndexMap; +pub type DemandMap = IndexMap; /// Demand for a given combination of commodity, region and time slice -type AllDemandMap = IndexMap<(CommodityID, RegionID, TimeSliceID), Flow>; - -/// Represents a set of markets which are invested in together. -#[derive(PartialEq, Debug, Clone, Eq, Hash)] -pub enum InvestmentSet { - /// Assets are selected for a single market using `select_assets_for_single_market` - Single((CommodityID, RegionID)), - /// Assets are selected for a group of markets which forms a cycle. - /// Experimental: handled by `select_assets_for_cycle` and guarded by the broken options - /// parameter. - Cycle(Vec<(CommodityID, RegionID)>), - /// Assets are selected for a layer of independent `InvestmentSet`s - Layer(Vec), -} - -impl InvestmentSet { - /// Recursively iterate over all markets contained in this `InvestmentSet`. - pub fn iter_markets<'a>( - &'a self, - ) -> Box + 'a> { - match self { - InvestmentSet::Single(market) => Box::new(std::iter::once(market)), - InvestmentSet::Cycle(markets) => Box::new(markets.iter()), - InvestmentSet::Layer(set) => Box::new(set.iter().flat_map(|s| s.iter_markets())), - } - } - - /// Selects assets for this investment set variant and passes through the shared - /// context needed by single-market, cycle, or layered selection. - /// - /// # Arguments - /// - /// * `model` – Simulation model supplying parameters, processes, and dispatch. - /// * `year` – Planning year being solved. - /// * `demand` – Net demand profiles available to all markets before selection. - /// * `existing_assets` – Assets already commissioned in the system. - /// * `prices` – Commodity price assumptions to use when valuing investments. - /// * `seen_markets` – Markets for which investments have already been settled. - /// * `previously_selected_assets` – Assets chosen in earlier investment sets. - /// * `writer` – Data sink used to log optimisation artefacts. - #[allow(clippy::too_many_arguments)] - fn select_assets( - &self, - model: &Model, - year: u32, - demand: &AllDemandMap, - existing_assets: &[AssetRef], - prices: &Prices, - seen_markets: &[(CommodityID, RegionID)], - previously_selected_assets: &[AssetRef], - writer: &mut DataWriter, - ) -> Result> { - match self { - InvestmentSet::Single((commodity_id, region_id)) => select_assets_for_single_market( - model, - commodity_id, - region_id, - year, - demand, - existing_assets, - prices, - writer, - ), - InvestmentSet::Cycle(markets) => { - debug!("Starting investment for cycle '{self}'"); - select_assets_for_cycle( - model, - markets, - year, - demand, - existing_assets, - prices, - seen_markets, - previously_selected_assets, - writer, - ) - .with_context(|| { - format!( - "Investments failed for market set {self} with cyclical dependencies. \ - Please note that the investment algorithm is currently experimental for \ - models with circular commodity dependencies and may not be able to find \ - a solution in all cases." - ) - }) - } - InvestmentSet::Layer(investment_sets) => { - debug!("Starting asset selection for layer '{self}'"); - let mut all_assets = Vec::new(); - for investment_set in investment_sets { - let assets = investment_set.select_assets( - model, - year, - demand, - existing_assets, - prices, - seen_markets, - previously_selected_assets, - writer, - )?; - all_assets.extend(assets); - } - debug!("Completed asset selection for layer '{self}'"); - Ok(all_assets) - } - } - } -} - -impl Display for InvestmentSet { - fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { - match self { - InvestmentSet::Single((commodity_id, region_id)) => { - write!(f, "{commodity_id}|{region_id}") - } - InvestmentSet::Cycle(markets) => { - write!( - f, - "({})", - markets.iter().map(|(c, r)| format!("{c}|{r}")).join(", ") - ) - } - InvestmentSet::Layer(ids) => { - write!(f, "[{}]", ids.iter().join(", ")) - } - } - } -} +pub type AllDemandMap = IndexMap<(CommodityID, RegionID, TimeSliceID), Flow>; /// Perform agent investment to determine capacity investment of new assets for next milestone year. /// @@ -196,10 +69,10 @@ pub fn perform_agent_investment( // markets that have been seen so far. let mut seen_markets = Vec::new(); - // Iterate over investment sets in the investment order for this year - for investment_set in investment_order { - // Select assets for this investment set - let selected_assets = investment_set.select_assets( + // Iterate over market sets in the investment order for this year + for market_set in investment_order { + // Select assets for this market set + let selected_assets = market_set.select_assets( model, year, &net_demand, @@ -211,7 +84,7 @@ pub fn perform_agent_investment( )?; // Update our list of seen markets - for market in investment_set.iter_markets() { + for market in market_set.iter_markets() { seen_markets.push(market.clone()); } @@ -219,7 +92,7 @@ pub fn perform_agent_investment( // **TODO**: this probably means there's no demand for the market, which we could // presumably preempt if selected_assets.is_empty() { - debug!("No assets selected for '{investment_set}'"); + debug!("No assets selected for '{market_set}'"); continue; } @@ -229,14 +102,14 @@ pub fn perform_agent_investment( // Perform dispatch optimisation with assets that have been selected so far // **TODO**: presumably we only need to do this for selected_assets, as assets added in // previous iterations should not change - debug!("Running post-investment dispatch for '{investment_set}'"); + debug!("Running post-investment dispatch for '{market_set}'"); // As upstream markets by definition will not yet have producers, we explicitly set // their prices using external values so that they don't appear free let solution = DispatchRun::new(model, &all_selected_assets, year) .with_market_balance_subset(&seen_markets) .with_input_prices(&prices.shadow) - .run(&format!("post {investment_set} investment"), writer)?; + .run(&format!("post {market_set} investment"), writer)?; // Update demand map with flows from newly added assets update_net_demand_map( @@ -249,210 +122,6 @@ pub fn perform_agent_investment( Ok(all_selected_assets) } -/// Select assets for a single market in a given year -/// -/// Returns a list of assets that are selected for investment for this market in this year. -#[allow(clippy::too_many_arguments)] -fn select_assets_for_single_market( - model: &Model, - commodity_id: &CommodityID, - region_id: &RegionID, - year: u32, - demand: &AllDemandMap, - existing_assets: &[AssetRef], - prices: &Prices, - writer: &mut DataWriter, -) -> Result> { - let commodity = &model.commodities[commodity_id]; - - let mut selected_assets = Vec::new(); - for (agent, commodity_portion) in - get_responsible_agents(model.agents.values(), commodity_id, region_id, year) - { - debug!( - "Running asset selection for agent '{}' in market '{}|{}'", - &agent.id, commodity_id, region_id - ); - - // Get demand portion for this market for this agent in this year - let demand_portion_for_market = get_demand_portion_for_market( - &model.time_slice_info, - demand, - commodity_id, - region_id, - commodity_portion, - ); - - // Existing and candidate assets from which to choose - let opt_assets = get_asset_options( - &model.time_slice_info, - existing_assets, - &demand_portion_for_market, - agent, - commodity, - region_id, - year, - ) - .collect::>(); - - // Calculate investment limits for candidate assets - let investment_limits = - calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); - - // Choose assets from among existing pool and candidates - let best_assets = select_best_assets( - model, - opt_assets, - investment_limits, - commodity, - agent, - region_id, - prices, - demand_portion_for_market, - year, - writer, - )?; - selected_assets.extend(best_assets); - } - - Ok(selected_assets) -} - -/// Iterates through the a pre-ordered set of markets forming a cycle, selecting assets for each -/// market in turn. -/// -/// Dispatch optimisation is performed after each market is visited to rebalance demand. -/// While dispatching, newly selected (`Ready`) assets are given flexible capacity (bounded by -/// `capacity_margin`) so small demand shifts caused by later markets can be absorbed. After all -/// markets have been visited once, the final set of assets is returned, applying any capacity -/// adjustments from the final full-system dispatch optimisation. -/// -/// Dispatch may fail at any point if new demands are encountered for previously visited markets, -/// and the `capacity_margin` is not sufficient to absorb the demand shift. At this point, the -/// simulation is terminated with an error prompting the user to increase the `capacity_margin`. -/// A longer-term solution (TODO) may be to trigger re-investment for the affected markets. Other -/// yet-to-implement features may also help to stabilise the cycle, such as capacity growth limits. -#[allow(clippy::too_many_arguments)] -fn select_assets_for_cycle( - model: &Model, - markets: &[(CommodityID, RegionID)], - year: u32, - demand: &AllDemandMap, - existing_assets: &[AssetRef], - prices: &Prices, - seen_markets: &[(CommodityID, RegionID)], - previously_selected_assets: &[AssetRef], - writer: &mut DataWriter, -) -> Result> { - // Precompute a joined string for logging - let markets_str = markets.iter().map(|(c, r)| format!("{c}|{r}")).join(", "); - - // Iterate over the markets to select assets - let mut current_demand = demand.clone(); - let mut assets_for_cycle = IndexMap::new(); - let mut last_solution = None; - for (idx, (commodity_id, region_id)) in markets.iter().enumerate() { - // Select assets for this market - let assets = select_assets_for_single_market( - model, - commodity_id, - region_id, - year, - ¤t_demand, - existing_assets, - prices, - writer, - )?; - assets_for_cycle.insert((commodity_id.clone(), region_id.clone()), assets); - - // Assemble full list of assets for dispatch (previously selected + all chosen so far) - let mut all_assets = previously_selected_assets.to_vec(); - let assets_for_cycle_flat: Vec<_> = assets_for_cycle - .values() - .flat_map(|v| v.iter().cloned()) - .collect(); - all_assets.extend_from_slice(&assets_for_cycle_flat); - - // We balance all previously seen markets plus all cycle markets up to and including this one - let mut markets_to_balance = seen_markets.to_vec(); - markets_to_balance.extend_from_slice(&markets[0..=idx]); - - // We allow all `Ready` state assets to have flexible capacity - let flexible_capacity_assets: Vec<_> = assets_for_cycle_flat - .iter() - .filter(|asset| matches!(asset.state(), AssetState::Ready { .. })) - .cloned() - .collect(); - - // Retrieve installable capacity limits for flexible capacity assets. - let mut agent_share_cache = HashMap::new(); - let capacity_limits = flexible_capacity_assets - .iter() - .filter_map(|asset| { - let agent_id = asset.agent_id().unwrap(); - let commodity_id = asset.primary_output_commodity().unwrap(); - let agent_share = *agent_share_cache - .entry((agent_id, commodity_id)) - .or_insert_with(|| { - model.agents[agent_id].commodity_portions[&(commodity_id.clone(), year)] - }); - asset - .max_installable_capacity(agent_share) - .map(|max_capacity| (asset.clone(), max_capacity)) - }) - .collect::>(); - - // Run dispatch - let solution = DispatchRun::new(model, &all_assets, year) - .with_market_balance_subset(&markets_to_balance) - .with_flexible_capacity_assets( - &flexible_capacity_assets, - Some(&capacity_limits), - // Gives newly selected cycle assets limited capacity wiggle-room; existing assets stay fixed. - model.parameters.capacity_margin, - ) - .run( - &format!("cycle ({markets_str}) post {commodity_id}|{region_id} investment"), - writer, - ) - .with_context(|| { - format!( - "Cycle balancing failed for cycle ({markets_str}), capacity_margin: {}. \ - Try increasing the capacity_margin.", - model.parameters.capacity_margin - ) - })?; - - // Calculate new net demand map with all assets selected so far - current_demand.clone_from(demand); - update_net_demand_map( - &mut current_demand, - &solution.create_flow_map(), - &assets_for_cycle_flat, - ); - last_solution = Some(solution); - } - - // Finally, update flexible capacity assets based on the final solution - let mut all_cycle_assets: Vec<_> = assets_for_cycle.into_values().flatten().collect(); - if let Some(solution) = last_solution { - let new_capacities: HashMap<_, _> = solution.iter_capacity().collect(); - for asset in &mut all_cycle_assets { - if let Some(new_capacity) = new_capacities.get(asset) { - debug!( - "Capacity of asset '{}' modified during cycle balancing ({} to {})", - asset.process_id(), - asset.total_capacity(), - new_capacity.total_capacity() - ); - asset.make_mut().set_capacity(*new_capacity); - } - } - } - - Ok(all_cycle_assets) -} - /// Flatten the preset commodity demands for a given year into a map of commodity, region and /// time slice to demand. /// @@ -500,7 +169,7 @@ fn flatten_preset_demands_for_year( /// /// TODO: this is a very flawed approach. The proper solution might be for agents to consider /// multiple commodities simultaneously, but that would require substantial work to implement. -fn update_net_demand_map(demand: &mut AllDemandMap, flows: &FlowMap, assets: &[AssetRef]) { +pub fn update_net_demand_map(demand: &mut AllDemandMap, flows: &FlowMap, assets: &[AssetRef]) { for ((asset, commodity_id, time_slice), flow) in flows { if assets.contains(asset) { let key = ( @@ -527,50 +196,6 @@ fn update_net_demand_map(demand: &mut AllDemandMap, flows: &FlowMap, assets: &[A } /// Get a portion of the demand profile for this market -fn get_demand_portion_for_market( - time_slice_info: &TimeSliceInfo, - demand: &AllDemandMap, - commodity_id: &CommodityID, - region_id: &RegionID, - commodity_portion: Dimensionless, -) -> DemandMap { - time_slice_info - .iter_ids() - .map(|time_slice| { - ( - time_slice.clone(), - commodity_portion - * *demand - .get(&(commodity_id.clone(), region_id.clone(), time_slice.clone())) - .unwrap_or(&Flow(0.0)), - ) - }) - .collect() -} - -/// Get the agents responsible for a given market in a given year along with the commodity -/// portion for which they are responsible -fn get_responsible_agents<'a, I>( - agents: I, - commodity_id: &'a CommodityID, - region_id: &'a RegionID, - year: u32, -) -> impl Iterator -where - I: Iterator, -{ - agents.filter_map(move |agent| { - if !agent.regions.contains(region_id) { - return None; - } - let portion = agent - .commodity_portions - .get(&(commodity_id.clone(), year))?; - - Some((agent, *portion)) - }) -} - /// Returns the minimum installed capacity required for `asset` to satisfy the demand that it can /// potentially serve, accounting for its activity constraints. /// @@ -584,7 +209,7 @@ where /// /// Selections whose maximum supply is zero are ignored. Such selections would otherwise imply an /// infinite capacity requirement and therefore provide no useful lower bound. -fn get_demand_limiting_capacity( +pub fn get_demand_limiting_capacity( time_slice_info: &TimeSliceInfo, asset: &Asset, commodity: &Commodity, @@ -654,57 +279,6 @@ fn get_demand_limiting_capacity( capacity } -/// Get options from existing and potential assets for the given parameters -fn get_asset_options<'a>( - time_slice_info: &'a TimeSliceInfo, - all_existing_assets: &'a [AssetRef], - demand: &'a DemandMap, - agent: &'a Agent, - commodity: &'a Commodity, - region_id: &'a RegionID, - year: u32, -) -> impl Iterator + 'a { - // Get existing assets which produce the commodity of interest - let existing_assets = all_existing_assets - .iter() - .filter_agent(&agent.id) - .filter_region(region_id) - .filter_primary_producers_of(&commodity.id) - .cloned(); - - // Get candidates assets which produce the commodity of interest - let candidate_assets = - get_candidate_assets(time_slice_info, demand, agent, region_id, commodity, year); - - chain(existing_assets, candidate_assets) -} - -/// Get candidate assets which produce a particular commodity for a given agent -fn get_candidate_assets<'a>( - time_slice_info: &'a TimeSliceInfo, - demand: &'a DemandMap, - agent: &'a Agent, - region_id: &'a RegionID, - commodity: &'a Commodity, - year: u32, -) -> impl Iterator + 'a { - agent - .iter_search_space(region_id, &commodity.id, year) - .map(move |process| { - let mut asset = - Asset::new_candidate(process.clone(), region_id.clone(), Capacity(0.0), year) - .unwrap(); - - // Set capacity based on demand - // This will serve as the upper limit when appraising the asset - let capacity = get_demand_limiting_capacity(time_slice_info, &asset, commodity, demand); - let asset_capacity = AssetCapacity::from_capacity(capacity, asset.unit_size()); - asset.set_capacity(asset_capacity); - - asset.into() - }) -} - /// Print debug message if there are multiple equally good outputs fn log_on_equal_appraisal_outputs( outputs: &[AppraisalOutput], @@ -742,37 +316,9 @@ fn log_on_equal_appraisal_outputs( } } -/// Calculate investment limits for an agent's candidate assets in a given year -/// -/// Investment limits are based on demand for the commodity (capacity cannot exceed that needed to -/// meet demand), and any annual addition limits specified by the process (scaled according to the -/// agent's portion of the commodity demand and the number of years elapsed since the previous -/// milestone year). -fn calculate_investment_limits_for_candidates( - opt_assets: &[AssetRef], - commodity_portion: Dimensionless, -) -> HashMap { - // Calculate limits for each candidate asset - opt_assets - .iter() - .filter(|asset| !asset.is_commissioned()) - .map(|asset| { - // Start off with the demand-limiting capacity (pre-calculated when creating candidate) - let mut cap = asset.capacity(); - - // Cap by the addition limits of the process, if specified - if let Some(limit_capacity) = asset.max_installable_capacity(commodity_portion) { - cap = cap.min(limit_capacity); - } - - (asset.clone(), cap) - }) - .collect() -} - /// Get the best assets for meeting demand for the given commodity #[allow(clippy::too_many_arguments)] -fn select_best_assets( +pub fn select_best_assets( model: &Model, mut opt_assets: Vec, investment_limits: HashMap, @@ -966,19 +512,14 @@ mod tests { use super::*; use crate::commodity::Commodity; use crate::fixture::{ - agent_id, asset, process, process_activity_limits_map, process_flows_map, - process_investment_constraints, process_parameter_map, region_id, svd_commodity, - time_slice, time_slice_info, time_slice_info2, - }; - use crate::process::{ - ActivityLimits, FlowType, Process, ProcessActivityLimitsMap, ProcessFlow, ProcessFlowsMap, - ProcessInvestmentConstraint, ProcessInvestmentConstraintsMap, ProcessParameterMap, + asset, process, process_activity_limits_map, process_flows_map, svd_commodity, time_slice, + time_slice_info, time_slice_info2, }; - use crate::region::RegionID; + use crate::process::{ActivityLimits, FlowType, Process, ProcessFlow}; use crate::time_slice::{TimeSliceID, TimeSliceInfo, TimeSliceSelection}; use crate::units::Dimensionless; - use crate::units::{ActivityPerCapacity, Capacity, Flow, FlowPerActivity, MoneyPerFlow}; - use indexmap::{IndexSet, indexmap}; + use crate::units::{Capacity, Flow, FlowPerActivity, MoneyPerFlow}; + use indexmap::indexmap; use rstest::rstest; use std::rc::Rc; @@ -1125,260 +666,4 @@ mod tests { assert_eq!(result, Capacity(20.0)); } - - #[rstest] - fn calculate_investment_limits_for_candidates_empty_list() { - // Test with empty list of assets - let opt_assets: Vec = vec![]; - let commodity_portion = Dimensionless(1.0); - - let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); - - assert!(result.is_empty()); - } - - #[rstest] - fn calculate_investment_limits_for_candidates_commissioned_assets_filtered( - process: Process, - region_id: RegionID, - agent_id: AgentID, - ) { - // Create a mix of commissioned and candidate assets - let process_rc = Rc::new(process); - let capacity = Capacity(10.0); - - // Create commissioned asset - should be filtered out - let commissioned_asset = Asset::new_commissioned( - agent_id.clone(), - process_rc.clone(), - region_id.clone(), - capacity, - 2015, - ) - .unwrap(); - - // Create candidate asset - should be included - let candidate_asset = - Asset::new_candidate(process_rc.clone(), region_id.clone(), capacity, 2015).unwrap(); - - let candidate_asset_ref = AssetRef::from(candidate_asset); - let opt_assets = vec![ - AssetRef::from(commissioned_asset), - candidate_asset_ref.clone(), - ]; - let commodity_portion = Dimensionless(1.0); - - let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); - - // Only the candidate asset should be in the result - assert_eq!(result.len(), 1); - assert!(result.contains_key(&candidate_asset_ref)); - } - - #[rstest] - fn calculate_investment_limits_for_candidates_no_investment_constraints( - process: Process, - region_id: RegionID, - ) { - // Create candidate asset without investment constraints - let process_rc = Rc::new(process); - let capacity = Capacity(15.0); - - let candidate_asset = Asset::new_candidate(process_rc, region_id, capacity, 2015).unwrap(); - - let opt_assets = vec![AssetRef::from(candidate_asset.clone())]; - let commodity_portion = Dimensionless(0.8); - - let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); - - // Should return the asset's original capacity since no constraints apply - assert_eq!(result.len(), 1); - let asset_ref = AssetRef::from(candidate_asset); - assert_eq!(result[&asset_ref], AssetCapacity::Continuous(capacity)); - } - - #[rstest] - // Asset capacity higher than constraint -> limited by constraint - #[case(Capacity(15.0), Capacity(10.0))] - // Asset capacity lower than constraint -> limited by asset capacity - #[case(Capacity(5.0), Capacity(5.0))] - fn calculate_investment_limits_for_candidates_with_constraints( - region_id: RegionID, - process_activity_limits_map: ProcessActivityLimitsMap, - process_flows_map: ProcessFlowsMap, - process_parameter_map: ProcessParameterMap, - #[case] asset_capacity: Capacity, - #[case] expected_limit: Capacity, - ) { - let region_ids: IndexSet = [region_id.clone()].into(); - - // Add investment constraint with addition limit - let constraint = ProcessInvestmentConstraint { - addition_limit: Some(Capacity(10.0)), - }; - let mut constraints = ProcessInvestmentConstraintsMap::new(); - constraints.insert((region_id.clone(), 2015), Rc::new(constraint)); - - let process = Process { - id: "constrained_process".into(), - description: "Process with constraints".into(), - years: 2010..=2020, - activity_limits: process_activity_limits_map, - flows: process_flows_map, - parameters: process_parameter_map, - regions: region_ids, - primary_output: None, - capacity_to_activity: ActivityPerCapacity(1.0), - investment_constraints: constraints, - unit_size: None, - }; - - let process_rc = Rc::new(process); - - let candidate_asset = - Asset::new_candidate(process_rc, region_id, asset_capacity, 2015).unwrap(); - - let opt_assets = vec![AssetRef::from(candidate_asset.clone())]; - let commodity_portion = Dimensionless(1.0); - - let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); - - // Should be limited by the minimum of asset capacity and constraint - assert_eq!(result.len(), 1); - let asset_ref = AssetRef::from(candidate_asset); - assert_eq!( - result[&asset_ref], - AssetCapacity::Continuous(expected_limit) - ); - } - - #[rstest] - fn calculate_investment_limits_for_candidates_multiple_assets( - region_id: RegionID, - process_activity_limits_map: ProcessActivityLimitsMap, - process_flows_map: ProcessFlowsMap, - process_parameter_map: ProcessParameterMap, - ) { - let region_ids: IndexSet = [region_id.clone()].into(); - - // Create first process with constraints - let constraint1 = ProcessInvestmentConstraint { - addition_limit: Some(Capacity(12.0)), - }; - let mut constraints1 = ProcessInvestmentConstraintsMap::new(); - constraints1.insert((region_id.clone(), 2015), Rc::new(constraint1)); - - let process1 = Process { - id: "process1".into(), - description: "First process".into(), - years: 2010..=2020, - activity_limits: process_activity_limits_map.clone(), - flows: process_flows_map.clone(), - parameters: process_parameter_map.clone(), - regions: region_ids.clone(), - primary_output: None, - capacity_to_activity: ActivityPerCapacity(1.0), - investment_constraints: constraints1, - unit_size: None, - }; - - // Create second process without constraints - let process2 = Process { - id: "process2".into(), - description: "Second process".into(), - years: 2010..=2020, - activity_limits: process_activity_limits_map, - flows: process_flows_map, - parameters: process_parameter_map, - regions: region_ids, - primary_output: None, - capacity_to_activity: ActivityPerCapacity(1.0), - investment_constraints: process_investment_constraints(), - unit_size: None, - }; - - let process1_rc = Rc::new(process1); - let process2_rc = Rc::new(process2); - - let candidate1 = - Asset::new_candidate(process1_rc, region_id.clone(), Capacity(20.0), 2015).unwrap(); - - let candidate2 = Asset::new_candidate(process2_rc, region_id, Capacity(8.0), 2015).unwrap(); - - let opt_assets = vec![ - AssetRef::from(candidate1.clone()), - AssetRef::from(candidate2.clone()), - ]; - let commodity_portion = Dimensionless(0.75); - - let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); - - // Should have both assets in result - assert_eq!(result.len(), 2); - - // First asset should be limited by constraint: 12.0 * 0.75 = 9.0 - let asset1_ref = AssetRef::from(candidate1); - assert_eq!( - result[&asset1_ref], - AssetCapacity::Continuous(Capacity(9.0)) - ); - - // Second asset should use its original capacity (no constraints) - let asset2_ref = AssetRef::from(candidate2); - assert_eq!( - result[&asset2_ref], - AssetCapacity::Continuous(Capacity(8.0)) - ); - } - - #[rstest] - fn calculate_investment_limits_for_candidates_discrete_capacity( - region_id: RegionID, - process_activity_limits_map: crate::process::ProcessActivityLimitsMap, - process_flows_map: crate::process::ProcessFlowsMap, - process_parameter_map: crate::process::ProcessParameterMap, - ) { - let region_ids: IndexSet = [region_id.clone()].into(); - - // Add investment constraint - let constraint = ProcessInvestmentConstraint { - addition_limit: Some(Capacity(35.0)), // Enough for 3.5 units at 10.0 each - }; - let mut constraints = ProcessInvestmentConstraintsMap::new(); - constraints.insert((region_id.clone(), 2015), Rc::new(constraint)); - - let process = Process { - id: "discrete_process".into(), - description: "Process with discrete units".into(), - years: 2010..=2020, - activity_limits: process_activity_limits_map, - flows: process_flows_map, - parameters: process_parameter_map, - regions: region_ids, - primary_output: None, - capacity_to_activity: ActivityPerCapacity(1.0), - investment_constraints: constraints, - unit_size: Some(Capacity(10.0)), // Discrete units of 10.0 capacity each - }; - - let process_rc = Rc::new(process); - let capacity = Capacity(50.0); // 5 units at 10.0 each - - let candidate_asset = Asset::new_candidate(process_rc, region_id, capacity, 2015).unwrap(); - - let opt_assets = vec![AssetRef::from(candidate_asset.clone())]; - let commodity_portion = Dimensionless(1.0); - - let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); - - // Should be limited by constraint and rounded down to whole units - // Constraint: 35.0, divided by unit size 10.0 = 3.5 -> floor to 3 units = 30.0 - assert_eq!(result.len(), 1); - let asset_ref = AssetRef::from(candidate_asset); - assert_eq!( - result[&asset_ref], - AssetCapacity::Discrete(3, Capacity(10.0)) - ); - assert_eq!(result[&asset_ref].total_capacity(), Capacity(30.0)); - } } diff --git a/src/simulation/market.rs b/src/simulation/market.rs new file mode 100644 index 000000000..f54e703aa --- /dev/null +++ b/src/simulation/market.rs @@ -0,0 +1,748 @@ +//! Code for creating sets of markets. +use super::optimisation::DispatchRun; +use crate::agent::Agent; +use crate::asset::{Asset, AssetCapacity, AssetIterator, AssetRef, AssetState}; +use crate::commodity::{Commodity, CommodityID}; +use crate::model::Model; +use crate::output::DataWriter; +use crate::region::RegionID; +use crate::simulation::investment::{ + AllDemandMap, DemandMap, get_demand_limiting_capacity, select_best_assets, + update_net_demand_map, +}; +use crate::simulation::prices::Prices; +use crate::time_slice::TimeSliceInfo; +use crate::units::{Capacity, Dimensionless, Flow}; +use anyhow::{Context, Result}; +use indexmap::IndexMap; +use itertools::{Itertools, chain}; +use log::debug; +use std::collections::HashMap; +use std::fmt::Display; + +/// Represents a set of markets which are invested in together. +#[derive(PartialEq, Debug, Clone, Eq, Hash)] +pub enum MarketSet { + /// Assets are selected for a single market using [`select_assets_for_single_market`] + Single((CommodityID, RegionID)), + /// Assets are selected for a group of markets which forms a cycle. + /// Experimental: handled by [`select_assets_for_cycle`]. May not work in every case. + Cycle(Vec<(CommodityID, RegionID)>), + /// Assets are selected for a layer of independent [`MarketSet`]s + Layer(Vec), +} + +impl MarketSet { + /// Recursively iterate over all markets contained in this `MarketSet`. + pub fn iter_markets<'a>( + &'a self, + ) -> Box + 'a> { + match self { + MarketSet::Single(market) => Box::new(std::iter::once(market)), + MarketSet::Cycle(markets) => Box::new(markets.iter()), + MarketSet::Layer(set) => Box::new(set.iter().flat_map(|s| s.iter_markets())), + } + } + + /// Selects assets for this market set variant and passes through the shared + /// context needed by single-market, cycle, or layered selection. + /// + /// # Arguments + /// + /// * `model` – Simulation model supplying parameters, processes, and dispatch. + /// * `year` – Planning year being solved. + /// * `demand` – Net demand profiles available to all markets before selection. + /// * `existing_assets` – Assets already commissioned in the system. + /// * `prices` – Commodity price assumptions to use when valuing investments. + /// * `seen_markets` – Markets for which investments have already been settled. + /// * `previously_selected_assets` – Assets chosen in earlier market sets. + /// * `writer` – Data sink used to log optimisation artefacts. + #[allow(clippy::too_many_arguments)] + pub fn select_assets( + &self, + model: &Model, + year: u32, + demand: &AllDemandMap, + existing_assets: &[AssetRef], + prices: &Prices, + seen_markets: &[(CommodityID, RegionID)], + previously_selected_assets: &[AssetRef], + writer: &mut DataWriter, + ) -> Result> { + match self { + MarketSet::Single((commodity_id, region_id)) => select_assets_for_single_market( + model, + commodity_id, + region_id, + year, + demand, + existing_assets, + prices, + writer, + ), + MarketSet::Cycle(markets) => { + debug!("Starting investment for cycle '{self}'"); + select_assets_for_cycle( + model, + markets, + year, + demand, + existing_assets, + prices, + seen_markets, + previously_selected_assets, + writer, + ) + .with_context(|| { + format!( + "Investments failed for market set {self} with cyclical dependencies. \ + Please note that the investment algorithm is currently experimental for \ + models with circular commodity dependencies and may not be able to find \ + a solution in all cases." + ) + }) + } + MarketSet::Layer(investment_sets) => { + debug!("Starting asset selection for layer '{self}'"); + let mut all_assets = Vec::new(); + for investment_set in investment_sets { + let assets = investment_set.select_assets( + model, + year, + demand, + existing_assets, + prices, + seen_markets, + previously_selected_assets, + writer, + )?; + all_assets.extend(assets); + } + debug!("Completed asset selection for layer '{self}'"); + Ok(all_assets) + } + } + } +} + +impl Display for MarketSet { + fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { + match self { + MarketSet::Single((commodity_id, region_id)) => { + write!(f, "{commodity_id}|{region_id}") + } + MarketSet::Cycle(markets) => { + write!( + f, + "({})", + markets.iter().map(|(c, r)| format!("{c}|{r}")).join(", ") + ) + } + MarketSet::Layer(ids) => { + write!(f, "[{}]", ids.iter().join(", ")) + } + } + } +} + +/// Calculate investment limits for an agent's candidate assets in a given year +/// +/// Investment limits are based on demand for the commodity (capacity cannot exceed that needed to +/// meet demand), and any annual addition limits specified by the process (scaled according to the +/// agent's portion of the commodity demand and the number of years elapsed since the previous +/// milestone year). +fn calculate_investment_limits_for_candidates( + opt_assets: &[AssetRef], + commodity_portion: Dimensionless, +) -> HashMap { + // Calculate limits for each candidate asset + opt_assets + .iter() + .filter(|asset| !asset.is_commissioned()) + .map(|asset| { + // Start off with the demand-limiting capacity (pre-calculated when creating candidate) + let mut cap = asset.capacity(); + + // Cap by the addition limits of the process, if specified + if let Some(limit_capacity) = asset.max_installable_capacity(commodity_portion) { + cap = cap.min(limit_capacity); + } + + (asset.clone(), cap) + }) + .collect() +} + +/// Get candidate assets which produce a particular commodity for a given agent +fn get_candidate_assets<'a>( + time_slice_info: &'a TimeSliceInfo, + demand: &'a DemandMap, + agent: &'a Agent, + region_id: &'a RegionID, + commodity: &'a Commodity, + year: u32, +) -> impl Iterator + 'a { + agent + .iter_search_space(region_id, &commodity.id, year) + .map(move |process| { + let mut asset = + Asset::new_candidate(process.clone(), region_id.clone(), Capacity(0.0), year) + .unwrap(); + + // Set capacity based on demand + // This will serve as the upper limit when appraising the asset + let capacity = get_demand_limiting_capacity(time_slice_info, &asset, commodity, demand); + let asset_capacity = AssetCapacity::from_capacity(capacity, asset.unit_size()); + asset.set_capacity(asset_capacity); + + asset.into() + }) +} + +/// Get options from existing and potential assets for the given parameters +fn get_asset_options<'a>( + time_slice_info: &'a TimeSliceInfo, + all_existing_assets: &'a [AssetRef], + demand: &'a DemandMap, + agent: &'a Agent, + commodity: &'a Commodity, + region_id: &'a RegionID, + year: u32, +) -> impl Iterator + 'a { + // Get existing assets which produce the commodity of interest + let existing_assets = all_existing_assets + .iter() + .filter_agent(&agent.id) + .filter_region(region_id) + .filter_primary_producers_of(&commodity.id) + .cloned(); + + // Get candidates assets which produce the commodity of interest + let candidate_assets = + get_candidate_assets(time_slice_info, demand, agent, region_id, commodity, year); + + chain(existing_assets, candidate_assets) +} + +/// Select assets for a single market in a given year +/// +/// Returns a list of assets that are selected for investment for this market in this year. +#[allow(clippy::too_many_arguments)] +fn select_assets_for_single_market( + model: &Model, + commodity_id: &CommodityID, + region_id: &RegionID, + year: u32, + demand: &AllDemandMap, + existing_assets: &[AssetRef], + prices: &Prices, + writer: &mut DataWriter, +) -> Result> { + let commodity = &model.commodities[commodity_id]; + + let mut selected_assets = Vec::new(); + for (agent, commodity_portion) in + get_responsible_agents(model.agents.values(), commodity_id, region_id, year) + { + debug!( + "Running asset selection for agent '{}' in market '{}|{}'", + &agent.id, commodity_id, region_id + ); + + // Get demand portion for this market for this agent in this year + let demand_portion_for_market = get_demand_portion_for_market( + &model.time_slice_info, + demand, + commodity_id, + region_id, + commodity_portion, + ); + + // Existing and candidate assets from which to choose + let opt_assets = get_asset_options( + &model.time_slice_info, + existing_assets, + &demand_portion_for_market, + agent, + commodity, + region_id, + year, + ) + .collect::>(); + + // Calculate investment limits for candidate assets + let investment_limits = + calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); + + // Choose assets from among existing pool and candidates + let best_assets = select_best_assets( + model, + opt_assets, + investment_limits, + commodity, + agent, + region_id, + prices, + demand_portion_for_market, + year, + writer, + )?; + selected_assets.extend(best_assets); + } + + Ok(selected_assets) +} + +/// Iterates through the a pre-ordered set of markets forming a cycle, selecting assets for each +/// market in turn. +/// +/// Dispatch optimisation is performed after each market is visited to rebalance demand. +/// While dispatching, newly selected (`Ready`) assets are given flexible capacity (bounded by +/// `capacity_margin`) so small demand shifts caused by later markets can be absorbed. After all +/// markets have been visited once, the final set of assets is returned, applying any capacity +/// adjustments from the final full-system dispatch optimisation. +/// +/// Dispatch may fail at any point if new demands are encountered for previously visited markets, +/// and the `capacity_margin` is not sufficient to absorb the demand shift. At this point, the +/// simulation is terminated with an error prompting the user to increase the `capacity_margin`. +/// A longer-term solution (TODO) may be to trigger re-investment for the affected markets. Other +/// yet-to-implement features may also help to stabilise the cycle, such as capacity growth limits. +#[allow(clippy::too_many_arguments)] +fn select_assets_for_cycle( + model: &Model, + markets: &[(CommodityID, RegionID)], + year: u32, + demand: &AllDemandMap, + existing_assets: &[AssetRef], + prices: &Prices, + seen_markets: &[(CommodityID, RegionID)], + previously_selected_assets: &[AssetRef], + writer: &mut DataWriter, +) -> Result> { + // Precompute a joined string for logging + let markets_str = markets.iter().map(|(c, r)| format!("{c}|{r}")).join(", "); + + // Iterate over the markets to select assets + let mut current_demand = demand.clone(); + let mut assets_for_cycle = IndexMap::new(); + let mut last_solution = None; + for (idx, (commodity_id, region_id)) in markets.iter().enumerate() { + // Select assets for this market + let assets = select_assets_for_single_market( + model, + commodity_id, + region_id, + year, + ¤t_demand, + existing_assets, + prices, + writer, + )?; + assets_for_cycle.insert((commodity_id.clone(), region_id.clone()), assets); + + // Assemble full list of assets for dispatch (previously selected + all chosen so far) + let mut all_assets = previously_selected_assets.to_vec(); + let assets_for_cycle_flat: Vec<_> = assets_for_cycle + .values() + .flat_map(|v| v.iter().cloned()) + .collect(); + all_assets.extend_from_slice(&assets_for_cycle_flat); + + // We balance all previously seen markets plus all cycle markets up to and including this one + let mut markets_to_balance = seen_markets.to_vec(); + markets_to_balance.extend_from_slice(&markets[0..=idx]); + + // We allow all `Ready` state assets to have flexible capacity + let flexible_capacity_assets: Vec<_> = assets_for_cycle_flat + .iter() + .filter(|asset| matches!(asset.state(), AssetState::Ready { .. })) + .cloned() + .collect(); + + // Retrieve installable capacity limits for flexible capacity assets. + let mut agent_share_cache = HashMap::new(); + let capacity_limits = flexible_capacity_assets + .iter() + .filter_map(|asset| { + let agent_id = asset.agent_id().unwrap(); + let commodity_id = asset.primary_output_commodity().unwrap(); + let agent_share = *agent_share_cache + .entry((agent_id, commodity_id)) + .or_insert_with(|| { + model.agents[agent_id].commodity_portions[&(commodity_id.clone(), year)] + }); + asset + .max_installable_capacity(agent_share) + .map(|max_capacity| (asset.clone(), max_capacity)) + }) + .collect::>(); + + // Run dispatch + let solution = DispatchRun::new(model, &all_assets, year) + .with_market_balance_subset(&markets_to_balance) + .with_flexible_capacity_assets( + &flexible_capacity_assets, + Some(&capacity_limits), + // Gives newly selected cycle assets limited capacity wiggle-room; existing assets stay fixed. + model.parameters.capacity_margin, + ) + .run( + &format!("cycle ({markets_str}) post {commodity_id}|{region_id} investment"), + writer, + ) + .with_context(|| { + format!( + "Cycle balancing failed for cycle ({markets_str}), capacity_margin: {}. \ + Try increasing the capacity_margin.", + model.parameters.capacity_margin + ) + })?; + + // Calculate new net demand map with all assets selected so far + current_demand.clone_from(demand); + update_net_demand_map( + &mut current_demand, + &solution.create_flow_map(), + &assets_for_cycle_flat, + ); + last_solution = Some(solution); + } + + // Finally, update flexible capacity assets based on the final solution + let mut all_cycle_assets: Vec<_> = assets_for_cycle.into_values().flatten().collect(); + if let Some(solution) = last_solution { + let new_capacities: HashMap<_, _> = solution.iter_capacity().collect(); + for asset in &mut all_cycle_assets { + if let Some(new_capacity) = new_capacities.get(asset) { + debug!( + "Capacity of asset '{}' modified during cycle balancing ({} to {})", + asset.process_id(), + asset.total_capacity(), + new_capacity.total_capacity() + ); + asset.make_mut().set_capacity(*new_capacity); + } + } + } + + Ok(all_cycle_assets) +} + +fn get_demand_portion_for_market( + time_slice_info: &TimeSliceInfo, + demand: &AllDemandMap, + commodity_id: &CommodityID, + region_id: &RegionID, + commodity_portion: Dimensionless, +) -> DemandMap { + time_slice_info + .iter_ids() + .map(|time_slice| { + ( + time_slice.clone(), + commodity_portion + * *demand + .get(&(commodity_id.clone(), region_id.clone(), time_slice.clone())) + .unwrap_or(&Flow(0.0)), + ) + }) + .collect() +} + +/// Get the agents responsible for a given market in a given year along with the commodity +/// portion for which they are responsible +fn get_responsible_agents<'a, I>( + agents: I, + commodity_id: &'a CommodityID, + region_id: &'a RegionID, + year: u32, +) -> impl Iterator +where + I: Iterator, +{ + agents.filter_map(move |agent| { + if !agent.regions.contains(region_id) { + return None; + } + let portion = agent + .commodity_portions + .get(&(commodity_id.clone(), year))?; + + Some((agent, *portion)) + }) +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::agent::AgentID; + use crate::fixture::{ + agent_id, process, process_activity_limits_map, process_flows_map, + process_investment_constraints, process_parameter_map, region_id, + }; + use crate::process::{ + Process, ProcessActivityLimitsMap, ProcessFlowsMap, ProcessInvestmentConstraint, + ProcessInvestmentConstraintsMap, ProcessParameterMap, + }; + use crate::region::RegionID; + use crate::units::Dimensionless; + use crate::units::{ActivityPerCapacity, Capacity}; + use indexmap::IndexSet; + use rstest::rstest; + use std::rc::Rc; + + #[rstest] + fn calculate_investment_limits_for_candidates_empty_list() { + // Test with empty list of assets + let opt_assets: Vec = vec![]; + let commodity_portion = Dimensionless(1.0); + + let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); + + assert!(result.is_empty()); + } + + #[rstest] + fn calculate_investment_limits_for_candidates_commissioned_assets_filtered( + process: Process, + region_id: RegionID, + agent_id: AgentID, + ) { + // Create a mix of commissioned and candidate assets + let process_rc = Rc::new(process); + let capacity = Capacity(10.0); + + // Create commissioned asset - should be filtered out + let commissioned_asset = Asset::new_commissioned( + agent_id.clone(), + process_rc.clone(), + region_id.clone(), + capacity, + 2015, + ) + .unwrap(); + + // Create candidate asset - should be included + let candidate_asset = + Asset::new_candidate(process_rc.clone(), region_id.clone(), capacity, 2015).unwrap(); + + let candidate_asset_ref = AssetRef::from(candidate_asset); + let opt_assets = vec![ + AssetRef::from(commissioned_asset), + candidate_asset_ref.clone(), + ]; + let commodity_portion = Dimensionless(1.0); + + let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); + + // Only the candidate asset should be in the result + assert_eq!(result.len(), 1); + assert!(result.contains_key(&candidate_asset_ref)); + } + + #[rstest] + fn calculate_investment_limits_for_candidates_no_investment_constraints( + process: Process, + region_id: RegionID, + ) { + // Create candidate asset without investment constraints + let process_rc = Rc::new(process); + let capacity = Capacity(15.0); + + let candidate_asset = Asset::new_candidate(process_rc, region_id, capacity, 2015).unwrap(); + + let opt_assets = vec![AssetRef::from(candidate_asset.clone())]; + let commodity_portion = Dimensionless(0.8); + + let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); + + // Should return the asset's original capacity since no constraints apply + assert_eq!(result.len(), 1); + let asset_ref = AssetRef::from(candidate_asset); + assert_eq!(result[&asset_ref], AssetCapacity::Continuous(capacity)); + } + + #[rstest] + // Asset capacity higher than constraint -> limited by constraint + #[case(Capacity(15.0), Capacity(10.0))] + // Asset capacity lower than constraint -> limited by asset capacity + #[case(Capacity(5.0), Capacity(5.0))] + fn calculate_investment_limits_for_candidates_with_constraints( + region_id: RegionID, + process_activity_limits_map: ProcessActivityLimitsMap, + process_flows_map: ProcessFlowsMap, + process_parameter_map: ProcessParameterMap, + #[case] asset_capacity: Capacity, + #[case] expected_limit: Capacity, + ) { + let region_ids: IndexSet = [region_id.clone()].into(); + + // Add investment constraint with addition limit + let constraint = ProcessInvestmentConstraint { + addition_limit: Some(Capacity(10.0)), + }; + let mut constraints = ProcessInvestmentConstraintsMap::new(); + constraints.insert((region_id.clone(), 2015), Rc::new(constraint)); + + let process = Process { + id: "constrained_process".into(), + description: "Process with constraints".into(), + years: 2010..=2020, + activity_limits: process_activity_limits_map, + flows: process_flows_map, + parameters: process_parameter_map, + regions: region_ids, + primary_output: None, + capacity_to_activity: ActivityPerCapacity(1.0), + investment_constraints: constraints, + unit_size: None, + }; + + let process_rc = Rc::new(process); + + let candidate_asset = + Asset::new_candidate(process_rc, region_id, asset_capacity, 2015).unwrap(); + + let opt_assets = vec![AssetRef::from(candidate_asset.clone())]; + let commodity_portion = Dimensionless(1.0); + + let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); + + // Should be limited by the minimum of asset capacity and constraint + assert_eq!(result.len(), 1); + let asset_ref = AssetRef::from(candidate_asset); + assert_eq!( + result[&asset_ref], + AssetCapacity::Continuous(expected_limit) + ); + } + + #[rstest] + fn calculate_investment_limits_for_candidates_multiple_assets( + region_id: RegionID, + process_activity_limits_map: ProcessActivityLimitsMap, + process_flows_map: ProcessFlowsMap, + process_parameter_map: ProcessParameterMap, + ) { + let region_ids: IndexSet = [region_id.clone()].into(); + + // Create first process with constraints + let constraint1 = ProcessInvestmentConstraint { + addition_limit: Some(Capacity(12.0)), + }; + let mut constraints1 = ProcessInvestmentConstraintsMap::new(); + constraints1.insert((region_id.clone(), 2015), Rc::new(constraint1)); + + let process1 = Process { + id: "process1".into(), + description: "First process".into(), + years: 2010..=2020, + activity_limits: process_activity_limits_map.clone(), + flows: process_flows_map.clone(), + parameters: process_parameter_map.clone(), + regions: region_ids.clone(), + primary_output: None, + capacity_to_activity: ActivityPerCapacity(1.0), + investment_constraints: constraints1, + unit_size: None, + }; + + // Create second process without constraints + let process2 = Process { + id: "process2".into(), + description: "Second process".into(), + years: 2010..=2020, + activity_limits: process_activity_limits_map, + flows: process_flows_map, + parameters: process_parameter_map, + regions: region_ids, + primary_output: None, + capacity_to_activity: ActivityPerCapacity(1.0), + investment_constraints: process_investment_constraints(), + unit_size: None, + }; + + let process1_rc = Rc::new(process1); + let process2_rc = Rc::new(process2); + + let candidate1 = + Asset::new_candidate(process1_rc, region_id.clone(), Capacity(20.0), 2015).unwrap(); + + let candidate2 = Asset::new_candidate(process2_rc, region_id, Capacity(8.0), 2015).unwrap(); + + let opt_assets = vec![ + AssetRef::from(candidate1.clone()), + AssetRef::from(candidate2.clone()), + ]; + let commodity_portion = Dimensionless(0.75); + + let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); + + // Should have both assets in result + assert_eq!(result.len(), 2); + + // First asset should be limited by constraint: 12.0 * 0.75 = 9.0 + let asset1_ref = AssetRef::from(candidate1); + assert_eq!( + result[&asset1_ref], + AssetCapacity::Continuous(Capacity(9.0)) + ); + + // Second asset should use its original capacity (no constraints) + let asset2_ref = AssetRef::from(candidate2); + assert_eq!( + result[&asset2_ref], + AssetCapacity::Continuous(Capacity(8.0)) + ); + } + #[rstest] + fn calculate_investment_limits_for_candidates_discrete_capacity( + region_id: RegionID, + process_activity_limits_map: crate::process::ProcessActivityLimitsMap, + process_flows_map: crate::process::ProcessFlowsMap, + process_parameter_map: crate::process::ProcessParameterMap, + ) { + let region_ids: IndexSet = [region_id.clone()].into(); + + // Add investment constraint + let constraint = ProcessInvestmentConstraint { + addition_limit: Some(Capacity(35.0)), // Enough for 3.5 units at 10.0 each + }; + let mut constraints = ProcessInvestmentConstraintsMap::new(); + constraints.insert((region_id.clone(), 2015), Rc::new(constraint)); + + let process = Process { + id: "discrete_process".into(), + description: "Process with discrete units".into(), + years: 2010..=2020, + activity_limits: process_activity_limits_map, + flows: process_flows_map, + parameters: process_parameter_map, + regions: region_ids, + primary_output: None, + capacity_to_activity: ActivityPerCapacity(1.0), + investment_constraints: constraints, + unit_size: Some(Capacity(10.0)), // Discrete units of 10.0 capacity each + }; + + let process_rc = Rc::new(process); + let capacity = Capacity(50.0); // 5 units at 10.0 each + + let candidate_asset = Asset::new_candidate(process_rc, region_id, capacity, 2015).unwrap(); + + let opt_assets = vec![AssetRef::from(candidate_asset.clone())]; + let commodity_portion = Dimensionless(1.0); + + let result = calculate_investment_limits_for_candidates(&opt_assets, commodity_portion); + + // Should be limited by constraint and rounded down to whole units + // Constraint: 35.0, divided by unit size 10.0 = 3.5 -> floor to 3 units = 30.0 + assert_eq!(result.len(), 1); + let asset_ref = AssetRef::from(candidate_asset); + assert_eq!( + result[&asset_ref], + AssetCapacity::Discrete(3, Capacity(10.0)) + ); + assert_eq!(result[&asset_ref].total_capacity(), Capacity(30.0)); + } +} diff --git a/src/simulation/prices.rs b/src/simulation/prices.rs index 750fc30fb..7445624cd 100644 --- a/src/simulation/prices.rs +++ b/src/simulation/prices.rs @@ -3,7 +3,7 @@ use crate::asset::AssetRef; use crate::commodity::{CommodityID, CommodityMap, PricingStrategy}; use crate::model::Model; use crate::region::RegionID; -use crate::simulation::investment::InvestmentSet; +use crate::simulation::market::MarketSet; use crate::simulation::optimisation::Solution; use crate::time_slice::{TimeSliceID, TimeSliceInfo, TimeSliceSelection}; use crate::units::{Activity, Dimensionless, Flow, MoneyPerActivity, MoneyPerFlow, UnitType, Year}; @@ -118,10 +118,10 @@ pub fn calculate_prices(model: &Model, solution: &Solution, year: u32) -> Result // Lazily computed only if at least one FullCost market is encountered. let mut annual_activities: Option> = None; - // Iterate over investment sets in reverse order - for investment_set in model.investment_order[&year].iter().rev() { - price_investment_set( - investment_set, + // Iterate over market sets in reverse order + for market_set in model.investment_order[&year].iter().rev() { + price_market_set( + market_set, model, solution, year, @@ -137,9 +137,9 @@ pub fn calculate_prices(model: &Model, solution: &Solution, year: u32) -> Result }) } -/// Calculate prices for the markets in an investment set, updating `market_prices`. -fn price_investment_set( - investment_set: &InvestmentSet, +/// Calculate prices for the markets in an market set, updating `market_prices`. +fn price_market_set( + market_set: &MarketSet, model: &Model, solution: &Solution, year: u32, @@ -147,8 +147,8 @@ fn price_investment_set( annual_activities: &mut Option>, market_prices: &mut PriceMap, ) { - match investment_set { - InvestmentSet::Single(market) => { + match market_set { + MarketSet::Single(market) => { price_markets( model, solution, @@ -159,7 +159,7 @@ fn price_investment_set( market_prices, ); } - InvestmentSet::Cycle(markets) => { + MarketSet::Cycle(markets) => { price_cycle( model, solution, @@ -170,9 +170,9 @@ fn price_investment_set( market_prices, ); } - InvestmentSet::Layer(investment_sets) => { - for set in investment_sets { - price_investment_set( + MarketSet::Layer(market_sets) => { + for set in market_sets { + price_market_set( set, model, solution,