Update share_hybrid

src/G1_oemof_create_model.py

@@ -449,18 +449,25 @@ def build(experiment, case_dict):
         logging.info(
             "Added constraint: Stability though actual generation of diesel generators and backup through batteries."
         )
+
+        if critical_constraint is True:
+            critical_demand_sinks = [sink_demand_ac_critical, sink_demand_ac_reducable]
+        else:
+            critical_demand_sinks = None
+
         constraints_custom.hybrid(
             model,
             case_dict,
             experiment=experiment,
             storage=storage,
             sink_demand=sink_demand_ac,
-            demand_ac_critical=demand_ac_critical,
+            critical_demand=critical_demand_sinks,
             genset=genset,
             pcc_consumption=pointofcoupling_consumption,
             source_shortage=source_shortage,
             el_bus_ac=bus_electricity_ac,
             el_bus_dc=bus_electricity_dc,
+
         )
     else:
         logging.warning(
@@ -473,9 +480,9 @@ def build(experiment, case_dict):
     if case_dict[CRITICAL_CONSTRAINT] is True:
         logging.info("Added constraint: Critical demand fulfilled at all timesteps")
 
-        if case_dict[STABILITY_CONSTRAINT] != False:
-            raise ValueError(
-                "At the moment you cannot use the stability constraint with critical demand"
+        if case_dict[STABILITY_CONSTRAINT] is True:
+            logging.warning(
+                "At the moment the stability constraint share_hybrid together with critical demand might be erroneous, use at your own risks!"
             )
 
         # list of assets which bring energy into the AC bus which could be used to fullfill critical demand

src/G2b_constraints_custom.py

@@ -255,37 +255,40 @@ def hybrid(
     source_shortage,
     el_bus_ac,
     el_bus_dc,
-    demand_ac_critical=None,
+    critical_demand=None,
 ):
 
     stability_limit = experiment[SHORTAGE_LIMIT]
-    print(case_dict[ALLOW_SHORTAGE])
-    print(case_dict[GENSET_FIXED_CAPACITY])
-    print(case_dict[STORAGE_FIXED_CAPACITY])
-    import ipdb
-
-    ipdb.set_trace()
 
     # ac_critical_demand = model.flow[el_bus_ac, demand_ac_critical, t]
 
     def stability_rule_capacity(model, t):
         expr = 0
         ## ------- Get demand at t ------- #
-        # TODO adapt this
+
         demand = model.flow[el_bus_ac, sink_demand, t]
+
+        if critical_demand is not None:
+            sink_demand_ac_critical, sink_demand_ac_reducable = critical_demand
+            demand += model.flow[el_bus_ac, sink_demand_ac_critical, t]
+            demand += model.flow[el_bus_ac, sink_demand_ac_reducable, t]
+
         expr += -stability_limit * demand
 
-        # Current state
+        # old state
         # -stability*demand_nc + stability*shortage + SUM_i genset_i + storage_out (can be negative or positive) >= 0
 
-        # what do we want?
-        # a) -stability*(demand_nc + demand_c) + stabililty*shortage + SUM_i genset_i + storage_out (can be negative or positive) >= 0
-        # b) -stability*(demand_c) + stability*shortage + SUM_i genset_i + storage_out (can be negative or positive) >= 0
-
+        # current state
+        # a) -stability*(demand_nc + demand_nc reducable  + demand_c) + stabililty*shortage2 + SUM_i genset_i + storage_out (can be negative or positive) >= 0
+        # shortage2 is shortage + max_allowed_reducable_demand - reducable_demand = shortage + non_supplied_reducable_demand
         ## ------- Get shortage at t------- #
         if case_dict[ALLOW_SHORTAGE] is True:
             # TODO this is considered
-            shortage = model.flow[source_shortage, el_bus_ac, t]
+            if critical_demand is not None:
+                max_allowed_demand_reduction = sink_demand_ac_reducable.inputs[el_bus_ac].max[t]
+                shortage = model.flow[source_shortage, el_bus_ac, t] - model.flow[el_bus_ac, sink_demand_ac_reducable, t] + max_allowed_demand_reduction
+            else:
+                shortage = model.flow[source_shortage, el_bus_ac, t]
             expr += stability_limit * shortage
 
         ## ------- Generation Diesel ------- #
@@ -335,19 +338,32 @@ def stability_rule_power(model, t):
         expr = 0
         ## ------- Get demand at t ------- #
         demand = model.flow[el_bus_ac, sink_demand, t]
+
+        if critical_demand is not None:
+            sink_demand_ac_critical, sink_demand_ac_reducable = critical_demand
+            demand += model.flow[el_bus_ac, sink_demand_ac_critical, t]
+            demand += model.flow[el_bus_ac, sink_demand_ac_reducable, t]
+
+
         expr += -stability_limit * demand
 
-        # Current state
+        # old state
         # -stability_limit*demand_nc + stability*shortage + SUM_i genset_i + storage_out (can be negative or positive) >= 0
 
-        # what do we want?
-        # a) -stability*(demand_nc + demand_c) + stabililty*shortage + SUM_i genset_i + storage_out (can be negative or positive) >= 0
-        # b) -stability*(demand_c) + stability*shortage + SUM_i genset_i + storage_out (can be negative or positive) >= 0
+        # current state
+        # a) -stability*(demand_nc + demand_nc reducable  + demand_c) + stabililty*shortage2 + SUM_i genset_i + storage_out (can be negative or positive) >= 0
+        # shortage2 is shortage + max_allowed_reducable_demand - reducable_demand = shortage + non_supplied_reducable_demand
 
         ## ------- Get shortage at t------- #
         if case_dict[ALLOW_SHORTAGE] is True:
             # TODO this is considered
-            shortage = model.flow[source_shortage, el_bus_ac, t]
+            if critical_demand is not None:
+                max_allowed_demand_reduction = sink_demand_ac_reducable.inputs[el_bus_ac].max[t]
+
+                shortage = model.flow[source_shortage, el_bus_ac, t] - model.flow[
+                    el_bus_ac, sink_demand_ac_reducable, t] + max_allowed_demand_reduction
+            else:
+                shortage = model.flow[source_shortage, el_bus_ac, t]
             expr += +stability_limit * shortage
 
         ## ------- Generation Diesel ------- #
@@ -746,8 +762,11 @@ def critical(
 
     def meet_critical_ac_demand_rule(model, t):
         """
-        wind + genset + inverter + pcc_consumption + shortage - rectifier - pcc_feedin - non-critical demand - critical demand >= 0
+        wind + genset + inverter + pcc_consumption + shortage - rectifier - pcc_feedin - non-critical demand - non-critical demand reducable - critical demand >= 0
         """
+
+        # TODO maybe implement wind + genset + inverter + pcc_consumption - rectifier - pcc_feedin  - supplied_critical demand >= 0
+        # because if currently critical demand is 0 then the equation is satisfied, but we don't want the critical demand to be 0
         ac_production = 0
         for ac_asset in ac_generation_assets:
             ac_production += model.flow[ac_asset, el_bus_ac, t]
@@ -762,7 +781,7 @@ def meet_critical_ac_demand_rule(model, t):
 
     def meet_critical_dc_demand_rule(model, t):
         """
-        PV + storage-discharge + shortage + rectifier - inverter - storage-charge - non-critical demand - critical demand >= 0
+        PV + storage-discharge + shortage + rectifier - inverter - storage-charge - non-critical demand - non-critical demand reducable - critical demand >= 0
         """
         dc_production = 0
         for dc_asset in dc_generation_assets:

src/G3_oemof_evaluate.py

@@ -258,7 +258,7 @@ def get_shortage(
     logging.debug("Evaluate flow: shortage")
 
     # Get flow
-    shortage = pd.Series([0 for i in e_flows_df.index], index=e_flows_df.index)
+    total_demand_reduction = pd.Series([0 for i in e_flows_df.index], index=e_flows_df.index)
 
     critical_constraint = case_dict.get(CRITICAL_CONSTRAINT, False)
 
@@ -278,9 +278,9 @@ def get_shortage(
             e_flows_df = join_e_flows_df(shortage_ac, DEMAND_SHORTAGE_AC, e_flows_df)
 
             if case_dict[EVALUATION_PERSPECTIVE] == AC_SYSTEM:
-                shortage += shortage_ac
+                total_demand_reduction += shortage_ac
             else:
-                shortage += shortage_ac / experiment[INVERTER_DC_AC_EFFICIENCY]
+                total_demand_reduction += shortage_ac / experiment[INVERTER_DC_AC_EFFICIENCY]
 
         if electricity_bus_dc != None:
 
@@ -297,9 +297,9 @@ def get_shortage(
             e_flows_df = join_e_flows_df(shortage_dc, DEMAND_SHORTAGE_DC, e_flows_df)
 
             if case_dict[EVALUATION_PERSPECTIVE] == AC_SYSTEM:
-                shortage += shortage_dc / experiment[RECTIFIER_AC_DC_EFFICIENCY]
+                total_demand_reduction += shortage_dc / experiment[RECTIFIER_AC_DC_EFFICIENCY]
             else:
-                shortage += shortage_dc
+                total_demand_reduction += shortage_dc
 
             if critical_constraint is True:
 
@@ -312,19 +312,21 @@ def get_shortage(
                 share_reducable_demand_not_supplied = max_allowed_demand_reduction - e_flows_df[DEMAND_NON_CRITICAL_REDUCABLE]
                 # update the demand reduction: the total demand reduction, or demand shortage, is the sum of what was
                 # provided by the shortage source and what could not be supplied by the reducable demand sink
-                shortage += share_reducable_demand_not_supplied
 
-                under_supply = shortage > max_allowed_demand_reduction
+                e_flows_df = join_e_flows_df(share_reducable_demand_not_supplied, "reducable_demand_not_supplied", e_flows_df)
+                total_demand_reduction += share_reducable_demand_not_supplied
+
+                under_supply = total_demand_reduction > max_allowed_demand_reduction
                 if under_supply.any():
                     ts = "\n".join([str(d) for d in under_supply.loc[under_supply == True].index.values])
                     logging.warning(f"The total demand reduction, or demand shortage, exceeds the allowed {100 *case_dict[MAX_SHORTAGE]}% of the non-critical demand:\n {ts}")
-                under_supply = shortage > e_flows_df[DEMAND_NON_CRITICAL]
+                under_supply = total_demand_reduction > e_flows_df[DEMAND_NON_CRITICAL]
                 if under_supply.any():
                     ts = "\n".join([str(d) for d in under_supply.loc[under_supply == True].index.values])
                     logging.error(f"A portion of the critical demand could not be provided :\n {ts}")
 
 
-        demand_supplied = e_flows_df[DEMAND] - shortage
+        demand_supplied = e_flows_df[DEMAND] - total_demand_reduction
 
 
         annual_value(
@@ -334,9 +336,9 @@ def get_shortage(
             case_dict,
         )
         annual_value(
-            TOTAL_DEMAND_SHORTAGE_ANNUAL_KWH, shortage, oemof_results, case_dict
+            TOTAL_DEMAND_SHORTAGE_ANNUAL_KWH, total_demand_reduction, oemof_results, case_dict
         )
-        e_flows_df = join_e_flows_df(shortage, DEMAND_SHORTAGE, e_flows_df)
+        e_flows_df = join_e_flows_df(total_demand_reduction, DEMAND_SHORTAGE, e_flows_df)
         e_flows_df = join_e_flows_df(demand_supplied, DEMAND_SUPPLIED, e_flows_df)
     else:
         oemof_results.update(

src/G3a_economic_evaluation.py

@@ -487,12 +487,13 @@ def expenditures_shortage(oemof_results, experiment):
     )
 
     if experiment[INCLUDE_SHORTAGE_PENALTY_COSTS_IN_LCOE] is True:
-        oemof_results.update(
-            {
-                ANNUITY: oemof_results[ANNUITY]
-                + oemof_results[EXPENDITURES_SHORTAGE_ANNUAL]
-            }
-        )
+        pass
+        # oemof_results.update(
+        #     {
+        #         ANNUITY: oemof_results[ANNUITY]
+        #         + oemof_results[EXPENDITURES_SHORTAGE_ANNUAL]
+        #     }
+        # )
     else:
         oemof_results.update(
             {

src/G4_output_functions.py

@@ -167,7 +167,7 @@ def save_mg_flows(experiment, case_dict, e_flows_df, filename):
         flows_connected_to_electricity_mg_bus = (
             [DEMAND_NON_CRITICAL, DEMAND_CRITICAL, DEMAND_NON_CRITICAL_REDUCABLE]
             + flows_connected_to_electricity_mg_bus[:2]
-            + [DEMAND_AC_CRITICAL, DEMAND_DC_CRITICAL]
+            + [DEMAND_AC_CRITICAL, DEMAND_DC_CRITICAL, "reducable_demand_not_supplied"]
             + flows_connected_to_electricity_mg_bus[2:]
         )