|
@@ -0,0 +1,194 @@
|
|
|
+import pandas as pd
|
|
|
+import numpy as np
|
|
|
+
|
|
|
+from matplotlib import pyplot as plt
|
|
|
+import matplotlib.dates as mdates
|
|
|
+
|
|
|
+import yaml
|
|
|
+
|
|
|
+from mix_simul.scenarios import Scenario
|
|
|
+
|
|
|
+with open("scenarios/rte.yml", 'r') as f:
|
|
|
+ rte = yaml.load(f, Loader=yaml.FullLoader)
|
|
|
+
|
|
|
+potential = pd.read_parquet("data/potential.parquet")
|
|
|
+potential = potential.loc['1985-01-01 00:00:00':'2015-01-01 00:00:00', :]
|
|
|
+potential.fillna(0, inplace=True)
|
|
|
+
|
|
|
+begin = "2012-01-01"
|
|
|
+end = "2015-01-01"
|
|
|
+flexibility = False
|
|
|
+
|
|
|
+potential = potential.loc[(
|
|
|
+ slice(f'{begin} 00:00:00', f'{end} 00:00:00'), 'FR'), :]
|
|
|
+
|
|
|
+# intermittent sources potential
|
|
|
+p = potential[["onshore", "offshore", "solar"]].to_xarray().to_array()
|
|
|
+p = np.insert(p, 3, 0.68, axis=0) # nuclear power-like
|
|
|
+
|
|
|
+times = potential.index.get_level_values(0)
|
|
|
+
|
|
|
+fig, axes = plt.subplots(nrows=6, ncols=2, sharex="col", sharey=True)
|
|
|
+w, h = fig.get_size_inches()
|
|
|
+fig.set_figwidth(w*1.5)
|
|
|
+fig.set_figheight(h*1.5)
|
|
|
+
|
|
|
+fig_storage, axes_storage = plt.subplots(nrows=6, ncols=2, sharex="col", sharey=True)
|
|
|
+fig_storage.set_figwidth(w*1.5)
|
|
|
+fig_storage.set_figheight(h*1.5)
|
|
|
+
|
|
|
+fig_dispatch, axes_dispatch = plt.subplots(nrows=6, ncols=2, sharex="col", sharey=True)
|
|
|
+fig_dispatch.set_figwidth(w*1.5)
|
|
|
+fig_dispatch.set_figheight(h*1.5)
|
|
|
+
|
|
|
+date_fmt = mdates.DateFormatter('%d/%m')
|
|
|
+
|
|
|
+
|
|
|
+# for step in np.linspace(start, stop, 2050-2022, True)[::-1]:
|
|
|
+row = 0
|
|
|
+for scenario in rte:
|
|
|
+ rte[scenario]["flexibility_power"] = 0
|
|
|
+ scenario_model = Scenario(**rte[scenario])
|
|
|
+ S, load, production, gap, storage, dp = scenario_model.run(times, p)
|
|
|
+
|
|
|
+ print(f"{scenario}:", S, gap.max(), np.quantile(gap, 0.95))
|
|
|
+ print(f"exports: {np.minimum(np.maximum(-gap, 0), 39).sum()/1000} TWh; imports: {np.minimum(np.maximum(gap, 0), 39).sum()/1000} TWh")
|
|
|
+ print(f"dispatchable: " + ", ".join([f"{dp[i].sum()/1000:.2f} TWh" for i in range(dp.shape[0])]))
|
|
|
+
|
|
|
+ potential['load'] = load
|
|
|
+ potential['production'] = production
|
|
|
+ potential['available'] = production-np.diff(storage.sum(axis=0), append=0)
|
|
|
+
|
|
|
+ for i in range(3):
|
|
|
+ potential[f"storage_{i}"] = np.diff(storage[i,:], append=0)#storage[i,:]/1000
|
|
|
+ potential[f"storage_{i}"] = storage[i,:]/1000
|
|
|
+
|
|
|
+ for i in range(dp.shape[0]):
|
|
|
+ potential[f"dispatch_{i}"] = dp[i,:]
|
|
|
+
|
|
|
+ potential["dispatch"] = dp.sum(axis=0)
|
|
|
+
|
|
|
+ data = [
|
|
|
+ potential.loc[(slice('2013-12-01 00:00:00',
|
|
|
+ '2014-01-01 00:00:00'), 'FR'), :],
|
|
|
+ potential.loc[(slice('2013-06-01 00:00:00',
|
|
|
+ '2013-07-01 00:00:00'), 'FR'), :]
|
|
|
+ ]
|
|
|
+
|
|
|
+ months = [
|
|
|
+ "Février",
|
|
|
+ "Juin"
|
|
|
+ ]
|
|
|
+
|
|
|
+ labels = [
|
|
|
+ "load (GW)",
|
|
|
+ "production (GW)",
|
|
|
+ "available power (production-storage) (GW)",
|
|
|
+ "power deficit"
|
|
|
+ ]
|
|
|
+
|
|
|
+ labels_storage = [
|
|
|
+ "Batteries (TWh)",
|
|
|
+ "STEP (TWh)",
|
|
|
+ "P2G (TWh)"
|
|
|
+ ]
|
|
|
+
|
|
|
+ labels_dispatch = [
|
|
|
+ "Hydro (GW)",
|
|
|
+ "Biomass (GW)",
|
|
|
+ "Thermal (GW)"
|
|
|
+ ]
|
|
|
+
|
|
|
+ for col in range(2):
|
|
|
+ ax = axes[row, col]
|
|
|
+ ax.plot(data[col].index.get_level_values(0), data[col]
|
|
|
+ ["load"], label="adjusted load (GW)", lw=1)
|
|
|
+ ax.plot(data[col].index.get_level_values(0), data[col]
|
|
|
+ ["production"], label="production (GW)", ls="dotted", lw=1)
|
|
|
+ ax.plot(data[col].index.get_level_values(0), data[col]["available"],
|
|
|
+ label="available power (production-d(storage)/dt) (GW)", lw=1)
|
|
|
+
|
|
|
+ ax.fill_between(
|
|
|
+ data[col].index.get_level_values(0),
|
|
|
+ data[col]["available"],
|
|
|
+ data[col]["load"],
|
|
|
+ where=data[col]["load"] > data[col]["available"],
|
|
|
+ color='red',
|
|
|
+ alpha=0.15
|
|
|
+ )
|
|
|
+
|
|
|
+ ax.xaxis.set_major_formatter(date_fmt)
|
|
|
+ ax.text(
|
|
|
+ 0.5, 0.87, f"Scénario {scenario} ({months[col]})", ha='center', transform=ax.transAxes)
|
|
|
+
|
|
|
+ ax.set_ylim(25, 225)
|
|
|
+
|
|
|
+ ax = axes_storage[row, col]
|
|
|
+ for i in np.arange(3):
|
|
|
+ if i == 2:
|
|
|
+ base = 0
|
|
|
+ else:
|
|
|
+ base = np.sum([data[col][f"storage_{j}"] for j in np.arange(i+1,3)], axis=0)
|
|
|
+
|
|
|
+ ax.fill_between(data[col].index.get_level_values(0), base, base+data[col]
|
|
|
+ [f"storage_{i}"], label=f"storage {i}", alpha=0.5)
|
|
|
+
|
|
|
+ ax.plot(data[col].index.get_level_values(0), base+data[col]
|
|
|
+ [f"storage_{i}"], label=f"storage {i}", lw=0.25)
|
|
|
+
|
|
|
+ ax.xaxis.set_major_formatter(date_fmt)
|
|
|
+ ax.text(
|
|
|
+ 0.5, 0.87, f"Scénario {scenario} ({months[col]})", ha='center', transform=ax.transAxes)
|
|
|
+
|
|
|
+ ax = axes_dispatch[row, col]
|
|
|
+ for i in range(dp.shape[0]):
|
|
|
+ if i == 0:
|
|
|
+ base = 0
|
|
|
+ else:
|
|
|
+ base = np.sum([data[col][f"dispatch_{j}"] for j in np.arange(i)], axis=0)
|
|
|
+
|
|
|
+ ax.fill_between(data[col].index.get_level_values(0), base, base+data[col]
|
|
|
+ [f"dispatch_{i}"], label=f"dispatch {i}", alpha=0.5)
|
|
|
+
|
|
|
+ ax.plot(data[col].index.get_level_values(0), base+data[col]
|
|
|
+ [f"dispatch_{i}"], label=f"dispatch {i}", lw=0.25)
|
|
|
+
|
|
|
+ ax.xaxis.set_major_formatter(date_fmt)
|
|
|
+ ax.text(
|
|
|
+ 0.5, 0.87, f"Scénario {scenario} ({months[col]})", ha='center', transform=ax.transAxes)
|
|
|
+
|
|
|
+ row += 1
|
|
|
+
|
|
|
+for label in axes[-1, 0].get_xmajorticklabels() + axes[-1, 1].get_xmajorticklabels():
|
|
|
+ label.set_rotation(30)
|
|
|
+ label.set_horizontalalignment("right")
|
|
|
+
|
|
|
+for label in axes_storage[-1, 0].get_xmajorticklabels() + axes_storage[-1, 1].get_xmajorticklabels():
|
|
|
+ label.set_rotation(30)
|
|
|
+ label.set_horizontalalignment("right")
|
|
|
+
|
|
|
+for label in axes_dispatch[-1, 0].get_xmajorticklabels() + axes_dispatch[-1, 1].get_xmajorticklabels():
|
|
|
+ label.set_rotation(30)
|
|
|
+ label.set_horizontalalignment("right")
|
|
|
+
|
|
|
+flex = "With" if flexibility else "Without"
|
|
|
+
|
|
|
+plt.subplots_adjust(wspace=0, hspace=0)
|
|
|
+fig.suptitle(f"Simulations based on {begin}--{end} weather data.\n{flex} consumption flexibility; no nuclear seasonality (unrealistic)")
|
|
|
+fig.text(1, 0, 'Lucas Gautheron', ha="right")
|
|
|
+fig.legend(labels, loc='lower right', bbox_to_anchor=(1, -0.1),
|
|
|
+ ncol=len(labels), bbox_transform=fig.transFigure)
|
|
|
+fig.savefig("output.png", bbox_inches="tight", dpi=200)
|
|
|
+
|
|
|
+fig_storage.suptitle(f"Simulations based on {begin}--{end} weather data.\n{flex} consumption flexibility; no nuclear seasonality (unrealistic)")
|
|
|
+fig_storage.text(1, 0, 'Lucas Gautheron', ha="right")
|
|
|
+fig_storage.legend(labels_storage, loc='lower right', bbox_to_anchor=(1, -0.1),
|
|
|
+ ncol=len(labels_storage), bbox_transform=fig_storage.transFigure)
|
|
|
+fig_storage.savefig("output_storage.png", bbox_inches="tight", dpi=200)
|
|
|
+
|
|
|
+fig_dispatch.suptitle(f"Simulations based on {begin}--{end} weather data.\n{flex} consumption flexibility; no nuclear seasonality (unrealistic)")
|
|
|
+fig_dispatch.text(1, 0, 'Lucas Gautheron', ha="right")
|
|
|
+fig_dispatch.legend(labels_dispatch, loc='lower right', bbox_to_anchor=(1, -0.1),
|
|
|
+ ncol=len(labels_dispatch), bbox_transform=fig_dispatch.transFigure)
|
|
|
+fig_dispatch.savefig("output_dispatch.png", bbox_inches="tight", dpi=200)
|
|
|
+plt.show()
|