Note
You can download this example as a Jupyter notebook or start it in interactive mode.
Power to Gas with Heat Coupling
This is an example for power to gas with optional coupling to heat sector (via boiler OR Combined-Heat-and-Power (CHP))
A location has an electric, gas and heat bus. The primary source is wind power, which can be converted to gas. The gas can be stored to convert into electricity or heat (with either a boiler or a CHP).
[1]:
import pypsa
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from pyomo.environ import Constraint
%matplotlib inline
Combined-Heat-and-Power (CHP) parameterisation
This setup follows http://www.ea-energianalyse.dk/reports/student-reports/integration_of_50_percent_wind%20power.pdf pages 35-6 which follows http://www.sciencedirect.com/science/article/pii/030142159390282K
[2]:
#ratio between max heat output and max electric output
nom_r = 1.
#backpressure limit
c_m = 0.75
#marginal loss for each additional generation of heat
c_v = 0.15
Graph for the case that max heat output equals max electric output
[3]:
fig, ax = plt.subplots(figsize=(9,5))
t = 0.01
ph = np.arange(0,1.0001,t)
ax.plot(ph,c_m*ph)
ax.set_xlabel("P_heat_out")
ax.set_ylabel("P_elec_out")
ax.grid(True)
ax.set_xlim([0,1.1])
ax.set_ylim([0,1.1])
ax.text(0.1, 0.7, "Allowed output", color="r")
ax.plot(ph, 1-c_v*ph)
for i in range(1,10):
k = 0.1*i
x = np.arange(0,k/(c_m+c_v),t)
ax.plot(x,k-c_v*x,color="g",alpha=0.5)
ax.text(0.05, 0.41, "iso-fuel-lines", color="g", rotation=-7)
ax.fill_between(ph, c_m*ph, 1-c_v*ph, facecolor="r", alpha=0.5)
fig.tight_layout()
Optimisation
[4]:
network = pypsa.Network()
network.set_snapshots(pd.date_range("2016-01-01 00:00","2016-01-01 03:00",freq="H"))
network.add("Bus", "0", carrier="AC")
network.add("Bus", "0 gas", carrier="gas")
network.add("Carrier", "wind")
network.add("Carrier", "gas", co2_emissions=0.2)
network.add("GlobalConstraint",
"co2_limit",
sense="<=",
constant=0.)
network.add("Generator",
"wind turbine",
bus="0",
carrier="wind",
p_nom_extendable=True,
p_max_pu=[0.,0.2,0.7,0.4],
capital_cost=1000)
network.add("Load",
"load",
bus="0",
p_set=5.)
network.add("Link",
"P2G",
bus0="0",
bus1="0 gas",
efficiency=0.6,
capital_cost=1000,
p_nom_extendable=True)
network.add("Link",
"generator",
bus0="0 gas",
bus1="0",
efficiency=0.468,
capital_cost=400,
p_nom_extendable=True)
network.add("Store",
"gas depot",
bus="0 gas",
e_cyclic=True,
e_nom_extendable=True)
Add heat sector
[5]:
network.add("Bus", "0 heat", carrier="heat")
network.add("Carrier", "heat")
network.add("Load",
"heat load",
bus="0 heat",
p_set=10.)
network.add("Link",
"boiler",
bus0="0 gas",
bus1="0 heat",
efficiency=0.9,
capital_cost=300,
p_nom_extendable=True)
network.add("Store",
"water tank",
bus="0 heat",
e_cyclic=True,
e_nom_extendable=True)
Add CHP constraints
[6]:
#Guarantees ISO fuel lines, i.e. fuel consumption p_b0 + p_g0 = constant along p_g1 + c_v p_b1 = constant
network.links.at["boiler","efficiency"] = network.links.at["generator","efficiency"]/c_v
def extra_functionality(network,snapshots):
#Guarantees heat output and electric output nominal powers are proportional
network.model.chp_nom = Constraint(rule=lambda model : network.links.at["generator","efficiency"]*nom_r*model.link_p_nom["generator"] == network.links.at["boiler","efficiency"]*model.link_p_nom["boiler"])
#Guarantees c_m p_b1 \leq p_g1
def backpressure(model,snapshot):
return c_m*network.links.at["boiler","efficiency"]*model.link_p["boiler",snapshot] <= network.links.at["generator","efficiency"]*model.link_p["generator",snapshot]
network.model.backpressure = Constraint(list(snapshots),rule=backpressure)
#Guarantees p_g1 +c_v p_b1 \leq p_g1_nom
def top_iso_fuel_line(model,snapshot):
return model.link_p["boiler",snapshot] + model.link_p["generator",snapshot] <= model.link_p_nom["generator"]
network.model.top_iso_fuel_line = Constraint(list(snapshots),rule=top_iso_fuel_line)
[7]:
network.lopf(network.snapshots, extra_functionality=extra_functionality)
network.objective
INFO:pypsa.opf:Performed preliminary steps
INFO:pypsa.opf:Building pyomo model using `kirchhoff` formulation
---------------------------------------------------------------------------
RuntimeError Traceback (most recent call last)
/tmp/ipykernel_636/3585628486.py in <module>
----> 1 network.lopf(network.snapshots, extra_functionality=extra_functionality)
2 network.objective
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/components.py in lopf(self, snapshots, pyomo, solver_name, solver_options, solver_logfile, formulation, keep_files, extra_functionality, multi_investment_periods, **kwargs)
646
647 if pyomo:
--> 648 return network_lopf(self, **args)
649 else:
650 return network_lopf_lowmem(self, **args)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/opf.py in network_lopf(network, snapshots, solver_name, solver_io, skip_pre, extra_functionality, multi_investment_periods, solver_logfile, solver_options, keep_files, formulation, ptdf_tolerance, free_memory, extra_postprocessing)
1652
1653 network_lopf_build_model(network, snapshots, skip_pre=skip_pre,
-> 1654 formulation=formulation, ptdf_tolerance=ptdf_tolerance)
1655
1656 if extra_functionality is not None:
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/opf.py in network_lopf_build_model(network, snapshots, skip_pre, formulation, ptdf_tolerance)
1470 define_sub_network_balance_constraints(network,snapshots)
1471
-> 1472 define_global_constraints(network,snapshots)
1473
1474 define_linear_objective(network, snapshots)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/opf.py in define_global_constraints(network, snapshots)
1123
1124 l_constraint(network.model, "global_constraints",
-> 1125 global_constraints, list(network.global_constraints.index))
1126
1127
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pypsa/opt.py in l_constraint(model, name, constraints, *args)
206 raise KeyError('`sense` must be one of "==","<=",">=","><"; got: {}'.format(sense))
207
--> 208 v._data[i] = _GeneralConstraintData(constr_expr, v)
209
210
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pyomo/core/base/constraint.py in __init__(self, expr, component)
291 self._expr = None
292 if expr is not None:
--> 293 self.set_value(expr)
294
295 def __getstate__(self):
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pyomo/core/base/constraint.py in set_value(self, expr)
540 "\n sum(model.costs) == model.income"
541 "\n (0, model.price[item], 50)"
--> 542 % (self.name, str(expr)))
543 raise ValueError(msg)
544 #
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pyomo/core/base/component.py in name(self)
275 def name(self):
276 """Get the fully qualifed component name."""
--> 277 return self.getname(fully_qualified=True)
278
279 # Adding a setter here to help users adapt to the new
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pyomo/core/base/component.py in getname(self, fully_qualified, name_buffer, relative_to)
922 return base + index_repr(idx)
923 #
--> 924 raise RuntimeError("Fatal error: cannot find the component data in "
925 "the owning component's _data dictionary.")
926
RuntimeError: Fatal error: cannot find the component data in the owning component's _data dictionary.
Inspection
[8]:
network.loads_t.p
[8]:
| snapshot |
|---|
| 2016-01-01 00:00:00 |
| 2016-01-01 01:00:00 |
| 2016-01-01 02:00:00 |
| 2016-01-01 03:00:00 |
[9]:
network.links.p_nom_opt
[9]:
P2G 0.0
generator 0.0
boiler 0.0
Name: p_nom_opt, dtype: float64
[10]:
#CHP is dimensioned by the heat demand met in three hours when no wind
4*10./3/network.links.at["boiler","efficiency"]
[10]:
4.273504273504273
[11]:
#elec is set by the heat demand
28.490028*0.15
[11]:
4.2735042
[12]:
network.links_t.p0
[12]:
| snapshot |
|---|
| 2016-01-01 00:00:00 |
| 2016-01-01 01:00:00 |
| 2016-01-01 02:00:00 |
| 2016-01-01 03:00:00 |
[13]:
network.links_t.p1
[13]:
| snapshot |
|---|
| 2016-01-01 00:00:00 |
| 2016-01-01 01:00:00 |
| 2016-01-01 02:00:00 |
| 2016-01-01 03:00:00 |
[14]:
pd.DataFrame({attr: network.stores_t[attr]["gas depot"] for attr in ["p","e"]})
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3360 try:
-> 3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
KeyError: 'gas depot'
The above exception was the direct cause of the following exception:
KeyError Traceback (most recent call last)
/tmp/ipykernel_636/1697226736.py in <module>
----> 1 pd.DataFrame({attr: network.stores_t[attr]["gas depot"] for attr in ["p","e"]})
/tmp/ipykernel_636/1697226736.py in <dictcomp>(.0)
----> 1 pd.DataFrame({attr: network.stores_t[attr]["gas depot"] for attr in ["p","e"]})
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/frame.py in __getitem__(self, key)
3456 if self.columns.nlevels > 1:
3457 return self._getitem_multilevel(key)
-> 3458 indexer = self.columns.get_loc(key)
3459 if is_integer(indexer):
3460 indexer = [indexer]
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
-> 3363 raise KeyError(key) from err
3364
3365 if is_scalar(key) and isna(key) and not self.hasnans:
KeyError: 'gas depot'
[15]:
pd.DataFrame({attr: network.stores_t[attr]["water tank"] for attr in ["p","e"]})
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3360 try:
-> 3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
KeyError: 'water tank'
The above exception was the direct cause of the following exception:
KeyError Traceback (most recent call last)
/tmp/ipykernel_636/1237326891.py in <module>
----> 1 pd.DataFrame({attr: network.stores_t[attr]["water tank"] for attr in ["p","e"]})
/tmp/ipykernel_636/1237326891.py in <dictcomp>(.0)
----> 1 pd.DataFrame({attr: network.stores_t[attr]["water tank"] for attr in ["p","e"]})
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/frame.py in __getitem__(self, key)
3456 if self.columns.nlevels > 1:
3457 return self._getitem_multilevel(key)
-> 3458 indexer = self.columns.get_loc(key)
3459 if is_integer(indexer):
3460 indexer = [indexer]
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
-> 3363 raise KeyError(key) from err
3364
3365 if is_scalar(key) and isna(key) and not self.hasnans:
KeyError: 'water tank'
[16]:
pd.DataFrame({attr: network.links_t[attr]["boiler"] for attr in ["p0","p1"]})
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3360 try:
-> 3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/_libs/index.pyx in pandas._libs.index.IndexEngine.get_loc()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
pandas/_libs/hashtable_class_helper.pxi in pandas._libs.hashtable.PyObjectHashTable.get_item()
KeyError: 'boiler'
The above exception was the direct cause of the following exception:
KeyError Traceback (most recent call last)
/tmp/ipykernel_636/3866892312.py in <module>
----> 1 pd.DataFrame({attr: network.links_t[attr]["boiler"] for attr in ["p0","p1"]})
/tmp/ipykernel_636/3866892312.py in <dictcomp>(.0)
----> 1 pd.DataFrame({attr: network.links_t[attr]["boiler"] for attr in ["p0","p1"]})
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/frame.py in __getitem__(self, key)
3456 if self.columns.nlevels > 1:
3457 return self._getitem_multilevel(key)
-> 3458 indexer = self.columns.get_loc(key)
3459 if is_integer(indexer):
3460 indexer = [indexer]
~/checkouts/readthedocs.org/user_builds/pypsa-docs-staging/envs/latest/lib/python3.7/site-packages/pandas/core/indexes/base.py in get_loc(self, key, method, tolerance)
3361 return self._engine.get_loc(casted_key)
3362 except KeyError as err:
-> 3363 raise KeyError(key) from err
3364
3365 if is_scalar(key) and isna(key) and not self.hasnans:
KeyError: 'boiler'
[17]:
network.stores.loc["gas depot"]
[17]:
attribute
bus 0 gas
type
carrier gas
e_nom 0.0
e_nom_extendable True
e_nom_min 0.0
e_nom_max inf
e_min_pu 0.0
e_max_pu 1.0
e_initial 0.0
e_initial_per_period False
e_cyclic True
e_cyclic_per_period True
p_set 0.0
q_set 0.0
sign 1.0
marginal_cost 0.0
capital_cost 0.0
standing_loss 0.0
build_year 0
lifetime inf
e_nom_opt 0.0
Name: gas depot, dtype: object
[18]:
network.generators.loc["wind turbine"]
[18]:
attribute
bus 0
control Slack
type
p_nom 0.0
p_nom_extendable True
p_nom_min 0.0
p_nom_max inf
p_min_pu 0.0
p_max_pu 1.0
p_set 0.0
q_set 0.0
sign 1.0
carrier wind
marginal_cost 0.0
build_year 0
lifetime inf
capital_cost 1000.0
efficiency 1.0
committable False
start_up_cost 0.0
shut_down_cost 0.0
min_up_time 0
min_down_time 0
up_time_before 1
down_time_before 0
ramp_limit_up NaN
ramp_limit_down NaN
ramp_limit_start_up 1.0
ramp_limit_shut_down 1.0
p_nom_opt 0.0
Name: wind turbine, dtype: object
[19]:
network.links.p_nom_opt
[19]:
P2G 0.0
generator 0.0
boiler 0.0
Name: p_nom_opt, dtype: float64
Calculate the overall efficiency of the CHP
[20]:
eta_elec = network.links.at["generator","efficiency"]
r = 1/c_m
#P_h = r*P_e
(1+r)/((1/eta_elec)*(1+c_v*r))
[20]:
0.9099999999999999