Add endogenous transport option

config/config.default.yaml

@@ -393,7 +393,12 @@ sector:
   bev_charge_rate: 0.011
   bev_avail_max: 0.95
   bev_avail_mean: 0.8
+  bev_charge_avail: 5
   v2g: true
+  endogenous_transport: false
+  EV_consumption_1car: 0.01 
+  ICE_consumption_1car: 0.033 
+  H2_consumption_1car: 0.02 
   land_transport_fuel_cell_share:
     2020: 0
     2025: 0
@@ -418,7 +423,8 @@ sector:
     2040: 0.3
     2045: 0.15
     2050: 0
-  transport_fuel_cell_efficiency: 0.5
+  transport_electric_vehicle_efficiency: 0.99
+  transport_fuel_cell_efficiency: 0.71 #0.5
   transport_internal_combustion_efficiency: 0.3
   agriculture_machinery_electric_share: 0
   agriculture_machinery_oil_share: 1
@@ -1058,4 +1064,5 @@ plotting:
     DC: "#8a1caf"
     DC-DC: "#8a1caf"
     DC link: "#8a1caf"
+    land transport demand: '#2596be'
     load: "#dd2e23"

doc/configtables/sector.csv

@@ -24,6 +24,10 @@ bev_charge_rate,MWh,float,The power consumption for one electric vehicle (EV) in
 bev_avail_max,--,float,The maximum share plugged-in availability for passenger electric vehicles.
 bev_avail_mean,--,float,The average share plugged-in availability for passenger electric vehicles.
 v2g,--,"{true, false}",Allows feed-in to grid from EV battery
+endogenous transport,--,"{true, false}",Allows to optimize land transport shares endogenously
+EV_consumption_1car,MWh_elec/hour,float,assuming 0.2 kWh/km https://github.com/PyPSA/pypsa-eur/blob/1fbe971ab8dab60d972d3a7b905b9cec7171c0ad/config/config.default.yaml#L382 and velocity of 50km/h
+ICE_consumption_1car,Mwh_oil/hour,float,"assuming 0.66 kWh_oil/km and 50 km/h (with the link efficiency 0.3, they are equivalent to 0.01 MWh_elec/hour"
+H2_consumption_1car,Mwh_H2/hour,float,"assuming 0.4 kWh_H2/km and 50km/h) (with efficiency 0.5, they are equivalent to 0.01 MWh_elec/hour"
 land_transport_fuel_cell _share,--,Dictionary with planning horizons as keys.,The share of vehicles that uses fuel cells in a given year
 land_transport_electric _share,--,Dictionary with planning horizons as keys.,The share of vehicles that uses electric vehicles (EV) in a given year
 land_transport_ice _share,--,Dictionary with planning horizons as keys.,The share of vehicles that uses internal combustion engines (ICE) in a given year. What is not EV or FCEV is oil-fuelled ICE.

doc/release_notes.rst

@@ -212,6 +212,8 @@ PyPSA-Eur 0.9.0 (5th January 2024)
   different directory using the new ``root_dir`` argument
   (https://github.com/PyPSA/pypsa-eur/pull/771).
 
+* Add option to optimize fuel type shares of land transport endogenously in configuration file (https://github.com/PyPSA/pypsa-eur/pull/734)
+
 * Rule ``purge`` now initiates a dialog to confirm if purge is desired
   (https://github.com/PyPSA/pypsa-eur/pull/745).
 

rules/build_electricity.smk

@@ -306,7 +306,7 @@ rule build_renewable_profiles:
         BENCHMARKS + "build_renewable_profiles_{technology}"
     threads: ATLITE_NPROCESSES
     resources:
-        mem_mb=ATLITE_NPROCESSES * 5000,
+        mem_mb=ATLITE_NPROCESSES * 9000,
     wildcard_constraints:
         technology="(?!hydro).*",  # Any technology other than hydro
     conda:
@@ -460,7 +460,7 @@ rule simplify_network:
         BENCHMARKS + "simplify_network/elec_s{simpl}"
     threads: 1
     resources:
-        mem_mb=12000,
+        mem_mb=24000,
     conda:
         "../envs/environment.yaml"
     script:

rules/build_sector.smk

@@ -424,7 +424,7 @@ rule build_industrial_production_per_country:
         + "industrial_production_per_country.csv",
     threads: 8
     resources:
-        mem_mb=1000,
+        mem_mb=5000,
     log:
         LOGS + "build_industrial_production_per_country.log",
     benchmark:

rules/retrieve.smk

@@ -378,4 +378,4 @@ if config["enable"]["retrieve"]:
         conda:
             "../envs/environment.yaml"
         script:
-            "../scripts/retrieve_monthly_fuel_prices.py"
+            "../scripts/retrieve_monthly_fuel_prices.py"

rules/solve_myopic.smk

@@ -23,6 +23,8 @@ rule add_existing_baseyear:
         existing_solar="data/existing_infrastructure/solar_capacity_IRENA.csv",
         existing_onwind="data/existing_infrastructure/onwind_capacity_IRENA.csv",
         existing_offwind="data/existing_infrastructure/offwind_capacity_IRENA.csv",
+        existing_transport=RESOURCES + "transport_data_s{simpl}_{clusters}.csv",
+        temp_air_total=RESOURCES + "temp_air_total_elec_s{simpl}_{clusters}.nc",
     output:
         RESULTS
         + "prenetworks-brownfield/elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}.nc",
@@ -108,6 +110,7 @@ rule solve_sector_network_myopic:
     resources:
         mem_mb=config["solving"]["mem"],
         walltime=config["solving"].get("walltime", "12:00:00"),
+        disk_mb=180000 #500000
     benchmark:
         (
             BENCHMARKS

scripts/add_brownfield.py

@@ -67,6 +67,16 @@ def add_brownfield(n, n_p, year):
                 & (c.df[f"{attr}_nom_opt"] < threshold)
             ],
         )
+        check_transport = ( (snakemake.config["sector"]["land_transport_electric_share"][year] is not None ) 
+        or (snakemake.config["sector"]["land_transport_fuel_cell_share"][year] is not None)
+        or (snakemake.config["sector"]["land_transport_ice_share"][year] is not None))
+
+        if not snakemake.config["sector"]["endogenous_transport"] or check_transport:
+            n_p.mremove(
+                c.name,
+                c.df.index[c.df.carrier.str.contains("land transport" or "V2G" or "EV battery storage" or "BEV charger")
+                ],
+            )
 
         # copy over assets but fix their capacity
         c.df[f"{attr}_nom"] = c.df[f"{attr}_nom_opt"]
@@ -119,6 +129,92 @@ def add_brownfield(n, n_p, year):
             n.links.loc[new_pipes, "p_nom"] = 0.0
             n.links.loc[new_pipes, "p_nom_min"] = 0.0
 
+def adjust_EVs(n, n_p, year):
+    # set p_min_pu and p_max_pu for solved network, so that only the EV link for the current time horizon
+    # is not constraint (constraining all land transport link with p_min_pu/p_max_pu while p_nom_extentable=True leads to infeasible by gurobi
+
+    for car in ["EV","fuel cell","oil"]:
+        cartype = "land transport " + car
+        if not n.links[(n.links.carrier==cartype) ].index.empty:
+            lifetime_EV = n.links.lifetime[n.links[(n.links.carrier==cartype) ].index[0]]
+            i = 0
+            while (year-lifetime_EV+i) < year:
+                if not n.links_t.p_min_pu[(n.links.filter(like=cartype +"-"+str(int(year-lifetime_EV+i)),axis=0)).index].empty:
+                    #p_set =  n_p.loads_t.p[n.loads[n.loads.carrier.str.contains('land transport demand')].index]
+                    #eff = n.links_t.efficiency[(n.links.filter(like="land transport EV-"+str(int(year-lifetime_EV+i)),axis=0)).index]
+                    #p_set = p_set.add_suffix(' EV-'+str(int(year-lifetime_EV+i)))
+                    #p_set = p_set.drop([col for col in p_set.columns if col in p_set.columns and col not in eff.columns], axis=1)
+                    #pnom = (p_set.divide(eff)).max()
+                    pu=n.links_t.p_min_pu[(n.links.filter(like=cartype +"-"+str(int(year-lifetime_EV+i)),axis=0)).index] #p_set.divide(eff)/pnom
+                    n.links_t.p_max_pu[(n.links.filter(like=cartype +"-"+str(int(year-lifetime_EV+i)),axis=0)).index] = pu.values
+                    #n.links_t.p_max_pu[(n.links.filter(like="land transpo-"+str(int(year-lifetime_EV+i)),axis=0)).index] = pu.values
+                i = i+1   
+
+    #Get existing capacity of EV batteries
+    ev_battery_storage = n.stores[n.stores.carrier.str.contains("Li ion|EV battery storage")]
+    sum_ev_battery = ev_battery_storage.groupby("bus")["e_nom"].sum()
+    print("EV battery storage",ev_battery_storage.groupby("bus")["e_nom"].sum())
+    sum_ev_battery.index = sum_ev_battery.index.str.replace(' EV battery','',regex=False)
+
+    #Set e_nom_max to remaining capacity in relation to maximum capacity of fully electric car fleet of 2015
+    ev_store = n.stores.carrier.str.contains('Li ion|EV battery storage') 
+    ev_store_ext = n.stores[ev_store].query("e_nom_extendable").index
+
+    #Caclulate existing EV land transport capacity 
+    ev_land_transport = n.links[n.links.carrier.str.contains("land transport EV")]
+    sum_ev_transport = ev_land_transport.groupby("bus0")["p_nom"].sum()
+    sum_ev_transport.index = sum_ev_transport.index.str.replace(' EV battery', '', regex=False)
+    print("sum ev transport",sum_ev_transport)
+
+    bev_availability = snakemake.config['sector']["bev_availability"]
+    for o in opts:
+
+        if "bevavail" not in o:
+            continue
+        oo = o.split("+")
+        bev_availability = float(oo[1])
+    ev_transport_to_store = sum_ev_transport/snakemake.config['sector']['EV_consumption_1car']*snakemake.config['sector']['bev_energy']*snakemake.config['sector']['bev_charge_avail']*bev_availability
+    print(sum_ev_battery)
+
+    if sum(sum_ev_battery):
+        enommax = n.stores.loc[ev_store_ext,'e_nom_max']
+        enommax.index = enommax.index.str.replace(' EV battery storage-'+str(year), '', regex=False) 
+        print(enommax)
+        #subtract existing battery capacity but add capacity if more EVs are used in land transport than is reflected in EV battery capacity
+        add_enommax = ev_transport_to_store-sum_ev_battery
+        add_enommax[add_enommax<0] = 0
+        newenommax = enommax-sum_ev_battery+add_enommax
+        newenommax[newenommax < 0] = 0
+        n.stores.loc[ev_store_ext,'e_nom_max'] = newenommax.values
+
+    BEV_charger = n.links[n.links.carrier.str.contains("BEV charger")]
+    sum_bev_charger = BEV_charger.groupby("bus0")["p_nom"].sum()
+    print("sum bev charger",sum_bev_charger)
+
+    sum_bev_charger.index = sum_bev_charger.index.str.replace(' low voltage', '', regex=False)
+
+    print("sum bev charger",sum_bev_charger)
+
+
+
+    ev_transport_to_bev = sum_ev_transport/snakemake.config['sector']['EV_consumption_1car']*snakemake.config['sector']['bev_charge_rate']*snakemake.config['sector']['bev_charge_avail']
+    print("ev transport to bev",ev_transport_to_bev)
+    bev_charger_ext = n.links[n.links.carrier.str.contains("BEV charger")].query("p_nom_extendable").index
+    pnommax = n.links.loc[bev_charger_ext,"p_nom_max"]
+    pnommax.index = pnommax.index.str.replace(' BEV charger-'+str(year), '', regex=False)
+    add_pnommax = ev_transport_to_bev-sum_bev_charger
+    print("add pnomax",add_pnommax)
+    add_pnommax[add_pnommax < 0] = 0
+    newpnommax = pnommax - sum_bev_charger + add_pnommax
+    print(pnommax, sum_bev_charger, add_pnommax)
+    newpnommax[newpnommax < 0] = 0
+    print(newpnommax)
+    print(len(newpnommax.values))
+    print(newpnommax.values.shape)
+    n.links.loc[bev_charger_ext,'p_nom_max'] = newpnommax.values
+    print(newpnommax)
+
+
 
 def disable_grid_expansion_if_LV_limit_hit(n):
     if not "lv_limit" in n.global_constraints.index:
@@ -177,5 +273,9 @@ def disable_grid_expansion_if_LV_limit_hit(n):
 
     disable_grid_expansion_if_LV_limit_hit(n)
 
+    opts = snakemake.wildcards.sector_opts.split("-")
+    if "T" in opts: # and snakemake.config["sector"]["endogenous_transport"]:
+        adjust_EVs(n, n_p, year)
+
     n.meta = dict(snakemake.config, **dict(wildcards=dict(snakemake.wildcards)))
     n.export_to_netcdf(snakemake.output[0])

scripts/add_electricity.py

@@ -127,6 +127,7 @@ def add_missing_carriers(n, carriers):
     Function to add missing carriers to the network without raising errors.
     """
     missing_carriers = set(carriers) - set(n.carriers.index)
+
     if len(missing_carriers) > 0:
         n.madd("Carrier", missing_carriers)
 
@@ -473,22 +474,38 @@ def attach_conventional_generators(
     # Define generators using modified ppl DataFrame
     caps = ppl.groupby("carrier").p_nom.sum().div(1e3).round(2)
     logger.info(f"Adding {len(ppl)} generators with capacities [GW] \n{caps}")
-
-    n.madd(
-        "Generator",
-        ppl.index,
-        carrier=ppl.carrier,
-        bus=ppl.bus,
-        p_nom_min=ppl.p_nom.where(ppl.carrier.isin(conventional_carriers), 0),
-        p_nom=ppl.p_nom.where(ppl.carrier.isin(conventional_carriers), 0),
-        p_nom_extendable=ppl.carrier.isin(extendable_carriers["Generator"]),
-        efficiency=ppl.efficiency,
-        marginal_cost=marginal_cost,
-        capital_cost=ppl.capital_cost,
-        build_year=ppl.datein.fillna(0).astype(int),
-        lifetime=(ppl.dateout - ppl.datein).fillna(np.inf),
-        **committable_attrs,
-    )
+    if snakemake.config["foresight"]== "myopic":
+        n.madd(
+            "Generator",
+            ppl.index,
+            carrier=ppl.carrier,
+            bus=ppl.bus,
+            p_nom_min=0, #ppl.p_nom.where(ppl.carrier.isin(conventional_carriers), 0),
+            p_nom=ppl.p_nom.where(ppl.carrier.isin(conventional_carriers), 0),
+            p_nom_extendable=ppl.carrier.isin(extendable_carriers["Generator"]),
+            efficiency=ppl.efficiency,
+            marginal_cost=marginal_cost,
+            capital_cost=ppl.capital_cost,
+            build_year=ppl.datein.fillna(0).astype(int),
+            lifetime=(ppl.dateout - ppl.datein).fillna(np.inf),
+            **committable_attrs,
+        )
+    else:
+        n.madd(
+            "Generator",
+            ppl.index,
+            carrier=ppl.carrier,
+            bus=ppl.bus,
+            p_nom_min=ppl.p_nom.where(ppl.carrier.isin(conventional_carriers), 0),
+            p_nom=ppl.p_nom.where(ppl.carrier.isin(conventional_carriers), 0),
+            p_nom_extendable=ppl.carrier.isin(extendable_carriers["Generator"]),
+            efficiency=ppl.efficiency,
+            marginal_cost=marginal_cost,
+            capital_cost=ppl.capital_cost,
+            build_year=ppl.datein.fillna(0).astype(int),
+            lifetime=(ppl.dateout - ppl.datein).fillna(np.inf),
+            **committable_attrs,
+        )
 
     for carrier in set(conventional_params) & set(carriers):
         # Generators with technology affected
@@ -878,6 +895,10 @@ def attach_line_rating(
     conventional_inputs = {
         k: v for k, v in snakemake.input.items() if k.startswith("conventional_")
     }
+    conventional_inputs2 = {
+        k: v for k, v in snakemake.input.items() if k.startswith("extendable_")
+    }
+    conventional_inputs.update(conventional_inputs2)
 
     if params.conventional["unit_commitment"]:
         unit_commitment = pd.read_csv(snakemake.input.unit_commitment, index_col=0)
@@ -958,6 +979,6 @@ def attach_line_rating(
         )
 
     sanitize_carriers(n, snakemake.config)
-
+   
     n.meta = snakemake.config
     n.export_to_netcdf(snakemake.output[0])