add option to take today's district heating share

config.default.yaml

@@ -141,6 +141,12 @@ existing_capacities:
 sector:
   central: true
   central_fraction: 0.6
+  district_heating_increase: true
+  dh_strength:
+    2020: 0   # starting at today's share
+    2030: 0.2
+    2040: 0.5
+    2050: 1    # maximum possible share defined in central fraction
   bev_dsm_restriction_value: 0.75  #Set to 0 for no restriction on BEV DSM
   bev_dsm_restriction_time: 7  #Time at which SOC of BEV has to be dsm_restriction_value
   transport_heating_deadband_upper: 20.

scripts/prepare_sector_network.py

@@ -1179,12 +1179,11 @@ def add_heat(n, costs):
 
     sectors = ["residential", "services"]
 
-    nodes = create_nodes_for_heat_sector()
+
+    nodes, dist_fraction, urban_fraction = create_nodes_for_heat_sector()
 
     #NB: must add costs of central heating afterwards (EUR 400 / kWpeak, 50a, 1% FOM from Fraunhofer ISE)
 
-    urban_fraction = options['central_fraction'] * pop_layout["urban"] / pop_layout[["urban", "rural"]].sum(axis=1)
-
     # exogenously reduce space heat demand
     if options["reduce_space_heat_exogenously"]:
         dE = get(options["reduce_space_heat_exogenously_factor"], investment_year)
@@ -1215,10 +1214,17 @@ def add_heat(n, costs):
         ## Add heat load
 
         for sector in sectors:
+            # heat demand weighting
             if "rural" in name:
                 factor = 1 - urban_fraction[nodes[name]]
-            elif "urban" in name:
+            elif "urban central" in name:
-                factor = urban_fraction[nodes[name]]
+                factor = dist_fraction[nodes[name]]
+            elif "urban decentral" in name:
+                factor = urban_fraction[nodes[name]] - \
+                    dist_fraction[nodes[name]]
+            else:
+                factor = None
+
             if sector in name:
                 heat_load = heat_demand[[sector + " water",sector + " space"]].groupby(level=1,axis=1).sum()[nodes[name]].multiply(factor)
 
@@ -1504,23 +1510,54 @@ def create_nodes_for_heat_sector():
     # urban are areas with high heating density
     # urban can be split into district heating (central) and individual heating (decentral)
 
+    ct_urban = pop_layout.urban.groupby(pop_layout["ct"]).sum()
+    pop_layout["urban_ct_fraction"] = pop_layout["urban"] / \
+            pop_layout["ct"].map(ct_urban.get)
+    # todays district heating share per country
+    dist_heat_share_ct = pd.read_csv(snakemake.input.dh_share, index_col=0,
+                                  usecols=[0,1]).dropna()/100
+    dist_heat_share = pop_layout.ct.map(dist_heat_share_ct["district heating share"])
+
+
     sectors = ["residential", "services"]
 
     nodes = {}
+    urban_fraction = pop_layout["urban"] / \
+                     (pop_layout[["urban", "rural"]].sum(axis=1))
+
     for sector in sectors:
         nodes[sector + " rural"] = pop_layout.index
-
-        if options["central"]:
-            # TODO: this looks hardcoded, move to config
-            urban_decentral_ct = pd.Index(["ES", "GR", "PT", "IT", "BG"])
-            nodes[sector + " urban decentral"] = pop_layout.index[pop_layout.ct.isin(urban_decentral_ct)]
-        else:
         nodes[sector + " urban decentral"] = pop_layout.index
 
-    # for central nodes, residential and services are aggregated
+    if options["central"] and not options['district_heating_increase']:
-    nodes["urban central"] = pop_layout.index.symmetric_difference(nodes["residential urban decentral"])
+        central_fraction = options['central_fraction']
+        dist_fraction = central_fraction * urban_fraction
+        nodes["urban central"] = dist_fraction.index
 
-    return nodes
+    if options['district_heating_increase']:  # take current district heating share
+        dist_fraction = dist_heat_share * \
+            pop_layout["urban_ct_fraction"] / pop_layout["fraction"]
+        nodes["urban central"] = dist_fraction.index
+        # if district heating share larger than urban fraction -> set urban
+        # fraction to district heating share
+        urban_fraction = pd.concat([urban_fraction, dist_fraction],
+                                   axis=1).max(axis=1)
+        diff = urban_fraction - dist_fraction
+        dist_fraction += diff * get(options["dh_strength"], investment_year)
+        print("************************************")
+        print(
+            "the current DH share compared to the maximum possible is increased \
+               \n by a factor of ",
+            get(options["dh_strength"], investment_year),
+            "resulting DH share: ",
+            dist_fraction)
+        print("**********************************")
+
+    else:
+        dist_fraction = urban_fraction * 0
+        nodes["urban central"] = dist_fraction.index
+
+    return nodes, dist_fraction, urban_fraction
 
 
 def add_biomass(n, costs):
@@ -1781,7 +1818,7 @@ def add_industry(n, costs):
 
     if options["oil_boilers"]:
 
-        nodes_heat = create_nodes_for_heat_sector()
+        nodes_heat = create_nodes_for_heat_sector()[0]
 
         for name in ["residential rural", "services rural", "residential urban decentral", "services urban decentral"]:
 
@@ -1973,17 +2010,18 @@ def limit_individual_line_extension(n, maxext):
     hvdc = n.links.index[n.links.carrier == 'DC']
     n.links.loc[hvdc, 'p_nom_max'] = n.links.loc[hvdc, 'p_nom'] + maxext
 
-
+#%%
 if __name__ == "__main__":
     if 'snakemake' not in globals():
         from helper import mock_snakemake
         snakemake = mock_snakemake(
             'prepare_sector_network',
             simpl='',
-            clusters=48,
+            opts="",
+            clusters="37",
             lv=1.0,
             sector_opts='Co2L0-168H-T-H-B-I-solar3-dist1',
-            planning_horizons=2020,
+            planning_horizons="2020",
         )
 
     logging.basicConfig(level=snakemake.config['logging_level'])
@@ -1998,7 +2036,7 @@ if __name__ == "__main__":
     n = pypsa.Network(snakemake.input.network, override_component_attrs=overrides)
 
     pop_layout = pd.read_csv(snakemake.input.clustered_pop_layout, index_col=0)
-    Nyears = n.snapshot_weightings.generators.sum() / 8760
+    Nyears = n.snapshot_weightings.sum() / 8760
 
     costs = prepare_costs(snakemake.input.costs,
                           snakemake.config['costs']['USD2013_to_EUR2013'],