diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml index ec0d5770c..099d5fc1d 100644 --- a/.github/workflows/CI.yml +++ b/.github/workflows/CI.yml @@ -22,5 +22,5 @@ jobs: version: ${{ matrix.julia-version }} - uses: julia-actions/julia-buildpkg@latest # - uses: mxschmitt/action-tmate@v3 # for interactive debugging - - run: julia --project=. -e 'using Pkg; Pkg.activate("test"); Pkg.rm("Xpress"); Pkg.activate("."); using TestEnv; TestEnv.activate(); cd("test"); include("runtests.jl")' + - run: julia --project=. -e 'using Pkg; Pkg.activate("test"); Pkg.rm("Xpress"); Pkg.activate("."); using TestEnv; TestEnv.activate(); Pkg.add(PackageSpec(name="GhpGhx", url="https://github.com/NREL/GhpGhx.jl.git")); using GhpGhx; cd("test"); include("runtests.jl")' shell: bash diff --git a/CHANGELOG.md b/CHANGELOG.md index 22cb5f36d..146ecfb98 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -106,6 +106,7 @@ Classify the change according to the following categories: ## v0.48.2 ### Added +- Add new optional parameter **max_ton** to GHP module to allow user to size GHP smaller than peak load - Battery residual value if choosing replacement strategy for degradation - Add new **ElectricStorage** parameters **max_duration_hours** and **min_duration_hours** to bound the energy duration of battery storage ### Changed diff --git a/src/REopt.jl b/src/REopt.jl index bf8819c93..8aa8ccfca 100644 --- a/src/REopt.jl +++ b/src/REopt.jl @@ -58,11 +58,12 @@ function __init__() end const EXISTING_BOILER_EFFICIENCY = 0.8 -const GAL_PER_M3 = 264.172 # [gal/m^3] -const KWH_PER_GAL_DIESEL = 40.7 # [kWh/gal_diesel] higher heating value of diesel -const KWH_PER_MMBTU = 293.07107 # [kWh/mmbtu] +const GAL_PER_M3 = 264.172 # [gal/m^3] +const KWH_PER_GAL_DIESEL = 40.7 # [kWh/gal_diesel] higher heating value of diesel +const KWH_PER_MMBTU = 293.07107 # [kWh/mmbtu] const KWH_THERMAL_PER_TONHOUR = 3.51685 -const TONNE_PER_LB = 1/2204.62 # [tonne/lb] +const TONNE_PER_LB = 1/2204.62 # [tonne/lb] +const TONNE_PER_MMBTU_HOUR = 0.012 # [tonne/mmbtu-hour] const FUEL_TYPES = ["natural_gas", "landfill_bio_gas", "propane", "diesel_oil"] const BIG_NUMBER = 1.0e10 #used for max size. TODO use this number elsewhere. const PRIME_MOVERS = ["recip_engine", "micro_turbine", "combustion_turbine", "fuel_cell"] #TODO replace `prime_movers` references in CHP code diff --git a/src/core/ghp.jl b/src/core/ghp.jl index d7f8f58a3..f9f151c21 100644 --- a/src/core/ghp.jl +++ b/src/core/ghp.jl @@ -19,11 +19,14 @@ struct with outer constructor: installed_cost_ghx_per_ft::Float64 = 14.0 installed_cost_building_hydronic_loop_per_sqft = 1.70 om_cost_per_sqft_year::Float64 = -0.51 - building_sqft::Float64 # Required input + building_sqft::Float64 # Required input space_heating_efficiency_thermal_factor::Float64 = NaN # Default depends on building and location - cooling_efficiency_thermal_factor::Float64 = NaN # Default depends on building and location + cooling_efficiency_thermal_factor::Float64 = NaN # Default depends on building and location ghpghx_response::Dict = Dict() can_serve_dhw::Bool = false + max_ton::Real # Maximum heat pump capacity size. Default at a big number + max_number_of_boreholes::Real # Maximum GHX size + load_served_by_ghp::String # "scaled" or "nonpeak" macrs_option_years::Int = 5 macrs_bonus_fraction::Float64 = 0.6 @@ -77,6 +80,9 @@ Base.@kwdef mutable struct GHP <: AbstractGHP can_serve_space_heating::Bool = true can_serve_process_heat::Bool = false can_supply_steam_turbine::Bool = false + max_ton::Real = BIG_NUMBER + max_number_of_boreholes::Real = BIG_NUMBER + load_served_by_ghp::String = "nonpeak" aux_heater_type::String = "electric" is_ghx_hybrid::Bool = false @@ -154,7 +160,7 @@ function GHP(response::Dict, d::Dict) end # incentives = IncentivesNoProdBased(**d_mod) - setup_installed_cost_curve!(ghp, response) + setup_installed_cost_curve!(d, ghp, response) setup_om_cost!(ghp) @@ -168,10 +174,10 @@ function GHP(response::Dict, d::Dict) end """ - setup_installed_cost_curve!(response::Dict, ghp::GHP) + setup_installed_cost_curve!(d::Dict, response::Dict, ghp::GHP) """ -function setup_installed_cost_curve!(ghp::GHP, response::Dict) +function setup_installed_cost_curve!(d::Dict, ghp::GHP, response::Dict) big_number = 1.0e10 # GHX and GHP sizing metrics for cost calculations total_ghx_ft = response["outputs"]["number_of_boreholes"] * response["outputs"]["length_boreholes_ft"] diff --git a/src/core/scenario.jl b/src/core/scenario.jl index 43169fda6..77c8cf26b 100644 --- a/src/core/scenario.jl +++ b/src/core/scenario.jl @@ -480,7 +480,7 @@ function Scenario(d::Dict; flex_hvac_from_json=false) space_heating_thermal_load_reduction_with_ghp_kw = zeros(8760 * settings.time_steps_per_hour) cooling_thermal_load_reduction_with_ghp_kw = zeros(8760 * settings.time_steps_per_hour) eval_ghp = false - get_ghpghx_from_input = false + get_ghpghx_from_input = false if haskey(d, "GHP") && haskey(d["GHP"],"building_sqft") eval_ghp = true if haskey(d["GHP"], "ghpghx_responses") && !isempty(d["GHP"]["ghpghx_responses"]) @@ -620,13 +620,101 @@ function Scenario(d::Dict; flex_hvac_from_json=false) end ghpghx_results = Dict() + try # Call GhpGhx.jl to size GHP and GHX @info "Starting GhpGhx.jl" # Call GhpGhx.jl to size GHP and GHX + # If user provides undersized GHP, calculate load to send to GhpGhx.jl, and load to send to REopt for backup + thermal_load_ton = ghpghx_inputs["heating_thermal_load_mmbtu_per_hr"] .* 1/TONNE_PER_MMBTU_HOUR + if haskey(ghpghx_inputs,"cooling_thermal_load_ton") + cooling_load_ton = ghpghx_inputs["cooling_thermal_load_ton"] + thermal_load_ton .+= cooling_load_ton + end + peak_thermal_load_ton = maximum(thermal_load_ton) + if haskey(d["GHP"],"max_ton") && peak_thermal_load_ton > d["GHP"]["max_ton"] + @info "User entered undersized GHP. Calculating load that can be served by user specified undersized GHP" + # When user specifies undersized GHP, calculate the load to be served by GHP and send the rest to REopt + if !haskey(d["GHP"], "load_served_by_ghp") + d["GHP"]["load_served_by_ghp"] = "nonpeak" + end + # If user choose to scale down total load (load_served_by_ghp="scaled"), calculate the ratio of the udersized GHP size and peak load + if d["GHP"]["load_served_by_ghp"] == "scaled" + @info "GHP served scaled down of total thermal load" + peak_ratio = d["GHP"]["max_ton"] / peak_thermal_load_ton + # Scale the total load profile down by the peak_ratio and use this scaled down load to rerun GhpGhx.jl + heating_load_mmbtu = ghpghx_inputs["heating_thermal_load_mmbtu_per_hr"] + heating_load_mmbtu = heating_load_mmbtu .* peak_ratio + ghpghx_inputs["heating_thermal_load_mmbtu_per_hr"] = heating_load_mmbtu + if haskey(ghpghx_inputs,"cooling_thermal_load_ton") + ghpghx_inputs["cooling_thermal_load_ton"] = cooling_load_ton .* peak_ratio + end + elseif d["GHP"]["load_served_by_ghp"] == "nonpeak" + @info "GHP serves all thermal load below peak" + heating_load_mmbtu = ghpghx_inputs["heating_thermal_load_mmbtu_per_hr"] + # if cooling load is included, cut down total thermal load and send as much heating load to GhpGhx.jl as possible + if haskey(ghpghx_inputs,"cooling_thermal_load_ton") + # If total thermal load (heating + cooling) is more than user-defined GHP size, + # first reduce heating load as much as possible while keeping cooling load the same + if peak_thermal_load_ton > d["GHP"]["max_ton"] + thermal_load_ton[thermal_load_ton .>= d["GHP"]["max_ton"]] .= d["GHP"]["max_ton"] + heating_load_ton = thermal_load_ton .- cooling_load_ton + # Make sure that the reduced heating load is not negative + heating_load_ton[heating_load_ton .< 0] .= 0 + # If the updated peak thermal load is still more than user-defined GHP size, + # reduce cooling load as well + updated_thermal_load_ton = heating_load_ton .+ cooling_load_ton + updated_peak_thermal_load_ton = maximum(updated_thermal_load_ton) + if updated_peak_thermal_load_ton > d["GHP"]["max_ton"] + updated_thermal_load_ton[updated_thermal_load_ton .>= d["GHP"]["max_ton"]] .= d["GHP"]["max_ton"] + cooling_load_ton = updated_thermal_load_ton .- heating_load_ton + ghpghx_inputs["cooling_thermal_load_ton"] = cooling_load_ton + end + heating_load_mmbtu = heating_load_ton .* TONNE_PER_MMBTU_HOUR + ghpghx_inputs["heating_thermal_load_mmbtu_per_hr"] = heating_load_mmbtu + end + # if cooling load is not included, cut down heating load only and send to GhpGhx.jl + else + heating_load_mmbtu[heating_load_mmbtu .>= d["GHP"]["max_ton"] * TONNE_PER_MMBTU_HOUR] .= d["GHP"]["max_ton"] * TONNE_PER_MMBTU_HOUR + ghpghx_inputs["heating_thermal_load_mmbtu_per_hr"] = heating_load_mmbtu + end + end + end results, inputs_params = GhpGhx.ghp_model(ghpghx_inputs) + # If max_number_of_boreholes is specified, check if number of boreholes sized by GhpGhx.jl greater than user-specified max_number_of_boreholes, + # and if max_number_of_boreholes is less, reduce thermal load served by GHP until max_number_of_boreholes = number of boreholses sized by GhpGhx.jl + if haskey(d["GHP"],"max_number_of_boreholes") + optimal_number_of_boreholes = GhpGhx.get_results_for_reopt(results, inputs_params)["number_of_boreholes"] + if optimal_number_of_boreholes > d["GHP"]["max_number_of_boreholes"] + @info "Max number of boreholes specified is less than number of boreholes sized in GhpGhx.jl, reducing thermal load served by GHP further" + max_iter = 10 + for iter = 1:max_iter + borehole_ratio = d["GHP"]["max_number_of_boreholes"] / optimal_number_of_boreholes + heating_load_mmbtu .*= borehole_ratio + if haskey(ghpghx_inputs,"cooling_thermal_load_ton") + cooling_load_ton .*= borehole_ratio + # if cooling load is not included, cut down heating load only and send to GhpGhx.jl + end + ghpghx_inputs["heating_thermal_load_mmbtu_per_hr"] = heating_load_mmbtu + ghpghx_inputs["cooling_thermal_load_ton"] = cooling_load_ton + + # Rerun GhpGhx.jl + results, inputs_params = GhpGhx.ghp_model(ghpghx_inputs) + optimal_number_of_boreholes = GhpGhx.get_results_for_reopt(results, inputs_params)["number_of_boreholes"] + # Solution is found if the new optimal number of boreholes sized by GhpGhx.jl = user-specified max number of boreholes, + # Otherwise, continue solving until reaching max iteration + if -0.5 < optimal_number_of_boreholes - d["GHP"]["max_number_of_boreholes"] < 0.5 + break + else + iter += 1 + end + end + end + end + # Create a dictionary of the results data needed for REopt ghpghx_results = GhpGhx.get_results_for_reopt(results, inputs_params) + # Return results from GhpGhx.jl without load scaling if user does not provide GHP size or if user entered GHP size is greater than GHP size output @info "GhpGhx.jl model solved" #with status $(results["status"])." catch e @info e @@ -644,9 +732,7 @@ function Scenario(d::Dict; flex_hvac_from_json=false) end append!(ghp_option_list, [GHP(ghpghx_response, ghp_inputs_removed_ghpghx_params)]) # Print out ghpghx_response for loading into a future run without running GhpGhx.jl again - #open("scenarios/ghpghx_response.json","w") do f - # JSON.print(f, ghpghx_response) - #end + # open("scenarios/ghpghx_response.json","w") do f end # If ghpghx_responses is included in inputs, do NOT run GhpGhx.jl model and use already-run ghpghx result as input to REopt elseif eval_ghp && get_ghpghx_from_input diff --git a/src/results/existing_boiler.jl b/src/results/existing_boiler.jl index 9d44f552c..82833b3ee 100644 --- a/src/results/existing_boiler.jl +++ b/src/results/existing_boiler.jl @@ -5,7 +5,7 @@ - `fuel_consumption_series_mmbtu_per_hour` # Fuel consumption series [MMBtu/hr] - `annual_fuel_consumption_mmbtu` # Fuel consumed in a year [MMBtu] - `thermal_production_series_mmbtu_per_hour` # Thermal energy production series [MMBtu/hr] -- `annual_thermal_production_mmbtu` # Thermal power production to TES (HotThermalStorage) series [MMBtu/hr] +- `annual_thermal_production_mmbtu` # Thermal power production in a year [MMBtu] - `thermal_to_storage_series_mmbtu_per_hour` # Thermal power production to TES (HotThermalStorage) series [MMBtu/hr] - `thermal_to_steamturbine_series_mmbtu_per_hour` # Thermal power production to SteamTurbine series [MMBtu/hr] - `thermal_to_load_series_mmbtu_per_hour` # Thermal power production to serve the heating load series [MMBtu/hr] diff --git a/src/results/ghp.jl b/src/results/ghp.jl index cab60c115..212c6252c 100644 --- a/src/results/ghp.jl +++ b/src/results/ghp.jl @@ -15,6 +15,9 @@ GHP results: - `thermal_to_space_heating_load_series_mmbtu_per_hour` - `thermal_to_dhw_load_series_mmbtu_per_hour` - `thermal_to_load_series_ton` +- `annual_thermal_production_mmbtu` # GHP's heating thermal power production in a year [MMBtu] +- `annual_thermal_production_tonhour` # GHP's cooling thermal power production in a year [ton] + """ function add_ghp_results(m::JuMP.AbstractModel, p::REoptInputs, d::Dict; _n="") @@ -25,6 +28,7 @@ function add_ghp_results(m::JuMP.AbstractModel, p::REoptInputs, d::Dict; _n="") # r["size_heat_pump_ton"] = 0.0 # r["size_wwhp_heating_pump_ton"] = 0.0 # r["size_wwhp_cooling_pump_ton"] = 0.0 + if ghp_option_chosen >= 1 r["ghpghx_chosen_outputs"] = p.s.ghp_option_list[ghp_option_chosen].ghpghx_response["outputs"] @@ -42,13 +46,22 @@ function add_ghp_results(m::JuMP.AbstractModel, p::REoptInputs, d::Dict; _n="") r["space_heating_thermal_load_reduction_with_ghp_mmbtu_per_hour"] = round.(value.(HeatingThermalReductionWithGHP) ./ KWH_PER_MMBTU, digits=3) @expression(m, CoolingThermalReductionWithGHP[ts in p.time_steps], sum(p.cooling_thermal_load_reduction_with_ghp_kw[g,ts] * m[Symbol("binGHP"*_n)][g] for g in p.ghp_options)) + + @expression(m, HeatingThermalLoadServedWithGHP[ts in p.time_steps], + sum(p.ghp_heating_thermal_load_served_kw[g,ts] * m[Symbol("binGHP"*_n)][g] for g in p.ghp_options)) + @expression(m, CoolingThermalLoadServedWithGHP[ts in p.time_steps], + sum(p.ghp_cooling_thermal_load_served_kw[g,ts] * m[Symbol("binGHP"*_n)][g] for g in p.ghp_options)) + r["cooling_thermal_load_reduction_with_ghp_ton"] = round.(value.(CoolingThermalReductionWithGHP) ./ KWH_THERMAL_PER_TONHOUR, digits=3) r["ghx_residual_value_present_value"] = value(m[:ResidualGHXCapCost]) r["avoided_capex_by_ghp_present_value"] = value(m[:AvoidedCapexByGHP]) - r["thermal_to_space_heating_load_series_mmbtu_per_hour"] = d["HeatingLoad"]["space_heating_thermal_load_series_mmbtu_per_hour"] .- r["space_heating_thermal_load_reduction_with_ghp_mmbtu_per_hour"] - r["thermal_to_load_series_ton"] = d["CoolingLoad"]["load_series_ton"] .- r["cooling_thermal_load_reduction_with_ghp_ton"] + r["thermal_to_space_heating_load_series_mmbtu_per_hour"] = round.(value.(HeatingThermalLoadServedWithGHP) ./ KWH_PER_MMBTU, digits=3) + r["thermal_to_load_series_ton"] = round.(value.(CoolingThermalLoadServedWithGHP) ./ KWH_THERMAL_PER_TONHOUR, digits=3) + r["annual_thermal_production_mmbtu"] = sum(r["thermal_to_space_heating_load_series_mmbtu_per_hour"]) + r["annual_thermal_production_tonhour"] = sum(r["thermal_to_load_series_ton"]) if p.s.ghp_option_list[ghp_option_chosen].can_serve_dhw r["thermal_to_dhw_load_series_mmbtu_per_hour"] = d["HeatingLoad"]["dhw_thermal_load_series_mmbtu_per_hour"] + r["annual_thermal_production_mmbtu"] = r["annual_thermal_production_mmbtu"] + sum(r["thermal_to_dhw_load_series_mmbtu_per_hour"]) else r["thermal_to_dhw_load_series_mmbtu_per_hour"] = zeros(length(p.time_steps)) end @@ -60,6 +73,8 @@ function add_ghp_results(m::JuMP.AbstractModel, p::REoptInputs, d::Dict; _n="") r["thermal_to_space_heating_load_series_mmbtu_per_hour"] = zeros(length(p.time_steps)) r["thermal_to_load_series_ton"] = zeros(length(p.time_steps)) r["thermal_to_dhw_load_series_mmbtu_per_hour"] = zeros(length(p.time_steps)) + r["annual_thermal_production_mmbtu"] = 0.0 + r["annual_thermal_production_tonhour"] = 0.0 end d["GHP"] = r nothing diff --git a/test/runtests.jl b/test/runtests.jl index b4d03c178..be573e9d1 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2189,6 +2189,8 @@ else # run HiGHS tests 6. Hybrid GHP capability functions as expected 7. Check GHP LCC calculation for URBANopt 8. Check GHX LCC calculation for URBANopt + 9. Allow User-defined max GHP size + 10. Allow User-defined max GHP size and max number of boreholes """ # Load base inputs @@ -2310,8 +2312,39 @@ else # run HiGHS tests @test ghx_lccc_initial - boreholes*boreholes_len*14 ≈ 0.0 atol = 0.01 # GHP size must be 0 @test ghp_size ≈ 0.0 atol = 0.01 - # LCCC should be around be around 52% of initial capital cost due to incentive and bonus + # LCCC should be around 52% of initial capital cost due to incentive and bonus @test ghx_lccc/ghx_lccc_initial ≈ 0.518 atol = 0.01 + + # User specified GHP size + input_presizedGHP = deepcopy(input_data) + input_presizedGHP["GHP"]["max_ton"] = 300 + input_presizedGHP["GHP"]["heatpump_capacity_sizing_factor_on_peak_load"] = 1.0 + delete!(input_presizedGHP["GHP"], "ghpghx_responses") + # Rerun REopts + s_presizedGHP = Scenario(input_presizedGHP) + inputs_presizedGHP = REoptInputs(s_presizedGHP) + m1 = Model(optimizer_with_attributes(HiGHS.Optimizer, "output_flag" => false, "log_to_console" => false, "mip_rel_gap" => 0.01)) + m2 = Model(optimizer_with_attributes(HiGHS.Optimizer, "output_flag" => false, "log_to_console" => false, "mip_rel_gap" => 0.01)) + results = run_reopt([m1,m2], inputs_presizedGHP) + # GHP output size should equal user-defined GHP size + output_GHP_size = sum(results["GHP"]["size_heat_pump_ton"]) + @test output_GHP_size ≈ 300.00 atol=0.1 + + # User specified max GHP and GHX sizes + input_presizedGHPGHX = deepcopy(input_presizedGHP) + input_presizedGHPGHX["GHP"]["max_number_of_boreholes"] = 400 + # Rerun REopts + s_presizedGHPGHX = Scenario(input_presizedGHPGHX) + inputs_presizedGHPGHX = REoptInputs(s_presizedGHPGHX) + m1 = Model(optimizer_with_attributes(HiGHS.Optimizer, "output_flag" => false, "log_to_console" => false, "mip_rel_gap" => 0.01)) + m2 = Model(optimizer_with_attributes(HiGHS.Optimizer, "output_flag" => false, "log_to_console" => false, "mip_rel_gap" => 0.01)) + results = run_reopt([m1,m2], inputs_presizedGHPGHX) + # GHP output size should equal user-defined GHP size + output_GHP_size = results["GHP"]["size_heat_pump_ton"] + output_GHX_size = results["GHP"]["ghpghx_chosen_outputs"]["number_of_boreholes"] + @test output_GHX_size ≈ 400.00 atol=0.5 + @test output_GHP_size < 300.00 + finalize(backend(m)) empty!(m) GC.gc() @@ -3669,7 +3702,7 @@ else # run HiGHS tests @test s.space_heating_load.loads_kw[end-24+1:end] == s.space_heating_load.loads_kw[1:24] @test s.cooling_load.loads_kw_thermal[end-24+1:end] == s.cooling_load.loads_kw_thermal[1:24] end - + @testset "After-tax savings and capital cost metric for alternative payback calculation" begin """ Check alignment between REopt simple_payback_years and a simple X/Y payback metric with