Headroom/reserve rate connections

h2integrate/converters/grid/grid.py

@@ -100,10 +100,26 @@ def setup(self):
             "electricity_sold",
             val=0.0,
             shape=n_timesteps,
-            units="kW",
+            units=self.commodity_rate_units,
             desc="Electricity sold to the grid",
         )
 
+        self.add_output(
+            "electricity_headroom_sold",
+            val=0.0,
+            shape=n_timesteps,
+            units=self.commodity_rate_units,
+            desc="Reserve capacity that could be sold to the grid",
+        )
+
+        self.add_output(
+            "electricity_headroom_out",
+            val=0.0,
+            shape=n_timesteps,
+            units=self.commodity_rate_units,
+            desc="Reserve capacity that could be bought from the grid",
+        )
+
         self.add_output(
             "electricity_unmet_demand",
             val=0.0,
@@ -148,6 +164,8 @@ def compute(self, inputs, outputs):
         max_production = (
             inputs["interconnection_size"] * len(outputs["electricity_out"]) * (self.dt / 3600)
         )
+        outputs["electricity_headroom_sold"] = interconnection_size - electricity_sold
+        outputs["electricity_headroom_out"] = interconnection_size - electricity_bought
         outputs["rated_electricity_production"] = inputs["interconnection_size"]
         outputs["total_electricity_produced"] = np.sum(outputs["electricity_out"]) * (
             self.dt / 3600

h2integrate/converters/grid/test/test_grid.py

@@ -118,6 +118,21 @@ def test_grid_performance_outputs(plant_config, subtests):
     with subtests.test(f"{commodity}_out length"):
         assert len(prob.get_val(f"comp.{commodity}_out", units=commodity_rate_units)) == n_timesteps
 
+    # Test that interconnect output headroom is greater than zero (plant oversized) and less than the rating
+    with subtests.test(f"0 < {commodity}_headroom_out < rated_{commodity}_production"):
+        assert np.all(prob.get_val(f"comp.{commodity}_headroom_out", units="MW") >= 0)
+        assert np.all(
+            prob.get_val(f"comp.{commodity}_headroom_out", units="MW")
+            <= prob.get_val(f"comp.rated_{commodity}_production", units="MW")
+        )
+
+    # Test that interconnect sales headroom is at the rating (nothing fed to grid)
+    with subtests.test(f"{commodity}_headroom_sold == rated_{commodity}_production"):
+        assert np.all(
+            prob.get_val(f"comp.{commodity}_headroom_sold", units="MW")
+            == prob.get_val(f"comp.rated_{commodity}_production", units="MW")
+        )
+
     # Test default values
     with subtests.test("operational_life default value"):
         assert prob.get_val("comp.operational_life", units="yr") == plant_life
@@ -220,6 +235,13 @@ def test_selling_electricity(plant_config, n_timesteps):
     actual_in = prob.get_val("grid.electricity_in")
     np.testing.assert_array_almost_equal(actual_in, electricity_in)
 
+    # The headroom should be the difference between electricity_in and the rating
+    headroom = prob.get_val("grid.electricity_headroom_sold")
+    np.testing.assert_array_almost_equal(
+        headroom,
+        tech_config["model_inputs"]["shared_parameters"]["interconnection_size"] - electricity_in,
+    )
+
 
 @pytest.mark.unit
 @pytest.mark.parametrize("n_timesteps", [10])
@@ -281,6 +303,13 @@ def test_varying_demand_profile(plant_config, n_timesteps):
     expected = np.clip(demand, 0, 100000)
     np.testing.assert_array_almost_equal(electricity_out, expected)
 
+    # The output headroom should be the difference between electricity_out and the rating
+    headroom = prob.get_val("grid.electricity_headroom_out")
+    np.testing.assert_array_almost_equal(
+        headroom,
+        tech_config["model_inputs"]["shared_parameters"]["interconnection_size"] - electricity_out,
+    )
+
 
 @pytest.mark.unit
 @pytest.mark.parametrize("n_timesteps", [10])
@@ -413,6 +442,10 @@ def test_grid_integration_dt_1800(subtests, tmp_path):
         annual_energy = h2i.prob.get_val("grid.annual_electricity_produced", units="kW*h/year")
         assert annual_energy == pytest.approx(expected_annual)
 
+    with subtests.test("headroom calculation accurately represents rating less demand"):
+        headroom = h2i.prob.get_val("grid.electricity_headroom_out", units="kW")
+        assert headroom == pytest.approx(tech_config["technologies"]["grid"]["model_inputs"]["shared_parameters"]["interconnection_size"] - demand)
+
 
 @pytest.mark.unit
 @pytest.mark.parametrize("n_timesteps", [8760])

h2integrate/converters/natural_gas/natural_gas_cc_ct.py

@@ -119,6 +119,13 @@ def setup(self):
             desc="Natural gas input energy",
         )
 
+        self.add_output(
+            f"{self.commodity}_headroom_out",  # unused but "accessible" portion of rated power
+            val=0.0,
+            shape=self.n_timesteps,
+            units=self.commodity_rate_units,
+        )
+
         self.add_output(
             "unmet_electricity_demand",
             val=0.0,
@@ -170,6 +177,10 @@ def compute(self, inputs, outputs):
 
         outputs["electricity_out"] = electricity_out
         outputs["natural_gas_consumed"] = natural_gas_consumed
+        outputs["electricity_headroom_out"] = np.minimum(  # we are limitied by either
+            natural_gas_available / heat_rate_mmbtu_per_mwh,  # the power available in the natural gas supply
+            system_capacity,  # or the rated power of the system
+        ) - electricity_out  # and subtracting out what we're using gives the available excess capacity
 
         outputs["rated_electricity_production"] = inputs["system_capacity"]
 

h2integrate/converters/natural_gas/test/test_natural_gas_models.py

@@ -167,6 +167,14 @@ def test_ngcc_performance_outputs(plant_config, ngcc_performance_params, subtest
     with subtests.test(f"{commodity}_out length"):
         assert len(prob.get_val(f"comp.{commodity}_out", units=commodity_rate_units)) == n_timesteps
 
+    # Test that headroom is greater than zero (plant oversized) and less than the rating
+    with subtests.test(f"0 < {commodity}_headroom_out < rated_{commodity}_production"):
+        assert np.all(prob.get_val(f"comp.{commodity}_headroom_out", units="MW") >= 0)
+        assert np.all(
+            prob.get_val(f"comp.{commodity}_headroom_out", units="MW")
+            <= prob.get_val(f"comp.rated_{commodity}_production", units="MW")
+        )
+
     # Test default values
     with subtests.test("operational_life default value"):
         assert prob.get_val("comp.operational_life", units="yr") == plant_life
@@ -210,6 +218,13 @@ def test_ngcc_performance(plant_config, ngcc_performance_params, subtests):
         # Check average output is 100 MW
         assert pytest.approx(np.mean(electricity_out), rel=1e-6) == 100.0
 
+    headroom_out = prob.get_val("electricity_headroom_out")
+
+    with subtests.test("NGCC Headroom Output"):
+        # Headroom should be capacity less expected output (here zero)
+        expected_headroom = ngcc_performance_params["system_capacity_mw"] - expected_output
+        assert np.allclose(headroom_out, expected_headroom, rtol=1.0e-6)
+
 
 @pytest.mark.regression
 def test_ngct_performance(plant_config, ngct_performance_params, subtests):
@@ -247,6 +262,13 @@ def test_ngct_performance(plant_config, ngct_performance_params, subtests):
         # Check average output is 50 MW
         assert pytest.approx(np.mean(electricity_out), rel=1e-6) == 50.0
 
+    headroom_out = prob.get_val("electricity_headroom_out")
+
+    with subtests.test("NGCT Headroom Output"):
+        # Headroom should be capacity less expected output
+        expected_headroom = ngct_performance_params["system_capacity_mw"] - expected_output
+        assert np.allclose(headroom_out, expected_headroom, rtol=1.0e-6)
+
 
 @pytest.mark.unit
 def test_ngcc_cost(plant_config, ngcc_cost_params, subtests):
@@ -397,3 +419,11 @@ def test_ngcc_performance_demand(plant_config, ngcc_performance_params, subtests
             pytest.approx(np.max(electricity_out), rel=1e-6)
             == ngcc_performance_params["system_capacity_mw"]
         )
+
+    headroom_out = prob.get_val("electricity_headroom_out")
+
+    with subtests.test("NGCC Electricity Headroom"):
+        # Headroom should be capacity less expected output (here zero)
+        expected_headroom = ngcc_performance_params["system_capacity_mw"] - expected_output
+        assert np.allclose(headroom_out, expected_headroom, rtol=1.0e-6)
+

h2integrate/storage/storage_performance_model.py

@@ -137,6 +137,13 @@ def setup(self):
 
         super().setup()
 
+        self.add_output(
+            f"{self.commodity}_headroom_out",
+            val=0.0,
+            shape=self.n_timesteps,
+            units=self.commodity_rate_units,
+        )  # this represents the excess capacity that could have been dumped
+
     def compute(self, inputs, outputs, discrete_inputs=[], discrete_outputs=[]):
         """Run the storage performance model."""
         self.current_soc = self.config.init_soc_fraction
@@ -150,3 +157,25 @@ def compute(self, inputs, outputs, discrete_inputs=[], discrete_outputs=[]):
         outputs = self.run_storage(
             charge_rate, discharge_rate, storage_capacity, inputs, outputs, discrete_inputs
         )
+
+        # add headroom calculation
+        headroom_discharge = (
+            (outputs["SOC"]/100.0 - self.config.min_soc_fraction) * storage_capacity / self.dt_hr
+        )  # i *could've* dumped the state of charge by this much
+
+        available_discharge = np.maximum(
+            0.0,  # at worst no discharge is available
+            np.minimum(
+                discharge_rate,  # the fundamental limit on discharge rate is the max
+                headroom_discharge,  # but that could be limited by the available capacity
+            ),
+        )  # this is the max i could discharge right now
+
+        outputs[f"{self.commodity}_headroom_out"] = (
+            available_discharge*self.config.discharge_efficiency  # i could dump this much power out total
+            - outputs[f"{self.commodity}_out"]  #  remove current discharge, ADD charge also (not sure if accounting is correct)
+            # # the below was my first attempt, to ignore charging, but I think the
+            # # current charging current should be treated as "available" and *should*
+            # # be accounted as reserve power
+            # - np.maximum(0.0, outputs[f"{self.commodity}_out"])  #  remove current discharge, throw away current charging?
+        )

h2integrate/storage/test/test_storage_performance_model.py

@@ -104,6 +104,32 @@ def test_generic_storage_with_simple_control_dmd_lessthan_charge_rate(plant_conf
             == performance_model_config["init_soc_fraction"]
         )
 
+    with subtests.test("Headroom is non-negative"):
+        assert np.all(prob.get_val("storage.hydrogen_headroom_out") >= 0.0)
+
+    with subtests.test("Headroom doesn't exceed rated plus output"):
+        # when rating-limited, headroom can't exceed the max discharge rate
+        # minus the output (positive if discharging, negative if charging)
+        assert np.all(
+            prob.get_val("storage.hydrogen_headroom_out")
+            <= prob.model.storage.config.max_discharge_rate
+            - prob.get_val("storage.hydrogen_out")
+        )
+
+    with subtests.test("Headroom doesn't exceed liquidatable capacity plus output"):
+        # when rating-limited, headroom can't exceed the max discharge rate
+        # minus the output (positive if discharging, negative if charging),
+        # which represents, respectively, rating that is not available as
+        # headroom and excess power capacity that can be diverted to generation
+        assert np.all(
+            prob.get_val("storage.hydrogen_headroom_out")
+            <= prob.model.storage.config.max_capacity*(
+                prob.get_val("storage.SOC")/100.0
+                - prob.model.storage.config.min_soc_fraction
+            )
+            - prob.get_val("storage.hydrogen_out")
+        )
+
     indx_soc_increase = np.argwhere(
         np.diff(prob.model.get_val("storage.SOC", units="unitless"), prepend=True) > 0
     ).flatten()