Merge pull request #404 from NREL/enforce-waste-heat-tracking

Enforce waste heat tracking

CHANGELOG.md

@@ -23,6 +23,7 @@ Classify the change according to the following categories:
     ### Deprecated
     ### Removed
 
+
 ## v0.47.0
 ### Added
 - Added inputs options and handling for ProcessHeatLoad for scaling annual or monthly fuel consumption values with reference hourly profiles, same as other loads
@@ -34,10 +35,17 @@ Classify the change according to the following categories:
 - Updated `src/core/doe_commercial_reference_building_loads.jl` to include **ProcessHeatLoad** for built-in load handling.
 - Refactored various functions to ensure **ProcessHeatLoad** is processed correctly in line with other heating loads.
 - When the URDB response `energyratestructure` has a "unit" value that is not "kWh", throw an error instead of averaging rates in each energy tier.
-- Changed default Financial **owner_tax_rate_fraction** and **offtaker_tax_rate_fraction** from 0.257 to 0.26 to align with API and user manual defaults. 
+- Refactored heating flow constraints to be in ./src/constraints/thermal_tech_constraints.jl instead of its previous separate locations in the storage and turbine constraints.
 ### Fixed
 - Updated the PV result **lifecycle_om_cost_after_tax** to account for the third-party factor for third-party ownership analyses.
 - Convert `max_electric_load_kw` to _Float64_ before passing to function `get_chp_defaults_prime_mover_size_class`
+- Fixed a bug in which excess heat from one heating technology resulted in waste heat from another technology.
+- Modified thermal waste heat constraints for heating technologies to avoid errors in waste heat results tracking.
+
+## v0.46.2
+### Changed
+- When the URDB response `energyratestructure` has a "unit" value that is not "kWh", throw an error instead of averaging rates in each energy tier.
+- Changed default Financial **owner_tax_rate_fraction** and **offtaker_tax_rate_fraction** from 0.257 to 0.26 to align with API and user manual defaults. 
 
 ## v0.46.1
 ### Changed

src/constraints/storage_constraints.jl

@@ -103,60 +103,6 @@ end
 
 function add_hot_thermal_storage_dispatch_constraints(m, p, b; _n="")
 
-    # # Constraint (4f)-1: (Hot) Thermal production sent to storage or grid must be less than technology's rated production
-	# # Constraint (4f)-1a: BoilerTechs
-	for t in p.techs.boiler
-		if !isempty(p.techs.steam_turbine) && (t in p.techs.can_supply_steam_turbine)
-            @constraint(m, [b in p.s.storage.types.hot, q in p.heating_loads, ts in p.time_steps],
-                    m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts] + m[Symbol("dvThermalToSteamTurbine"*_n)][t,q,ts]  <=
-                    m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
-                    )
-        else
-            @constraint(m, [b in p.s.storage.types.hot, q in p.heating_loads, ts in p.time_steps],
-                    m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts]  <=
-                    m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
-                    )
-        end
-    end
-
-    if !isempty(p.techs.electric_heater)
-        for t in p.techs.electric_heater
-            if !isempty(p.techs.steam_turbine) && (t in p.techs.can_supply_steam_turbine)
-                @constraint(m, [b in p.s.storage.types.hot, q in p.heating_loads, ts in p.time_steps],
-                        m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts] + m[Symbol("dvThermalToSteamTurbine"*_n)][t,q,ts]  <=
-                        m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
-                        )
-            else
-                @constraint(m, [b in p.s.storage.types.hot, q in p.heating_loads, ts in p.time_steps],
-                        m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts]  <=
-                        m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
-                        )
-            end
-        end
-    end
-
-    # Constraint (4f)-1b: SteamTurbineTechs
-	if !isempty(p.techs.steam_turbine)
-		@constraint(m, SteamTurbineTechProductionFlowCon[b in p.s.storage.types.hot, t in p.techs.steam_turbine, q in p.heating_loads, ts in p.time_steps],
-			m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts] <=  m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
-			)
-	end
-
-    # # Constraint (4g): CHP Thermal production sent to storage or grid must be less than technology's rated production
-	if !isempty(p.techs.chp)
-		if !isempty(p.techs.steam_turbine) && p.s.chp.can_supply_steam_turbine
-            @constraint(m, CHPTechProductionFlowCon[b in p.s.storage.types.hot, t in p.techs.chp, q in p.heating_loads, ts in p.time_steps],
-                    m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts] + m[Symbol("dvProductionToWaste"*_n)][t,q,ts] + m[Symbol("dvThermalToSteamTurbine"*_n)][t,q,ts] <=
-                    m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
-                    )
-        else
-            @constraint(m, CHPTechProductionFlowCon[b in p.s.storage.types.hot, t in p.techs.chp, q in p.heating_loads, ts in p.time_steps],
-                    m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts] + m[Symbol("dvProductionToWaste"*_n)][t,q,ts]  <=
-                    m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
-                    )
-        end
-	end
-
     # Constraint (4j)-1: Reconcile state-of-charge for (hot) thermal storage
 	@constraint(m, [b in p.s.storage.types.hot, ts in p.time_steps],
     m[Symbol("dvStoredEnergy"*_n)][b,ts] == m[Symbol("dvStoredEnergy"*_n)][b,ts-1] + (1/p.s.settings.time_steps_per_hour) * (

src/constraints/thermal_tech_constraints.jl

@@ -24,6 +24,43 @@ function add_boiler_tech_constraints(m, p; _n="")
 end
 
 function add_heating_tech_constraints(m, p; _n="")
+    # Constraint (7_heating_flow): Flows to Steam turbine, waste, and turbine must be less than or equal to total production
+    if !isempty(p.techs.steam_turbine)
+        if !isempty(p.s.storage.types.hot)
+            @constraint(m, [t in p.techs.can_supply_steam_turbine, q in p.heating_loads, ts in p.time_steps],
+                sum(m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts] for b in p.s.storage.types.hot) + m[Symbol("dvThermalToSteamTurbine"*_n)][t,q,ts] + m[Symbol("dvProductionToWaste"*_n)][t,q,ts]  <=
+                m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
+            )
+            if !isempty(setdiff(union(p.techs.heating, p.techs.chp),p.techs.can_supply_steam_turbine))
+                @constraint(m, [t in setdiff(union(p.techs.heating,p.techs.chp),p.techs.can_supply_steam_turbine), q in p.heating_loads, ts in p.time_steps],
+                    sum(m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts] for b in p.s.storage.types.hot) + m[Symbol("dvProductionToWaste"*_n)][t,q,ts]  <=  
+                    m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
+                )
+            end
+        else
+            @constraint(m, [t in p.techs.can_supply_steam_turbine, q in p.heating_loads, ts in p.time_steps],
+                m[Symbol("dvThermalToSteamTurbine"*_n)][t,q,ts] + m[Symbol("dvProductionToWaste"*_n)][t,q,ts]  <=
+                m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
+            )
+            if !isempty(setdiff(union(p.techs.heating, p.techs.chp),p.techs.can_supply_steam_turbine))
+                @constraint(m, [t in setdiff(union(p.techs.heating,p.techs.chp),p.techs.can_supply_steam_turbine), q in p.heating_loads, ts in p.time_steps],
+                    m[Symbol("dvProductionToWaste"*_n)][t,q,ts]  <=  m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
+                )
+            end
+        end
+    else
+        if !isempty(p.s.storage.types.hot)
+            @constraint(m, [t in union(p.techs.heating, p.techs.chp), q in p.heating_loads, ts in p.time_steps],
+                sum(m[Symbol("dvHeatToStorage"*_n)][b,t,q,ts] for b in p.s.storage.types.hot) + m[Symbol("dvProductionToWaste"*_n)][t,q,ts]  <=
+                m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
+            )
+        else
+            @constraint(m, [t in union(p.techs.heating, p.techs.chp), q in p.heating_loads, ts in p.time_steps],
+                m[Symbol("dvProductionToWaste"*_n)][t,q,ts]  <=  m[Symbol("dvHeatingProduction"*_n)][t,q,ts]
+            )
+        end
+    end
+
     # Constraint (7_heating_prod_size): Production limit based on size for non-electricity-producing heating techs
     if !isempty(setdiff(p.techs.heating, union(p.techs.elec, p.techs.ghp)))
         @constraint(m, [t in setdiff(p.techs.heating, union(p.techs.elec, p.techs.ghp)), ts in p.time_steps],
@@ -35,19 +72,23 @@ function add_heating_tech_constraints(m, p; _n="")
         if !(t in p.techs.can_serve_space_heating)
             for ts in p.time_steps
                 fix(m[Symbol("dvHeatingProduction"*_n)][t,"SpaceHeating",ts], 0.0, force=true)
+                fix(m[Symbol("dvProductionToWaste"*_n)][t,"SpaceHeating",ts], 0.0, force=true)
             end
         end
         if !(t in p.techs.can_serve_dhw)
             for ts in p.time_steps
                 fix(m[Symbol("dvHeatingProduction"*_n)][t,"DomesticHotWater",ts], 0.0, force=true)
+                fix(m[Symbol("dvProductionToWaste"*_n)][t,"DomesticHotWater",ts], 0.0, force=true)
             end
         end
         if !(t in p.techs.can_serve_process_heat)
             for ts in p.time_steps
                 fix(m[Symbol("dvHeatingProduction"*_n)][t,"ProcessHeat",ts], 0.0, force=true)
+                fix(m[Symbol("dvProductionToWaste"*_n)][t,"ProcessHeat",ts], 0.0, force=true)
             end
         end
     end
+    # Enfore
 end
 
 function no_existing_boiler_production(m, p; _n="")

src/core/reopt.jl

@@ -636,7 +636,11 @@ function add_variables!(m::JuMP.AbstractModel, p::REoptInputs)
 	end
 
     if !isempty(p.techs.steam_turbine)
-        @variable(m, dvThermalToSteamTurbine[p.techs.can_supply_steam_turbine, p.heating_loads, p.time_steps] >= 0)
+		if !isempty(p.techs.can_supply_steam_turbine)
+	        @variable(m, dvThermalToSteamTurbine[p.techs.can_supply_steam_turbine, p.heating_loads, p.time_steps] >= 0)
+		else
+			throw(@error("Steam turbine is present, but set p.techs.can_supply_steam_turbine is empty."))
+		end
     end
 
 	if !isempty(p.s.electric_utility.outage_durations) # add dvUnserved Load if there is at least one outage

test/runtests.jl

@@ -659,10 +659,10 @@ else  # run HiGHS tests
                 data_sizing = JSON.parsefile("./scenarios/chp_sizing.json")
                 s = Scenario(data_sizing)
                 inputs = REoptInputs(s)
-                m = Model(optimizer_with_attributes(HiGHS.Optimizer, "output_flag" => false, "log_to_console" => false, "mip_rel_gap" => 0.01))
+                m = Model(optimizer_with_attributes(HiGHS.Optimizer, "output_flag" => false, "log_to_console" => false, "mip_rel_gap" => 0.01, "presolve" => "on"))
                 results = run_reopt(m, inputs)
 
-                @test round(results["CHP"]["size_kw"], digits=0) ≈ 330.0 atol=20.0
+                @test round(results["CHP"]["size_kw"], digits=0) ≈ 400.0 atol=50.0
                 @test round(results["Financial"]["lcc"], digits=0) ≈ 1.3476e7 rtol=1.0e-2
             end
 
@@ -695,6 +695,7 @@ else  # run HiGHS tests
                     [0, init_capex_chp_expected * data_cost_curve["CHP"]["federal_itc_fraction"]])
 
                 #PV
+                data_cost_curve["PV"] = Dict()
                 data_cost_curve["PV"]["min_kw"] = 1500
                 data_cost_curve["PV"]["max_kw"] = 1500
                 data_cost_curve["PV"]["installed_cost_per_kw"] = 1600
@@ -823,6 +824,13 @@ else  # run HiGHS tests
                 @test results["CHP"]["annual_thermal_production_mmbtu"] ≈ 149136.6 rtol=1e-5
                 @test results["ElectricTariff"]["lifecycle_demand_cost_after_tax"] ≈ 5212.7 rtol=1e-5
             end
+
+            @testset "CHP to Waste Heat" begin
+                m = Model(optimizer_with_attributes(HiGHS.Optimizer, "output_flag" => false, "log_to_console" => false, "presolve" => "on"))
+                d = JSON.parsefile("./scenarios/chp_waste.json")
+                results = run_reopt(m, d)
+                @test sum(results["CHP"]["thermal_curtailed_series_mmbtu_per_hour"]) ≈ 4174.455 atol=1e-3
+            end
         end
 
         @testset "FlexibleHVAC" begin

test/scenarios/chp_sizing.json

@@ -58,12 +58,5 @@
                 "duration_hours": 0
             }
         ]
-    },
-    "PV": {
-        "max_kw": 0.0
-      },
-    "ElectricStorage": {
-        "max_kw": 0.0,
-        "max_kwh": 0.0
-    }    
+    }
 }

test/scenarios/chp_waste.json

@@ -0,0 +1,176 @@
+{
+    "user_uuid": "1d2ef71e-fd93-4c4a-b5c3-1485a87f772e",
+    "webtool_uuid": "73518b4b-2836-4a44-9852-a4889c33b76a",
+    "Settings": {
+        "solver_name": "HiGHS",
+        "off_grid_flag": false,
+        "include_climate_in_objective": false,
+        "include_health_in_objective": false
+    },
+    "Meta": {
+        "address": "800 Hoffman Rd Watertown WI 53094 USA"
+    },
+    "Site": {
+        "latitude": 43.1752702,
+        "longitude": -88.7344012,
+        "include_exported_renewable_electricity_in_total": true,
+        "include_exported_elec_emissions_in_total": true,
+        "land_acres": 16.33,
+        "roof_squarefeet": 0
+    },
+    "ElectricLoad": {
+        "monthly_totals_kwh": [
+            217080.0,
+            209585.0,
+            225612.0,
+            225612.0,
+            248200.0,
+            242739.0,
+            240361.0,
+            279924.0,
+            254856.0,
+            217403.0,
+            177556.0,
+            222634.0
+        ],
+        "doe_reference_name": "FlatLoad"
+    },
+    "ElectricTariff": {
+        "monthly_energy_rates": [
+            0.098,
+            0.11,
+            0.119,
+            0.107,
+            0.112,
+            0.116,
+            0.124,
+            0.116,
+            0.118,
+            0.116,
+            0.107,
+            0.101
+        ],
+        "monthly_demand_rates": [
+            0.0,
+            0.0,
+            0.0,
+            0.0,
+            0.0,
+            0.0,
+            0.0,
+            0.0,
+            0.0,
+            0.0,
+            0.0,
+            0.0
+        ]
+    },
+    "Generator": {
+        "max_kw": 0,
+        "existing_kw": 0
+    },
+    "CHP": {
+        "can_net_meter": true,
+        "prime_mover": "recip_engine",
+        "size_class": 2,
+        "fuel_type": "natural_gas",
+        "electric_efficiency_full_load": 0.34,
+        "electric_efficiency_half_load": 0.34,
+        "macrs_option_years": 5,
+        "macrs_bonus_fraction": 0.6,
+        "can_wholesale": false,
+        "fuel_cost_per_mmbtu": [
+            6.26,
+            8.36,
+            8.35,
+            4.61,
+            4.66,
+            5.22,
+            5.83,
+            5.37,
+            4.69,
+            4.65,
+            5.07,
+            6.04
+        ]
+    },
+    "DomesticHotWaterLoad": {
+        "monthly_mmbtu": [
+            135.097,
+            117.915,
+            96.579,
+            57.283,
+            15.684,
+            14.794,
+            19.994,
+            28.147,
+            37.784,
+            46.926,
+            147.828,
+            162.636
+        ],
+        "addressable_load_fraction": [
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0
+        ],
+        "doe_reference_name": "FlatLoad"
+    },
+    "SpaceHeatingLoad": {
+        "monthly_mmbtu": [
+            821.003,
+            716.585,
+            586.921,
+            348.117,
+            95.316,
+            89.906,
+            121.506,
+            171.053,
+            229.616,
+            285.174,
+            898.372,
+            988.364
+        ],
+        "addressable_load_fraction": [
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0,
+            1.0
+        ],
+        "doe_reference_name": "FlatLoad"
+    },
+    "ExistingBoiler": {
+        "fuel_type": "natural_gas",
+        "fuel_cost_per_mmbtu": [
+            6.26,
+            8.36,
+            8.35,
+            4.61,
+            4.66,
+            5.22,
+            5.83,
+            5.37,
+            4.69,
+            4.65,
+            5.07,
+            6.04
+        ]
+    }
+}