Skip to content

Commit

Permalink
Merge pull request #1182 from GbotemiB/linopy-sec
Browse files Browse the repository at this point in the history 
Linopy Integration for Sec
  • Loading branch information
@ekatef
ekatef authored Nov 21, 2024
2 parents 6770144 + b0d9cb7 commit 92c931b
Showing 1 changed file with 70 additions and 70 deletions.
140 changes: 70 additions & 70 deletions scripts/solve_network.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -462,7 +462,6 @@ def add_battery_constraints(n):
1 * charger_size - efficiency * discharger_size = 0
"""
nodes = n.buses.index[n.buses.carrier == "battery"]

if nodes.empty:
return
vars_link = n.model["Link-p_nom"]
Expand Down Expand Up @@ -609,14 +608,19 @@ def add_h2_network_cap(n, cap):
h2_network = n.links.loc[n.links.carrier == "H2 pipeline"]
if h2_network.index.empty:
return
h2_network_cap = get_var(n, "Link", "p_nom")
subset_index = h2_network.index.intersection(h2_network_cap.index)
lhs = linexpr(
(h2_network.loc[subset_index, "length"], h2_network_cap[subset_index])
h2_network_cap = n.model["Link-p_nom"]
h2_network_cap_index = h2_network_cap.indexes["Link-ext"]
subset_index = h2_network.index.intersection(h2_network_cap_index)
diff_index = h2_network_cap_index.difference(subset_index)
if len(diff_index) > 0:
logger.warning(
f"Impossible to set a limit for H2 pipelines extension for the following links: {diff_index}"
)
lhs = (
h2_network_cap.loc[subset_index] * h2_network.loc[subset_index, "length"]
).sum()
# lhs = linexpr((1, h2_network_cap[h2_network.index])).sum()
rhs = cap * 1000
define_constraints(n, lhs, "<=", rhs, "h2_network_cap")
n.model.add_constraints(lhs <= rhs, name="h2_network_cap")


def H2_export_yearly_constraint(n):
Expand All @@ -637,9 +641,10 @@ def H2_export_yearly_constraint(n):
index=n.snapshots,
columns=res_index,
)
res = join_exprs(
linexpr((weightings, get_var(n, "Generator", "p")[res_index]))
) # single line sum
capacity_variable = n.model["Generator-p"]

# single line sum
res = (weightings * capacity_variable.loc[res_index]).sum()

load_ind = n.loads[n.loads.carrier == "AC"].index.intersection(
n.loads_t.p_set.columns
Expand Down Expand Up @@ -667,7 +672,7 @@ def H2_export_yearly_constraint(n):
else:
rhs = h2_export * (1 / 0.7)

con = define_constraints(n, lhs, ">=", rhs, "H2ExportConstraint", "RESproduction")
n.model.add_constraints(lhs >= rhs, name="H2ExportConstraint-RESproduction")


def monthly_constraints(n, n_ref):
Expand All @@ -690,15 +695,17 @@ def monthly_constraints(n, n_ref):
index=n.snapshots,
columns=res_index,
)
capacity_variable = n.model["Generator-p"]

res = linexpr((weightings, get_var(n, "Generator", "p")[res_index])).sum(
axis=1
) # single line sum
# single line sum
res = (weightings * capacity_variable[res_index]).sum(axis=1)
res = res.groupby(res.index.month).sum()

electrolysis = get_var(n, "Link", "p")[
link_p = n.model["Link-p"]
electrolysis = link_p.loc[
n.links.index[n.links.index.str.contains("H2 Electrolysis")]
]

weightings_electrolysis = pd.DataFrame(
np.outer(
n.snapshot_weightings["generators"], [1.0] * len(electrolysis.columns)
Expand All @@ -707,7 +714,7 @@ def monthly_constraints(n, n_ref):
columns=electrolysis.columns,
)

elec_input = linexpr((-allowed_excess * weightings_electrolysis, electrolysis)).sum(
elec_input = ((-allowed_excess * weightings_electrolysis) * electrolysis).sum(
axis=1
)

Expand All @@ -730,16 +737,16 @@ def monthly_constraints(n, n_ref):
for i in range(len(res.index)):
lhs = res.iloc[i] + "\n" + elec_input.iloc[i]
rhs = res_ref.iloc[i] + elec_input_ref.iloc[i]
con = define_constraints(
n, lhs, ">=", rhs, f"RESconstraints_{i}", f"REStarget_{i}"
n.model.add_constraints(
lhs >= rhs, name=f"RESconstraints_{i}-REStarget_{i}"
)

else:
for i in range(len(res.index)):
lhs = res.iloc[i] + "\n" + elec_input.iloc[i]

con = define_constraints(
n, lhs, ">=", 0.0, f"RESconstraints_{i}", f"REStarget_{i}"
n.model.add_constraints(
lhs >= 0.0, name=f"RESconstraints_{i}-REStarget_{i}"
)
# else:
# logger.info("ignoring H2 export constraint as wildcard is set to 0")
Expand All @@ -760,57 +767,46 @@ def add_chp_constraints(n):
electric = n.links.index[electric_bool]
heat = n.links.index[heat_bool]

electric_ext = n.links.index[electric_bool & n.links.p_nom_extendable]
heat_ext = n.links.index[heat_bool & n.links.p_nom_extendable]
electric_ext = n.links[electric_bool].query("p_nom_extendable").index
heat_ext = n.links[heat_bool].query("p_nom_extendable").index

electric_fix = n.links.index[electric_bool & ~n.links.p_nom_extendable]
heat_fix = n.links.index[heat_bool & ~n.links.p_nom_extendable]
electric_fix = n.links[electric_bool].query("~p_nom_extendable").index
heat_fix = n.links[heat_bool].query("~p_nom_extendable").index

link_p = get_var(n, "Link", "p")
p = n.model["Link-p"] # dimension: [time, link]

# output ratio between heat and electricity and top_iso_fuel_line for extendable
if not electric_ext.empty:
link_p_nom = get_var(n, "Link", "p_nom")

# ratio of output heat to electricity set by p_nom_ratio
lhs = linexpr(
(
n.links.loc[electric_ext, "efficiency"]
* n.links.loc[electric_ext, "p_nom_ratio"],
link_p_nom[electric_ext],
),
(-n.links.loc[heat_ext, "efficiency"].values, link_p_nom[heat_ext].values),
)

define_constraints(n, lhs, "=", 0, "chplink", "fix_p_nom_ratio")
p_nom = n.model["Link-p_nom"]

# top_iso_fuel_line for extendable
lhs = linexpr(
(1, link_p[heat_ext]),
(1, link_p[electric_ext].values),
(-1, link_p_nom[electric_ext].values),
lhs = (
p_nom.loc[electric_ext]
* (n.links.p_nom_ratio * n.links.efficiency)[electric_ext].values
- p_nom.loc[heat_ext] * n.links.efficiency[heat_ext].values
)
n.model.add_constraints(lhs == 0, name="chplink-fix_p_nom_ratio")

define_constraints(n, lhs, "<=", 0, "chplink", "top_iso_fuel_line_ext")
rename = {"Link-ext": "Link"}
lhs = (
p.loc[:, electric_ext]
+ p.loc[:, heat_ext]
- p_nom.rename(rename).loc[electric_ext]
)
n.model.add_constraints(lhs <= 0, name="chplink-top_iso_fuel_line_ext")

# top_iso_fuel_line for fixed
if not electric_fix.empty:
# top_iso_fuel_line for fixed
lhs = linexpr((1, link_p[heat_fix]), (1, link_p[electric_fix].values))

rhs = n.links.loc[electric_fix, "p_nom"].values

define_constraints(n, lhs, "<=", rhs, "chplink", "top_iso_fuel_line_fix")
lhs = p.loc[:, electric_fix] + p.loc[:, heat_fix]
rhs = n.links.p_nom[electric_fix]
n.model.add_constraints(lhs <= rhs, name="chplink-top_iso_fuel_line_fix")

# back-pressure
if not electric.empty:
# backpressure
lhs = linexpr(
(
n.links.loc[electric, "c_b"].values * n.links.loc[heat, "efficiency"],
link_p[heat],
),
(-n.links.loc[electric, "efficiency"].values, link_p[electric].values),
lhs = (
p.loc[:, heat] * (n.links.efficiency[heat] * n.links.c_b[electric].values)
- p.loc[:, electric] * n.links.efficiency[electric]
)

define_constraints(n, lhs, "<=", 0, "chplink", "backpressure")
n.model.add_constraints(lhs <= rhs, name="chplink-backpressure")


def add_co2_sequestration_limit(n, sns):
Expand All @@ -819,25 +815,24 @@ def add_co2_sequestration_limit(n, sns):
if co2_stores.empty:
return

vars_final_co2_stored = get_var(n, "Store", "e").loc[sns[-1], co2_stores]
vars_final_co2_stored = n.model["Store-e"].loc[sns[-1], co2_stores]

lhs = linexpr((1, vars_final_co2_stored)).sum()
lhs = (1 * vars_final_co2_stored).sum()
rhs = (
n.config["sector"].get("co2_sequestration_potential", 5) * 1e6
) # TODO change 200 limit (Europe)

name = "co2_sequestration_limit"
define_constraints(
n, lhs, "<=", rhs, "GlobalConstraint", "mu", axes=pd.Index([name]), spec=name
)

n.model.add_constraints(lhs <= rhs, name=f"GlobalConstraint-{name}")


def set_h2_colors(n):
blue_h2 = get_var(n, "Link", "p")[
blue_h2 = n.model["Link-p"].loc[
n.links.index[n.links.index.str.contains("blue H2")]
]

pink_h2 = get_var(n, "Link", "p")[
pink_h2 = n.model["Link-p"].loc[
n.links.index[n.links.index.str.contains("pink H2")]
]

Expand Down Expand Up @@ -869,16 +864,16 @@ def set_h2_colors(n):
columns=pink_h2.columns,
)

total_blue = linexpr((weightings_blue, blue_h2)).sum().sum()
total_blue = (weightings_blue * blue_h2).sum().sum()

total_pink = linexpr((weightings_pink, pink_h2)).sum().sum()
total_pink = (weightings_pink * pink_h2).sum().sum()

rhs_blue = load_h2 * snakemake.config["sector"]["hydrogen"]["blue_share"]
rhs_pink = load_h2 * snakemake.config["sector"]["hydrogen"]["pink_share"]

define_constraints(n, total_blue, "=", rhs_blue, "blue_h2_share")
n.model.add_constraints(total_blue == rhs_blue, name="blue_h2_share")

define_constraints(n, total_pink, "=", rhs_pink, "pink_h2_share")
n.model.add_constraints(total_pink == rhs_pink, name="pink_h2_share")


def add_existing(n):
Expand Down Expand Up @@ -939,8 +934,13 @@ def extra_functionality(n, snapshots):
for o in opts:
if "EQ" in o:
add_EQ_constraints(n, o)

add_battery_constraints(n)

if snakemake.config["sector"]["chp"]:
logger.info("setting CHP constraints")
add_chp_constraints(n)

if (
snakemake.config["policy_config"]["hydrogen"]["temporal_matching"]
== "h2_yearly_matching"
Expand Down

0 comments on commit 92c931b

Please sign in to comment.