diff --git a/CHANGELOG.md b/CHANGELOG.md index a2e0d1944..8a9c2292c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -26,6 +26,11 @@ Classify the change according to the following categories: ##### Removed ### Patches +## heuristic-dispatch-option +### Minor udpates +#### Added +- Add **ElectricStorage** inputs field **dispatch_options** with heuristic options + ## fixed-soc ### Minor udpates #### Added diff --git a/julia_src/Manifest.toml b/julia_src/Manifest.toml index dc825f89c..ceb0887ae 100644 --- a/julia_src/Manifest.toml +++ b/julia_src/Manifest.toml @@ -948,9 +948,11 @@ version = "1.11.0" [[deps.REopt]] deps = ["ArchGDAL", "CSV", "CoolProp", "DataFrames", "Dates", "DelimitedFiles", "HTTP", "JLD", "JSON", "JuMP", "LinDistFlow", "LinearAlgebra", "Logging", "MathOptInterface", "Requires", "Roots", "Statistics", "TestEnv"] -git-tree-sha1 = "e807b62ab249fca7b59ac6ab7bb71725f504cc5f" +git-tree-sha1 = "ec4dad03ea876407e69a384cecc62c29be858405" +repo-rev = "dispatch-options" +repo-url = "https://github.com/NatLabRockies/REopt.jl.git" uuid = "d36ad4e8-d74a-4f7a-ace1-eaea049febf6" -version = "0.59.2" +version = "0.59.3" [[deps.Random]] deps = ["SHA"] diff --git a/julia_src/http.jl b/julia_src/http.jl index 9eab6396d..c89172a49 100644 --- a/julia_src/http.jl +++ b/julia_src/http.jl @@ -8,6 +8,7 @@ DotEnv.load!() const test_nrel_developer_api_key = ENV["NREL_DEVELOPER_API_KEY"] ENV["NREL_DEVELOPER_EMAIL"] = "reopt@nlr.gov" +include("mpc.jl") include("os_solvers.jl") @@ -64,6 +65,7 @@ function reopt(req::HTTP.Request) ENV["NREL_DEVELOPER_API_KEY"] = test_nrel_developer_api_key delete!(d, "api_key") end + settings = d["Settings"] solver_name = get(settings, "solver_name", "HiGHS") if solver_name == "Xpress" && !(xpress_installed=="True") @@ -71,6 +73,34 @@ function reopt(req::HTTP.Request) @warn "Changing solver_name from Xpress to $solver_name because Xpress is not installed. Next time Specify Settings.solver_name = 'HiGHS' or 'Cbc' or 'SCIP'" end + + # ---- API-only battery heuristic dispatch strategy: "daily_foresight_optimized" ---- + # When ElectricStorage.dispatch_strategy == "daily_foresight_optimized", if needed, first run REopt to get optimal sizing of PV and battery. + # Then run the MPC rolling-horizon loop to get a SOC profile (skip MPC dispatch if optimal battery size is 0). + # Then set ElectricStorage.fixed_soc_series_fraction = MPC SOC and fix PV and BESS sizing before running the main REopt optimization. + electric_storage = get(d, "ElectricStorage", Dict()) + if get(electric_storage, "dispatch_strategy", nothing) == "daily_foresight_optimized" + try + @info "Running MPC to obtain daily foresight optimized battery dispatch profile." + mpc_results = get_mpc_results!(d; solver_name=solver_name) + # TODO: Cache sizing run results and avoid a second call to REopt? Are those the same results? + if get(mpc_results, "skip_mpc", false) == true + @info "Cannot execute daily_foresight_optimized battery dispatch because optimal battery size is 0. Setting dispatch strategy to 'optimized'." + d["ElectricStorage"]["dispatch_strategy"] = "optimized" + else + soc = mpc_results["ElectricStorage"]["soc_series_fraction"] + d["ElectricStorage"]["fixed_soc_series_fraction"] = soc + d["ElectricStorage"]["dispatch_strategy"] = "custom_soc" + end + catch e + @error "MPC pre-solve failed" exception=(e, catch_backtrace()) + return HTTP.Response(500, JSON.json(Dict( + "error" => "MPC pre-solve failed: " * sprint(showerror, e), + "reopt_version" => string(pkgversion(reoptjl)), + ))) + end + end + timeout_seconds = pop!(settings, "timeout_seconds") optimality_tolerance = pop!(settings, "optimality_tolerance") solver_attributes = SolverAttributes(timeout_seconds, optimality_tolerance) @@ -810,6 +840,63 @@ function job_no_xpress(req::HTTP.Request) return HTTP.Response(500, JSON.json(error_response)) end +""" + mpc(req::HTTP.Request) + +HTTP endpoint for rolling-horizon Model Predictive Control (MPC) dispatch optimization. + +This endpoint performs a full-year rolling-horizon MPC dispatch for PV + ElectricStorage systems, +optimizing daily dispatch using a 24-hour look-ahead window. Runs the MPC dispatch loop via +`get_mpc_results!` and returns dispatch results and cost metrics as JSON. + +Arguments: + req::HTTP.Request: REopt inputs dictionary + +Returns JSON dictionary containing: + - MPC: Metadata (time_steps_per_hour, horizon_time_steps) + - PV: Size and dispatch series (to load, storage, grid, curtailed) + - ElectricStorage: Sizes and state-of-charge series + - ElectricUtility: Grid dispatch series and emissions + - ElectricLoad: Load profile used + - ElectricTariff: Energy and demand costs, peak demands by month/ratchet + - status: "optimal" + - reopt_version: Version of REopt.jl used + +""" +function mpc(req::HTTP.Request) + d = JSON.parse(String(req.body)) + error_response = Dict() + results = Dict() + try + if !isempty(get(d, "api_key", "")) + ENV["NREL_DEVELOPER_API_KEY"] = pop!(d, "api_key") + else + ENV["NREL_DEVELOPER_API_KEY"] = test_nrel_developer_api_key + delete!(d, "api_key") + end + + solver_name = get(get(d, "Settings", Dict()), "solver_name", "HiGHS") + if solver_name == "Xpress" && !(xpress_installed=="True") + solver_name = "HiGHS" + @warn "Changing solver_name from Xpress to $solver_name because Xpress is not installed. Next time + Specify Settings.solver_name = 'HiGHS' or 'Cbc' or 'SCIP'" + end + + results = get_mpc_results!(d; solver_name=solver_name) + catch e + @error "MPC failed" exception=(e, catch_backtrace()) + error_response["error"] = sprint(showerror, e) + error_response["reopt_version"] = string(pkgversion(reoptjl)) + end + GC.gc() + if isempty(error_response) + @info "MPC ran successfully." + return HTTP.Response(200, JSON.json(results)) + else + return HTTP.Response(500, JSON.json(error_response)) + end +end + # define REST endpoints to dispatch to "service" functions const ROUTER = HTTP.Router() @@ -820,6 +907,7 @@ else end HTTP.register!(ROUTER, "POST", "/reopt", reopt) HTTP.register!(ROUTER, "POST", "/erp", erp) +HTTP.register!(ROUTER, "POST", "/mpc", mpc) HTTP.register!(ROUTER, "POST", "/ghpghx", ghpghx) HTTP.register!(ROUTER, "GET", "/chp_defaults", chp_defaults) HTTP.register!(ROUTER, "GET", "/avert_emissions_profile", avert_emissions_profile) diff --git a/julia_src/mpc.jl b/julia_src/mpc.jl new file mode 100644 index 000000000..e9f8b5c1d --- /dev/null +++ b/julia_src/mpc.jl @@ -0,0 +1,604 @@ +# REopt®, Copyright (c) Alliance for Sustainable Energy, LLC. See also https://github.com/NatLabRockies/REopt_API/blob/master/LICENSE. +# ============================================================================ +# MPC (Model Predictive Control) endpoint +# ---------------------------------------------------------------------------- +# Rolling-horizon dispatch with one day look-ahead using REopt.run_mpc. Assumes: +# * PV and ElectricStorage are the only available technologies +# * Only perfect forecast scenarios are modeled (no forecast errors) +# * Sizing: Perform a sizing run first if PV or ElectricStorage sizes are not provided (storage sizes must be > 0) +# * PV: Use PVWatts if `PV.production_factor_series` is not provided +# * ElectricLoad: Use commercial reference profiles if `ElectricLoad.loads_kw` is not provided +# * Code wraps with Jan 1 data to determine Dec 31 dispatch (leap-year inputs are not supported) +# * MPC settings are currently hard coded (e.g., forecast horizon, control horizon, optimization horizon) +# ============================================================================ + +""" + get_month_transition_timesteps(time_steps_per_hour) + +Return an array of length 12 specifying the index marking the start of each month in a non-leap year +""" +function get_month_transition_timesteps(time_steps_per_hour::Int) + # hours in each month for a non-leap year + hours_per_month = [744, 672, 744, 720, 744, 720, 744, 744, 720, 744, 720, 744] + starts = Vector{Int}(undef, 12) + starts[1] = 1 + for m in 2:12 + starts[m] = starts[m-1] + hours_per_month[m-1] * time_steps_per_hour + end + return starts +end + +""" + slice_data(arr, idx, end_idx) + +Returns arr[idx:end_idx] but wrap around to the start of arr if end_idx > length(arr). +""" +function slice_data(arr::AbstractVector, idx::Int, end_idx::Int) + n = length(arr) + if end_idx <= n + return arr[idx:end_idx] + else + wrap_len = end_idx - n + return vcat(arr[idx:n], arr[1:wrap_len]) + end +end + + +""" + generate_pv_production_factors(d, time_steps_per_hour) + +Generate a PV production factor series using PVWatts by calling REopt.get_production_factor. +This is called by get_mpc_results! only when the user does not provide a custom production_factor_series and REopt is not called for sizing. +""" +function generate_pv_production_factors(d::Dict, time_steps_per_hour::Int) + site = get(d, "Site", Dict()) + lat = Float64(site["latitude"]) + lon = Float64(site["longitude"]) + + @info "MPC: PV.production_factor_series not provided; generating using PVWatts through REopt.jl (lat=$(lat), lon=$(lon))." + + pv = get(d, "PV", Dict()) + pv_pf_kwargs = (:array_type, :tilt, :module_type, :losses, :azimuth, :gcr, + :radius, :name, :location, :dc_ac_ratio, :inv_eff) + kwargs = Dict{Symbol,Any}(Symbol(k) => v for (k, v) in pv if Symbol(k) in pv_pf_kwargs) + pv_struct = reoptjl.PV(; latitude = lat, kwargs...) + pv_production_factor_series = reoptjl.get_production_factor(pv_struct, lat, lon; + time_steps_per_hour = time_steps_per_hour) + return Vector{Float64}(pv_production_factor_series) +end + +""" + build_mpc_response(status; skip_mpc=false, messages=Dict(), result_dict=Dict()) + +Build a consistent MPC response envelope with status, version info, and messages. +Ensures all response paths (success, error, skip) have uniform structure. +""" +function build_mpc_response(status::String; skip_mpc=false, messages=Dict(), result_dict=Dict()) + response = Dict( + "status" => status, + "reopt_version" => string(pkgversion(reoptjl)), + "Messages" => messages, + "skip_mpc" => skip_mpc, + ) + # Merge result data if provided (for success case) + return merge(response, result_dict) +end + +""" + get_technology_sizes!(d::Dict, model_inputs::reoptjl.REoptInputs, solver_settings::Dict) + +Determine PV and ElectricStorage sizes for the MPC loop. If min_kw != max_kw and/or min_kwh != max_kwh, +call REopt to size technologies. User input battery sizes must also be greater than zero. + d is mutated in-place to update PV and ElectricStorage sizes. + model_inputs has already been processed to remove inputs not used in REopt.jl +""" +function get_technology_sizes!(d::Dict, model_inputs::reoptjl.REoptInputs, solver_settings::Dict) + # pv and batt are updated in-place from the dictionary `d` and used in the final REopt run + pv = get!(d, "PV", Dict()) + batt = get!(d, "ElectricStorage", Dict()) + + function is_fixed(dct, lo_key, hi_key) + lo = get(dct, lo_key, nothing) + hi = get(dct, hi_key, nothing) + return lo !== nothing && hi !== nothing && Float64(lo) == Float64(hi) + end + + # Check storage sizes are greater than zero + if Float64(get(batt, "max_kw", 1.0)) <= 0.0 || Float64(get(batt, "max_kwh", 1.0)) <= 0.0 + error("ElectricStorage max_kw and max_kwh must both be greater than zero " * + "to run the daily_foresight_optimized dispatch option.") + end + + # Check if both PV and BESS sizes fixed + pv_fixed = is_fixed(pv, "min_kw", "max_kw") + batt_fixed = is_fixed(batt, "min_kw", "max_kw") && + is_fixed(batt, "min_kwh", "max_kwh") + + if pv_fixed && batt_fixed + pv_kw = Float64(pv["min_kw"]) + batt_kw = Float64(batt["min_kw"]) + batt_kwh = Float64(batt["min_kwh"]) + return (pv_kw = pv_kw, batt_kw = batt_kw, batt_kwh = batt_kwh, skip_mpc = false, pv_production_factor_series = nothing) + end + + @info "MPC: PV and/or ElectricStorage sizes are not specified — running REopt sizing first." + # TODO: Should we "remove tiers" here for sizing or allow for optimizing with tiers? + + m = get_solver_model(get_solver_model_type(solver_settings["solver_name"]), solver_settings["solver_attributes"]) + + # model_inputs = reoptjl.REoptInputs(sizing_post) + sizing_results = reoptjl.run_reopt(m, model_inputs) + + if get(sizing_results, "status", "") != "optimal" + status = get(sizing_results, "status", "unknown") + msgs = get(sizing_results, "Messages", Dict()) + errs = get(msgs, "errors", []) + warns = get(msgs, "warnings", []) + error("MPC sizing pre-step did not solve (status = $(status)). " * + "REopt errors: $(errs). REopt warnings: $(warns).") + end + + pv_kw = Float64(get(get(sizing_results, "PV", Dict()), "size_kw", 0.0)) + batt_kw = Float64(get(get(sizing_results, "ElectricStorage", Dict()), "size_kw", 0.0)) + batt_kwh = Float64(get(get(sizing_results, "ElectricStorage", Dict()), "size_kwh", 0.0)) + pv_production_factor_series = get(get(sizing_results, "PV", Dict()), "production_factor_series", nothing) + + # Skip the MPC loop if no battery is sized + if batt_kw <= 0.0 || batt_kwh <= 0.0 + return (pv_kw = pv_kw, batt_kw = 0.0, batt_kwh = 0.0, skip_mpc = true, pv_production_factor_series = pv_production_factor_series) + end + + # Fix inputs for final REopt run in http.jl + pv["min_kw"] = pv_kw + pv["max_kw"] = pv_kw + batt["min_kw"] = batt_kw + batt["max_kw"] = batt_kw + batt["min_kwh"] = batt_kwh + batt["max_kwh"] = batt_kwh + + @info "MPC: REopt sizing solved with PV = $(pv_kw) kW and battery = $(batt_kw) kW / $(batt_kwh) kWh." + return (; pv_kw, batt_kw, batt_kwh, skip_mpc = false, pv_production_factor_series) +end + +function get_mpc_results!(d::Dict; solver_name::String="HiGHS")::Dict + """ + Run a full-year rolling-horizon MPC dispatch for PV + ElectricStorage by + calling `REopt.run_mpc` once per timestep with a 24-hour look-ahead. + + Inputs: + d::Dict, REopt inputs dictionary (will be modified in-place to store computed fixed sizes and PV production factors.) + + Returns JSON dictionary containing: + - status: "optimal", "error", or "skipped" + - reopt_version: Version of REopt.jl used + - Messages: Dict with optional errors, warnings, or info + - skip_mpc: Boolean indicating if MPC was skipped + + For "optimal" status, also includes: + - MPC: Metadata (time_steps_per_hour, horizon_time_steps) + - PV: Size and dispatch series (to load, storage, grid, curtailed) + - ElectricStorage: Sizes and state-of-charge series + - ElectricUtility: Grid dispatch series and emissions + - ElectricLoad: Load profile used + - ElectricTariff: Separate cost components (total_energy_cost, total_export_benefit, + total_tou_demand_cost, total_non_tou_monthly_demand_cost), a combined total_electricity_bill, + per-timestep energy/export series, and peak demands by month/ratchet + + """ + + ## Validation on allowable inputs for MPC ## + # Error if any techs other than PV and ElectricStorage are provided + mpc_allowed_keys = Set(["PV", "ElectricStorage", "ElectricLoad", "ElectricTariff", "ElectricUtility", "Site", "Settings", "Financial"]) + unsupported_keys = setdiff(keys(d), mpc_allowed_keys) + if !isempty(unsupported_keys) + error("When using MPC (daily_foresight_optimized dispatch), only PV and ElectricStorage are supported technologies. " * + "Unsupported inputs found: $(join(unsupported_keys, ", ")).") + end + # PV can be provided as a Dict (single system) or an array of Dicts. + # MPC only supports a single PV, so error on multiple PVs and normalize a + # one-element array down to a Dict so downstream code can treat d["PV"] as a Dict. + # TODO: Handle multiple PVs + if haskey(d, "PV") && isa(d["PV"], AbstractArray) + if length(d["PV"]) > 1 + error("MPC: Multiple PV systems are not supported.") + elseif length(d["PV"]) == 1 + d["PV"] = d["PV"][1] + else + delete!(d, "PV") # empty PV array -> treat as no PV + end + end + + # Error if unsupported CO2/renewable-fraction constraints are set + _site_input = get(d, "Site", Dict()) + if !isnothing(get(_site_input, "CO2_emissions_reduction_min_fraction", nothing)) + error("MPC: Site.CO2_emissions_reduction_min_fraction is not supported in MPC runs.") + end + if get(_site_input, "include_grid_renewable_fraction_in_RE_constraints", false) == true + error("MPC: Site.include_grid_renewable_fraction_in_RE_constraints is not supported in MPC runs.") + end + if get(_site_input, "include_exported_elec_emissions_in_total", true) == false + error("MPC: Site.include_exported_elec_emissions_in_total = false is not supported in MPC runs.") + end + if get(_site_input, "include_exported_renewable_electricity_in_total", true) == false + error("MPC: Site.include_exported_renewable_electricity_in_total = false is not supported in MPC runs.") + end + + # TODO: Add warnings for REopt inputs and scenarios that are not modeled in MPC (e.g., coincident peak charges, demand lookback, etc.) + @warn "Using MPC to determine dispatch. MPC does not model: tiered electricity rates; rates will be flattened to the first tier." + + # TODO: Error if rate tariff contains lookbacks or coincident peak charges. + + # TODO: Test with outage inputs before enabling this warning. + # # Warning for outage inputs (MPC does not model outages) + # _utility_input = get(d, "ElectricUtility", Dict()) + # if any(k -> haskey(_utility_input, k), ("outage_start_time_step", "outage_start_time_steps", "outage_durations")) + # @warn "MPC: Outage inputs detected (outage_start_time_step, outage_start_time_steps, outage_durations). " * + # "MPC does not model outages; these inputs will be ignored." + # end + + ## Set up MPC inputs ## + + # TODO: MPC horizons and timeout are currently hard coded + settings = get!(d, "Settings", Dict()) + time_steps_per_hour = Int(get(settings, "time_steps_per_hour", 1)) + length_of_data = 8760 * time_steps_per_hour + horizon = 24 * time_steps_per_hour + per_iter_timeout_s = 30.0 + + # Update sizing_post to remove solver settings and dispatch inputs, to be able to validate inputs using REoptInputs + sizing_post = deepcopy(d) + settings = get(sizing_post, "Settings", Dict()) + solver_settings=Dict() + delete!(settings, "run_bau") # Remove run_bau from sizing run + solver_settings["timeout_seconds"] = pop!(settings, "timeout_seconds", 600) # only gets used in sizing run. + solver_settings["optimality_tolerance"] = pop!(settings, "optimality_tolerance", 0.001) # Update to a higher value if solve time becomes an issue + solver_settings["solver_attributes"] = SolverAttributes(solver_settings["timeout_seconds"], solver_settings["optimality_tolerance"]) + solver_settings["solver_name"] = solver_name + # Delete inputs specific to the heuristic battery dispatch run + if haskey(sizing_post, "ElectricStorage") + delete!(sizing_post["ElectricStorage"], "dispatch_strategy") + delete!(sizing_post["ElectricStorage"], "fixed_soc_series_fraction") + end + + # Process and validate inputs using REoptInputs + model_inputs = nothing + try + model_inputs = reoptjl.REoptInputs(sizing_post) + + # REoptInputs returns an error Dict (rather than throwing) when input validation fails. + # Surface those messages instead of falling through to get_technology_sizes!, which expects + # a REoptInputs and would otherwise raise a confusing MethodError. + if isa(model_inputs, Dict) + @error "REopt input validation failed during MPC pre-solve." messages=get(model_inputs, "Messages", Dict()) + return build_mpc_response( + "error", + messages = get(model_inputs, "Messages", Dict("errors" => ["REopt input validation failed."])) + ) + else + @info "Successfully processed REopt inputs." + end + catch e + @error "Something went wrong during REopt inputs processing!" exception=(e, catch_backtrace()) + return build_mpc_response( + "error", + messages = Dict("errors" => [sprint(showerror, e)]) + ) + end + + + + # MPC requires fixed PV and battery sizes. If not provided, call REopt first in a sizing run. + technology_sizes = get_technology_sizes!(d, model_inputs, solver_settings) + s = model_inputs.s # Access the processed Scenario struct + + # Skip MPC if no battery is optimally sized + if technology_sizes.skip_mpc + return build_mpc_response( + "skipped", + skip_mpc = true, + messages = Dict("info" => ["No battery was optimally sized in REopt pre-step; MPC dispatch not needed."]) + ) + end + + # Sizes from user or initial REopt run + pv_kw = technology_sizes.pv_kw + batt_kw = technology_sizes.batt_kw + batt_kwh = technology_sizes.batt_kwh + + # Note: REoptInputs does not provide PV production factors if user doesn't specify custom values + # PV production is NOT levelized in MPC but IS levelized in REopt. This will result in a slight mis-match. + if !isempty(s.pvs) # Get prod factors if PV considered. + if !isnothing(s.pvs[1].production_factor_series) + pv_prod_factor = Float64.(s.pvs[1].production_factor_series) + elseif technology_sizes.pv_production_factor_series !== nothing + pv_prod_factor = Float64.(technology_sizes.pv_production_factor_series) + else + pv_prod_factor = generate_pv_production_factors(d, time_steps_per_hour) + end + # Avoid another PVWatts call in the final REopt run. + d["PV"]["production_factor_series"] = pv_prod_factor + else + pv_prod_factor = zeros(Float64, length_of_data) + end + + loads_kw = Float64.(s.electric_load.loads_kw) + + # Extract tariff inputs relevant to MPC (use first tier only if tiered rates) + # TODO: Are all of these relevant? Any missing inputs? + # TODO: Implement lookback (demand_lookback_months, demand_lookback_percent, demand_lookback_range) Ignoring coincident peak charges for now + energy_rates = Float64.(s.electric_tariff.energy_rates[:, 1]) + monthly_demand_rates = isempty(s.electric_tariff.monthly_demand_rates) ? + zeros(Float64, 12) : Float64.(s.electric_tariff.monthly_demand_rates[:, 1]) + tou_demand_rates = isempty(s.electric_tariff.tou_demand_rates) ? Float64[] : Float64.(s.electric_tariff.tou_demand_rates[:, 1]) + tou_demand_ratchet_time_steps = [Int.(v) for v in s.electric_tariff.tou_demand_ratchet_time_steps] + # TODO: Monthly demand has never been tested + n_tou_ratchets = length(tou_demand_rates) # Number of TOU ratchets + tou_previous_peak_demands = zeros(Float64, n_tou_ratchets) # Tracks past TOU peak demand per ratchet + monthly_previous_peak_demands = zeros(Float64, 12) # Tracks past monthly peak demand + + println("monthly_demand_rates: ", monthly_demand_rates) + + + # Extract storage efficiency and SOC defaults from processed inputs + rect_eff = Float64(s.storage.attr["ElectricStorage"].rectifier_efficiency_fraction) + inv_eff = Float64(s.storage.attr["ElectricStorage"].inverter_efficiency_fraction) + int_eff = Float64(s.storage.attr["ElectricStorage"].internal_efficiency_fraction) + charge_eff = rect_eff * sqrt(int_eff) + discharge_eff = inv_eff * sqrt(int_eff) + soc_0 = Float64(s.storage.attr["ElectricStorage"].soc_init_fraction) + soc_min = Float64(s.storage.attr["ElectricStorage"].soc_min_fraction) + + # Extract emissions defaults (or use user input if provided) + co2_grid_emissions_series = Float64.(s.electric_utility.emissions_factor_series_lb_CO2_per_kwh) + + # --- Export / net metering setup (mirror the sizing-run scenario) --- + # NEM is enabled in MPC when ElectricUtility.net_metering_limit_kw > 0. + nm_limit_kw = Float64(s.electric_utility.net_metering_limit_kw) + + # WHL (net billing) rate: export_rates[:WHL] in the processed scenario is a negative "cost"; + # MPCElectricTariff expects a positive wholesale_rate (it negates internally). + whl_rate_full = (:WHL in s.electric_tariff.export_bins) ? + -1.0 .* Float64.(s.electric_tariff.export_rates[:WHL]) : nothing + + # PV export capability comes from the processed PV struct. + pv_can_net_meter = !isempty(s.pvs) ? Bool(s.pvs[1].can_net_meter) : false + pv_can_wholesale = !isempty(s.pvs) ? Bool(s.pvs[1].can_wholesale) : false + + # Battery export capability comes from the processed ElectricStorage struct + _batt_attr = s.storage.attr["ElectricStorage"] + batt_can_net_meter = Bool(_batt_attr.can_net_meter) + batt_can_wholesale = Bool(_batt_attr.can_wholesale) + + # NEM export is credited (approximately) at the retail energy rate; used for cost reporting below. + nem_active = nm_limit_kw > 0 && (pv_can_net_meter || batt_can_net_meter) + + month_starts = get_month_transition_timesteps(time_steps_per_hour) + + # ts_to_month = 8760 array specifying which month each timestep falls in (1-12) + ts_to_month = Vector{Int}(undef, length_of_data) + for m in 1:12 + s_idx = month_starts[m] + e = m < 12 ? month_starts[m+1] - 1 : length_of_data + ts_to_month[s_idx:e] .= m + end + + # ts_to_ratchet = 8760 array specifying which ratchet each timestep falls in + ts_to_ratchet = zeros(Int, length_of_data) + for (t, ratchet_ts) in enumerate(tou_demand_ratchet_time_steps), g in ratchet_ts + if 1 <= g <= length_of_data + ts_to_ratchet[g] = t + end + end + + # Saved dispatch series (first timestep of each MPC loop) + dispatch_series = Dict( + "PV" => Dict( + "electric_to_load_series_kw" => Float64[], + "electric_to_storage_series_kw" => Float64[], + "electric_to_grid_series_kw" => Float64[], + "electric_curtailed_series_kw" => Float64[], + ), + "ElectricStorage" => Dict( + "storage_to_load_series_kw" => Float64[], + "storage_to_grid_series_kw" => Float64[], + "soc_series_fraction" => Float64[], + ), + "ElectricUtility" => Dict( + "electric_to_load_series_kw" => Float64[], + "electric_to_storage_series_kw" => Float64[], + "emissions_series_lb_CO2" => Float64[], + ), + "ElectricLoad" => Dict( + "load_series_kw" => Float64[], + ), + ) + energy_cost_series = Float64[] # grid purchase (energy) charges per timestep + export_benefit_series = Float64[] # NEM/WHL export credits per timestep (positive = revenue) + total_energy_cost = 0.0 + total_export_benefit = 0.0 + soc_init_frac = soc_0 + + # Build MPC post + function build_mpc_post(current_horizon_pv, current_horizon_load, current_horizon_energy_rates, + current_horizon_emissions, current_horizon_tou_ts, current_horizon_monthly_ts, + tou_previous_peak_demands, monthly_previous_peak_demands, soc_init_frac, + current_horizon_whl_rate) + tariff = Dict( + "energy_rates" => current_horizon_energy_rates, + "tou_demand_rates" => tou_demand_rates, + "tou_demand_time_steps" => current_horizon_tou_ts, + "tou_previous_peak_demands" => tou_previous_peak_demands, + "monthly_demand_rates" => monthly_demand_rates, + "time_steps_monthly" => current_horizon_monthly_ts, + "monthly_previous_peak_demands" => monthly_previous_peak_demands, + ) + # WHL (net billing) export: MPCElectricTariff reads a positive `wholesale_rate` and + # builds the :WHL export bin when it is provided. + if current_horizon_whl_rate !== nothing + tariff["wholesale_rate"] = current_horizon_whl_rate + end + return Dict( + "PV" => Dict( + "size_kw" => pv_kw, + "production_factor_series" => current_horizon_pv, + "can_net_meter" => pv_can_net_meter, + "can_wholesale" => pv_can_wholesale, + ), + "ElectricStorage" => Dict( + "size_kw" => batt_kw, + "size_kwh" => batt_kwh, + "charge_efficiency" => charge_eff, + "discharge_efficiency" => discharge_eff, + "soc_init_fraction" => soc_init_frac, + "soc_min_fraction" => soc_min, + "can_net_meter" => batt_can_net_meter, + "can_wholesale" => batt_can_wholesale, + ), + "ElectricLoad" => Dict( + "loads_kw" => current_horizon_load, + ), + "ElectricTariff" => tariff, + # net_metering_limit_kw drives the NEM export bin in MPCElectricTariff (NEM on if > 0) + "ElectricUtility" => Dict( + "net_metering_limit_kw" => nm_limit_kw, + "emissions_factor_series_lb_CO2_per_kwh" => current_horizon_emissions, + ), + ) + end + + @info "MPC: starting rolling-horizon optimization ($(length_of_data) iterations, horizon = $(horizon) timesteps)" + for idx in 1:length_of_data + end_ts = idx + horizon - 1 + + current_horizon_pv = slice_data(pv_prod_factor, idx, end_ts) + current_horizon_load = slice_data(loads_kw, idx, end_ts) + current_horizon_energy_rates = slice_data(energy_rates, idx, end_ts) + current_horizon_emissions = slice_data(co2_grid_emissions_series, idx, end_ts) + current_horizon_whl_rate = whl_rate_full === nothing ? nothing : slice_data(whl_rate_full, idx, end_ts) + + # List of length n_tou_ratchets, specifies which ts of the current horizon are in each TOU ratchet + # by placing values 1 to horizon into the corresponding element of the array based on ratchet number + current_horizon_tou_ts = [Int[] for _ in 1:n_tou_ratchets] + + # 12 element list, each element for one month of the year. Specifies which timesteps of the current horizon are + # in each month by placing values 1 - horizon into the corresponding element of the array based on month number + current_horizon_monthly_ts = [Int[] for _ in 1:12] + for k in 1:horizon + g = idx + k - 1 + if g > length_of_data + g -= length_of_data + end + push!(current_horizon_monthly_ts[ts_to_month[g]], k) + ratchet = ts_to_ratchet[g] + if ratchet > 0 + push!(current_horizon_tou_ts[ratchet], k) + end + end + + post = build_mpc_post(current_horizon_pv, current_horizon_load, current_horizon_energy_rates, + current_horizon_emissions, current_horizon_tou_ts, current_horizon_monthly_ts, + tou_previous_peak_demands, monthly_previous_peak_demands, soc_init_frac, + current_horizon_whl_rate + ) + + model = get_solver_model(get_solver_model_type(solver_name), + SolverAttributes(per_iter_timeout_s, solver_settings["optimality_tolerance"])) + result = reoptjl.run_mpc(model, post) + + # Assume perfect forecast; save first timestep of results as the executed state + pv_res = result["PV"] + batt_res = result["ElectricStorage"] + util_res = result["ElectricUtility"] + + pv_to_load = pv_res["electric_to_load_series_kw"][1] + pv_to_batt = pv_res["electric_to_storage_series_kw"][1] + pv_to_grid = haskey(pv_res, "electric_to_grid_series_kw") ? pv_res["electric_to_grid_series_kw"][1] : 0.0 + pv_curtailed = haskey(pv_res, "electric_curtailed_series_kw") ? pv_res["electric_curtailed_series_kw"][1] : 0.0 + batt_to_load = batt_res["storage_to_load_series_kw"][1] + batt_to_grid = haskey(batt_res, "storage_to_grid_series_kw") ? batt_res["storage_to_grid_series_kw"][1] : 0.0 + batt_soc = batt_res["soc_series_fraction"][1] + util_to_load = util_res["electric_to_load_series_kw"][1] + util_to_batt = util_res["electric_to_storage_series_kw"][1] + grid_power = max(util_to_load + util_to_batt, 0.0) + + push!(dispatch_series["PV"]["electric_to_load_series_kw"], pv_to_load) + push!(dispatch_series["PV"]["electric_to_storage_series_kw"], pv_to_batt) + push!(dispatch_series["PV"]["electric_to_grid_series_kw"], pv_to_grid) + push!(dispatch_series["PV"]["electric_curtailed_series_kw"], pv_curtailed) + push!(dispatch_series["ElectricStorage"]["storage_to_load_series_kw"], batt_to_load) + push!(dispatch_series["ElectricStorage"]["storage_to_grid_series_kw"], batt_to_grid) + push!(dispatch_series["ElectricStorage"]["soc_series_fraction"], batt_soc) + push!(dispatch_series["ElectricUtility"]["electric_to_load_series_kw"], util_to_load) + push!(dispatch_series["ElectricUtility"]["electric_to_storage_series_kw"], util_to_batt) + push!(dispatch_series["ElectricUtility"]["emissions_series_lb_CO2"], + co2_grid_emissions_series[idx] * grid_power / time_steps_per_hour) + push!(dispatch_series["ElectricLoad"]["load_series_kw"], loads_kw[idx]) + + # Running electricity costs (reported as separate components; see results below) + # Energy (grid purchase) charge for this timestep + step_energy_charge = grid_power * energy_rates[idx] / time_steps_per_hour + + # Export credit for this timestep. NEM credits at the retail energy rate; otherwise WHL credits + # at the wholesale rate. NOTE: when both NEM and WHL are available the model chooses per horizon + # and the executed-timestep bin split is not returned, so this can misprice such timesteps. + step_export_kw = pv_to_grid + batt_to_grid + export_rate_idx = nem_active ? energy_rates[idx] : + (whl_rate_full !== nothing ? whl_rate_full[idx] : 0.0) + step_export_benefit = step_export_kw * export_rate_idx / time_steps_per_hour + + push!(energy_cost_series, step_energy_charge) + push!(export_benefit_series, step_export_benefit) + total_energy_cost += step_energy_charge + total_export_benefit += step_export_benefit + + soc_init_frac = batt_soc + + # Update monthly and TOU peak demand as max(current ts grid_power, previous max) + current_month = ts_to_month[idx] + monthly_previous_peak_demands[current_month] = max(grid_power, monthly_previous_peak_demands[current_month]) + + if n_tou_ratchets > 0 + current_ratchet = ts_to_ratchet[idx] + if current_ratchet > 0 + tou_previous_peak_demands[current_ratchet] = max(grid_power, tou_previous_peak_demands[current_ratchet]) + end + end + + end + + @info "MPC looping completed." + + # Calculate final demand costs + monthly_demand_cost_total = sum(monthly_previous_peak_demands .* monthly_demand_rates) + tou_demand_cost_total = n_tou_ratchets > 0 ? + sum(tou_previous_peak_demands .* tou_demand_rates) : 0.0 + + return build_mpc_response( + "optimal", + result_dict = Dict( + "MPC" => Dict( + "time_steps_per_hour" => time_steps_per_hour, + "horizon_time_steps" => horizon, + ), + "PV" => merge(Dict("size_kw" => pv_kw), dispatch_series["PV"]), + "ElectricStorage" => merge(Dict("size_kw" => batt_kw, "size_kwh" => batt_kwh), dispatch_series["ElectricStorage"]), + "ElectricUtility" => dispatch_series["ElectricUtility"], + "ElectricLoad" => dispatch_series["ElectricLoad"], + "ElectricTariff" => Dict( + # --- Cost components (before tax); combine for the total bill below --- + "total_energy_cost" => total_energy_cost, # grid purchases (energy) only + "total_export_benefit" => total_export_benefit, # NEM/WHL credits (positive = revenue) + "total_tou_demand_cost" => tou_demand_cost_total, + "total_non_tou_monthly_demand_cost" => monthly_demand_cost_total, + # --- Total electricity bill = energy charge - export benefit + demand charges --- + "total_electricity_bill" => total_energy_cost - total_export_benefit + + tou_demand_cost_total + monthly_demand_cost_total, + # --- Per-timestep series (energy/exports only; demand is a peak-based charge) --- + "energy_cost_series_per_timestep" => energy_cost_series, + "export_benefit_series_per_timestep" => export_benefit_series, + "tou_peaks_by_ratchet_kw" => tou_previous_peak_demands, + "monthly_peaks_kw" => monthly_previous_peak_demands, + ), + ) + ) +end diff --git a/reoptjl/custom_table_config.py b/reoptjl/custom_table_config.py index 4cd3ec479..4d52f13cc 100644 --- a/reoptjl/custom_table_config.py +++ b/reoptjl/custom_table_config.py @@ -804,6 +804,12 @@ "bau_value" : lambda df: safe_get(df, "outputs.ElectricStorage.storage_to_load_series_kw_bau"), "scenario_value": lambda df: safe_get(df, "outputs.ElectricStorage.storage_to_load_series_kw") }, + { + "label" : "Battery Exported to Grid (kWh/yr)", + "key" : "battery_exported_to_grid", + "bau_value" : lambda df: safe_get(df, "outputs.ElectricStorage.storage_to_grid_series_kw_bau"), + "scenario_value": lambda df: safe_get(df, "outputs.ElectricStorage.storage_to_grid_series_kw") + }, { "label" : "Generator Serving Load (kWh/yr)", "key" : "generator_serving_load", diff --git a/reoptjl/migrations/0117_merge_20260530_1640.py b/reoptjl/migrations/0117_merge_20260530_1640.py new file mode 100644 index 000000000..0302ff8e6 --- /dev/null +++ b/reoptjl/migrations/0117_merge_20260530_1640.py @@ -0,0 +1,14 @@ +# Generated by Django 4.2.26 on 2026-05-30 16:40 + +from django.db import migrations + + +class Migration(migrations.Migration): + + dependencies = [ + ('reoptjl', '0116_alter_apimeta_api_key_alter_apimeta_job_type'), + ('reoptjl', '0116_rename_state_of_health_electricstorageoutputs_state_of_health_series_fraction_and_more'), + ] + + operations = [ + ] diff --git a/reoptjl/migrations/0118_merge_20260629_2347.py b/reoptjl/migrations/0118_merge_20260629_2347.py new file mode 100644 index 000000000..e55228f28 --- /dev/null +++ b/reoptjl/migrations/0118_merge_20260629_2347.py @@ -0,0 +1,14 @@ +# Generated by Django 4.2.26 on 2026-06-29 23:47 + +from django.db import migrations + + +class Migration(migrations.Migration): + + dependencies = [ + ('reoptjl', '0117_merge_20260530_1640'), + ('reoptjl', '0117_merge_20260601_2056'), + ] + + operations = [ + ] diff --git a/reoptjl/migrations/0119_alter_electricstorageinputs_fixed_soc_series_fraction_tolerance.py b/reoptjl/migrations/0119_alter_electricstorageinputs_fixed_soc_series_fraction_tolerance.py new file mode 100644 index 000000000..2378b63f3 --- /dev/null +++ b/reoptjl/migrations/0119_alter_electricstorageinputs_fixed_soc_series_fraction_tolerance.py @@ -0,0 +1,19 @@ +# Generated by Django 4.2.26 on 2026-07-09 18:09 + +import django.core.validators +from django.db import migrations, models + + +class Migration(migrations.Migration): + + dependencies = [ + ('reoptjl', '0118_merge_20260629_2347'), + ] + + operations = [ + migrations.AlterField( + model_name='electricstorageinputs', + name='fixed_soc_series_fraction_tolerance', + field=models.FloatField(blank=True, help_text='Absolute tolerance on fixed_soc_series_fraction to avoid infeasible solutions when fixed_soc_series_fraction is provided.', null=True, validators=[django.core.validators.MinValueValidator(0), django.core.validators.MaxValueValidator(1)]), + ), + ] diff --git a/reoptjl/migrations/0120_electricstorageinputs_dispatch_strategy_and_more.py b/reoptjl/migrations/0120_electricstorageinputs_dispatch_strategy_and_more.py new file mode 100644 index 000000000..486adb62d --- /dev/null +++ b/reoptjl/migrations/0120_electricstorageinputs_dispatch_strategy_and_more.py @@ -0,0 +1,24 @@ +# Generated by Django 4.2.26 on 2026-07-09 20:32 + +import django.core.validators +from django.db import migrations, models + + +class Migration(migrations.Migration): + + dependencies = [ + ('reoptjl', '0119_alter_electricstorageinputs_fixed_soc_series_fraction_tolerance'), + ] + + operations = [ + migrations.AddField( + model_name='electricstorageinputs', + name='dispatch_strategy', + field=models.TextField(blank=True, choices=[('optimized', 'Optimized'), ('peak_shaving_look_ahead', 'Peak Shaving Look Ahead'), ('peak_shaving_look_behind', 'Peak Shaving Look Behind'), ('self_consumption', 'Self Consumption'), ('backup', 'Backup'), ('custom_soc', 'Custom Soc'), ('daily_foresight_optimized', 'Daily Foresight Optimized')], default='optimized', help_text='Electric storage dispatch strategy can be one of: optimized, peak_shaving_look_ahead, peak_shaving_look_behind, self_consumption, backup, custom_soc, daily_foresight_optimized', null=True), + ), + migrations.AlterField( + model_name='electricstorageinputs', + name='soc_min_fraction', + field=models.FloatField(blank=True, help_text='Minimum allowable battery state of charge as fraction of energy capacity.', validators=[django.core.validators.MinValueValidator(0), django.core.validators.MaxValueValidator(1.0)]), + ), + ] diff --git a/reoptjl/migrations/0121_electricstorageinputs_can_export_beyond_nem_limit_and_more.py b/reoptjl/migrations/0121_electricstorageinputs_can_export_beyond_nem_limit_and_more.py new file mode 100644 index 000000000..d1c39e257 --- /dev/null +++ b/reoptjl/migrations/0121_electricstorageinputs_can_export_beyond_nem_limit_and_more.py @@ -0,0 +1,34 @@ +# Generated by Django 4.2.26 on 2026-07-15 19:47 + +import django.contrib.postgres.fields +from django.db import migrations, models + + +class Migration(migrations.Migration): + + dependencies = [ + ('reoptjl', '0120_electricstorageinputs_dispatch_strategy_and_more'), + ] + + operations = [ + migrations.AddField( + model_name='electricstorageinputs', + name='can_export_beyond_nem_limit', + field=models.BooleanField(blank=True, default=False, help_text='True/False for if technology can export energy beyond the annual site load (and be compensated for that energy at the export_rate_beyond_net_metering_limit).'), + ), + migrations.AddField( + model_name='electricstorageinputs', + name='can_net_meter', + field=models.BooleanField(blank=True, default=False, help_text='True/False for if technology has option to participate in net metering agreement with utility. Note that a technology can only participate in either net metering or wholesale rates (not both).'), + ), + migrations.AddField( + model_name='electricstorageinputs', + name='can_wholesale', + field=models.BooleanField(blank=True, default=False, help_text='True/False for if technology has option to export energy that is compensated at the wholesale_rate. Note that a technology can only participate in either net metering or wholesale rates (not both).'), + ), + migrations.AddField( + model_name='electricstorageoutputs', + name='storage_to_grid_series_kw', + field=django.contrib.postgres.fields.ArrayField(base_field=models.FloatField(blank=True, null=True), blank=True, default=list, size=None), + ), + ] diff --git a/reoptjl/models.py b/reoptjl/models.py index cee7bfaa1..b3d0e14f0 100644 --- a/reoptjl/models.py +++ b/reoptjl/models.py @@ -3668,6 +3668,16 @@ class ElectricStorageInputs(BaseModel, models.Model): primary_key=True ) + ELECTRICSTORAGE_DISPATCH_STRATEGY = models.TextChoices('ELECTRICSTORAGE_DISPATCH_STRATEGY', ( + "optimized", + "peak_shaving_look_ahead", + "peak_shaving_look_behind", + "self_consumption", + "backup", + "custom_soc", + "daily_foresight_optimized" + )) + min_kw = models.FloatField( default=0, validators=[ @@ -3731,8 +3741,14 @@ class ElectricStorageInputs(BaseModel, models.Model): blank=True, help_text="Battery rectifier efficiency" ) + dispatch_strategy = models.TextField( + default=ELECTRICSTORAGE_DISPATCH_STRATEGY.optimized, + choices=ELECTRICSTORAGE_DISPATCH_STRATEGY.choices, + null=True, + blank=True, + help_text="Electric storage dispatch strategy can be one of: optimized, peak_shaving_look_ahead, peak_shaving_look_behind, self_consumption, backup, custom_soc, daily_foresight_optimized" + ) soc_min_fraction = models.FloatField( - default=0.2, validators=[ MinValueValidator(0), MaxValueValidator(1.0) @@ -3777,6 +3793,24 @@ class ElectricStorageInputs(BaseModel, models.Model): blank=True, help_text="Flag to set whether the battery can be charged from the grid, or just onsite generation." ) + can_net_meter = models.BooleanField( + default=False, + blank=True, + help_text=("True/False for if technology has option to participate in net metering agreement with utility. " + "Note that a technology can only participate in either net metering or wholesale rates (not both).") + ) + can_wholesale = models.BooleanField( + default=False, + blank=True, + help_text=("True/False for if technology has option to export energy that is compensated at the wholesale_rate. " + "Note that a technology can only participate in either net metering or wholesale rates (not both).") + ) + can_export_beyond_nem_limit = models.BooleanField( + default=False, + blank=True, + help_text=("True/False for if technology can export energy beyond the annual site load (and be compensated for " + "that energy at the export_rate_beyond_net_metering_limit).") + ) installed_cost_per_kw = models.FloatField( default=968.0, validators=[ @@ -3962,6 +3996,10 @@ class ElectricStorageOutputs(BaseModel, models.Model): models.FloatField(null=True, blank=True), blank=True, default=list ) + storage_to_grid_series_kw = ArrayField( + models.FloatField(null=True, blank=True), + blank=True, default=list + ) initial_capital_cost = models.FloatField(null=True, blank=True) maintenance_cost = models.FloatField(null=True, blank=True) state_of_health_series_fraction = ArrayField( diff --git a/reoptjl/test/posts/all_inputs_test.json b/reoptjl/test/posts/all_inputs_test.json index d53def1ad..3cc26956b 100644 --- a/reoptjl/test/posts/all_inputs_test.json +++ b/reoptjl/test/posts/all_inputs_test.json @@ -167,6 +167,7 @@ "internal_efficiency_fraction": 0.975, "inverter_efficiency_fraction": 0.96, "rectifier_efficiency_fraction": 0.96, + "dispatch_strategy": "optimized", "soc_min_fraction": 0.2, "soc_min_applies_during_outages": true, "soc_init_fraction": 0.5, diff --git a/reoptjl/validators.py b/reoptjl/validators.py index 99cdc456a..76dbd852a 100644 --- a/reoptjl/validators.py +++ b/reoptjl/validators.py @@ -363,6 +363,12 @@ def update_pv_defaults_offgrid(self, pvmodel): else: self.models["ElectricStorage"].soc_init_fraction = 1.0 + if self.models["ElectricStorage"].__getattribute__("soc_min_fraction") == None: + if self.models["ElectricStorage"].dispatch_strategy=="backup": + self.models["ElectricStorage"].soc_min_fraction = 0.8 + else: + self.models["ElectricStorage"].soc_min_fraction = 0.2 + if self.models["ElectricStorage"].__getattribute__("can_grid_charge") == None: if self.models["Settings"].off_grid_flag==False: self.models["ElectricStorage"].can_grid_charge = True