diff --git a/gtep/driver_esr.py b/gtep/driver_esr.py index 365117df..e3c34496 100644 --- a/gtep/driver_esr.py +++ b/gtep/driver_esr.py @@ -22,7 +22,16 @@ # Add data data_path = "./data/5bus" -data_object = ExpansionPlanningData() +data_object = ExpansionPlanningData( + stages=2, + num_reps=2, + num_commit=6, + num_dispatch=4, + duration_representative_period=6, + save_period_structure_file=False, + period_structure_json_file=None, + # period_structure_json_file="period_structure_from_gtep.json", +) data_object.load_prescient(data_path) # [ESR WIP: Add bus and cost data files to be used on the @@ -52,15 +61,8 @@ # Populate and create GTEP model mod_object = ExpansionPlanningModel( - stages=2, data=data_object, cost_data=data_processing_object, - num_reps=2, - len_reps=1, - num_commit=6, - num_dispatch=4, - # [ESR: in min by default, for now] - duration_dispatch=15, ) mod_object.config["include_investment"] = True diff --git a/gtep/gtep_data.py b/gtep/gtep_data.py index 251fcbc6..d6a1ea70 100644 --- a/gtep/gtep_data.py +++ b/gtep/gtep_data.py @@ -35,23 +35,33 @@ def __init__( len_reps=1, num_commit=24, num_dispatch=1, - duration_dispatch=60, + duration_representative_period=24, + save_period_structure_file=False, + period_structure_json_file=None, ): """Initialize generation & expansion planning data object. - :param stages: integer number of investment periods - :param num_reps: integer number of representative periods per investment period - :param len_reps: (for now integer) length of each representative period (in hours) - :param num_commit: integer number of commitment periods per representative period - :param num_dispatch: integer number of dispatch periods per commitment period - :param duration_dispatch: (for now integer) duration of each dispatch period (in minutes) + :param: stages: integer number of investment periods + :param: num_reps: integer number of representative periods per investment period + :param: num_commit: integer number of commitment periods per representative period + :param: num_dispatch: integer number of dispatch periods per commitment period + :param: duration_representative_period: duration of each representative period + (in hours) + :param: save_period_structure_file: (optional) If True, saves the generated + period structure as a JSON file in the data directory. Default is False. + :param: period_structure_json_file: (optional) Path to a JSON file in the data + directory specifying the period structure. Overrides scalar/list arguments + if provided. Default is None. + """ self.stages = stages self.num_reps = num_reps self.len_reps = len_reps self.num_commit = num_commit self.num_dispatch = num_dispatch - self.duration_dispatch = duration_dispatch + self.duration_representative_period = duration_representative_period + self.save_period_structure_file = save_period_structure_file + self.period_structure_json_file = period_structure_json_file def load_prescient( self, diff --git a/gtep/gtep_model.py b/gtep/gtep_model.py index 321ee19e..88ecd61d 100644 --- a/gtep/gtep_model.py +++ b/gtep/gtep_model.py @@ -27,6 +27,7 @@ import json import numpy as np import re +import os import pyomo.environ as pyo from pyomo.environ import units as u @@ -52,6 +53,13 @@ import gtep.model_library.storage as stor import gtep.model_library.transmission as transm +from gtep.utils import ( + _set_period_structure_dict, + check_period_structure_consistency, +) + +curr_dir = os.path.dirname(os.path.abspath(__file__)) + # Define what a USD is for pyomo units purposes. This will be set to a # base year and we will do NPV calculations based on automatic Pyomo # unit transformations. @@ -99,17 +107,64 @@ def __init__( :return: Pyomo model for full GTEP """ + self.config = _get_model_config() + self.timer = TicTocTimer() + self.stages = data.stages self.formulation = formulation self.data = data self.cost_data = cost_data - self.num_reps = data.num_reps - self.len_reps = data.len_reps - self.num_commit = data.num_commit - self.num_dispatch = data.num_dispatch - self.duration_dispatch = data.duration_dispatch - self.config = _get_model_config() - self.timer = TicTocTimer() + self.save_period_structure_file = data.save_period_structure_file + self.period_structure_json_file = data.period_structure_json_file + + # Calculate commitment and dispatch period durations using + # provided scalars for the duration of representative period + # and number of commitment and dispatch periods. Note: We are + # assuming that all periods within their parent are of equal + # length. + duration_commitment = data.duration_representative_period / data.num_commit + duration_dispatch = pyo.value( + pyo.units.convert( + (duration_commitment / data.num_dispatch) * u.hours, + to_units=u.minutes, + ) + ) + + # Generate the period structure dictionary from provided + # scalars, or load it from a .json file if specified. In + # either case, the period structure is returned as a + # dictionary. + period_dict = _set_period_structure_dict( + data.num_reps, + data.num_commit, + data.num_dispatch, + data.duration_representative_period, + duration_commitment, + duration_dispatch, + self.save_period_structure_file, + self.period_structure_json_file, + ) + + # Assign period structure attributes from the dictionary + self.num_reps = period_dict.get("number_representative", data.num_reps) + self.num_commit = period_dict["number_commitment"] + self.num_dispatch = period_dict["number_dispatch"] + self.duration_representative_period = period_dict[ + "duration_representative_period" + ] + self.duration_commitment = period_dict["duration_commitment"] + self.duration_dispatch = period_dict["duration_dispatch"] + + # Run a consistency check on commitment and dispatch + # durations. + check_period_structure_consistency( + self.num_reps, + self.num_commit, + self.num_dispatch, + self.duration_representative_period, + self.duration_commitment, + self.duration_dispatch, + ) _add_common_configs(self.config) _add_investment_configs(self.config) @@ -151,15 +206,18 @@ def create_model(self): m, self.stages, rep_per=[i for i in range(1, self.num_reps + 1)] ) - comps.add_model_parameters( + comps.add_model_parameters(m) + + create_stages( m, + self.stages, self.num_commit, self.num_dispatch, + self.duration_representative_period, + self.duration_commitment, self.duration_dispatch, ) - create_stages(m, self.stages) - obj_comp.create_objective_function(m) self.model = m @@ -202,8 +260,24 @@ def report_large_coefficients(self, outfile, magnitude_cutoff=1e5): with open(outfile, "w") as fil: json.dump(really_bad_var_coef_list, fil) + return ( + num_commit_dict, + num_dispatch_dict, + duration_rep_dict, + duration_commit_dict, + duration_dispatch_dict, + ) + -def create_stages(m, stages): +def create_stages( + m, + stages, + num_commit, + num_dispatch, + duration_representative_period, + duration_commitment, + duration_dispatch, +): """This method constructs the block structure for the Generation and Transmission Expansion Planning (GTEP) model. It creates investment, representative period, and commitment blocks for each @@ -247,6 +321,13 @@ def create_stages(m, stages): # representative_period. for representative_period in b_inv.representativePeriods: b_rep = b_inv.representativePeriod[representative_period] + + b_rep.representativePeriodLength = pyo.Param( + initialize=duration_representative_period[representative_period], + within=pyo.PositiveReals, + units=u.hr, + ) + b_rep.representative_date = m.data.representative_dates[ representative_period - 1 ] @@ -258,15 +339,12 @@ def create_stages(m, stages): # b_rep.day = int(broken_date[2]) b_rep.currentPeriod = representative_period - # [ESR NOTE: Include commitment blocks regardless of the - # value of include_commitment. When include_commitment is - # False, generators are assumed to be always online, and + # Include commitment blocks regardless of the value of + # include_commitment. When False, generators are on and # operational costs are determined solely by dispatch - # decisions. No redispatch for now.] - # if m.config["include_commitment"] or m.config["include_redispatch"]: - b_rep.commitmentPeriods = pyo.RangeSet( - m.numCommitmentPeriods[representative_period] - ) + # decisions. + n_commit = num_commit[representative_period] + b_rep.commitmentPeriods = pyo.RangeSet(n_commit) b_rep.commitmentPeriod = pyo.Block(b_rep.commitmentPeriods) # --.--.--.--.--.--.----.--.--.--.--.--.----.--.--.--.--.--.-- @@ -274,12 +352,20 @@ def create_stages(m, stages): # constraints for commitment_period in b_rep.commitmentPeriods: b_comm = b_rep.commitmentPeriod[commitment_period] + + b_comm.commitmentPeriodLength = pyo.Param( + initialize=duration_commitment[representative_period][ + commitment_period + ], + within=pyo.PositiveReals, + units=u.hr, + ) + b_comm.commitmentPeriod = commitment_period if m.config["include_redispatch"]: - b_comm.dispatchPeriods = pyo.RangeSet( - m.numDispatchPeriods[b_rep.currentPeriod] - ) + n_dispatch = num_dispatch[representative_period][commitment_period] + b_comm.dispatchPeriods = pyo.RangeSet(n_dispatch) b_comm.dispatchPeriod = pyo.Block(b_comm.dispatchPeriods) # [TODO: update properties for this time period!] @@ -299,21 +385,22 @@ def create_stages(m, stages): # [TODO: This feels REALLY inelegant and # bad. Check a better way of declaring these.] - for period in b_comm.dispatchPeriods: - b_comm.dispatchPeriod[period].periodLength = pyo.Param( - initialize=1, + for dispatch_period in b_comm.dispatchPeriods: + b_disp = b_comm.dispatchPeriod[dispatch_period] + b_disp.dispatchPeriodLength = pyo.Param( + initialize=duration_dispatch[representative_period][ + commitment_period + ][dispatch_period], within=pyo.PositiveReals, - # units=u.minutes, + units=u.minutes, ) disp.add_dispatch_variables( - b_comm.dispatchPeriod[period], - period, - m.dispatchPeriodLength, - ) - disp.add_dispatch_constraints( - b_comm.dispatchPeriod[period], period + b_disp, + dispatch_period, + b_disp.dispatchPeriodLength, ) + disp.add_dispatch_constraints(b_disp, dispatch_period) # =.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=.=. @@ -338,8 +425,8 @@ def create_stages(m, stages): rep_period.add_representative_period_variables(b_rep, representative_period) if m.config["include_commitment"]: - # These logical constraint ensure the state - # disjuncts stay consistent. + # These logical constraints ensure the state disjuncts + # stay consistent. rep_period.add_representative_period_logical_constraints( b_rep, representative_period ) diff --git a/gtep/gtep_solution.py b/gtep/gtep_solution.py index d7e5f773..e22b164d 100644 --- a/gtep/gtep_solution.py +++ b/gtep/gtep_solution.py @@ -69,7 +69,6 @@ def load_from_model(self, gtep_model): self.formulation = gtep_model.formulation # None (???) self.data = gtep_model.data # ModelData object self.num_reps = gtep_model.num_reps # int - self.len_reps = gtep_model.len_reps # int self.num_commit = gtep_model.num_commit # int self.num_dispatch = gtep_model.num_dispatch # int diff --git a/gtep/model_library/commitment.py b/gtep/model_library/commitment.py index a8be863b..fae75b41 100644 --- a/gtep/model_library/commitment.py +++ b/gtep/model_library/commitment.py @@ -28,9 +28,6 @@ def add_commitment_parameters(b, commitment_period, investmentStage): m = b.model() - b.commitmentPeriodLength = pyo.Param( - within=pyo.PositiveReals, default=1, units=u.hr - ) b.carbonTax = pyo.Param(default=0) # [ESR: Corrected to be in the commitment block "b", not in main @@ -174,7 +171,9 @@ def renewable_curtailment_cost(b): for com_per in b.representativePeriod[rep_per].commitmentPeriods: renewableCurtailmentRep += ( m.weights[rep_per] - * m.commitmentPeriodLength + * b.representativePeriod[rep_per] + .commitmentPeriod[com_per] + .commitmentPeriodLength * b.representativePeriod[rep_per] .commitmentPeriod[com_per] .renewableCurtailmentCommitment # in MW diff --git a/gtep/model_library/components.py b/gtep/model_library/components.py index 7d6bfdf1..4e06fd7a 100644 --- a/gtep/model_library/components.py +++ b/gtep/model_library/components.py @@ -106,7 +106,7 @@ def add_model_sets(m, stages, rep_per=["a", "b"], com_per=2, dis_per=2): ) -def add_model_parameters(m, num_commit, num_dispatch, duration_dispatch): +def add_model_parameters(m): """Creates and labels all the parameters in the GTEP model. This method ties input data directly to the model. @@ -118,33 +118,6 @@ def add_model_parameters(m, num_commit, num_dispatch, duration_dispatch): # configuration arg and not hardcoded values.] m.years = [2025, 2030, 2035] - # Add parameters related to the representative periods for the - # different stages - m.representativePeriodLength = pyo.Param( - m.representativePeriods, within=pyo.PositiveReals, default=24, units=u.hr - ) - m.numCommitmentPeriods = pyo.Param( - m.representativePeriods, - within=pyo.PositiveIntegers, - default=2, - initialize=num_commit, - ) - m.numDispatchPeriods = pyo.Param( - m.representativePeriods, - within=pyo.PositiveIntegers, - default=2, - initialize=num_dispatch, - ) - m.commitmentPeriodLength = pyo.Param( - within=pyo.PositiveReals, default=1, units=u.hr - ) - - # [TODO: Index by dispatch period? Certainly index by - # commitment period.] - m.dispatchPeriodLength = pyo.Param( - within=pyo.PositiveReals, initialize=duration_dispatch, units=u.minutes - ) - # Add power-related parameters m.thermalCapacity = pyo.Param( m.thermalGenerators, diff --git a/gtep/model_library/gen.py b/gtep/model_library/gen.py index ae0530c0..dbdde32a 100644 --- a/gtep/model_library/gen.py +++ b/gtep/model_library/gen.py @@ -281,9 +281,7 @@ def ramp_up_limits(disj, dispatchPeriod, generator): return ( b.dispatchPeriod[dispatchPeriod].thermalGeneration[generator] - b.dispatchPeriod[dispatchPeriod - 1].thermalGeneration[generator] - <= m.rampUpRates[generator] - * b.dispatchPeriod[dispatchPeriod].periodLength - * m.thermalCapacity[generator] + <= m.rampUpRates[generator] * m.thermalCapacity[generator] ) elif dispatchPeriod == 1 and commitment_period != 1: return ( @@ -291,9 +289,7 @@ def ramp_up_limits(disj, dispatchPeriod, generator): - r_p.commitmentPeriod[commitment_period - 1] .dispatchPeriod[b.dispatchPeriods.last()] .thermalGeneration[generator] - <= m.rampUpRates[generator] - * b.dispatchPeriod[dispatchPeriod].periodLength - * m.thermalCapacity[generator] + <= m.rampUpRates[generator] * m.thermalCapacity[generator] ) else: return pyo.Constraint.Skip @@ -309,7 +305,6 @@ def ramp_down_limits(disj, dispatchPeriod, generator): b.dispatchPeriod[dispatchPeriod - 1].thermalGeneration[generator] - b.dispatchPeriod[dispatchPeriod].thermalGeneration[generator] <= m.rampDownRates[generator] # in MW/min - * b.dispatchPeriod[dispatchPeriod].periodLength # in min * m.thermalCapacity[generator] # in MW ) elif dispatchPeriod == 1 and commitment_period != 1: @@ -319,7 +314,6 @@ def ramp_down_limits(disj, dispatchPeriod, generator): .thermalGeneration[generator] - b.dispatchPeriod[dispatchPeriod].thermalGeneration[generator] <= m.rampDownRates[generator] # in MW/min - * b.dispatchPeriod[dispatchPeriod].periodLength # in min * m.thermalCapacity[generator] # in MW ) else: @@ -369,12 +363,7 @@ def ramp_up_limits(disj, dispatchPeriod, generator): - b.dispatchPeriod[dispatchPeriod - 1].thermalGeneration[generator] <= max( pyo.value(m.thermalMin[generator]), - # [ESR: Make sure the time units are consistent - # here since we are only taking the value] pyo.value(m.rampUpRates[generator]) - * pyo.value( - b.dispatchPeriod[dispatchPeriod].periodLength - ) # in minutes * pyo.value(m.thermalCapacity[generator]), ) * u.MW @@ -387,12 +376,7 @@ def ramp_up_limits(disj, dispatchPeriod, generator): .thermalGeneration[generator] <= max( pyo.value(m.thermalMin[generator]), - # [ESR: Make sure the time units are consistent - # here since we are only taking the value] pyo.value(m.rampUpRates[generator]) - * pyo.value( - b.dispatchPeriod[dispatchPeriod].periodLength - ) # in minutes * pyo.value(m.thermalCapacity[generator]), ) * u.MW @@ -436,12 +420,7 @@ def ramp_down_limits(disj, dispatchPeriod, generator): - b.dispatchPeriod[dispatchPeriod].thermalGeneration[generator] <= max( pyo.value(m.thermalMin[generator]), - # [ESR: Make sure the time units are consistent - # here since we are taking the value only] pyo.value(m.rampDownRates[generator]) - * pyo.value( - b.dispatchPeriod[dispatchPeriod].periodLength - ) # in minutes * pyo.value(m.thermalCapacity[generator]), ) * u.MW @@ -454,12 +433,7 @@ def ramp_down_limits(disj, dispatchPeriod, generator): - b.dispatchPeriod[dispatchPeriod].thermalGeneration[generator] <= max( pyo.value(m.thermalMin[generator]), - # [ESR: Make sure the time units are consistent - # here since we are taking the value only] pyo.value(m.rampDownRates[generator]) - * pyo.value( - b.dispatchPeriod[dispatchPeriod].periodLength - ) # in minutes * pyo.value(m.thermalCapacity[generator]), ) * u.MW diff --git a/gtep/tests/unit/test_gtep_model.py b/gtep/tests/unit/test_gtep_model.py index 58cf77d7..e255506b 100644 --- a/gtep/tests/unit/test_gtep_model.py +++ b/gtep/tests/unit/test_gtep_model.py @@ -11,6 +11,8 @@ # for full copyright and license information. ################################################################################# +import json +from pathlib import Path from os.path import abspath, join, dirname import pytest @@ -23,6 +25,8 @@ _component_data_handlers, assert_units_equivalent, ) +from pyomo.common.tempfiles import TempfileManager + from gtep.gtep_model import ExpansionPlanningModel from gtep.gtep_data import ExpansionPlanningData from gtep.gtep_data_processing import DataProcessing @@ -50,20 +54,22 @@ def patch_unit_handlers(): def read_debug_model( stages=1, num_reps=3, - len_reps=24, num_commit=24, num_dispatch=4, - duration_dispatch=15, + duration_representative_period=24, + save_period_structure_file=False, + period_structure_json_file=None, ): curr_dir = dirname(abspath(__file__)) debug_data_path = abspath(join(curr_dir, "..", "..", "data", "5bus")) dataObject = ExpansionPlanningData( stages=stages, num_reps=num_reps, - len_reps=len_reps, num_commit=num_commit, num_dispatch=num_dispatch, - duration_dispatch=duration_dispatch, + duration_representative_period=duration_representative_period, + save_period_structure_file=save_period_structure_file, + period_structure_json_file=period_structure_json_file, ) dataObject.load_prescient(debug_data_path) return dataObject @@ -72,19 +78,21 @@ def read_debug_model( def prepare_model_and_cost_data( stages=1, num_reps=3, - len_reps=24, num_commit=24, num_dispatch=4, - duration_dispatch=15, + duration_representative_period=24, + save_period_structure_file=False, + period_structure_json_file=None, ): # Prepare model and cost data dataObject = read_debug_model( stages, num_reps, - len_reps, num_commit, num_dispatch, - duration_dispatch, + duration_representative_period, + save_period_structure_file, + period_structure_json_file, ) curr_dir = dirname(abspath(__file__)) data_path = abspath(join(curr_dir, "..", "..", "data", "costs")) @@ -114,38 +122,57 @@ def prepare_model_and_cost_data( @pytest.mark.usefixtures("patch_unit_handlers") class TestGTEP(unittest.TestCase): def test_model_init(self): - # Test that the ExpansionPlanningModel object can read a default dataset and init - # properly with default values, including building a Pyomo.ConcreteModel object - data_object = read_debug_model( - stages=1, - num_reps=3, - len_reps=24, - num_commit=24, - num_dispatch=4, - duration_dispatch=60, - ) + # Test that the ExpansionPlanningModel object can read a + # default dataset and init properly with default values, + # including building a Pyomo.ConcreteModel object + data_object = read_debug_model() modObject = ExpansionPlanningModel(data=data_object) self.assertIsInstance(modObject, ExpansionPlanningModel) + modObject.create_model() + self.assertIsInstance(modObject.model, ConcreteModel) self.assertEqual(modObject.stages, 1) self.assertEqual(modObject.formulation, None) self.assertIsInstance(modObject.model.md, ModelData) self.assertEqual(modObject.num_reps, 3) - self.assertEqual(modObject.len_reps, 24) - self.assertEqual(modObject.num_commit, 24) - self.assertEqual(modObject.num_dispatch, 4) - self.assertEqual(modObject.duration_dispatch, 60) + self.assertEqual(modObject.num_commit, {1: 24, 2: 24, 3: 24}) + self.assertEqual( + modObject.num_dispatch, + { + 1: {i: 4 for i in range(1, 25)}, + 2: {i: 4 for i in range(1, 25)}, + 3: {i: 4 for i in range(1, 25)}, + }, + ) + self.assertEqual( + modObject.duration_representative_period, {1: 24, 2: 24, 3: 24} + ) + self.assertEqual( + modObject.duration_commitment, + { + 1: {i: 1 for i in range(1, 25)}, + 2: {i: 1 for i in range(1, 25)}, + 3: {i: 1 for i in range(1, 25)}, + }, + ) + self.assertEqual( + modObject.duration_dispatch, + { + 1: {com: {disp: 15 for disp in range(1, 5)} for com in range(1, 25)}, + 2: {com: {disp: 15 for disp in range(1, 5)} for com in range(1, 25)}, + 3: {com: {disp: 15 for disp in range(1, 5)} for com in range(1, 25)}, + }, + ) - # Test that the ExpansionPlanningModel object can read a default dataset and init - # properly with non-default values + # Test that the ExpansionPlanningModel object can read a + # default dataset and init properly with non-default values data_object = read_debug_model( stages=2, num_reps=4, - len_reps=16, num_commit=12, - num_dispatch=12, - duration_dispatch=30, + num_dispatch=24, + duration_representative_period=24, ) modObject = ExpansionPlanningModel( data=data_object, @@ -157,10 +184,37 @@ def test_model_init(self): self.assertEqual(modObject.formulation, None) self.assertIsInstance(modObject.model.md, ModelData) self.assertEqual(modObject.num_reps, 4) - self.assertEqual(modObject.len_reps, 16) - self.assertEqual(modObject.num_commit, 12) - self.assertEqual(modObject.num_dispatch, 12) - self.assertEqual(modObject.duration_dispatch, 30) + self.assertEqual(modObject.num_commit, {1: 12, 2: 12, 3: 12, 4: 12}) + self.assertEqual( + modObject.num_dispatch, + { + 1: {i: 24 for i in range(1, 13)}, + 2: {i: 24 for i in range(1, 13)}, + 3: {i: 24 for i in range(1, 13)}, + 4: {i: 24 for i in range(1, 13)}, + }, + ) + self.assertEqual( + modObject.duration_representative_period, {1: 24, 2: 24, 3: 24, 4: 24} + ) + self.assertEqual( + modObject.duration_commitment, + { + 1: {i: 2 for i in range(1, 13)}, + 2: {i: 2 for i in range(1, 13)}, + 3: {i: 2 for i in range(1, 13)}, + 4: {i: 2 for i in range(1, 13)}, + }, + ) + self.assertEqual( + modObject.duration_dispatch, + { + 1: {com: {disp: 5 for disp in range(1, 25)} for com in range(1, 13)}, + 2: {com: {disp: 5 for disp in range(1, 25)} for com in range(1, 13)}, + 3: {com: {disp: 5 for disp in range(1, 25)} for com in range(1, 13)}, + 4: {com: {disp: 5 for disp in range(1, 25)} for com in range(1, 13)}, + }, + ) # We have expansion blocks and they are where and what we think they are expansion_blocks = modObject.model.component("investmentStage") @@ -196,9 +250,9 @@ def test_model_unit_consistency(self): data_object = read_debug_model( stages=2, num_reps=2, - len_reps=2, num_commit=2, num_dispatch=2, + duration_representative_period=2, ) modObject = ExpansionPlanningModel( data=data_object, @@ -210,54 +264,36 @@ def test_model_unit_consistency(self): assert_units_consistent(m) # Check that subset of model components have expected units + m_inv = m.investmentStage[1] + m_rep_period = m_inv.representativePeriod[1] + m_commit = m_rep_period.commitmentPeriod[1] + m_disp = m_commit.dispatchPeriod[1] + + assert_units_equivalent(m_rep_period.representativePeriodLength, u.h) + assert_units_equivalent(m_commit.commitmentPeriodLength, u.h) + assert_units_equivalent(m_disp.dispatchPeriodLength, u.min) assert_units_equivalent(m.renewable_capacity_enforcement[1, "10_PV"].expr, u.MW) - assert_units_equivalent( - m.investmentStage[1].renewable_curtailment_cost.expr, u.USD - ) - assert_units_equivalent( - m.investmentStage[1] - .representativePeriod[1] - .commitmentPeriod[1] - .dispatchPeriod[1] - .flow_balance["bus1"] - .expr, - u.MW, - ) - assert_units_equivalent(m.commitmentPeriodLength, u.h) - assert_units_equivalent(m.dispatchPeriodLength, u.min) + assert_units_equivalent(m_inv.renewable_curtailment_cost.expr, u.USD) + assert_units_equivalent(m_disp.flow_balance["bus1"].expr, u.MW) assert_units_equivalent(m.rampUpRates, u.dimensionless) assert_units_equivalent(m.varCost, u.USD / u.h / u.MW) + assert_units_equivalent(m_disp.spinningReserve, u.MW) assert_units_equivalent( - m.investmentStage[1] - .representativePeriod[1] - .commitmentPeriod[1] - .dispatchPeriod[1] - .spinningReserve, + m_commit.genOn["3_CT"].operating_limit_min[1].expr, u.MW, ) assert_units_equivalent( - m.investmentStage[1] - .representativePeriod[1] - .commitmentPeriod[1] - .genOn["3_CT"] - .operating_limit_min[1] - .expr, - u.MW, - ) - assert_units_equivalent( - m.investmentStage[1] - .representativePeriod[1] - .commitmentPeriod[1] - .genOn["4_STEAM"] - .max_spinning_reserve[1, "4_STEAM"] - .expr, + m_commit.genOn["4_STEAM"].max_spinning_reserve[1, "4_STEAM"].expr, u.MW, ) def test_solve_bigm(self): # Solve the debug model as is data_object = read_debug_model( - num_reps=1, len_reps=1, num_commit=1, num_dispatch=1 + num_reps=1, + num_commit=1, + num_dispatch=1, + duration_representative_period=1, ) modObject = ExpansionPlanningModel(data=data_object) modObject.create_model() @@ -275,9 +311,10 @@ def test_solve_bigm(self): modObject.results = opt.solve(modObject.model) - # previous successful objective values: 9207.95, 6078.86, 531860.15, 531883.43 + # Previous successful objective values: 9207.95, 6078.86, + # 531860.15, 531883.43 self.assertAlmostEqual( - value(modObject.model.total_cost_objective_rule), 531883.43, places=1 + value(modObject.model.total_cost_objective_rule), 2127462.53, places=1 ) assert_units_equivalent(modObject.model.total_cost_objective_rule.expr, u.USD) @@ -285,9 +322,9 @@ def test_no_investment(self): # Solve the debug model with no investment data_object = read_debug_model( num_reps=1, - len_reps=1, num_commit=1, num_dispatch=1, + duration_representative_period=1, ) modObject = ExpansionPlanningModel( data=data_object, @@ -312,7 +349,7 @@ def test_no_investment(self): # previous successful objective values: 531860.15, 531883.43 self.assertAlmostEqual( - value(modObject.model.total_cost_objective_rule), 531883.43, places=1 + value(modObject.model.total_cost_objective_rule), 2127462.53, places=1 ) assert_units_equivalent(modObject.model.total_cost_objective_rule.expr, u.USD) @@ -323,10 +360,9 @@ def test_with_cost_data_and_commitment(self): dataObject, dataProcessingObject = prepare_model_and_cost_data( stages=2, num_reps=2, - len_reps=1, num_commit=6, num_dispatch=4, - duration_dispatch=15, + duration_representative_period=6, ) # Populate and create GTEP model @@ -371,10 +407,9 @@ def test_with_cost_data_and_no_commitment(self): dataObject, dataProcessingObject = prepare_model_and_cost_data( stages=2, num_reps=2, - len_reps=1, num_commit=6, num_dispatch=4, - duration_dispatch=15, + duration_representative_period=6, ) # Populate and create GTEP model @@ -411,3 +446,129 @@ def test_with_cost_data_and_no_commitment(self): ) assert_units_equivalent(modObject.model.total_cost_objective_rule.expr, u.USD) + + def test_period_structure_from_scalars(self): + # Test period structure dictionary created using the provided + # scalars + dataObject, dataProcessingObject = prepare_model_and_cost_data( + num_reps=2, + num_commit=3, + num_dispatch=4, + duration_representative_period=3, + ) + + modObject = ExpansionPlanningModel( + data=dataObject, cost_data=dataProcessingObject + ) + + # Assert that all values are as expected + self.assertEqual(modObject.num_commit[1], 3) + self.assertEqual(modObject.num_dispatch[2][3], 4) + self.assertEqual(modObject.duration_representative_period[1], 3) + self.assertEqual(modObject.duration_commitment[1][2], 1) + self.assertEqual(modObject.duration_dispatch[2][3][4], 15) + + # Test custom period structure dictionary with irregular + # values. + period_dict = { + "number_representative": 2, + "number_commitment": {1: 2, 2: 3}, + "number_dispatch": {1: {1: 3, 2: 2}, 2: {1: 2, 2: 3, 3: 2}}, + "duration_representative_period": {1: 24, 2: 18}, + "duration_commitment": {1: {1: 12, 2: 12}, 2: {1: 6, 2: 6, 3: 6}}, + "duration_dispatch": { + 1: {1: {1: 360, 2: 180, 3: 180}, 2: {1: 360, 2: 360}}, + 2: { + 1: {1: 180, 2: 180}, + 2: {1: 120, 2: 120, 3: 120}, + 3: {1: 180, 2: 180}, + }, + }, + } + + with TempfileManager.new_context() as tempfile: + temp_dir = Path(tempfile.mkdtemp()) + json_path = temp_dir / "test_custom_period_structure.json" + with open(json_path, "w") as f: + json.dump(period_dict, f, indent=2) + + # Test that the model correctly reads and assigns the + # custom period structure values. For this, we instantiate + # the model using the temp .json file. + dataObject, dataProcessingObject = prepare_model_and_cost_data( + period_structure_json_file=str(json_path), + ) + + modObject = ExpansionPlanningModel( + data=dataObject, cost_data=dataProcessingObject + ) + + # Assert that we have the correct reading of the structure + self.assertEqual(modObject.num_reps, 2) + self.assertEqual(modObject.num_commit[2], 3) + self.assertEqual(modObject.num_dispatch[2][2], 3) + self.assertEqual(modObject.duration_commitment[2][2], 6) + self.assertEqual(modObject.duration_dispatch[2][2][3], 120) + self.assertEqual(modObject.duration_representative_period[2], 18) + self.assertEqual(modObject.duration_dispatch[1][1][2], 180) + self.assertEqual(modObject.duration_commitment[1][2], 12) + + def test_period_structure_consistency_errors(self): + + # Check that a consistency error is raised if + # commitment/dispatch period durations do not sum to the + # representative/commitment period duration when loading from + # a .json file. The test first creates and writes the .json + # file, then checks for the error. + period_dict = { + "number_representative": 2, + "number_commitment": {"1": 6, "2": 6}, + "number_dispatch": { + "1": {"1": 4, "2": 4, "3": 4, "4": 4, "5": 4, "6": 4}, + "2": {"1": 4, "2": 4, "3": 4, "4": 4, "5": 4, "6": 4}, + }, + "duration_representative_period": {"1": 6, "2": 6}, + "duration_commitment": { + "1": {"1": 1.0, "2": 1.0, "3": 1.0, "4": 1.0, "5": 1.0, "6": 3.0}, + "2": {"1": 1.0, "2": 1.0, "3": 2.0, "4": 1.0, "5": 1.0, "6": 1.0}, + }, + "duration_dispatch": { + "1": { + "1": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + "2": {"1": 25.0, "2": 15.0, "3": 15.0, "4": 15.0}, + "3": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + "4": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + "5": {"1": 15.0, "2": 15.0, "3": 5.0, "4": 15.0}, + "6": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + }, + "2": { + "1": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + "2": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + "3": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + "4": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + "5": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + "6": {"1": 15.0, "2": 15.0, "3": 15.0, "4": 15.0}, + }, + }, + } + + with TempfileManager.new_context() as tempfile: + temp_dir = Path(tempfile.mkdtemp()) + json_path = temp_dir / "test_consistency_errors_period_structure.json" + with open(json_path, "w") as f: + json.dump(period_dict, f, indent=2) + + # Instantiate the model using the temp .json file. + dataObject, dataProcessingObject = prepare_model_and_cost_data( + period_structure_json_file=str(json_path), + ) + + # Assert that the sum of commitment durations (20hr) does + # not match the representative period duration (18hr) + with self.assertRaises(ValueError) as cm: + ExpansionPlanningModel(data=dataObject, cost_data=dataProcessingObject) + + self.assertIn( + "Period structure consistency check failed:\n\nERROR: Found (2) mismatches for commitment period duration:", + str(cm.exception), + ) diff --git a/gtep/tests/unit/test_validation.py b/gtep/tests/unit/test_validation.py index 9cda178c..339e5d26 100644 --- a/gtep/tests/unit/test_validation.py +++ b/gtep/tests/unit/test_validation.py @@ -44,9 +44,9 @@ def get_solution_object(): data_object = ExpansionPlanningData( stages=2, num_reps=2, - len_reps=1, num_commit=6, num_dispatch=4, + duration_representative_period=6, ) data_object.load_prescient(str(input_data_source)) diff --git a/gtep/utils.py b/gtep/utils.py new file mode 100644 index 00000000..586ac989 --- /dev/null +++ b/gtep/utils.py @@ -0,0 +1,248 @@ +################################################################################# +# The Institute for the Design of Advanced Energy Systems Integrated Platform +# Framework (IDAES IP) was produced under the DOE Institute for the +# Design of Advanced Energy Systems (IDAES). +# +# Copyright (c) 2018-2026 by the software owners: The Regents of the +# University of California, through Lawrence Berkeley National Laboratory, +# National Technology & Engineering Solutions of Sandia, LLC, Carnegie Mellon +# University, West Virginia University Research Corporation, et al. +# All rights reserved. Please see the files COPYRIGHT.md and LICENSE.md +# for full copyright and license information. +################################################################################# + +import os +import json + +import pyomo.environ as pyo +from pyomo.environ import units as u + +curr_dir = os.path.dirname(os.path.abspath(__file__)) + + +def generate_period_structure_skeleton( + num_reps, + num_commit, + num_dispatch, + rep_duration=24, + com_duration=1, + disp_duration=15, +): + """This method generates a skeleton of the period structure + dictionary for user editing. The dictionary is nested as follows: + + { + "number_representative": , + "number_commitment": {rep: }, + "number_dispatch": {rep: {com: }}, + "duration_representative_period": {rep: }, + "duration_commitment": {rep: {com: }}, + "duration_dispatch": {rep: {com: {disp: }}} + } + + """ + + period_dict = { + "number_representative": num_reps, + "number_commitment": {rep: num_commit for rep in range(1, num_reps + 1)}, + "number_dispatch": { + rep: {com: num_dispatch for com in range(1, num_commit + 1)} + for rep in range(1, num_reps + 1) + }, + "duration_representative_period": { + rep: rep_duration for rep in range(1, num_reps + 1) + }, + "duration_commitment": { + rep: {com: com_duration for com in range(1, num_commit + 1)} + for rep in range(1, num_reps + 1) + }, + "duration_dispatch": { + rep: { + com: {disp: disp_duration for disp in range(1, num_dispatch + 1)} + for com in range(1, num_commit + 1) + } + for rep in range(1, num_reps + 1) + }, + } + + return period_dict + + +def save_period_structure_json(period_structure, filename): + """This method saves a period structure dictionary to a .json file + file. + + """ + + with open(filename, "w") as f: + json.dump(period_structure, f, indent=2) + + +def generate_period_structure_utils( + num_reps, + num_commit, + num_dispatch, + rep_duration=24, + com_duration=1, + disp_duration=15, + filename=None, +): + """This method generates a period structure skeleton and + optionally saves it as a .json file for user editing. + + :return: period structure dictionary; if a filename is provided, + also saves it as a .json file. + + """ + + period_dict = generate_period_structure_skeleton( + num_reps, num_commit, num_dispatch, rep_duration, com_duration, disp_duration + ) + + if filename: + save_period_structure_json(period_dict, filename) + + return period_dict + + +def _set_period_structure_dict( + num_reps, + num_commit, + num_dispatch, + duration_representative_period, + duration_commitment, + duration_dispatch, + save_period_structure_file, + period_structure_json_file, +): + """This method returns a period structure dictionary with keys for + the number and duration of representative, commitment, and + dispatch periods. + + If a JSON file is specified, it loads the period structure from + that file. Otherwise, generates the period structure from the + provided scalar arguments. Optionally saves the generated + structure to a JSON file. + + Returns: + + :dict: period structure dictionary + + """ + + # If a .json file with period structure data is provided, use + # it, otherwise, expand to a dictionary using the provided + # scalars. + + if period_structure_json_file is not None: + # Use provided .json file + json_path = os.path.abspath( + os.path.join(curr_dir, "data", period_structure_json_file) + ) + with open(json_path, "r") as f: + period_dict = json.load(f) + + # Helper function to recursively convert string keys to + # integers + def convert_keys_to_int(obj): + if isinstance(obj, dict): + return { + ( + int(k) if isinstance(k, str) and k.isdigit() else k + ): convert_keys_to_int(v) + for k, v in obj.items() + } + else: + return obj + + period_dict = convert_keys_to_int(period_dict) + + else: + # .json file not provided; expand period structure + # dictionary from scalar arguments. Optionally save the + # expanded dictionary as a .json file with a default name + # under the data directory. + filename = ( + os.path.abspath( + os.path.join(curr_dir, "data", "period_structure_from_gtep.json") + ) + if save_period_structure_file + else None + ) + + period_dict = generate_period_structure_utils( + num_reps, + num_commit, + num_dispatch, + duration_representative_period, + duration_commitment, + duration_dispatch, + filename=filename, + ) + if save_period_structure_file: + print( + f"\nINFO: Period structure dictionary generated from scalar period arguments has been written to '{filename}'.\n" + ) + + return period_dict + + +def check_period_structure_consistency( + num_reps, + num_commit, + num_dispatch, + duration_representative_period, + duration_commitment, + duration_dispatch, +): + """This method checks that the sum of commitment and dispatch + durations equals the representative and commitment period + duration. It raises ValueError with details if mismatches are + found. + + """ + + commitment_errors = [] + dispatch_errors = [] + for rep in range(1, num_reps + 1): + # Consistency check (1): Sum commitment durations (in hours) + commitment_sum_hr = sum( + duration_commitment[rep][com] for com in range(1, num_commit[rep] + 1) + ) + rep_period_hr = duration_representative_period[rep] + if abs(commitment_sum_hr - rep_period_hr) > 1e-6: + commitment_errors.append( + f" - Representative period {rep}: sum of commitment durations ({commitment_sum_hr} hr) != representative period duration ({rep_period_hr} hr)" + ) + + for com in range(1, num_commit[rep] + 1): + # Consistency check (2): Sum dispatch durations (in + # minutes) and convert to hours + dispatch_sum_hr = pyo.units.convert( + sum( + duration_dispatch[rep][com][disp] + for disp in range(1, num_dispatch[rep][com] + 1) + ) + * u.minutes, + to_units=u.hours, + ) + commitment_hr = duration_commitment[rep][com] + if abs(pyo.value(dispatch_sum_hr) - commitment_hr) > 1e-6: + dispatch_errors.append( + f" - Representative period {rep}, commitment period {com}: sum of dispatch durations ({pyo.value(dispatch_sum_hr)} hr) != commitment period duration ({commitment_hr} hr)" + ) + + # Raise an error if any mismatches were found + if commitment_errors or dispatch_errors: + msg = ["Period structure consistency check failed:\n"] + if commitment_errors: + msg.append( + f"ERROR: Found ({len(commitment_errors)}) mismatches for commitment period duration:\n" + + "\n".join(commitment_errors) + ) + if dispatch_errors: + msg.append( + f"ERROR: Found ({len(dispatch_errors)}) mismatches for dispatch period duration:\n" + + "\n".join(dispatch_errors) + ) + raise ValueError("\n".join(msg))