Skip to content

Commit

Permalink
add global indexing capabilities to min/max activity/capacity constra… (
Browse files Browse the repository at this point in the history 
TemoaProject#52)

* add global indexing capabilities to min/max activity/capacity constraints

* fix bug in tech_groups definition
  • Loading branch information
@kathjordan @jeff-ws
kathjordan authored and jeff-ws committed Feb 22, 2024
1 parent 2814590 commit 294d715
Show file tree
Hide file tree
Showing 2 changed files with 155 additions and 57 deletions.
15 changes: 8 additions & 7 deletions temoa_model/temoa_model.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -61,7 +61,8 @@ def temoa_create_model(name="Temoa"):
# RegionalIndices is the set of all the possible combinations of interregional
# exhanges plus original region indices. If tech_exchange is empty, RegionalIndices =regions.
M.RegionalIndices = Set(initialize=CreateRegionalIndices)

M.RegionalGlobalIndices = Set(initialize=RegionalGlobalInitializedIndices)

# Define technology-related sets
M.tech_resource = Set()
M.tech_production = Set()
Expand All @@ -77,7 +78,7 @@ def temoa_create_model(name="Temoa"):
M.tech_flex = Set(within=M.tech_all)
M.tech_exchange = Set(within=M.tech_all)
M.groups = Set(dimen=1) # Define groups for technologies
M.tech_groups = Set(within=M.RegionalIndices * M.groups * M.tech_all)
M.tech_groups = Set(within=M.RegionalGlobalIndices * M.groups * M.tech_all)
M.tech_annual = Set(within=M.tech_all) # Define techs with constant output
M.tech_variable = Set(within=M.tech_all) # Define techs for use with TechInputSplitAverage constraint, where techs have variable annual output but the user wishes to constrain them annually
M.tech_retirement = Set(within=M.tech_all) # Define techs for which economic retirement is an option
Expand Down Expand Up @@ -210,7 +211,7 @@ def temoa_create_model(name="Temoa"):
)

# Define parameters associated with user-defined constraints
M.RegionalGlobalIndices = Set(initialize=RegionalGlobalInitializedIndices)

M.MinCapacity = Param(M.RegionalIndices, M.time_optimize, M.tech_all)
M.MaxCapacity = Param(M.RegionalIndices, M.time_optimize, M.tech_all)
M.MinNewCapacity = Param(M.RegionalIndices, M.time_optimize, M.tech_all)
Expand All @@ -227,10 +228,10 @@ def temoa_create_model(name="Temoa"):
M.EmissionLimit = Param(M.RegionalGlobalIndices, M.time_optimize, M.commodity_emissions)
M.EmissionActivity_reitvo = Set(dimen=6, initialize=EmissionActivityIndices)
M.EmissionActivity = Param(M.EmissionActivity_reitvo)
M.MinActivityGroup = Param(M.RegionalIndices, M.time_optimize, M.groups)
M.MaxActivityGroup = Param(M.RegionalIndices, M.time_optimize, M.groups)
M.MinCapacityGroup = Param(M.RegionalIndices, M.time_optimize, M.groups)
M.MaxCapacityGroup = Param(M.RegionalIndices, M.time_optimize, M.groups)
M.MinActivityGroup = Param(M.RegionalGlobalIndices, M.time_optimize, M.groups)
M.MaxActivityGroup = Param(M.RegionalGlobalIndices, M.time_optimize, M.groups)
M.MinCapacityGroup = Param(M.RegionalGlobalIndices, M.time_optimize, M.groups)
M.MaxCapacityGroup = Param(M.RegionalGlobalIndices, M.time_optimize, M.groups)
M.MinNewCapacityGroup = Param(M.RegionalIndices, M.time_optimize, M.groups)
M.MaxNewCapacityGroup = Param(M.RegionalIndices, M.time_optimize, M.groups)
M.MinCapShare_rptg = Set(dimen=4, initialize=MinCapShareIndices)
Expand Down
197 changes: 147 additions & 50 deletions temoa_model/temoa_rules.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1890,28 +1890,57 @@ def MinActivityGroup_Constraint(M, r, p, g):
refers to the :code:`MinActivityGroup` parameter.
"""

activity_p = sum(
M.V_FlowOut[r, p, s, d, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r and _g == g and (r, p, S_t) in M.processVintages
for S_v in M.processVintages[r, p, S_t]
for S_i in M.processInputs[r, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r, p, S_t, S_v, S_i]
for s in M.time_season
for d in M.time_of_day
if (S_t not in M.tech_annual) and ((r, p, S_t) in M.processVintages.keys())
)
if r == 'global':
reg = M.regions
elif '+' in r:
reg = r.split('+')
else:
reg = [r]

activity_p = 0
activity_p_annual = 0
for r_i in reg:
if r == 'global':
activity_p += sum(
M.V_FlowOut[r_i, p, s, d, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r and _g == g and (r_i, p, S_t) in M.processVintages
for S_v in M.processVintages[r_i, p, S_t]
for S_i in M.processInputs[r_i, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r_i, p, S_t, S_v, S_i]
for s in M.time_season
for d in M.time_of_day
if (S_t not in M.tech_annual) and ((r_i, p, S_t) in M.processVintages.keys())
)

activity_p_annual = sum(
M.V_FlowOutAnnual[r, p, S_i, S_t, S_v, S_o] * M.MinGenGroupWeight[r, S_t, g]
for r in M.RegionalIndices
for S_t in M.tech_groups if (S_t in M.tech_annual) and ((r, p, S_t) in M.processVintages.keys())
for S_v in M.processVintages[r, p, S_t]
for S_i in M.processInputs[r, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r, p, S_t, S_v, S_i]
if (S_t in M.tech_annual) and ((r, p, S_t) in M.processVintages.keys())
)
activity_p_annual += sum(
M.V_FlowOutAnnual[r_i, p, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r and _g == g and (r_i, p, S_t) in M.processVintages
for S_v in M.processVintages[r_i, p, S_t]
for S_i in M.processInputs[r_i, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r_i, p, S_t, S_v, S_i]
if (S_t in M.tech_annual) and ((r_i, p, S_t) in M.processVintages.keys())
)

else:
activity_p += sum(
M.V_FlowOut[r_i, p, s, d, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r_i and _g == g and (r_i, p, S_t) in M.processVintages
for S_v in M.processVintages[r_i, p, S_t]
for S_i in M.processInputs[r_i, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r_i, p, S_t, S_v, S_i]
for s in M.time_season
for d in M.time_of_day
if (S_t not in M.tech_annual) and ((r_i, p, S_t) in M.processVintages.keys())
)

activity_p_annual += sum(
M.V_FlowOutAnnual[r_i, p, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r_i and _g == g and (r_i, p, S_t) in M.processVintages
for S_v in M.processVintages[r_i, p, S_t]
for S_i in M.processInputs[r_i, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r_i, p, S_t, S_v, S_i]
if (S_t in M.tech_annual) and ((r_i, p, S_t) in M.processVintages.keys())
)
min_act = value(M.MinActivityGroup[r, p, g])
expr = activity_p + activity_p_annual >= min_act
return expr
Expand All @@ -1930,25 +1959,57 @@ def MaxActivityGroup_Constraint(M, r, p, g):
refers to the :code:`MaxActivityGroup` parameter.
"""

activity_p = sum(
M.V_FlowOut[r, p, s, d, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r and _g == g and (r, p, S_t) in M.processVintages
for S_v in M.processVintages[r, p, S_t]
for S_i in M.processInputs[r, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r, p, S_t, S_v, S_i]
for s in M.time_season
for d in M.time_of_day
if (S_t not in M.tech_annual) and ((r, p, S_t) in M.processVintages.keys())
)
if r == 'global':
reg = M.regions
elif '+' in r:
reg = r.split('+')
else:
reg = [r]

activity_p_annual = sum(
M.V_FlowOutAnnual[r, p, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r and _g == g and (r, p, S_t) in M.processVintages
for S_v in M.processVintages[r, p, S_t]
for S_i in M.processInputs[r, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r, p, S_t, S_v, S_i]
if (S_t in M.tech_annual) and ((r, p, S_t) in M.processVintages.keys())
)
activity_p = 0
activity_p_annual = 0
for r_i in reg:
if r == 'global':
activity_p += sum(
M.V_FlowOut[r_i, p, s, d, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r and _g == g and (r_i, p, S_t) in M.processVintages
for S_v in M.processVintages[r_i, p, S_t]
for S_i in M.processInputs[r_i, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r_i, p, S_t, S_v, S_i]
for s in M.time_season
for d in M.time_of_day
if (S_t not in M.tech_annual) and ((r_i, p, S_t) in M.processVintages.keys())
)

activity_p_annual += sum(
M.V_FlowOutAnnual[r_i, p, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r and _g == g and (r_i, p, S_t) in M.processVintages
for S_v in M.processVintages[r_i, p, S_t]
for S_i in M.processInputs[r_i, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r_i, p, S_t, S_v, S_i]
if (S_t in M.tech_annual) and ((r_i, p, S_t) in M.processVintages.keys())
)

else:
activity_p += sum(
M.V_FlowOut[r_i, p, s, d, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r_i and _g == g and (r_i, p, S_t) in M.processVintages
for S_v in M.processVintages[r_i, p, S_t]
for S_i in M.processInputs[r_i, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r_i, p, S_t, S_v, S_i]
for s in M.time_season
for d in M.time_of_day
if (S_t not in M.tech_annual) and ((r_i, p, S_t) in M.processVintages.keys())
)

activity_p_annual += sum(
M.V_FlowOutAnnual[r_i, p, S_i, S_t, S_v, S_o]
for _r, _g, S_t in M.tech_groups if _r == r_i and _g == g and (r_i, p, S_t) in M.processVintages
for S_v in M.processVintages[r_i, p, S_t]
for S_i in M.processInputs[r_i, p, S_t, S_v]
for S_o in M.ProcessOutputsByInput[r_i, p, S_t, S_v, S_i]
if (S_t in M.tech_annual) and ((r_i, p, S_t) in M.processVintages.keys())
)

max_act = value(M.MaxActivityGroup[r, p, g])
expr = activity_p + activity_p_annual <= max_act
Expand Down Expand Up @@ -2031,14 +2092,31 @@ def MaxCapacityGroup_Constraint(M, r, p, g):
Similar to the :code:`MaxCapacity` constraint, but works on a group of technologies.
"""
max_cap = value(M.MaxCapacityGroup[r, p, g])
aggcap = sum(
M.V_CapacityAvailableByPeriodAndTech[r, p, t]
for _r, _g, t in M.tech_groups if _r == r and _g == g and (r, p, t) in M.V_CapacityAvailableByPeriodAndTech.keys()
)
expr = aggcap <= max_cap
if r == 'global':
reg = M.regions
elif '+' in r:
reg = r.split('+')
else:
reg = [r]
max_capgroup = value(M.MaxCapacityGroup[r, p, g])

cap = 0
for r_i in reg:
if r == 'global':
cap += sum(
M.V_CapacityAvailableByPeriodAndTech[r_i, p, t]
for _r, _g, t in M.tech_groups
if _r == r and _g == g and (r_i,p,t) in M.V_CapacityAvailableByPeriodAndTech.keys()
)
else:
cap += sum(
M.V_CapacityAvailableByPeriodAndTech[r_i, p, t]
for _r, _g, t in M.tech_groups
if _r == r_i and _g == g and (r_i,p,t) in M.V_CapacityAvailableByPeriodAndTech.keys()
)

expr = cap <= max_capgroup
return expr

def MinNewCapacity_Constraint(M, r, p, t):
r"""
The MinNewCapacity constraint sets a limit on the minimum newly installed capacity of a
Expand Down Expand Up @@ -2077,14 +2155,33 @@ def MinCapacityGroup_Constraint(M, r, p, g):
Similar to the :code:`MinCapacity` constraint, but works on a group of technologies.
"""
min_cap = value(M.MinCapacityGroup[r, p, g])
aggcap = sum(
M.V_CapacityAvailableByPeriodAndTech[r, p, t]
for _r, _g, t in M.tech_groups if _r == r and _g == g and (r, p, t) in M.V_CapacityAvailableByPeriodAndTech.keys()
)
expr = aggcap >= min_cap
if r == 'global':
reg = M.regions
elif '+' in r:
reg = r.split('+')
else:
reg = [r]
min_capgroup = value(M.MinCapacityGroup[r, p, g])

cap = 0
for r_i in reg:
if r == 'global':
cap += sum(
M.V_CapacityAvailableByPeriodAndTech[r_i, p, t]
for _r, _g, t in M.tech_groups
if _r == r and _g == g and (r_i,p,t) in M.V_CapacityAvailableByPeriodAndTech.keys()
)
else:
cap += sum(
M.V_CapacityAvailableByPeriodAndTech[r_i, p, t]
for _r, _g, t in M.tech_groups
if _r == r_i and _g == g and (r_i,p,t) in M.V_CapacityAvailableByPeriodAndTech.keys()
)

expr = cap >= min_capgroup
return expr


def MinNewCapacityGroup_Constraint(M, r, p, g):
r"""
Similar to the :code:`MinNewCapacity` constraint, but works on a group of technologies.
Expand Down

0 comments on commit 294d715

Please sign in to comment.