diff --git a/.git-blame-ignore-revs b/.git-blame-ignore-revs index 536e3289aa..c228d15181 100644 --- a/.git-blame-ignore-revs +++ b/.git-blame-ignore-revs @@ -1,3 +1,5 @@ # .git-blame-ignore-revs # Standardize code formatting across project (#673) ee3f08756584ba16a57bb701492270a7bf129b4d +# Update code formatting +730f91df23447e94177c3a9c3d4e553cb502e2bf \ No newline at end of file diff --git a/.github/workflows/format_suggestions.yml b/.github/workflows/format_suggestions.yml index 05e574dd68..dbd307846d 100644 --- a/.github/workflows/format_suggestions.yml +++ b/.github/workflows/format_suggestions.yml @@ -7,3 +7,8 @@ jobs: runs-on: ubuntu-latest steps: - uses: julia-actions/julia-format@v2 + continue-on-error: true + - name: Check on failures + if: steps.julia-format.outcome != 'success' + run: echo "There are formatting errors. Please check the logs above." + shell: bash \ No newline at end of file diff --git a/CHANGELOG.md b/CHANGELOG.md index 60e91ac0f2..8b4478add3 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -43,13 +43,18 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 non-zero in multi-stage GenX (#666) - Added condition number scaling added to objective function (#667) - Added versioned doc-pages for v0.3.6 and v0.4.0 + +- Added a warning message in write_costs_multistage mentioning th approximate value of costs currently. + +### Fixed - Add constraint in mga to compute total capacity in each zone from a given technology type (#681) - New settings parameter MGAAnnualGeneration to switch between different MGA formulations (#681) - Add validation for `Can_Retire` column in multi-stage GenX since the current implementation does not allow a resource to switch from can_retire = 0 to can_retire = 1 between stages. (#683) - Add tutorials for running GenX (#637 and #685) - -### Fixed +- Add writing of multistage stats during optimization with foresight (#687) +- Fix docstring in operational_reserves.jl (#690) +- Fix docstring in energy_share_requirement.jl (#692) - Set MUST_RUN=1 for RealSystemExample/small_hydro plants (#517). Previously these plants had no resource flag set, and so they did not contribute to the power balance. As these plants are now useful, the objective in these cases is slightly lower. @@ -70,6 +75,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Fixed outputting capital recovery cost to 0 if the remaining number of years is 0 (#666) - Updated the docstring for the initialize_cost_to_go function and adjusted the formula for the discount factor to reflect the code implementation (#672). - Fix write_multi_stage_cost.jl: add discount with OPEX multipliers to cUnmetPolicyPenalty (#679) +- Fix DF calculation in DDP to make it more generic for variable length stages (#680) - Fix write_power_balance.jl: add additional two columns ("VRE_Storage_Discharge" and "VRE_Storage_Charge") for VRE_STOR ### Changed diff --git a/docs/src/Model_Reference/write_outputs.md b/docs/src/Model_Reference/write_outputs.md index 1c50d00dcb..4bb02f563f 100644 --- a/docs/src/Model_Reference/write_outputs.md +++ b/docs/src/Model_Reference/write_outputs.md @@ -142,13 +142,20 @@ GenX.write_vre_stor_discharge ``` ## Write Multi-stage files +```@autodocs +Modules = [GenX] +Pages = ["write_multi_stage_outputs.jl"] +``` ```@docs GenX.write_multi_stage_costs GenX.write_multi_stage_stats GenX.write_multi_stage_settings GenX.write_multi_stage_network_expansion +GenX.write_multi_stage_capacities_discharge GenX.write_multi_stage_capacities_charge GenX.write_multi_stage_capacities_energy +GenX.create_multi_stage_stats_file +GenX.update_multi_stage_stats_file ``` ## Write maintenance files diff --git a/src/additional_tools/method_of_morris.jl b/src/additional_tools/method_of_morris.jl index d4628b55fa..ed5cd5ddfa 100644 --- a/src/additional_tools/method_of_morris.jl +++ b/src/additional_tools/method_of_morris.jl @@ -62,12 +62,12 @@ function calculate_spread(matrix) end function sample_matrices(p_range, - p_steps, - rng; - num_trajectory, - total_num_trajectory, - len_design_mat, - groups) + p_steps, + rng; + num_trajectory, + total_num_trajectory, + len_design_mat, + groups) matrix_array = [] println(num_trajectory) println(total_num_trajectory) @@ -85,13 +85,13 @@ function sample_matrices(p_range, end function my_gsa(f, - p_steps, - num_trajectory, - total_num_trajectory, - p_range::AbstractVector, - len_design_mat, - groups, - random) + p_steps, + num_trajectory, + total_num_trajectory, + p_range::AbstractVector, + len_design_mat, + groups, + random) rng = Random.default_rng() if !random Random.seed!(SEED) @@ -113,7 +113,8 @@ function my_gsa(f, pairs(eachcol(L[!, ((i - 1) * len_design_mat + 1):(i * len_design_mat)])))))) distinct_trajectories[i] = length(Matrix(DataFrame(unique(last, - pairs(eachcol(L[!, ((i - 1) * len_design_mat + 1):(i * len_design_mat)])))))[1, + pairs(eachcol(L[!, ((i - 1) * len_design_mat + 1):(i * len_design_mat)])))))[ + 1, :]) end end @@ -189,12 +190,12 @@ function my_gsa(f, effects) end function morris(EP::Model, - path::AbstractString, - setup::Dict, - inputs::Dict, - outpath::AbstractString, - OPTIMIZER; - random = true) + path::AbstractString, + setup::Dict, + inputs::Dict, + outpath::AbstractString, + OPTIMIZER; + random = true) # Reading the input parameters Morris_range = load_dataframe(joinpath(path, "Method_of_morris_range.csv")) @@ -215,11 +216,12 @@ function morris(EP::Model, column_f = isdefined(GenX, col_sym) ? getfield(GenX, col_sym) : r -> getproperty(r, col_sym) sigma = [sigma; - [column_f.(gen) .* (1 .+ - Morris_range[Morris_range[!, :Parameter] .== column, :Lower_bound] ./ 100) column_f.(gen) .* - (1 .+ - Morris_range[Morris_range[!, :Parameter] .== column, - :Upper_bound] ./ 100)]] + [column_f.(gen) .* (1 .+ + Morris_range[Morris_range[!, :Parameter] .== column, :Lower_bound] ./ + 100) column_f.(gen) .* + (1 .+ + Morris_range[Morris_range[!, :Parameter] .== column, + :Upper_bound] ./ 100)]] end sigma = sigma[2:end, :] diff --git a/src/case_runners/case_runner.jl b/src/case_runners/case_runner.jl index 0fd002b206..1c65298b5e 100644 --- a/src/case_runners/case_runner.jl +++ b/src/case_runners/case_runner.jl @@ -168,12 +168,7 @@ function run_genx_case_multistage!(case::AbstractString, mysetup::Dict, optimize ### Solve model println("Solving Model") - # Step 3) Run DDP Algorithm - ## Solve Model - model_dict, mystats_d, inputs_dict = run_ddp(model_dict, mysetup, inputs_dict) - - # Step 4) Write final outputs from each stage - + # Prepare folder for results outpath = get_default_output_folder(case) if mysetup["OverwriteResults"] == 1 @@ -188,6 +183,11 @@ function run_genx_case_multistage!(case::AbstractString, mysetup::Dict, optimize mkdir(outpath) end + # Step 3) Run DDP Algorithm + ## Solve Model + model_dict, mystats_d, inputs_dict = run_ddp(outpath, model_dict, mysetup, inputs_dict) + + # Step 4) Write final outputs from each stage for p in 1:mysetup["MultiStageSettingsDict"]["NumStages"] outpath_cur = joinpath(outpath, "results_p$p") write_outputs(model_dict[p], outpath_cur, mysetup, inputs_dict[p]) diff --git a/src/configure_settings/configure_settings.jl b/src/configure_settings/configure_settings.jl index 2054b2a50a..bd130e151f 100644 --- a/src/configure_settings/configure_settings.jl +++ b/src/configure_settings/configure_settings.jl @@ -82,7 +82,8 @@ function validate_settings!(settings::Dict{Any, Any}) if haskey(settings, "Reserves") Base.depwarn("""The Reserves setting has been deprecated. Please use the - OperationalReserves setting instead.""", :validate_settings!, force = true) + OperationalReserves setting instead.""", + :validate_settings!, force = true) settings["OperationalReserves"] = settings["Reserves"] delete!(settings, "Reserves") end diff --git a/src/load_inputs/load_dataframe.jl b/src/load_inputs/load_dataframe.jl index 64aa112968..bd212b75ef 100644 --- a/src/load_inputs/load_dataframe.jl +++ b/src/load_inputs/load_dataframe.jl @@ -64,7 +64,8 @@ function load_dataframe(dir::AbstractString, basenames::Vector{String})::DataFra target = look_for_file_with_alternate_case(dir, base) # admonish if target != FILENOTFOUND - Base.depwarn("""The filename '$target' is deprecated. '$best_basename' is preferred.""", + Base.depwarn( + """The filename '$target' is deprecated. '$best_basename' is preferred.""", :load_dataframe, force = true) return load_dataframe_from_file(joinpath(dir, target)) @@ -132,8 +133,8 @@ function load_dataframe_from_file(path)::DataFrame end function find_matrix_columns_in_dataframe(df::DataFrame, - columnprefix::AbstractString; - prefixseparator = '_')::Vector{Int} + columnprefix::AbstractString; + prefixseparator = '_')::Vector{Int} all_columns = names(df) # 2 is the length of the '_' connector plus one for indexing @@ -167,8 +168,8 @@ ESR_1, other_thing, ESR_3, ESR_2, ``` """ function extract_matrix_from_dataframe(df::DataFrame, - columnprefix::AbstractString; - prefixseparator = '_') + columnprefix::AbstractString; + prefixseparator = '_') all_columns = names(df) columnnumbers = find_matrix_columns_in_dataframe(df, columnprefix, @@ -193,12 +194,13 @@ function extract_matrix_from_dataframe(df::DataFrame, end function extract_matrix_from_resources(rs::Vector{T}, - columnprefix::AbstractString, - default = 0.0) where {T <: AbstractResource} + columnprefix::AbstractString, + default = 0.0) where {T <: AbstractResource} # attributes starting with columnprefix with a numeric suffix attributes_n = [attr - for attr in string.(attributes(rs[1])) if startswith(attr, columnprefix)] + for attr in string.(attributes(rs[1])) + if startswith(attr, columnprefix)] # sort the attributes by the numeric suffix sort!(attributes_n, by = x -> parse(Int, split(x, "_")[end])) diff --git a/src/load_inputs/load_demand_data.jl b/src/load_inputs/load_demand_data.jl index 4c0e8a0319..52c5bd7bf2 100644 --- a/src/load_inputs/load_demand_data.jl +++ b/src/load_inputs/load_demand_data.jl @@ -7,12 +7,13 @@ function get_demand_dataframe(path) DEMAND_COLUMN_PREFIX_DEPRECATED()[1:(end - 1)], prefixseparator = 'z') old_column_symbols = Symbol.(DEMAND_COLUMN_PREFIX_DEPRECATED() * string(i) - for i in old_columns) + for i in old_columns) if length(old_column_symbols) > 0 pref_prefix = DEMAND_COLUMN_PREFIX() dep_prefix = DEMAND_COLUMN_PREFIX_DEPRECATED() @info "$dep_prefix is deprecated. Use $pref_prefix." - new_column_symbols = Symbol.(DEMAND_COLUMN_PREFIX() * string(i) for i in old_columns) + new_column_symbols = Symbol.(DEMAND_COLUMN_PREFIX() * string(i) + for i in old_columns) rename!(df, Dict(old_column_symbols .=> new_column_symbols)) end return df diff --git a/src/load_inputs/load_inputs.jl b/src/load_inputs/load_inputs.jl index 1b8705ec4e..e17a390ced 100644 --- a/src/load_inputs/load_inputs.jl +++ b/src/load_inputs/load_inputs.jl @@ -117,8 +117,8 @@ Returns: - String: The directory path based on the setup parameters. """ function get_systemfiles_path(setup::Dict, - TDR_directory::AbstractString, - path::AbstractString) + TDR_directory::AbstractString, + path::AbstractString) if setup["TimeDomainReduction"] == 1 && time_domain_reduced_files_exist(TDR_directory) return TDR_directory else diff --git a/src/load_inputs/load_multistage_data.jl b/src/load_inputs/load_multistage_data.jl index 95395c726f..2a51aaf2e6 100644 --- a/src/load_inputs/load_multistage_data.jl +++ b/src/load_inputs/load_multistage_data.jl @@ -34,7 +34,7 @@ function scale_multistage_data!(multistage_in::DataFrame, scale_factor::Float64) :min_retired_cap_charge_dc_mw, :min_retired_cap_charge_ac_mw, :min_retired_cap_discharge_dc_mw, - :min_retired_cap_discharge_ac_mw, + :min_retired_cap_discharge_ac_mw ] scale_columns!(multistage_in, columns_to_scale, scale_factor) return nothing diff --git a/src/load_inputs/load_operational_reserves.jl b/src/load_inputs/load_operational_reserves.jl index 7aad7a74de..6b6d67cb78 100644 --- a/src/load_inputs/load_operational_reserves.jl +++ b/src/load_inputs/load_operational_reserves.jl @@ -19,7 +19,8 @@ function load_operational_reserves!(setup::Dict, path::AbstractString, inputs::D end for col in columns if col in all_columns - Base.depwarn("The column name $col in file $filename is deprecated; prefer $best", + Base.depwarn( + "The column name $col in file $filename is deprecated; prefer $best", :load_operational_reserves, force = true) return float(df[firstrow, col]) diff --git a/src/load_inputs/load_resources_data.jl b/src/load_inputs/load_resources_data.jl index 979d9e7c21..31c3d9d1c3 100644 --- a/src/load_inputs/load_resources_data.jl +++ b/src/load_inputs/load_resources_data.jl @@ -35,7 +35,8 @@ function _get_policyfile_info() min_cap_filenames = ["Resource_minimum_capacity_requirement.csv"] max_cap_filenames = ["Resource_maximum_capacity_requirement.csv"] - policyfile_info = (esr = (filenames = esr_filenames, + policyfile_info = ( + esr = (filenames = esr_filenames, setup_param = "EnergyShareRequirement"), cap_res = (filenames = cap_res_filenames, setup_param = "CapacityReserveMargin"), min_cap = (filenames = min_cap_filenames, setup_param = "MinCapReq"), @@ -100,7 +101,7 @@ function scale_resources_data!(resource_in::DataFrame, scale_factor::Float64) :min_retired_charge_cap_mw, # to GW :min_retired_energy_cap_mw, # to GW :start_cost_per_mw, # to $M/GW - :ccs_disposal_cost_per_metric_ton, :hydrogen_mwh_per_tonne, # to GWh/t + :ccs_disposal_cost_per_metric_ton, :hydrogen_mwh_per_tonne # to GWh/t ] scale_columns!(resource_in, columns_to_scale, scale_factor) @@ -183,8 +184,8 @@ Scales in-place the columns in `columns_to_scale` of a dataframe `df` by a `scal """ function scale_columns!(df::DataFrame, - columns_to_scale::Vector{Symbol}, - scale_factor::Float64) + columns_to_scale::Vector{Symbol}, + scale_factor::Float64) for column in columns_to_scale if string(column) in names(df) df[!, column] /= scale_factor @@ -300,8 +301,8 @@ Construct the array of resources from multiple files of different types located """ function create_resource_array(resource_folder::AbstractString, - resources_info::NamedTuple, - scale_factor::Float64 = 1.0) + resources_info::NamedTuple, + scale_factor::Float64 = 1.0) resource_id_offset = 0 resources = [] # loop over available types and load all resources in resource_folder @@ -530,7 +531,8 @@ function validate_policy_files(resource_policies_path::AbstractString, setup::Di policyfile_info = _get_policyfile_info() for (filenames, setup_param) in values(policyfile_info) if setup[setup_param] == 1 && - any(!isfile(joinpath(resource_policies_path, filename)) for filename in filenames) + any(!isfile(joinpath(resource_policies_path, filename)) + for filename in filenames) msg = string(setup_param, " is set to 1 in settings but the required file(s) ", filenames, @@ -579,7 +581,7 @@ function validate_policy_dataframe!(filename::AbstractString, policy_in::DataFra end # Check that all policy columns have names with format "[policy_name]_[tagnum]" if !all(any([occursin(Regex("$(y)") * r"_\d", col) for y in accepted_cols]) - for col in cols) + for col in cols) msg = "Columns in policy file $filename must have names with format \"[policy_name]_[tagnum]\", case insensitive. (e.g., ESR_1, Min_Cap_1, Max_Cap_2, etc.)." error(msg) end @@ -598,8 +600,8 @@ Adds a set of new attributes (names and corresponding values) to a resource. The """ function add_attributes_to_resource!(resource::AbstractResource, - new_symbols::Vector{Symbol}, - new_values::T) where {T <: DataFrameRow} + new_symbols::Vector{Symbol}, + new_values::T) where {T <: DataFrameRow} # loop over new attributes for (sym, value) in zip(new_symbols, new_values) # add attribute to resource @@ -643,8 +645,8 @@ Loads a single policy file and adds the columns as new attributes to resources i - `filename::AbstractString`: The name of the policy file. """ function add_policy_to_resources!(resources::Vector{<:AbstractResource}, - path::AbstractString, - filename::AbstractString) + path::AbstractString, + filename::AbstractString) policy_in = load_dataframe(path) # check if policy file has any attributes, validate column names validate_policy_dataframe!(filename, policy_in) @@ -663,7 +665,7 @@ Reads policy files and adds policies-related attributes to resources in the mode - `resources_path::AbstractString`: The path to the resources folder. """ function add_policies_to_resources!(resources::Vector{<:AbstractResource}, - resource_policy_path::AbstractString) + resource_policy_path::AbstractString) # get filename for each type of policy available in GenX policies_info = _get_policyfile_info() # loop over policy files @@ -690,7 +692,7 @@ Reads module dataframe and adds columns as new attributes to the resources in th - `module_in::DataFrame`: The dataframe with the columns to add to the resources. """ function add_module_to_resources!(resources::Vector{<:AbstractResource}, - module_in::DataFrame) + module_in::DataFrame) # add module columns to resources as new attributes add_df_to_resources!(resources, module_in) return nothing @@ -707,8 +709,8 @@ Reads module dataframes, loops over files and adds columns as new attributes to - `resources_path::AbstractString`: The path to the resources folder. """ function add_modules_to_resources!(resources::Vector{<:AbstractResource}, - setup::Dict, - resources_path::AbstractString) + setup::Dict, + resources_path::AbstractString) scale_factor = setup["ParameterScale"] == 1 ? ModelScalingFactor : 1.0 modules = Vector{DataFrame}() @@ -769,8 +771,8 @@ Reads piecewise fuel usage data from the vector of generators, create a PWFU_dat - `inputs::Dict`: The dictionary containing the input data """ function process_piecewisefuelusage!(setup::Dict, - gen::Vector{<:AbstractResource}, - inputs::Dict) + gen::Vector{<:AbstractResource}, + inputs::Dict) inputs["PWFU_Num_Segments"] = 0 inputs["THERM_COMMIT_PWFU"] = Int64[] @@ -846,7 +848,8 @@ function process_piecewisefuelusage!(setup::Dict, # create a PWFU_data that contain processed intercept and slope (i.e., heat rate) intercept_cols = [Symbol("pwfu_intercept_", i) for i in 1:num_segments] intercept_df = DataFrame(intercept_mat, Symbol.(intercept_cols)) - slope_cols = Symbol.(filter(colname -> startswith(string(colname), + slope_cols = Symbol.(filter( + colname -> startswith(string(colname), "pwfu_heat_rate_mmbtu_per_mwh"), collect(attributes(thermal_gen[1])))) sort!(slope_cols, by = x -> parse(Int, split(string(x), "_")[end])) @@ -962,9 +965,9 @@ Adds resources to the `inputs` `Dict` with the key "RESOURCES" together with sev """ function add_resources_to_input_data!(inputs::Dict, - setup::Dict, - case_path::AbstractString, - gen::Vector{<:AbstractResource}) + setup::Dict, + case_path::AbstractString, + gen::Vector{<:AbstractResource}) # Number of resources G = length(gen) @@ -1297,9 +1300,9 @@ Raises: DeprecationWarning: If the `Generators_data.csv` file is found, a deprecation warning is issued, together with an error message. """ function load_resources_data!(inputs::Dict, - setup::Dict, - case_path::AbstractString, - resources_path::AbstractString) + setup::Dict, + case_path::AbstractString, + resources_path::AbstractString) if isfile(joinpath(case_path, "Generators_data.csv")) msg = "The `Generators_data.csv` file was deprecated in release v0.4. " * "Please use the new interface for generators creation, and see the documentation for additional details." @@ -1332,9 +1335,9 @@ end Function for reading input parameters related to multi fuels """ function load_multi_fuels_data!(inputs::Dict, - gen::Vector{<:AbstractResource}, - setup::Dict, - path::AbstractString) + gen::Vector{<:AbstractResource}, + setup::Dict, + path::AbstractString) inputs["NUM_FUELS"] = num_fuels.(gen) # Number of fuels that this resource can use max_fuels = maximum(inputs["NUM_FUELS"]) inputs["FUEL_COLS"] = [Symbol(string("Fuel", f)) for f in 1:max_fuels] diff --git a/src/model/core/co2.jl b/src/model/core/co2.jl index 95c3dd27de..7b03d7920a 100644 --- a/src/model/core/co2.jl +++ b/src/model/core/co2.jl @@ -95,10 +95,10 @@ function co2!(EP::Model, inputs::Dict) else sum((1 - biomass(gen[y]) - co2_capture_fraction(gen[y])) * EP[:vMulFuels][y, i, t] * fuel_CO2[fuel_cols(gen[y], tag = i)] - for i in 1:max_fuels) + + for i in 1:max_fuels) + sum((1 - biomass(gen[y]) - co2_capture_fraction_startup(gen[y])) * EP[:vMulStartFuels][y, i, t] * fuel_CO2[fuel_cols(gen[y], tag = i)] - for i in 1:max_fuels) + for i in 1:max_fuels) end) # CO2 captured from power plants in "Generators_data.csv" @@ -116,7 +116,7 @@ function co2!(EP::Model, inputs::Dict) @expression(EP, eEmissionsCaptureByPlantYear[y in CCS], sum(omega[t] * eEmissionsCaptureByPlant[y, t] - for t in 1:T)) + for t in 1:T)) # add CO2 sequestration cost to objective function # when scale factor is on tCO2/MWh = > kt CO2/GWh @expression(EP, ePlantCCO2Sequestration[y in CCS], @@ -125,7 +125,7 @@ function co2!(EP::Model, inputs::Dict) @expression(EP, eZonalCCO2Sequestration[z = 1:Z], sum(ePlantCCO2Sequestration[y] - for y in intersect(resources_in_zone_by_rid(gen, z), CCS))) + for y in intersect(resources_in_zone_by_rid(gen, z), CCS))) @expression(EP, eTotaleCCO2Sequestration, sum(eZonalCCO2Sequestration[z] for z in 1:Z)) diff --git a/src/model/core/discharge/discharge.jl b/src/model/core/discharge/discharge.jl index d67881c942..17bc2187dd 100644 --- a/src/model/core/discharge/discharge.jl +++ b/src/model/core/discharge/discharge.jl @@ -42,9 +42,9 @@ function discharge!(EP::Model, inputs::Dict, setup::Dict) if setup["EnergyShareRequirement"] >= 1 @expression(EP, eESRDischarge[ESR = 1:inputs["nESR"]], +sum(inputs["omega"][t] * esr(gen[y], tag = ESR) * EP[:vP][y, t] - for y in ids_with_policy(gen, esr, tag = ESR), t in 1:T) + for y in ids_with_policy(gen, esr, tag = ESR), t in 1:T) -sum(inputs["dfESR"][z, ESR] * inputs["omega"][t] * inputs["pD"][t, z] - for t in 1:T, z in findall(x -> x > 0, inputs["dfESR"][:, ESR]))) + for t in 1:T, z in findall(x -> x > 0, inputs["dfESR"][:, ESR]))) add_similar_to_expression!(EP[:eESR], eESRDischarge) end end diff --git a/src/model/core/fuel.jl b/src/model/core/fuel.jl index 30c4f73ef2..e7665c65a7 100644 --- a/src/model/core/fuel.jl +++ b/src/model/core/fuel.jl @@ -114,11 +114,11 @@ function fuel!(EP::Model, inputs::Dict, setup::Dict) COFIRE_MAX = [findall(g -> max_cofire_cols(g, tag = i) < 1, gen[MULTI_FUELS]) for i in 1:max_fuels] COFIRE_MAX_START = [findall(g -> max_cofire_start_cols(g, tag = i) < 1, - gen[MULTI_FUELS]) for i in 1:max_fuels] + gen[MULTI_FUELS]) for i in 1:max_fuels] COFIRE_MIN = [findall(g -> min_cofire_cols(g, tag = i) > 0, gen[MULTI_FUELS]) for i in 1:max_fuels] COFIRE_MIN_START = [findall(g -> min_cofire_start_cols(g, tag = i) > 0, - gen[MULTI_FUELS]) for i in 1:max_fuels] + gen[MULTI_FUELS]) for i in 1:max_fuels] @variable(EP, vMulFuels[y in MULTI_FUELS, i = 1:max_fuels, t = 1:T]>=0) @variable(EP, vMulStartFuels[y in MULTI_FUELS, i = 1:max_fuels, t = 1:T]>=0) @@ -226,14 +226,15 @@ function fuel!(EP::Model, inputs::Dict, setup::Dict) if !isempty(MULTI_FUELS) @expression(EP, eFuelConsumption_multi[f in 1:NUM_FUEL, t in 1:T], sum((EP[:vMulFuels][y, i, t] + EP[:vMulStartFuels][y, i, t]) #i: fuel id - for i in 1:max_fuels, - y in intersect(resource_id.(gen[fuel_cols.(gen, tag = i) .== string(fuels[f])]), + for i in 1:max_fuels, + y in intersect( + resource_id.(gen[fuel_cols.(gen, tag = i) .== string(fuels[f])]), MULTI_FUELS))) end @expression(EP, eFuelConsumption_single[f in 1:NUM_FUEL, t in 1:T], sum(EP[:vFuel][y, t] + EP[:eStartFuel][y, t] - for y in intersect(resources_with_fuel(gen, fuels[f]), SINGLE_FUEL))) + for y in intersect(resources_with_fuel(gen, fuels[f]), SINGLE_FUEL))) @expression(EP, eFuelConsumption[f in 1:NUM_FUEL, t in 1:T], if !isempty(MULTI_FUELS) @@ -249,13 +250,15 @@ function fuel!(EP::Model, inputs::Dict, setup::Dict) ### only apply constraint to generators with fuel type other than None @constraint(EP, - cFuelCalculation_single[y in intersect(SINGLE_FUEL, setdiff(HAS_FUEL, THERM_COMMIT)), + cFuelCalculation_single[ + y in intersect(SINGLE_FUEL, setdiff(HAS_FUEL, THERM_COMMIT)), t = 1:T], EP[:vFuel][y, t] - EP[:vP][y, t] * heat_rate_mmbtu_per_mwh(gen[y])==0) if !isempty(MULTI_FUELS) @constraint(EP, - cFuelCalculation_multi[y in intersect(MULTI_FUELS, + cFuelCalculation_multi[ + y in intersect(MULTI_FUELS, setdiff(HAS_FUEL, THERM_COMMIT)), t = 1:T], sum(EP[:vMulFuels][y, i, t] / heat_rates[i][y] for i in 1:max_fuels) - @@ -278,19 +281,21 @@ function fuel!(EP::Model, inputs::Dict, setup::Dict) PiecewiseFuelUsage[y in THERM_COMMIT_PWFU, t = 1:T, seg in segs], EP[:vFuel][y, t]>=(EP[:vP][y, t] * segment_slope(y, seg) + - EP[:vCOMMIT][y, t] * segment_intercept(y, seg))) + EP[:vCOMMIT][y, t] * segment_intercept(y, seg))) end # constraint for fuel consumption at a constant heat rate @constraint(EP, - FuelCalculationCommit_single[y in intersect(setdiff(THERM_COMMIT, + FuelCalculationCommit_single[ + y in intersect(setdiff(THERM_COMMIT, THERM_COMMIT_PWFU), SINGLE_FUEL), t = 1:T], EP[:vFuel][y, t] - EP[:vP][y, t] * heat_rate_mmbtu_per_mwh(gen[y])==0) if !isempty(MULTI_FUELS) @constraint(EP, - FuelCalculationCommit_multi[y in intersect(setdiff(THERM_COMMIT, + FuelCalculationCommit_multi[ + y in intersect(setdiff(THERM_COMMIT, THERM_COMMIT_PWFU), MULTI_FUELS), t = 1:T], diff --git a/src/model/core/non_served_energy.jl b/src/model/core/non_served_energy.jl index 0686a92cab..458ac6348a 100644 --- a/src/model/core/non_served_energy.jl +++ b/src/model/core/non_served_energy.jl @@ -93,7 +93,7 @@ function non_served_energy!(EP::Model, inputs::Dict, setup::Dict) @expression(EP, eCapResMarBalanceNSE[res = 1:inputs["NCapacityReserveMargin"], t = 1:T], sum(EP[:vNSE][s, t, z] - for s in 2:SEG, z in findall(x -> x != 0, inputs["dfCapRes"][:, res]))) + for s in 2:SEG, z in findall(x -> x != 0, inputs["dfCapRes"][:, res]))) add_similar_to_expression!(EP[:eCapResMarBalance], eCapResMarBalanceNSE) end end diff --git a/src/model/core/operational_reserves.jl b/src/model/core/operational_reserves.jl index c77c8363a4..9616cf96b0 100644 --- a/src/model/core/operational_reserves.jl +++ b/src/model/core/operational_reserves.jl @@ -46,9 +46,9 @@ where $\epsilon^{contingency}$ is static contingency requirement in MWs. Option 2 (dynamic capacity-based contingency) is expressed by the following constraints: ```math \begin{aligned} - &Contingency \geq \Omega^{size}_{y,z} \times \alpha^{Contingency,Aux}_{y,z} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ - &\alpha^{Contingency,Aux}_{y,z} \leq \Delta^{\text{total}}_{y,z} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ - &\alpha^{Contingency,Aux}_{y,z} \geq M_y \times \Delta^{\text{total}}_{y,z} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ + & Contingency \geq \Omega^{size}_{y,z} \times \alpha^{Contingency,Aux}_{y,z} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ + & \alpha^{Contingency,Aux}_{y,z} \leq \Delta^{\text{total}}_{y,z} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ + & \Delta^{\text{total}}_{y,z} \leq M_y \times \alpha^{Contingency,Aux}_{y,z} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ \end{aligned} ``` @@ -60,7 +60,7 @@ Option 3 (dynamic commitment-based contingency) is expressed by the following se \begin{aligned} & Contingency \geq \Omega^{size}_{y,z} \times Contingency\_Aux_{y,z,t} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ & Contingency\_Aux_{y,z,t} \leq \nu_{y,z,t} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ - & Contingency\_Aux_{y,z,t} \geq M_y \times \nu_{y,z,t} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ + & \nu_{y,z,t} \leq M_y \times Contingency\_Aux_{y,z,t} & \forall y \in \mathcal{UC}, z \in \mathcal{Z}\\ \end{aligned} ``` @@ -227,7 +227,8 @@ function operational_reserves_core!(EP::Model, inputs::Dict, setup::Dict) systemwide_hourly_demand = sum(pDemand, dims = 2) function must_run_vre_generation(t) - sum(pP_Max(y, t) * EP[:eTotalCap][y] + sum( + pP_Max(y, t) * EP[:eTotalCap][y] for y in intersect(inputs["VRE"], inputs["MUST_RUN"]); init = 0) end diff --git a/src/model/core/transmission/dcopf_transmission.jl b/src/model/core/transmission/dcopf_transmission.jl index c833c532f2..de2dcfd5bf 100644 --- a/src/model/core/transmission/dcopf_transmission.jl +++ b/src/model/core/transmission/dcopf_transmission.jl @@ -41,7 +41,7 @@ function dcopf_transmission!(EP::Model, inputs::Dict, setup::Dict) cPOWER_FLOW_OPF[l = 1:L, t = 1:T], EP[:vFLOW][l, t]==inputs["pDC_OPF_coeff"][l] * - sum(inputs["pNet_Map"][l, z] * vANGLE[z, t] for z in 1:Z)) + sum(inputs["pNet_Map"][l, z] * vANGLE[z, t] for z in 1:Z)) # Bus angle limits (except slack bus) @constraints(EP, diff --git a/src/model/core/transmission/investment_transmission.jl b/src/model/core/transmission/investment_transmission.jl index 804588e22c..a61b628f91 100644 --- a/src/model/core/transmission/investment_transmission.jl +++ b/src/model/core/transmission/investment_transmission.jl @@ -73,7 +73,7 @@ function investment_transmission!(EP::Model, inputs::Dict, setup::Dict) @expression(EP, eTotalCNetworkExp, sum(vNEW_TRANS_CAP[l] * inputs["pC_Line_Reinforcement"][l] - for l in EXPANSION_LINES)) + for l in EXPANSION_LINES)) if MultiStage == 1 # OPEX multiplier to count multiple years between two model stages diff --git a/src/model/core/transmission/transmission.jl b/src/model/core/transmission/transmission.jl index 342b2d7610..ca5c637159 100644 --- a/src/model/core/transmission/transmission.jl +++ b/src/model/core/transmission/transmission.jl @@ -312,7 +312,7 @@ function transmission!(EP::Model, inputs::Dict, setup::Dict) @expression(EP, eESRTran[ESR = 1:inputs["nESR"]], sum(inputs["dfESR"][z, ESR] * sum(inputs["omega"][t] * EP[:eLosses_By_Zone][z, t] for t in 1:T) - for z in findall(x -> x > 0, inputs["dfESR"][:, ESR]))) + for z in findall(x -> x > 0, inputs["dfESR"][:, ESR]))) add_similar_to_expression!(EP[:eESR], -eESRTran) end end diff --git a/src/model/expression_manipulation.jl b/src/model/expression_manipulation.jl index 7e783eae55..33b0b8e8e2 100644 --- a/src/model/expression_manipulation.jl +++ b/src/model/expression_manipulation.jl @@ -26,8 +26,8 @@ This can lead to errors later if a method can only operate on expressions. We don't currently have a method to do this with non-contiguous indexing. """ function create_empty_expression!(EP::Model, - exprname::Symbol, - dims::NTuple{N, Int64}) where {N} + exprname::Symbol, + dims::NTuple{N, Int64}) where {N} temp = Array{AffExpr}(undef, dims) fill_with_zeros!(temp) EP[exprname] = temp @@ -67,7 +67,7 @@ In the future we could expand this to non AffExpr, using GenericAffExpr e.g. if we wanted to use Float32 instead of Float64 """ function fill_with_const!(arr::AbstractArray{GenericAffExpr{C, T}, dims}, - con::Real) where {C, T, dims} + con::Real) where {C, T, dims} for i in eachindex(arr) arr[i] = AffExpr(con) end @@ -80,7 +80,7 @@ end ###### ###### ###### ###### ###### ###### # function extract_time_series_to_expression(var::Matrix{VariableRef}, - set::AbstractVector{Int}) + set::AbstractVector{Int}) TIME_DIM = 2 time_range = 1:size(var)[TIME_DIM] @@ -90,13 +90,14 @@ function extract_time_series_to_expression(var::Matrix{VariableRef}, return expr end -function extract_time_series_to_expression(var::JuMP.Containers.DenseAxisArray{ - VariableRef, - 2, - Tuple{X, Base.OneTo{Int64}}, - Y, - }, - set::AbstractVector{Int}) where {X, Y} +function extract_time_series_to_expression( + var::JuMP.Containers.DenseAxisArray{ + VariableRef, + 2, + Tuple{X, Base.OneTo{Int64}}, + Y + }, + set::AbstractVector{Int}) where {X, Y} TIME_DIM = 2 time_range = var.axes[TIME_DIM] @@ -125,7 +126,7 @@ This will work on JuMP DenseContainers which do not have linear indexing from 1: However, the accessed parts of both arrays must have the same dimensions. """ function add_similar_to_expression!(expr1::AbstractArray{GenericAffExpr{C, T}, dim1}, - expr2::AbstractArray{V, dim2}) where {C, T, V, dim1, dim2} + expr2::AbstractArray{V, dim2}) where {C, T, V, dim1, dim2} # This is defined for Arrays of different dimensions # despite the fact it will definitely throw an error # because the error will tell the user / developer @@ -155,7 +156,7 @@ Add an entry of type `V` to an array of expressions, in-place. This will work on JuMP DenseContainers which do not have linear indexing from 1:length(arr). """ function add_term_to_expression!(expr1::AbstractArray{GenericAffExpr{C, T}, dims}, - expr2::V) where {C, T, V, dims} + expr2::V) where {C, T, V, dims} for i in eachindex(expr1) add_to_expression!(expr1[i], expr2) end @@ -173,7 +174,7 @@ Check that two arrays have the same dimensions. If not, return an error message which includes the dimensions of both arrays. """ function check_sizes_match(expr1::AbstractArray{C, dim1}, - expr2::AbstractArray{T, dim2}) where {C, T, dim1, dim2} + expr2::AbstractArray{T, dim2}) where {C, T, dim1, dim2} # After testing, this appears to be just as fast as a method for Array{GenericAffExpr{C,T}, dims} or Array{AffExpr, dims} if size(expr1) != size(expr2) error(" diff --git a/src/model/policies/cap_reserve_margin.jl b/src/model/policies/cap_reserve_margin.jl index 74052fabd4..bdbd077ce9 100755 --- a/src/model/policies/cap_reserve_margin.jl +++ b/src/model/policies/cap_reserve_margin.jl @@ -82,5 +82,5 @@ function cap_reserve_margin!(EP::Model, inputs::Dict, setup::Dict) EP[:eCapResMarBalance][res, t] >=sum(inputs["pD"][t, z] * (1 + inputs["dfCapRes"][z, res]) - for z in findall(x -> x != 0, inputs["dfCapRes"][:, res]))) + for z in findall(x -> x != 0, inputs["dfCapRes"][:, res]))) end diff --git a/src/model/policies/co2_cap.jl b/src/model/policies/co2_cap.jl index d14b69a265..4c8badbfe8 100644 --- a/src/model/policies/co2_cap.jl +++ b/src/model/policies/co2_cap.jl @@ -92,33 +92,33 @@ function co2_cap!(EP::Model, inputs::Dict, setup::Dict) if setup["CO2Cap"] == 1 @constraint(EP, cCO2Emissions_systemwide[cap = 1:inputs["NCO2Cap"]], sum(inputs["omega"][t] * EP[:eEmissionsByZone][z, t] - for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]), t in 1:T) - + for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]), t in 1:T) - vCO2Cap_slack[cap]<= sum(inputs["dfMaxCO2"][z, cap] - for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]))) + for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]))) ## (fulfilled) demand + Rate-based: Emissions constraint in terms of rate (tons/MWh) elseif setup["CO2Cap"] == 2 ##This part moved to non_served_energy.jl @constraint(EP, cCO2Emissions_systemwide[cap = 1:inputs["NCO2Cap"]], sum(inputs["omega"][t] * EP[:eEmissionsByZone][z, t] - for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]), t in 1:T) - + for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]), t in 1:T) - vCO2Cap_slack[cap]<= sum(inputs["dfMaxCO2Rate"][z, cap] * sum(inputs["omega"][t] * (inputs["pD"][t, z] - sum(EP[:vNSE][s, t, z] for s in 1:SEG)) - for t in 1:T) - for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap])) + + for t in 1:T) + for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap])) + sum(inputs["dfMaxCO2Rate"][z, cap] * setup["StorageLosses"] * EP[:eELOSSByZone][z] - for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]))) + for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]))) ## Generation + Rate-based: Emissions constraint in terms of rate (tons/MWh) elseif (setup["CO2Cap"] == 3) @constraint(EP, cCO2Emissions_systemwide[cap = 1:inputs["NCO2Cap"]], sum(inputs["omega"][t] * EP[:eEmissionsByZone][z, t] - for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]), t in 1:T) - + for z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]), t in 1:T) - vCO2Cap_slack[cap]<= sum(inputs["dfMaxCO2Rate"][z, cap] * inputs["omega"][t] * EP[:eGenerationByZone][z, t] - for t in 1:T, z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]))) + for t in 1:T, z in findall(x -> x == 1, inputs["dfCO2CapZones"][:, cap]))) end end diff --git a/src/model/policies/energy_share_requirement.jl b/src/model/policies/energy_share_requirement.jl index a4e3225c08..4ed3728777 100644 --- a/src/model/policies/energy_share_requirement.jl +++ b/src/model/policies/energy_share_requirement.jl @@ -4,7 +4,7 @@ This function establishes constraints that can be flexibily applied to define al These policies usually require that the annual MWh generation from a subset of qualifying generators has to be higher than a pre-specified percentage of demand from qualifying zones. The implementation allows for user to define one or multiple RPS/CES style minimum energy share constraints, where each constraint can cover different combination of model zones to mimic real-world policy implementation (e.g. multiple state policies, multiple RPS tiers or overlapping RPS and CES policies). -The number of energy share requirement constraints is specified by the user by the value of the GenX settings parameter ```EnergyShareRequirement``` (this value should be an integer >=0). +Including an energy share requirement constraint is specified by the user by the value of the GenX settings parameter ```EnergyShareRequirement``` (this value should either 0 or 1). For each constraint $p \in \mathcal{P}^{ESR}$, we define a subset of zones $z \in \mathcal{Z}^{ESR}_{p} \subset \mathcal{Z}$ that are eligible for trading renewable/clean energy credits to meet the corresponding renewable/clean energy requirement. For each energy share requirement constraint $p \in \mathcal{P}^{ESR}$, we specify the share of total demand in each eligible model zone, diff --git a/src/model/resources/curtailable_variable_renewable/curtailable_variable_renewable.jl b/src/model/resources/curtailable_variable_renewable/curtailable_variable_renewable.jl index 9b790789d6..779cce88e8 100644 --- a/src/model/resources/curtailable_variable_renewable/curtailable_variable_renewable.jl +++ b/src/model/resources/curtailable_variable_renewable/curtailable_variable_renewable.jl @@ -68,8 +68,8 @@ function curtailable_variable_renewable!(EP::Model, inputs::Dict, setup::Dict) # Note: inequality constraint allows curtailment of output below maximum level. @constraint(EP, [t = 1:T], - EP[:vP][y, - t]<=sum(inputs["pP_Max"][yy, t] * EP[:eTotalCap][yy] for yy in VRE_BINS)) + EP[:vP][y,t]<=sum(inputs["pP_Max"][yy, t] * EP[:eTotalCap][yy] + for yy in VRE_BINS)) end end @@ -80,7 +80,7 @@ function curtailable_variable_renewable!(EP::Model, inputs::Dict, setup::Dict) ##CO2 Polcy Module VRE Generation by zone @expression(EP, eGenerationByVRE[z = 1:Z, t = 1:T], # the unit is GW sum(EP[:vP][y, t] - for y in intersect(inputs["VRE"], resources_in_zone_by_rid(gen, z)))) + for y in intersect(inputs["VRE"], resources_in_zone_by_rid(gen, z)))) add_similar_to_expression!(EP[:eGenerationByZone], eGenerationByVRE) end diff --git a/src/model/resources/flexible_demand/flexible_demand.jl b/src/model/resources/flexible_demand/flexible_demand.jl index 37c459f58c..56e1487c7b 100644 --- a/src/model/resources/flexible_demand/flexible_demand.jl +++ b/src/model/resources/flexible_demand/flexible_demand.jl @@ -67,7 +67,7 @@ function flexible_demand!(EP::Model, inputs::Dict, setup::Dict) ## Power Balance Expressions ## @expression(EP, ePowerBalanceDemandFlex[t = 1:T, z = 1:Z], sum(-EP[:vP][y, t] + EP[:vCHARGE_FLEX][y, t] - for y in intersect(FLEX, resources_in_zone_by_rid(gen, z)))) + for y in intersect(FLEX, resources_in_zone_by_rid(gen, z)))) add_similar_to_expression!(EP[:ePowerBalance], ePowerBalanceDemandFlex) # Capacity Reserves Margin policy @@ -127,14 +127,12 @@ function flexible_demand!(EP::Model, inputs::Dict, setup::Dict) @constraint(EP, [t in 1:T], # cFlexibleDemandDelay: Constraints looks forward over next n hours, where n = max_flexible_demand_delay sum(EP[:vP][y, e] - for e in hoursafter(hours_per_subperiod, t, 1:max_flex_demand_delay))>=EP[:vS_FLEX][y, - t]) + for e in hoursafter(hours_per_subperiod, t, 1:max_flex_demand_delay))>=EP[:vS_FLEX][y,t]) @constraint(EP, [t in 1:T], # cFlexibleDemandAdvance: Constraint looks forward over next n hours, where n = max_flexible_demand_advance sum(EP[:vCHARGE_FLEX][y, e] - for e in hoursafter(hours_per_subperiod, t, 1:max_flex_demand_advance))>=-EP[:vS_FLEX][y, - t]) + for e in hoursafter(hours_per_subperiod, t, 1:max_flex_demand_advance))>=-EP[:vS_FLEX][y,t]) end end diff --git a/src/model/resources/hydro/hydro_inter_period_linkage.jl b/src/model/resources/hydro/hydro_inter_period_linkage.jl index 5fdab6287d..c1812173b4 100644 --- a/src/model/resources/hydro/hydro_inter_period_linkage.jl +++ b/src/model/resources/hydro/hydro_inter_period_linkage.jl @@ -81,17 +81,22 @@ function hydro_inter_period_linkage!(EP::Model, inputs::Dict) cHydroReservoirLongDurationStorageStart[w = 1:REP_PERIOD, y in STOR_HYDRO_LONG_DURATION], EP[:vS_HYDRO][y, - hours_per_subperiod * (w - 1) + 1]==(EP[:vS_HYDRO][y, hours_per_subperiod * w] - vdSOC_HYDRO[y, w]) - - (1 / efficiency_down(gen[y]) * EP[:vP][y, hours_per_subperiod * (w - 1) + 1]) - - EP[:vSPILL][y, hours_per_subperiod * (w - 1) + 1] + - inputs["pP_Max"][y, hours_per_subperiod * (w - 1) + 1] * EP[:eTotalCap][y]) + hours_per_subperiod * (w - 1) + 1]==(EP[:vS_HYDRO][y, hours_per_subperiod * w] - + vdSOC_HYDRO[y, w]) - + (1 / efficiency_down(gen[y]) * EP[:vP][ + y, hours_per_subperiod * (w - 1) + 1]) - + EP[:vSPILL][ + y, hours_per_subperiod * (w - 1) + 1] + + inputs["pP_Max"][ + y, hours_per_subperiod * (w - 1) + 1] * EP[:eTotalCap][y]) # Storage at beginning of period w = storage at beginning of period w-1 + storage built up in period w (after n representative periods) ## Multiply storage build up term from prior period with corresponding weight @constraint(EP, cHydroReservoirLongDurationStorage[y in STOR_HYDRO_LONG_DURATION, r in MODELED_PERIODS_INDEX], vSOC_HYDROw[y, - mod1(r + 1, NPeriods)]==vSOC_HYDROw[y, r] + vdSOC_HYDRO[y, dfPeriodMap[r, :Rep_Period_Index]]) + mod1(r + 1, NPeriods)]==vSOC_HYDROw[y, r] + + vdSOC_HYDRO[y, dfPeriodMap[r, :Rep_Period_Index]]) # Storage at beginning of each modeled period cannot exceed installed energy capacity @constraint(EP, @@ -104,7 +109,6 @@ function hydro_inter_period_linkage!(EP::Model, inputs::Dict) @constraint(EP, cHydroReservoirLongDurationStorageSub[y in STOR_HYDRO_LONG_DURATION, r in REP_PERIODS_INDEX], - vSOC_HYDROw[y, - r]==EP[:vS_HYDRO][y, hours_per_subperiod * dfPeriodMap[r, :Rep_Period_Index]] - - vdSOC_HYDRO[y, dfPeriodMap[r, :Rep_Period_Index]]) + vSOC_HYDROw[y,r]==EP[:vS_HYDRO][y, hours_per_subperiod * dfPeriodMap[r, :Rep_Period_Index]] - + vdSOC_HYDRO[y, dfPeriodMap[r, :Rep_Period_Index]]) end diff --git a/src/model/resources/hydro/hydro_res.jl b/src/model/resources/hydro/hydro_res.jl index ce9b2c69f5..f4b90df76f 100644 --- a/src/model/resources/hydro/hydro_res.jl +++ b/src/model/resources/hydro/hydro_res.jl @@ -127,8 +127,8 @@ function hydro_res!(EP::Model, inputs::Dict, setup::Dict) cHydroReservoirStart[y in CONSTRAINTSET, t in START_SUBPERIODS], EP[:vS_HYDRO][y, t]==EP[:vS_HYDRO][y, hoursbefore(p, t, 1)] - - (1 / efficiency_down(gen[y]) * EP[:vP][y, t]) - vSPILL[y, t] + - inputs["pP_Max"][y, t] * EP[:eTotalCap][y]) + (1 / efficiency_down(gen[y]) * EP[:vP][y, t]) - vSPILL[y, t] + + inputs["pP_Max"][y, t] * EP[:eTotalCap][y]) ### Constraints commmon to all reservoir hydro (y in set HYDRO_RES) ### @constraints(EP, diff --git a/src/model/resources/hydrogen/electrolyzer.jl b/src/model/resources/hydrogen/electrolyzer.jl index bfd21d505d..08945a62b7 100644 --- a/src/model/resources/hydrogen/electrolyzer.jl +++ b/src/model/resources/hydrogen/electrolyzer.jl @@ -101,7 +101,7 @@ function electrolyzer!(EP::Model, inputs::Dict, setup::Dict) @expression(EP, ePowerBalanceElectrolyzers[t in 1:T, z in 1:Z], sum(EP[:vUSE][y, t] - for y in intersect(ELECTROLYZERS, resources_in_zone_by_rid(gen, z)))) + for y in intersect(ELECTROLYZERS, resources_in_zone_by_rid(gen, z)))) # Electrolyzers consume electricity so their vUSE is subtracted from power balance EP[:ePowerBalance] -= ePowerBalanceElectrolyzers @@ -145,7 +145,7 @@ function electrolyzer!(EP::Model, inputs::Dict, setup::Dict) @constraint(EP, cHydrogenMin[y in ELECTROLYZERS], sum(inputs["omega"][t] * EP[:vUSE][y, t] / hydrogen_mwh_per_tonne(gen[y]) - for t in 1:T)>=electrolyzer_min_kt(gen[y]) * kt_to_t) + for t in 1:T)>=electrolyzer_min_kt(gen[y]) * kt_to_t) ### Remove vP (electrolyzers do not produce power so vP = 0 for all periods) @constraints(EP, begin @@ -161,14 +161,9 @@ function electrolyzer!(EP::Model, inputs::Dict, setup::Dict) QUALIFIED_SUPPLY = ids_with(gen, qualified_hydrogen_supply) @constraint(EP, cHourlyMatching[z in HYDROGEN_ZONES, t in 1:T], sum(EP[:vP][y, t] - for y in intersect(resources_in_zone_by_rid(gen, z), QUALIFIED_SUPPLY))>=sum(EP[:vUSE][y, - t] for y in intersect(resources_in_zone_by_rid(gen, - z), - ELECTROLYZERS)) + sum(EP[:vCHARGE][y, - t] for y in intersect(resources_in_zone_by_rid(gen, - z), - QUALIFIED_SUPPLY, - STORAGE))) + for y in intersect(resources_in_zone_by_rid(gen, z), QUALIFIED_SUPPLY))>=sum(EP[:vUSE][y,t] + for y in intersect(resources_in_zone_by_rid(gen,z), ELECTROLYZERS)) + sum(EP[:vCHARGE][y,t] + for y in intersect(resources_in_zone_by_rid(gen,z), QUALIFIED_SUPPLY, STORAGE))) end ### Energy Share Requirement Policy ### @@ -179,7 +174,7 @@ function electrolyzer!(EP::Model, inputs::Dict, setup::Dict) @expression(EP, eElectrolyzerESR[ESR in 1:inputs["nESR"]], sum(inputs["omega"][t] * EP[:vUSE][y, t] - for y in intersect(ELECTROLYZERS, ids_with_policy(gen, esr, tag = ESR)), + for y in intersect(ELECTROLYZERS, ids_with_policy(gen, esr, tag = ESR)), t in 1:T)) EP[:eESR] -= eElectrolyzerESR end diff --git a/src/model/resources/maintenance.jl b/src/model/resources/maintenance.jl index 37d9c0ce82..1a9bda5f04 100644 --- a/src/model/resources/maintenance.jl +++ b/src/model/resources/maintenance.jl @@ -59,9 +59,9 @@ end maintenance_begin_hours: collection of hours in which maintenance is allowed to start """ function controlling_maintenance_start_hours(p::Int, - t::Int, - maintenance_duration::Int, - maintenance_begin_hours) + t::Int, + maintenance_duration::Int, + maintenance_begin_hours) controlled_hours = hoursbefore(p, t, 0:(maintenance_duration - 1)) return intersect(controlled_hours, maintenance_begin_hours) end @@ -102,16 +102,16 @@ end Adds constraints which act on the vCOMMIT-like variable. """ function maintenance_formulation!(EP::Model, - inputs::Dict, - resource_component::AbstractString, - r_id::Int, - maint_begin_cadence::Int, - maint_dur::Int, - maint_freq_years::Int, - cap::Float64, - vcommit::Symbol, - ecap::Symbol, - integer_operational_unit_commitment::Bool) + inputs::Dict, + resource_component::AbstractString, + r_id::Int, + maint_begin_cadence::Int, + maint_dur::Int, + maint_freq_years::Int, + cap::Float64, + vcommit::Symbol, + ecap::Symbol, + integer_operational_unit_commitment::Bool) T = 1:inputs["T"] hours_per_subperiod = inputs["hours_per_subperiod"] diff --git a/src/model/resources/resources.jl b/src/model/resources/resources.jl index 1e12375064..8eb07c98e8 100644 --- a/src/model/resources/resources.jl +++ b/src/model/resources/resources.jl @@ -350,7 +350,7 @@ function ids_with_positive(rs::Vector{T}, name::Symbol) where {T <: AbstractReso end function ids_with_positive(rs::Vector{T}, - name::AbstractString) where {T <: AbstractResource} + name::AbstractString) where {T <: AbstractResource} return ids_with_positive(rs, Symbol(lowercase(name))) end @@ -429,8 +429,8 @@ julia> existing_cap_mw(gen[21]) ``` """ function ids_with(rs::Vector{T}, - f::Function, - default = default_zero) where {T <: AbstractResource} + f::Function, + default = default_zero) where {T <: AbstractResource} return findall(r -> f(r) != default, rs) end @@ -460,16 +460,16 @@ julia> existing_cap_mw(gen[21]) ``` """ function ids_with(rs::Vector{T}, - name::Symbol, - default = default_zero) where {T <: AbstractResource} + name::Symbol, + default = default_zero) where {T <: AbstractResource} # if the getter function exists in GenX then use it, otherwise get the attribute directly f = isdefined(GenX, name) ? getfield(GenX, name) : r -> getproperty(r, name) return ids_with(rs, f, default) end function ids_with(rs::Vector{T}, - name::AbstractString, - default = default_zero) where {T <: AbstractResource} + name::AbstractString, + default = default_zero) where {T <: AbstractResource} return ids_with(rs, Symbol(lowercase(name)), default) end @@ -487,8 +487,8 @@ Function for finding resources in a vector `rs` where the policy specified by `f - `ids (Vector{Int64})`: The vector of resource ids with a positive value for policy `f` and tag `tag`. """ function ids_with_policy(rs::Vector{T}, - f::Function; - tag::Int64) where {T <: AbstractResource} + f::Function; + tag::Int64) where {T <: AbstractResource} return findall(r -> f(r, tag = tag) > 0, rs) end @@ -506,8 +506,8 @@ Function for finding resources in a vector `rs` where the policy specified by `n - `ids (Vector{Int64})`: The vector of resource ids with a positive value for policy `name` and tag `tag`. """ function ids_with_policy(rs::Vector{T}, - name::Symbol; - tag::Int64) where {T <: AbstractResource} + name::Symbol; + tag::Int64) where {T <: AbstractResource} # if the getter function exists in GenX then use it, otherwise get the attribute directly if isdefined(GenX, name) f = getfield(GenX, name) @@ -517,8 +517,8 @@ function ids_with_policy(rs::Vector{T}, end function ids_with_policy(rs::Vector{T}, - name::AbstractString; - tag::Int64) where {T <: AbstractResource} + name::AbstractString; + tag::Int64) where {T <: AbstractResource} return ids_with_policy(rs, Symbol(lowercase(name)), tag = tag) end @@ -744,9 +744,7 @@ function no_unit_commitment(rs::Vector{T}) where {T <: AbstractResource} end # Operational Reserves -function ids_with_regulation_reserve_requirements(rs::Vector{ - T, -}) where {T <: AbstractResource} +function ids_with_regulation_reserve_requirements(rs::Vector{T}) where {T <: AbstractResource} findall(r -> reg_max(r) > 0, rs) end function ids_with_spinning_reserve_requirements(rs::Vector{T}) where {T <: AbstractResource} @@ -826,7 +824,8 @@ vre(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, Vre), rs) Returns the indices of all electrolyzer resources in the vector `rs`. """ -electrolyzer(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, +electrolyzer(rs::Vector{T}) where {T <: AbstractResource} = findall( + r -> isa(r, Electrolyzer), rs) electrolyzer_min_kt(r::Electrolyzer) = r.electrolyzer_min_kt @@ -869,7 +868,8 @@ self_discharge(r::VreStorage) = r.self_disch Returns the indices of all co-located solar resources in the vector `rs`. """ -solar(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, VreStorage) && +solar(rs::Vector{T}) where {T <: AbstractResource} = findall( + r -> isa(r, VreStorage) && r.solar != 0, rs) @@ -878,7 +878,8 @@ solar(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, VreStor Returns the indices of all co-located wind resources in the vector `rs`. """ -wind(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, VreStorage) && +wind(rs::Vector{T}) where {T <: AbstractResource} = findall( + r -> isa(r, VreStorage) && r.wind != 0, rs) @@ -886,7 +887,8 @@ wind(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, VreStora storage_dc_discharge(rs::Vector{T}) where T <: AbstractResource Returns the indices of all co-located storage resources in the vector `rs` that discharge DC. """ -storage_dc_discharge(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, +storage_dc_discharge(rs::Vector{T}) where {T <: AbstractResource} = findall( + r -> isa(r, VreStorage) && r.stor_dc_discharge >= 1, rs) function storage_sym_dc_discharge(rs::Vector{T}) where {T <: AbstractResource} @@ -900,7 +902,8 @@ end storage_dc_charge(rs::Vector{T}) where T <: AbstractResource Returns the indices of all co-located storage resources in the vector `rs` that charge DC. """ -storage_dc_charge(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, +storage_dc_charge(rs::Vector{T}) where {T <: AbstractResource} = findall( + r -> isa(r, VreStorage) && r.stor_dc_charge >= 1, rs) function storage_sym_dc_charge(rs::Vector{T}) where {T <: AbstractResource} @@ -914,7 +917,8 @@ end storage_ac_discharge(rs::Vector{T}) where T <: AbstractResource Returns the indices of all co-located storage resources in the vector `rs` that discharge AC. """ -storage_ac_discharge(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, +storage_ac_discharge(rs::Vector{T}) where {T <: AbstractResource} = findall( + r -> isa(r, VreStorage) && r.stor_ac_discharge >= 1, rs) function storage_sym_ac_discharge(rs::Vector{T}) where {T <: AbstractResource} @@ -928,7 +932,8 @@ end storage_ac_charge(rs::Vector{T}) where T <: AbstractResource Returns the indices of all co-located storage resources in the vector `rs` that charge AC. """ -storage_ac_charge(rs::Vector{T}) where {T <: AbstractResource} = findall(r -> isa(r, +storage_ac_charge(rs::Vector{T}) where {T <: AbstractResource} = findall( + r -> isa(r, VreStorage) && r.stor_ac_charge >= 1, rs) function storage_sym_ac_charge(rs::Vector{T}) where {T <: AbstractResource} @@ -1085,7 +1090,7 @@ Find R_ID's of resources with retrofit cluster id `cluster_id`. - `Vector{Int64}`: The vector of resource ids in the retrofit cluster. """ function resources_in_retrofit_cluster_by_rid(rs::Vector{<:AbstractResource}, - cluster_id::String) + cluster_id::String) return resource_id.(rs[retrofit_id.(rs) .== cluster_id]) end @@ -1158,7 +1163,7 @@ Check if all retrofit options in the retrofit cluster of the retrofit resource ` - `Bool`: True if all retrofit options contribute to min retirement, otherwise false. """ function has_all_options_contributing(retrofit_res::AbstractResource, - rs::Vector{T}) where {T <: AbstractResource} + rs::Vector{T}) where {T <: AbstractResource} retro_id = retrofit_id(retrofit_res) return isempty(intersect(resources_in_retrofit_cluster_by_rid(rs, retro_id), ids_retrofit_options(rs), @@ -1198,7 +1203,7 @@ Check if all retrofit options in the retrofit cluster of the retrofit resource ` - `Bool`: True if all retrofit options do not contribute to min retirement, otherwise false. """ function has_all_options_not_contributing(retrofit_res::AbstractResource, - rs::Vector{T}) where {T <: AbstractResource} + rs::Vector{T}) where {T <: AbstractResource} retro_id = retrofit_id(retrofit_res) return isempty(intersect(resources_in_retrofit_cluster_by_rid(rs, retro_id), ids_retrofit_options(rs), diff --git a/src/model/resources/retrofits/retrofits.jl b/src/model/resources/retrofits/retrofits.jl index 6920dffd79..164cf59025 100644 --- a/src/model/resources/retrofits/retrofits.jl +++ b/src/model/resources/retrofits/retrofits.jl @@ -27,23 +27,31 @@ function retrofit(EP::Model, inputs::Dict) RETROFIT_IDS = inputs["RETROFIT_IDS"] # Set of unique IDs for retrofit resources @expression(EP, eRetrofittedCapByRetroId[id in RETROFIT_IDS], - sum(cap_size(gen[y]) * EP[:vRETROFITCAP][y] for y in intersect(RETROFIT_CAP, + sum( + cap_size(gen[y]) * EP[:vRETROFITCAP][y] + for y in intersect(RETROFIT_CAP, COMMIT, resources_in_retrofit_cluster_by_rid(gen, id)); init = 0) - +sum(EP[:vRETROFITCAP][y] for y in setdiff(intersect(RETROFIT_CAP, + +sum( + EP[:vRETROFITCAP][y] + for y in setdiff( + intersect(RETROFIT_CAP, resources_in_retrofit_cluster_by_rid(gen, id)), COMMIT); init = 0)) @expression(EP, eRetrofitCapByRetroId[id in RETROFIT_IDS], - sum(cap_size(gen[y]) * EP[:vCAP][y] * (1 / retrofit_efficiency(gen[y])) + sum( + cap_size(gen[y]) * EP[:vCAP][y] * (1 / retrofit_efficiency(gen[y])) for y in intersect(RETROFIT_OPTIONS, COMMIT, resources_in_retrofit_cluster_by_rid(gen, id)); init = 0) - +sum(EP[:vCAP][y] * (1 / retrofit_efficiency(gen[y])) - for y in setdiff(intersect(RETROFIT_OPTIONS, + +sum( + EP[:vCAP][y] * (1 / retrofit_efficiency(gen[y])) + for y in setdiff( + intersect(RETROFIT_OPTIONS, resources_in_retrofit_cluster_by_rid(gen, id)), COMMIT); init = 0)) diff --git a/src/model/resources/storage/long_duration_storage.jl b/src/model/resources/storage/long_duration_storage.jl index 1708332784..3730d3dff1 100644 --- a/src/model/resources/storage/long_duration_storage.jl +++ b/src/model/resources/storage/long_duration_storage.jl @@ -106,17 +106,21 @@ function long_duration_storage!(EP::Model, inputs::Dict, setup::Dict) cSoCBalLongDurationStorageStart[w = 1:REP_PERIOD, y in STOR_LONG_DURATION], EP[:vS][y, hours_per_subperiod * (w - 1) + 1]==(1 - self_discharge(gen[y])) * - (EP[:vS][y, hours_per_subperiod * w] - vdSOC[y, w]) - - - (1 / efficiency_down(gen[y]) * EP[:vP][y, hours_per_subperiod * (w - 1) + 1]) + - (efficiency_up(gen[y]) * EP[:vCHARGE][y, hours_per_subperiod * (w - 1) + 1])) + (EP[:vS][y, hours_per_subperiod * w] - + vdSOC[y, w]) + - + (1 / efficiency_down(gen[y]) * EP[:vP][ + y, hours_per_subperiod * (w - 1) + 1]) + + (efficiency_up(gen[y]) * EP[:vCHARGE][ + y, hours_per_subperiod * (w - 1) + 1])) # Storage at beginning of period w = storage at beginning of period w-1 + storage built up in period w (after n representative periods) ## Multiply storage build up term from prior period with corresponding weight @constraint(EP, cSoCBalLongDurationStorage[y in STOR_LONG_DURATION, r in MODELED_PERIODS_INDEX], vSOCw[y, - mod1(r + 1, NPeriods)]==vSOCw[y, r] + vdSOC[y, dfPeriodMap[r, :Rep_Period_Index]]) + mod1(r + 1, NPeriods)]==vSOCw[y, r] + + vdSOC[y, dfPeriodMap[r, :Rep_Period_Index]]) # Storage at beginning of each modeled period cannot exceed installed energy capacity @constraint(EP, @@ -130,7 +134,7 @@ function long_duration_storage!(EP::Model, inputs::Dict, setup::Dict) cSoCBalLongDurationStorageSub[y in STOR_LONG_DURATION, r in REP_PERIODS_INDEX], vSOCw[y, r]==EP[:vS][y, hours_per_subperiod * dfPeriodMap[r, :Rep_Period_Index]] - - vdSOC[y, dfPeriodMap[r, :Rep_Period_Index]]) + vdSOC[y, dfPeriodMap[r, :Rep_Period_Index]]) # Capacity Reserve Margin policy if CapacityReserveMargin > 0 @@ -144,12 +148,15 @@ function long_duration_storage!(EP::Model, inputs::Dict, setup::Dict) cVSoCBalLongDurationStorageStart[w = 1:REP_PERIOD, y in STOR_LONG_DURATION], EP[:vCAPRES_socinreserve][y, hours_per_subperiod * (w - 1) + 1]==(1 - self_discharge(gen[y])) * - (EP[:vCAPRES_socinreserve][y, hours_per_subperiod * w] - vCAPRES_dsoc[y, w]) - + - (1 / efficiency_down(gen[y]) * - EP[:vCAPRES_discharge][y, hours_per_subperiod * (w - 1) + 1]) - - (efficiency_up(gen[y]) * - EP[:vCAPRES_charge][y, hours_per_subperiod * (w - 1) + 1])) + (EP[:vCAPRES_socinreserve][ + y, hours_per_subperiod * w] - vCAPRES_dsoc[y, w]) + + + (1 / efficiency_down(gen[y]) * + EP[:vCAPRES_discharge][ + y, hours_per_subperiod * (w - 1) + 1]) - + (efficiency_up(gen[y]) * + EP[:vCAPRES_charge][ + y, hours_per_subperiod * (w - 1) + 1])) # Storage held in reserve at beginning of period w = storage at beginning of period w-1 + storage built up in period w (after n representative periods) ## Multiply storage build up term from prior period with corresponding weight @@ -157,15 +164,16 @@ function long_duration_storage!(EP::Model, inputs::Dict, setup::Dict) cVSoCBalLongDurationStorage[y in STOR_LONG_DURATION, r in MODELED_PERIODS_INDEX], vCAPRES_socw[y, - mod1(r + 1, NPeriods)]==vCAPRES_socw[y, r] + vCAPRES_dsoc[y, dfPeriodMap[r, :Rep_Period_Index]]) + mod1(r + 1, NPeriods)]==vCAPRES_socw[y, r] + + vCAPRES_dsoc[y, dfPeriodMap[r, :Rep_Period_Index]]) # Initial reserve storage level for representative periods must also adhere to sub-period storage inventory balance # Initial storage = Final storage - change in storage inventory across representative period @constraint(EP, cVSoCBalLongDurationStorageSub[y in STOR_LONG_DURATION, r in REP_PERIODS_INDEX], - vCAPRES_socw[y, - r]==EP[:vCAPRES_socinreserve][y, - hours_per_subperiod * dfPeriodMap[r, :Rep_Period_Index]] - vCAPRES_dsoc[y, dfPeriodMap[r, :Rep_Period_Index]]) + vCAPRES_socw[y,r]==EP[:vCAPRES_socinreserve][y, + hours_per_subperiod * dfPeriodMap[r, :Rep_Period_Index]] - + vCAPRES_dsoc[y, dfPeriodMap[r, :Rep_Period_Index]]) # energy held in reserve at the beginning of each modeled period acts as a lower bound on the total energy held in storage @constraint(EP, diff --git a/src/model/resources/storage/storage.jl b/src/model/resources/storage/storage.jl index ed6c0630ce..c4f9e6d90d 100644 --- a/src/model/resources/storage/storage.jl +++ b/src/model/resources/storage/storage.jl @@ -167,8 +167,8 @@ function storage!(EP::Model, inputs::Dict, setup::Dict) @expression(EP, eESRStor[ESR = 1:inputs["nESR"]], sum(inputs["dfESR"][z, ESR] * sum(EP[:eELOSS][y] - for y in intersect(resources_in_zone_by_rid(gen, z), STOR_ALL)) - for z in findall(x -> x > 0, inputs["dfESR"][:, ESR]))) + for y in intersect(resources_in_zone_by_rid(gen, z), STOR_ALL)) + for z in findall(x -> x > 0, inputs["dfESR"][:, ESR]))) add_similar_to_expression!(EP[:eESR], -eESRStor) end end @@ -178,14 +178,14 @@ function storage!(EP::Model, inputs::Dict, setup::Dict) @expression(EP, eCapResMarBalanceStor[res = 1:inputs["NCapacityReserveMargin"], t = 1:T], sum(derating_factor(gen[y], tag = res) * (EP[:vP][y, t] - EP[:vCHARGE][y, t]) - for y in STOR_ALL)) + for y in STOR_ALL)) if StorageVirtualDischarge > 0 @expression(EP, eCapResMarBalanceStorVirtual[res = 1:inputs["NCapacityReserveMargin"], t = 1:T], sum(derating_factor(gen[y], tag = res) * (EP[:vCAPRES_discharge][y, t] - EP[:vCAPRES_charge][y, t]) - for y in STOR_ALL)) + for y in STOR_ALL)) add_similar_to_expression!(eCapResMarBalanceStor, eCapResMarBalanceStorVirtual) end add_similar_to_expression!(EP[:eCapResMarBalance], eCapResMarBalanceStor) diff --git a/src/model/resources/storage/storage_all.jl b/src/model/resources/storage/storage_all.jl index cc002c78f5..84000f190c 100644 --- a/src/model/resources/storage/storage_all.jl +++ b/src/model/resources/storage/storage_all.jl @@ -49,9 +49,10 @@ function storage_all!(EP::Model, inputs::Dict, setup::Dict) # Energy losses related to technologies (increase in effective demand) @expression(EP, eELOSS[y in STOR_ALL], - sum(inputs["omega"][t] * EP[:vCHARGE][y, t] for t in 1:T)-sum(inputs["omega"][t] * - EP[:vP][y, t] - for t in 1:T)) + sum(inputs["omega"][t] * EP[:vCHARGE][y, t] + for t in 1:T)-sum(inputs["omega"][t] * + EP[:vP][y, t] + for t in 1:T)) ## Objective Function Expressions ## @@ -92,7 +93,7 @@ function storage_all!(EP::Model, inputs::Dict, setup::Dict) # Term to represent net dispatch from storage in any period @expression(EP, ePowerBalanceStor[t = 1:T, z = 1:Z], sum(EP[:vP][y, t] - EP[:vCHARGE][y, t] - for y in intersect(resources_in_zone_by_rid(gen, z), STOR_ALL))) + for y in intersect(resources_in_zone_by_rid(gen, z), STOR_ALL))) add_similar_to_expression!(EP[:ePowerBalance], ePowerBalanceStor) ### Constraints ### diff --git a/src/model/resources/thermal/thermal.jl b/src/model/resources/thermal/thermal.jl index ef5f9df385..2735a094a7 100644 --- a/src/model/resources/thermal/thermal.jl +++ b/src/model/resources/thermal/thermal.jl @@ -23,14 +23,15 @@ function thermal!(EP::Model, inputs::Dict, setup::Dict) ##CO2 Polcy Module Thermal Generation by zone @expression(EP, eGenerationByThermAll[z = 1:Z, t = 1:T], # the unit is GW sum(EP[:vP][y, t] - for y in intersect(inputs["THERM_ALL"], resources_in_zone_by_rid(gen, z)))) + for y in intersect(inputs["THERM_ALL"], resources_in_zone_by_rid(gen, z)))) add_similar_to_expression!(EP[:eGenerationByZone], eGenerationByThermAll) # Capacity Reserves Margin policy if setup["CapacityReserveMargin"] > 0 ncapres = inputs["NCapacityReserveMargin"] @expression(EP, eCapResMarBalanceThermal[capres in 1:ncapres, t in 1:T], - sum(derating_factor(gen[y], tag = capres) * EP[:eTotalCap][y] for y in THERM_ALL)) + sum(derating_factor(gen[y], tag = capres) * EP[:eTotalCap][y] + for y in THERM_ALL)) add_similar_to_expression!(EP[:eCapResMarBalance], eCapResMarBalanceThermal) MAINT = ids_with_maintenance(gen) diff --git a/src/model/resources/thermal/thermal_commit.jl b/src/model/resources/thermal/thermal_commit.jl index 347bdfbed5..ca7a74c7ba 100644 --- a/src/model/resources/thermal/thermal_commit.jl +++ b/src/model/resources/thermal/thermal_commit.jl @@ -154,7 +154,8 @@ function thermal_commit!(EP::Model, inputs::Dict, setup::Dict) ## Power Balance Expressions ## @expression(EP, ePowerBalanceThermCommit[t = 1:T, z = 1:Z], - sum(EP[:vP][y, t] for y in intersect(THERM_COMMIT, resources_in_zone_by_rid(gen, z)))) + sum(EP[:vP][y, t] + for y in intersect(THERM_COMMIT, resources_in_zone_by_rid(gen, z)))) add_similar_to_expression!(EP[:ePowerBalance], ePowerBalanceThermCommit) ### Constraints ### @@ -189,7 +190,8 @@ function thermal_commit!(EP::Model, inputs::Dict, setup::Dict) (EP[:vCOMMIT][y, t] - EP[:vSTART][y, t]) + min(inputs["pP_Max"][y, t], - max(min_power(gen[y]), ramp_up_fraction(gen[y]))) * cap_size(gen[y]) * EP[:vSTART][y, t] + max(min_power(gen[y]), ramp_up_fraction(gen[y]))) * + cap_size(gen[y]) * EP[:vSTART][y, t] - min_power(gen[y]) * cap_size(gen[y]) * EP[:vSHUT][y, t]) @@ -203,7 +205,8 @@ function thermal_commit!(EP::Model, inputs::Dict, setup::Dict) min_power(gen[y]) * cap_size(gen[y]) * EP[:vSTART][y, t] + min(inputs["pP_Max"][y, t], - max(min_power(gen[y]), ramp_down_fraction(gen[y]))) * cap_size(gen[y]) * EP[:vSHUT][y, t]) + max(min_power(gen[y]), ramp_down_fraction(gen[y]))) * + cap_size(gen[y]) * EP[:vSHUT][y, t]) ### Minimum and maximum power output constraints (Constraints #7-8) if setup["OperationalReserves"] == 1 @@ -376,7 +379,7 @@ end Eliminates the contribution of a plant to the capacity reserve margin while it is down for maintenance. """ function thermal_maintenance_capacity_reserve_margin_adjustment!(EP::Model, - inputs::Dict) + inputs::Dict) gen = inputs["RESOURCES"] T = inputs["T"] # Number of time steps (hours) @@ -387,18 +390,18 @@ function thermal_maintenance_capacity_reserve_margin_adjustment!(EP::Model, maint_adj = @expression(EP, [capres in 1:ncapres, t in 1:T], sum(thermal_maintenance_capacity_reserve_margin_adjustment(EP, - inputs, - y, - capres, - t) for y in applicable_resources)) + inputs, + y, + capres, + t) for y in applicable_resources)) add_similar_to_expression!(EP[:eCapResMarBalance], maint_adj) end function thermal_maintenance_capacity_reserve_margin_adjustment(EP::Model, - inputs::Dict, - y::Int, - capres::Int, - t) + inputs::Dict, + y::Int, + capres::Int, + t) gen = inputs["RESOURCES"] resource_component = resource_name(gen[y]) capresfactor = derating_factor(gen[y], tag = capres) diff --git a/src/model/resources/thermal/thermal_no_commit.jl b/src/model/resources/thermal/thermal_no_commit.jl index 1a75eb0980..dd1253e0bd 100644 --- a/src/model/resources/thermal/thermal_no_commit.jl +++ b/src/model/resources/thermal/thermal_no_commit.jl @@ -58,7 +58,7 @@ function thermal_no_commit!(EP::Model, inputs::Dict, setup::Dict) ## Power Balance Expressions ## @expression(EP, ePowerBalanceThermNoCommit[t = 1:T, z = 1:Z], sum(EP[:vP][y, t] - for y in intersect(THERM_NO_COMMIT, resources_in_zone_by_rid(gen, z)))) + for y in intersect(THERM_NO_COMMIT, resources_in_zone_by_rid(gen, z)))) add_similar_to_expression!(EP[:ePowerBalance], ePowerBalanceThermNoCommit) ### Constraints ### diff --git a/src/model/resources/vre_stor/vre_stor.jl b/src/model/resources/vre_stor/vre_stor.jl index 911fdca66b..0e1ca1709e 100644 --- a/src/model/resources/vre_stor/vre_stor.jl +++ b/src/model/resources/vre_stor/vre_stor.jl @@ -131,7 +131,7 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) for t in 1:T, z in 1:Z if !isempty(resources_in_zone_by_rid(gen_VRE_STOR, z)) ePowerBalance_VRE_STOR[t, z] += sum(EP[:vP][y, t] - for y in resources_in_zone_by_rid(gen_VRE_STOR, + for y in resources_in_zone_by_rid(gen_VRE_STOR, z)) end end @@ -173,17 +173,17 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) @expression(EP, eESRVREStor[ESR = 1:inputs["nESR"]], sum(inputs["omega"][t] * esr_vrestor(gen[y], tag = ESR) * EP[:vP_SOLAR][y, t] * by_rid(y, :etainverter) - for y in intersect(SOLAR, ids_with_policy(gen, esr_vrestor, tag = ESR)), + for y in intersect(SOLAR, ids_with_policy(gen, esr_vrestor, tag = ESR)), t in 1:T) +sum(inputs["omega"][t] * esr_vrestor(gen[y], tag = ESR) * EP[:vP_WIND][y, t] - for y in intersect(WIND, ids_with_policy(gen, esr_vrestor, tag = ESR)), + for y in intersect(WIND, ids_with_policy(gen, esr_vrestor, tag = ESR)), t in 1:T)) EP[:eESR] += eESRVREStor if IncludeLossesInESR == 1 @expression(EP, eESRVREStorLosses[ESR = 1:inputs["nESR"]], sum(inputs["dfESR"][z, ESR] * sum(EP[:eELOSS_VRE_STOR][y] - for y in intersect(STOR, resources_in_zone_by_rid(gen_VRE_STOR, z))) - for z in findall(x -> x > 0, inputs["dfESR"][:, ESR]))) + for y in intersect(STOR, resources_in_zone_by_rid(gen_VRE_STOR, z))) + for z in findall(x -> x > 0, inputs["dfESR"][:, ESR]))) EP[:eESR] -= eESRVREStorLosses end end @@ -191,19 +191,21 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) # Minimum Capacity Requirement if MinCapReq == 1 @expression(EP, eMinCapResSolar[mincap = 1:inputs["NumberOfMinCapReqs"]], - sum(by_rid(y, :etainverter) * EP[:eTotalCap_SOLAR][y] for y in intersect(SOLAR, + sum(by_rid(y, :etainverter) * EP[:eTotalCap_SOLAR][y] + for y in intersect(SOLAR, ids_with_policy(gen_VRE_STOR, min_cap_solar, tag = mincap)))) EP[:eMinCapRes] += eMinCapResSolar @expression(EP, eMinCapResWind[mincap = 1:inputs["NumberOfMinCapReqs"]], - sum(EP[:eTotalCap_WIND][y] for y in intersect(WIND, + sum(EP[:eTotalCap_WIND][y] + for y in intersect(WIND, ids_with_policy(gen_VRE_STOR, min_cap_wind, tag = mincap)))) EP[:eMinCapRes] += eMinCapResWind if !isempty(inputs["VS_ASYM_AC_DISCHARGE"]) @expression(EP, eMinCapResACDis[mincap = 1:inputs["NumberOfMinCapReqs"]], sum(EP[:eTotalCapDischarge_AC][y] - for y in intersect(inputs["VS_ASYM_AC_DISCHARGE"], + for y in intersect(inputs["VS_ASYM_AC_DISCHARGE"], ids_with_policy(gen_VRE_STOR, min_cap_stor, tag = mincap)))) EP[:eMinCapRes] += eMinCapResACDis end @@ -211,7 +213,7 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) if !isempty(inputs["VS_ASYM_DC_DISCHARGE"]) @expression(EP, eMinCapResDCDis[mincap = 1:inputs["NumberOfMinCapReqs"]], sum(EP[:eTotalCapDischarge_DC][y] - for y in intersect(inputs["VS_ASYM_DC_DISCHARGE"], + for y in intersect(inputs["VS_ASYM_DC_DISCHARGE"], ids_with_policy(gen_VRE_STOR, min_cap_stor, tag = mincap)))) EP[:eMinCapRes] += eMinCapResDCDis end @@ -219,7 +221,7 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) if !isempty(inputs["VS_SYM_AC"]) @expression(EP, eMinCapResACStor[mincap = 1:inputs["NumberOfMinCapReqs"]], sum(by_rid(y, :power_to_energy_ac) * EP[:eTotalCap_STOR][y] - for y in intersect(inputs["VS_SYM_AC"], + for y in intersect(inputs["VS_SYM_AC"], ids_with_policy(gen_VRE_STOR, min_cap_stor, tag = mincap)))) EP[:eMinCapRes] += eMinCapResACStor end @@ -227,7 +229,7 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) if !isempty(inputs["VS_SYM_DC"]) @expression(EP, eMinCapResDCStor[mincap = 1:inputs["NumberOfMinCapReqs"]], sum(by_rid(y, :power_to_energy_dc) * EP[:eTotalCap_STOR][y] - for y in intersect(inputs["VS_SYM_DC"], + for y in intersect(inputs["VS_SYM_DC"], ids_with_policy(gen_VRE_STOR, min_cap_stor, tag = mincap)))) EP[:eMinCapRes] += eMinCapResDCStor end @@ -236,19 +238,21 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) # Maximum Capacity Requirement if MaxCapReq == 1 @expression(EP, eMaxCapResSolar[maxcap = 1:inputs["NumberOfMaxCapReqs"]], - sum(by_rid(y, :etainverter) * EP[:eTotalCap_SOLAR][y] for y in intersect(SOLAR, + sum(by_rid(y, :etainverter) * EP[:eTotalCap_SOLAR][y] + for y in intersect(SOLAR, ids_with_policy(gen_VRE_STOR, max_cap_solar, tag = maxcap)))) EP[:eMaxCapRes] += eMaxCapResSolar @expression(EP, eMaxCapResWind[maxcap = 1:inputs["NumberOfMaxCapReqs"]], - sum(EP[:eTotalCap_WIND][y] for y in intersect(WIND, + sum(EP[:eTotalCap_WIND][y] + for y in intersect(WIND, ids_with_policy(gen_VRE_STOR, max_cap_wind, tag = maxcap)))) EP[:eMaxCapRes] += eMaxCapResWind if !isempty(inputs["VS_ASYM_AC_DISCHARGE"]) @expression(EP, eMaxCapResACDis[maxcap = 1:inputs["NumberOfMaxCapReqs"]], sum(EP[:eTotalCapDischarge_AC][y] - for y in intersect(inputs["VS_ASYM_AC_DISCHARGE"], + for y in intersect(inputs["VS_ASYM_AC_DISCHARGE"], ids_with_policy(gen_VRE_STOR, max_cap_stor, tag = maxcap)))) EP[:eMaxCapRes] += eMaxCapResACDis end @@ -256,7 +260,7 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) if !isempty(inputs["VS_ASYM_DC_DISCHARGE"]) @expression(EP, eMaxCapResDCDis[maxcap = 1:inputs["NumberOfMaxCapReqs"]], sum(EP[:eTotalCapDischarge_DC][y] - for y in intersect(inputs["VS_ASYM_DC_DISCHARGE"], + for y in intersect(inputs["VS_ASYM_DC_DISCHARGE"], ids_with_policy(gen_VRE_STOR, max_cap_stor, tag = maxcap)))) EP[:eMaxCapRes] += eMaxCapResDCDis end @@ -264,7 +268,7 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) if !isempty(inputs["VS_SYM_AC"]) @expression(EP, eMaxCapResACStor[maxcap = 1:inputs["NumberOfMaxCapReqs"]], sum(by_rid(y, :power_to_energy_ac) * EP[:eTotalCap_STOR][y] - for y in intersect(inputs["VS_SYM_AC"], + for y in intersect(inputs["VS_SYM_AC"], ids_with_policy(gen_VRE_STOR, max_cap_stor, tag = maxcap)))) EP[:eMaxCapRes] += eMaxCapResACStor end @@ -272,7 +276,7 @@ function vre_stor!(EP::Model, inputs::Dict, setup::Dict) if !isempty(inputs["VS_SYM_DC"]) @expression(EP, eMaxCapResDCStor[maxcap = 1:inputs["NumberOfMaxCapReqs"]], sum(by_rid(y, :power_to_energy_dc) * EP[:eTotalCap_STOR][y] - for y in intersect(inputs["VS_SYM_DC"], + for y in intersect(inputs["VS_SYM_DC"], ids_with_policy(gen_VRE_STOR, max_cap_stor, tag = maxcap)))) EP[:eMaxCapRes] += eMaxCapResDCStor end @@ -1124,7 +1128,8 @@ function stor_vre_stor!(EP::Model, inputs::Dict, setup::Dict) # SoC expressions @expression(EP, eSoCBalStart_VRE_STOR[y in CONSTRAINTSET, t in START_SUBPERIODS], vS_VRE_STOR[y, - t + hours_per_subperiod - 1]-self_discharge(gen[y]) * vS_VRE_STOR[y, t + hours_per_subperiod - 1]) + t + hours_per_subperiod - 1]-self_discharge(gen[y]) * + vS_VRE_STOR[y, t + hours_per_subperiod - 1]) @expression(EP, eSoCBalInterior_VRE_STOR[y in STOR, t in INTERIOR_SUBPERIODS], vS_VRE_STOR[y, t - 1]-self_discharge(gen[y]) * vS_VRE_STOR[y, t - 1]) # Expression for energy losses related to technologies (increase in effective demand) @@ -1179,7 +1184,7 @@ function stor_vre_stor!(EP::Model, inputs::Dict, setup::Dict) for y in AC_DISCHARGE EP[:eELOSS_VRE_STOR][y] -= sum(inputs["omega"][t] * vP_AC_DISCHARGE[y, t] - for t in 1:T) + for t in 1:T) for t in 1:T EP[:eInvACBalance][y, t] += vP_AC_DISCHARGE[y, t] end @@ -1208,7 +1213,7 @@ function stor_vre_stor!(EP::Model, inputs::Dict, setup::Dict) for z in 1:Z, t in 1:T if !isempty(resources_in_zone_by_rid(gen_VRE_STOR, z)) EP[:ePowerBalance_VRE_STOR][t, z] -= sum(vCHARGE_VRE_STOR[y, t] - for y in intersect(resources_in_zone_by_rid(gen_VRE_STOR, + for y in intersect(resources_in_zone_by_rid(gen_VRE_STOR, z), STOR)) end @@ -1219,7 +1224,7 @@ function stor_vre_stor!(EP::Model, inputs::Dict, setup::Dict) # From CO2 Policy module @expression(EP, eELOSSByZone_VRE_STOR[z = 1:Z], sum(EP[:eELOSS_VRE_STOR][y] - for y in intersect(resources_in_zone_by_rid(gen_VRE_STOR, z), STOR))) + for y in intersect(resources_in_zone_by_rid(gen_VRE_STOR, z), STOR))) add_similar_to_expression!(EP[:eELOSSByZone], eELOSSByZone_VRE_STOR) ### CONSTRAINTS ### @@ -1378,7 +1383,8 @@ function lds_vre_stor!(EP::Model, inputs::Dict) cVreStorSoCBalLongDurationStorage[y in VS_LDS, r in MODELED_PERIODS_INDEX], EP[:vSOCw_VRE_STOR][y, mod1(r + 1, NPeriods)]==EP[:vSOCw_VRE_STOR][y, r] + - EP[:vdSOC_VRE_STOR][y, dfPeriodMap[r, :Rep_Period_Index]]) + EP[:vdSOC_VRE_STOR][ + y, dfPeriodMap[r, :Rep_Period_Index]]) # Constraint 3: Storage at beginning of each modeled period cannot exceed installed energy capacity @constraint(EP, @@ -1389,10 +1395,10 @@ function lds_vre_stor!(EP::Model, inputs::Dict) # Initial storage = Final storage - change in storage inventory across representative period @constraint(EP, cVreStorSoCBalLongDurationStorageSub[y in VS_LDS, r in REP_PERIODS_INDEX], - EP[:vSOCw_VRE_STOR][y, - r]==EP[:vS_VRE_STOR][y, hours_per_subperiod * dfPeriodMap[r, :Rep_Period_Index]] - - - EP[:vdSOC_VRE_STOR][y, dfPeriodMap[r, :Rep_Period_Index]]) + EP[:vSOCw_VRE_STOR][y,r]==EP[:vS_VRE_STOR][ + y, hours_per_subperiod * dfPeriodMap[r, :Rep_Period_Index]] + - + EP[:vdSOC_VRE_STOR][y, dfPeriodMap[r, :Rep_Period_Index]]) end @doc raw""" @@ -1546,9 +1552,11 @@ function investment_charge_vre_stor!(EP::Model, inputs::Dict, setup::Dict) by_rid(rid, sym) = by_rid_res(rid, sym, gen_VRE_STOR) if !isempty(VS_ASYM_DC_DISCHARGE) - MAX_DC_DISCHARGE = intersect(ids_with_nonneg(gen_VRE_STOR, max_cap_discharge_dc_mw), + MAX_DC_DISCHARGE = intersect( + ids_with_nonneg(gen_VRE_STOR, max_cap_discharge_dc_mw), VS_ASYM_DC_DISCHARGE) - MIN_DC_DISCHARGE = intersect(ids_with_positive(gen_VRE_STOR, + MIN_DC_DISCHARGE = intersect( + ids_with_positive(gen_VRE_STOR, min_cap_discharge_dc_mw), VS_ASYM_DC_DISCHARGE) @@ -1742,9 +1750,11 @@ function investment_charge_vre_stor!(EP::Model, inputs::Dict, setup::Dict) end if !isempty(VS_ASYM_AC_DISCHARGE) - MAX_AC_DISCHARGE = intersect(ids_with_nonneg(gen_VRE_STOR, max_cap_discharge_ac_mw), + MAX_AC_DISCHARGE = intersect( + ids_with_nonneg(gen_VRE_STOR, max_cap_discharge_ac_mw), VS_ASYM_AC_DISCHARGE) - MIN_AC_DISCHARGE = intersect(ids_with_positive(gen_VRE_STOR, + MIN_AC_DISCHARGE = intersect( + ids_with_positive(gen_VRE_STOR, min_cap_discharge_ac_mw), VS_ASYM_AC_DISCHARGE) @@ -2186,7 +2196,7 @@ function vre_stor_capres!(EP::Model, inputs::Dict, setup::Dict) eCapResMarBalanceStor_VRE_STOR[res = 1:inputs["NCapacityReserveMargin"], t = 1:T], (sum(derating_factor(gen[y], tag = res) * by_rid(y, :etainverter) * inputs["pP_Max_Solar"][y, t] * EP[:eTotalCap_SOLAR][y] - for y in inputs["VS_SOLAR"]) + for y in inputs["VS_SOLAR"]) + sum(derating_factor(gen[y], tag = res) * inputs["pP_Max_Wind"][y, t] * EP[:eTotalCap_WIND][y] for y in inputs["VS_WIND"]) @@ -2195,26 +2205,29 @@ function vre_stor_capres!(EP::Model, inputs::Dict, setup::Dict) (EP[:vP_DC_DISCHARGE][y, t]) for y in DC_DISCHARGE) + sum(derating_factor(gen[y], tag = res) * (EP[:vP_AC_DISCHARGE][y, t]) - for y in AC_DISCHARGE) + for y in AC_DISCHARGE) - sum(derating_factor(gen[y], tag = res) * (EP[:vP_DC_CHARGE][y, t]) / - by_rid(y, :etainverter) for y in DC_CHARGE) + by_rid(y, :etainverter) + for y in DC_CHARGE) -sum(derating_factor(gen[y], tag = res) * (EP[:vP_AC_CHARGE][y, t]) - for y in AC_CHARGE))) + for y in AC_CHARGE))) if StorageVirtualDischarge > 0 @expression(EP, - eCapResMarBalanceStor_VRE_STOR_Virtual[res = 1:inputs["NCapacityReserveMargin"], + eCapResMarBalanceStor_VRE_STOR_Virtual[ + res = 1:inputs["NCapacityReserveMargin"], t = 1:T], (sum(derating_factor(gen[y], tag = res) * by_rid(y, :etainverter) * (vCAPRES_DC_DISCHARGE[y, t]) for y in DC_DISCHARGE) + sum(derating_factor(gen[y], tag = res) * (vCAPRES_AC_DISCHARGE[y, t]) - for y in AC_DISCHARGE) + for y in AC_DISCHARGE) - sum(derating_factor(gen[y], tag = res) * (vCAPRES_DC_CHARGE[y, t]) / - by_rid(y, :etainverter) for y in DC_CHARGE) + by_rid(y, :etainverter) + for y in DC_CHARGE) -sum(derating_factor(gen[y], tag = res) * (vCAPRES_AC_CHARGE[y, t]) - for y in AC_CHARGE))) + for y in AC_CHARGE))) add_similar_to_expression!(eCapResMarBalanceStor_VRE_STOR, eCapResMarBalanceStor_VRE_STOR_Virtual) end @@ -2304,7 +2317,8 @@ function vre_stor_capres!(EP::Model, inputs::Dict, setup::Dict) AC_CHARGE_CONSTRAINTSET = intersect(AC_CHARGE, VS_LDS) for w in 1:REP_PERIOD for y in DC_DISCHARGE_CONSTRAINTSET - eVreStorVSoCBalLongDurationStorageStart[y, w] += EP[:vCAPRES_DC_DISCHARGE][y, + eVreStorVSoCBalLongDurationStorageStart[y, w] += EP[:vCAPRES_DC_DISCHARGE][ + y, hours_per_subperiod * (w - 1) + 1] / by_rid(y, :eff_down_dc) end for y in DC_CHARGE_CONSTRAINTSET @@ -2313,7 +2327,8 @@ function vre_stor_capres!(EP::Model, inputs::Dict, setup::Dict) hours_per_subperiod * (w - 1) + 1] end for y in AC_DISCHARGE_CONSTRAINTSET - eVreStorVSoCBalLongDurationStorageStart[y, w] += EP[:vCAPRES_AC_DISCHARGE][y, + eVreStorVSoCBalLongDurationStorageStart[y, w] += EP[:vCAPRES_AC_DISCHARGE][ + y, hours_per_subperiod * (w - 1) + 1] / by_rid(y, :eff_down_ac) end for y in AC_CHARGE_CONSTRAINTSET @@ -2340,15 +2355,16 @@ function vre_stor_capres!(EP::Model, inputs::Dict, setup::Dict) cVreStorVSoCBalLongDurationStorage[y in VS_LDS, r in MODELED_PERIODS_INDEX], vCAPCONTRSTOR_VSOCw_VRE_STOR[y, mod1(r + 1, NPeriods)]==vCAPCONTRSTOR_VSOCw_VRE_STOR[y, r] + - vCAPCONTRSTOR_VdSOC_VRE_STOR[y, dfPeriodMap[r, :Rep_Period_Index]]) + vCAPCONTRSTOR_VdSOC_VRE_STOR[ + y, dfPeriodMap[r, :Rep_Period_Index]]) # Constraint 3: Initial reserve storage level for representative periods must also adhere to sub-period storage inventory balance # Initial storage = Final storage - change in storage inventory across representative period @constraint(EP, cVreStorVSoCBalLongDurationStorageSub[y in VS_LDS, r in REP_PERIODS_INDEX], - vCAPCONTRSTOR_VSOCw_VRE_STOR[y, - r]==EP[:vCAPRES_VS_VRE_STOR][y, - hours_per_subperiod * dfPeriodMap[r, :Rep_Period_Index]] - vCAPCONTRSTOR_VdSOC_VRE_STOR[y, dfPeriodMap[r, :Rep_Period_Index]]) + vCAPCONTRSTOR_VSOCw_VRE_STOR[y,r]==EP[:vCAPRES_VS_VRE_STOR][y, + hours_per_subperiod * dfPeriodMap[r, :Rep_Period_Index]] - + vCAPCONTRSTOR_VdSOC_VRE_STOR[y, dfPeriodMap[r, :Rep_Period_Index]]) # Constraint 4: Energy held in reserve at the beginning of each modeled period acts as a lower bound on the total energy held in storage @constraint(EP, @@ -2780,22 +2796,25 @@ function vre_stor_operational_reserves!(EP::Model, inputs::Dict, setup::Dict) eRegReqVreStor[t = 1:T], inputs["pReg_Req_VRE"] * sum(inputs["pP_Max_Solar"][y, t] * EP[:eTotalCap_SOLAR][y] * - by_rid(y, :etainverter) for y in SOLAR_REG) + by_rid(y, :etainverter) + for y in SOLAR_REG) +inputs["pReg_Req_VRE"] * sum(inputs["pP_Max_Wind"][y, t] * EP[:eTotalCap_WIND][y] for y in WIND_REG)) @expression(EP, eRsvReqVreStor[t = 1:T], inputs["pRsv_Req_VRE"] * sum(inputs["pP_Max_Solar"][y, t] * EP[:eTotalCap_SOLAR][y] * - by_rid(y, :etainverter) for y in SOLAR_RSV) + by_rid(y, :etainverter) + for y in SOLAR_RSV) +inputs["pRsv_Req_VRE"] * sum(inputs["pP_Max_Wind"][y, t] * EP[:eTotalCap_WIND][y] for y in WIND_RSV)) if !isempty(VRE_STOR_REG) @constraint(EP, cRegVreStor[t = 1:T], - sum(EP[:vREG][y, t] for y in inputs["REG"])>=EP[:eRegReq][t] + - eRegReqVreStor[t]) + sum(EP[:vREG][y, t] + for y in inputs["REG"])>=EP[:eRegReq][t] + + eRegReqVreStor[t]) end if !isempty(VRE_STOR_RSV) @constraint(EP, diff --git a/src/multi_stage/configure_multi_stage_inputs.jl b/src/multi_stage/configure_multi_stage_inputs.jl index d2f024d9a8..f74d379ed3 100644 --- a/src/multi_stage/configure_multi_stage_inputs.jl +++ b/src/multi_stage/configure_multi_stage_inputs.jl @@ -81,8 +81,8 @@ inputs: returns: dictionary containing updated model inputs, to be used in the generate\_model() method. """ function configure_multi_stage_inputs(inputs_d::Dict, - settings_d::Dict, - NetworkExpansion::Int64) + settings_d::Dict, + NetworkExpansion::Int64) gen = inputs_d["RESOURCES"] # Parameter inputs when multi-year discounting is activated @@ -121,31 +121,38 @@ function configure_multi_stage_inputs(inputs_d::Dict, # Conduct 1. and 2. for any co-located VRE-STOR resources if !isempty(inputs_d["VRE_STOR"]) gen_VRE_STOR = gen.VreStorage - gen_VRE_STOR.inv_cost_inverter_per_mwyr = compute_overnight_capital_cost(settings_d, + gen_VRE_STOR.inv_cost_inverter_per_mwyr = compute_overnight_capital_cost( + settings_d, inv_cost_inverter_per_mwyr.(gen_VRE_STOR), capital_recovery_period_dc.(gen_VRE_STOR), tech_wacc_dc.(gen_VRE_STOR)) - gen_VRE_STOR.inv_cost_solar_per_mwyr = compute_overnight_capital_cost(settings_d, + gen_VRE_STOR.inv_cost_solar_per_mwyr = compute_overnight_capital_cost( + settings_d, inv_cost_solar_per_mwyr.(gen_VRE_STOR), capital_recovery_period_solar.(gen_VRE_STOR), tech_wacc_solar.(gen_VRE_STOR)) - gen_VRE_STOR.inv_cost_wind_per_mwyr = compute_overnight_capital_cost(settings_d, + gen_VRE_STOR.inv_cost_wind_per_mwyr = compute_overnight_capital_cost( + settings_d, inv_cost_wind_per_mwyr.(gen_VRE_STOR), capital_recovery_period_wind.(gen_VRE_STOR), tech_wacc_wind.(gen_VRE_STOR)) - gen_VRE_STOR.inv_cost_discharge_dc_per_mwyr = compute_overnight_capital_cost(settings_d, + gen_VRE_STOR.inv_cost_discharge_dc_per_mwyr = compute_overnight_capital_cost( + settings_d, inv_cost_discharge_dc_per_mwyr.(gen_VRE_STOR), capital_recovery_period_discharge_dc.(gen_VRE_STOR), tech_wacc_discharge_dc.(gen_VRE_STOR)) - gen_VRE_STOR.inv_cost_charge_dc_per_mwyr = compute_overnight_capital_cost(settings_d, + gen_VRE_STOR.inv_cost_charge_dc_per_mwyr = compute_overnight_capital_cost( + settings_d, inv_cost_charge_dc_per_mwyr.(gen_VRE_STOR), capital_recovery_period_charge_dc.(gen_VRE_STOR), tech_wacc_charge_dc.(gen_VRE_STOR)) - gen_VRE_STOR.inv_cost_discharge_ac_per_mwyr = compute_overnight_capital_cost(settings_d, + gen_VRE_STOR.inv_cost_discharge_ac_per_mwyr = compute_overnight_capital_cost( + settings_d, inv_cost_discharge_ac_per_mwyr.(gen_VRE_STOR), capital_recovery_period_discharge_ac.(gen_VRE_STOR), tech_wacc_discharge_ac.(gen_VRE_STOR)) - gen_VRE_STOR.inv_cost_charge_ac_per_mwyr = compute_overnight_capital_cost(settings_d, + gen_VRE_STOR.inv_cost_charge_ac_per_mwyr = compute_overnight_capital_cost( + settings_d, inv_cost_charge_ac_per_mwyr.(gen_VRE_STOR), capital_recovery_period_charge_ac.(gen_VRE_STOR), tech_wacc_charge_ac.(gen_VRE_STOR)) diff --git a/src/multi_stage/dual_dynamic_programming.jl b/src/multi_stage/dual_dynamic_programming.jl index 9703784473..05ecd5154e 100644 --- a/src/multi_stage/dual_dynamic_programming.jl +++ b/src/multi_stage/dual_dynamic_programming.jl @@ -131,7 +131,7 @@ returns: * stats\_d – Dictionary which contains the run time, upper bound, and lower bound of each DDP iteration. * inputs\_d – Dictionary of inputs for each model stage, generated by the load\_inputs() method, modified by this method. """ -function run_ddp(models_d::Dict, setup::Dict, inputs_d::Dict) +function run_ddp(outpath::AbstractString, models_d::Dict, setup::Dict, inputs_d::Dict) settings_d = setup["MultiStageSettingsDict"] num_stages = settings_d["NumStages"] # Total number of investment planning stages EPSILON = settings_d["ConvergenceTolerance"] # Tolerance @@ -142,10 +142,13 @@ function run_ddp(models_d::Dict, setup::Dict, inputs_d::Dict) ic = 0 # Iteration Counter results_d = Dict() # Dictionary to store the results to return - stats_d = Dict() # Dictionary to store the statistics (total time, upper bound, and lower bound for each iteration) times_a = [] # Array to store the total time of each iteration upper_bounds_a = [] # Array to store the upper bound of each iteration lower_bounds_a = [] # Array to store the lower bound of each iteration + stats_d = Dict() # Dictionary to store the statistics (total time, upper bound, and lower bound for each iteration) + stats_d["TIMES"] = times_a + stats_d["UPPER_BOUNDS"] = upper_bounds_a + stats_d["LOWER_BOUNDS"] = lower_bounds_a # Step a.i) Initialize cost-to-go function for t = 1:num_stages for t in 1:num_stages @@ -181,10 +184,6 @@ function run_ddp(models_d::Dict, setup::Dict, inputs_d::Dict) println(string("Lower Bound = ", z_lower)) println("***********") - stats_d["TIMES"] = times_a - stats_d["UPPER_BOUNDS"] = upper_bounds_a - statd_d["LOWER_BOUNDS"] = lower_bounds_a - return models_d, stats_d, inputs_d end @@ -232,10 +231,6 @@ function run_ddp(models_d::Dict, setup::Dict, inputs_d::Dict) println(string("Upper Bound = ", z_upper)) println(string("Lower Bound = ", z_lower)) println("***********") - - stats_d["TIMES"] = times_a - stats_d["UPPER_BOUNDS"] = upper_bounds_a - stats_d["LOWER_BOUNDS"] = lower_bounds_a return models_d, stats_d, inputs_d end ### @@ -253,6 +248,7 @@ function run_ddp(models_d::Dict, setup::Dict, inputs_d::Dict) end append!(upper_bounds_a, z_upper) # Store current iteration upper bound + update_multi_stage_stats_file(outpath, ic, z_upper, z_lower, NaN, new_row = true) # Step f) Backward pass for t = num_stages:2 for t in num_stages:-1:2 @@ -274,10 +270,13 @@ function run_ddp(models_d::Dict, setup::Dict, inputs_d::Dict) # Step g) Recalculate lower bound and go back to c) z_lower = objective_value(models_d[1]) append!(lower_bounds_a, z_lower) # Store current iteration lower bound + update_multi_stage_stats_file(outpath, ic, z_upper, z_lower, NaN) # Step h) Store the total time of the current iteration (in seconds) ddp_iteration_time = time() - ddp_prev_time append!(times_a, ddp_iteration_time) + update_multi_stage_stats_file(outpath, ic, z_upper, z_lower, ddp_iteration_time) + ddp_prev_time = time() end @@ -313,41 +312,9 @@ function run_ddp(models_d::Dict, setup::Dict, inputs_d::Dict) end ##### END of final forward pass - stats_d["TIMES"] = times_a - stats_d["UPPER_BOUNDS"] = upper_bounds_a - stats_d["LOWER_BOUNDS"] = lower_bounds_a - return models_d, stats_d, inputs_d end -@doc raw""" - write_multi_stage_outputs(stats_d::Dict, outpath::String, settings_d::Dict) - -This function calls various methods which write multi-stage modeling outputs as .csv files. - -inputs: - - * stats\_d – Dictionary which contains the run time, upper bound, and lower bound of each DDP iteration. - * outpath – String which represents the path to the Results directory. - * settings\_d - Dictionary containing settings configured in the GenX settings genx\_settings.yml file as well as the multi-stage settings file multi\_stage\_settings.yml. -""" -function write_multi_stage_outputs(stats_d::Dict, - outpath::String, - settings_d::Dict, - inputs_dict::Dict) - multi_stage_settings_d = settings_d["MultiStageSettingsDict"] - - write_multi_stage_capacities_discharge(outpath, multi_stage_settings_d) - write_multi_stage_capacities_charge(outpath, multi_stage_settings_d) - write_multi_stage_capacities_energy(outpath, multi_stage_settings_d) - if settings_d["NetworkExpansion"] == 1 - write_multi_stage_network_expansion(outpath, multi_stage_settings_d) - end - write_multi_stage_costs(outpath, multi_stage_settings_d, inputs_dict) - multi_stage_settings_d["Myopic"] == 0 && write_multi_stage_stats(outpath, stats_d) - write_multi_stage_settings(outpath, settings_d) -end - @doc raw""" fix_initial_investments(EP_prev::Model, EP_cur::Model, start_cap_d::Dict) @@ -362,9 +329,9 @@ inputs: returns: JuMP model with updated linking constraints. """ function fix_initial_investments(EP_prev::Model, - EP_cur::Model, - start_cap_d::Dict, - inputs_d::Dict) + EP_cur::Model, + start_cap_d::Dict, + inputs_d::Dict) ALL_CAP = union(inputs_d["RET_CAP"], inputs_d["NEW_CAP"]) # Set of all resources subject to inter-stage capacity tracking # start_cap_d dictionary contains the starting capacity expression name (e) as a key, @@ -401,9 +368,9 @@ inputs: returns: JuMP model with updated linking constraints. """ function fix_capacity_tracking(EP_prev::Model, - EP_cur::Model, - cap_track_d::Dict, - cur_stage::Int) + EP_cur::Model, + cap_track_d::Dict, + cur_stage::Int) # cap_track_d dictionary contains the endogenous retirement tracking array variable name (v) as a key, # and the associated linking constraint name (c) as a value @@ -517,9 +484,9 @@ inputs: returns: JuMP expression representing a sum of Benders cuts for linking capacity investment variables to be added to the cost-to-go function. """ function generate_cut_component_track(EP_cur::Model, - EP_next::Model, - var_name::Symbol, - constr_name::Symbol) + EP_next::Model, + var_name::Symbol, + constr_name::Symbol) next_dual_value = Float64[] cur_inv_value = Float64[] cur_inv_var = [] @@ -560,9 +527,9 @@ inputs: returns: JuMP expression representing a sum of Benders cuts for linking capacity investment variables to be added to the cost-to-go function. """ function generate_cut_component_inv(EP_cur::Model, - EP_next::Model, - expr_name::Symbol, - constr_name::Symbol) + EP_next::Model, + expr_name::Symbol, + constr_name::Symbol) next_dual_value = Float64[] cur_inv_value = Float64[] cur_inv_var = [] @@ -593,7 +560,7 @@ The updated objective function $OBJ^{*}$ returned by this method takes the form: where $OBJ$ is the original objective function. $OBJ$ is scaled by two terms. The first is a discount factor (applied only in the non-myopic case), which discounts costs associated with the model stage $p$ to year-0 dollars: ```math \begin{aligned} - DF = \frac{1}{(1+WACC)^{L_{p}*(p-1)}} + DF = \frac{1}{(1+WACC)^{\sum^{(p-1)}_{k=0}L_{k}}} \end{aligned} ``` where $WACC$ is the weighted average cost of capital, and $L_{p}$ is the length of each stage in years (both set in multi\_stage\_settings.yml) @@ -617,6 +584,10 @@ returns: JuMP model with updated objective function. """ function initialize_cost_to_go(settings_d::Dict, EP::Model, inputs::Dict) cur_stage = settings_d["CurStage"] # Current DDP Investment Planning Stage + cum_years = 0 + for stage_count in 1:(cur_stage - 1) + cum_years += settings_d["StageLengths"][stage_count] + end stage_len = settings_d["StageLengths"][cur_stage] wacc = settings_d["WACC"] # Interest Rate and also the discount rate unless specified other wise myopic = settings_d["Myopic"] == 1 # 1 if myopic (only one forward pass), 0 if full DDP @@ -629,7 +600,7 @@ function initialize_cost_to_go(settings_d::Dict, EP::Model, inputs::Dict) ### No discount factor or OPEX multiplier applied in myopic case as costs are left annualized. @objective(EP, Min, EP[:eObj]) else - DF = 1 / (1 + wacc)^(stage_len * (cur_stage - 1)) # Discount factor applied all to costs in each stage ### + DF = 1 / (1 + wacc)^(cum_years) # Discount factor applied all to costs in each stage ### # Initialize the cost-to-go variable @variable(EP, vALPHA>=0) @objective(EP, Min, DF * OPEXMULT * EP[:eObj]+vALPHA) diff --git a/src/multi_stage/endogenous_retirement.jl b/src/multi_stage/endogenous_retirement.jl index b88ac93e1d..29c549d719 100644 --- a/src/multi_stage/endogenous_retirement.jl +++ b/src/multi_stage/endogenous_retirement.jl @@ -23,9 +23,9 @@ function get_retirement_stage(cur_stage::Int, lifetime::Int, stage_lens::Array{I end function update_cumulative_min_ret!(inputs_d::Dict, - t::Int, - Resource_Set::String, - RetCap::Symbol) + t::Int, + Resource_Set::String, + RetCap::Symbol) gen_name = "RESOURCES" CumRetCap = Symbol("cum_" * String(RetCap)) # if the getter function exists in GenX then use it, otherwise get the attribute directly @@ -170,10 +170,10 @@ In other words, it is the largest index $r \in \{1, ..., (p-1)\}$ such that: ``` """ function endogenous_retirement_discharge!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (Discharge) Module") gen = inputs["RESOURCES"] @@ -211,7 +211,7 @@ function endogenous_retirement_discharge!(EP::Model, @expression(EP, eNewCapTrack[y in RET_CAP], sum(EP[:vCAPTRACK][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrack[y in RET_CAP], if y in COMMIT cum_min_retired_cap_mw(gen[y]) / cap_size(gen[y]) @@ -258,10 +258,10 @@ function endogenous_retirement_discharge!(EP::Model, end function endogenous_retirement_charge!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (Charge) Module") gen = inputs["RESOURCES"] @@ -293,7 +293,7 @@ function endogenous_retirement_charge!(EP::Model, @expression(EP, eNewCapTrackCharge[y in RET_CAP_CHARGE], sum(EP[:vCAPTRACKCHARGE][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackCharge[y in RET_CAP_CHARGE], cum_min_retired_charge_cap_mw(gen[y])) @@ -324,10 +324,10 @@ function endogenous_retirement_charge!(EP::Model, end function endogenous_retirement_energy!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (Energy) Module") gen = inputs["RESOURCES"] @@ -359,7 +359,7 @@ function endogenous_retirement_energy!(EP::Model, @expression(EP, eNewCapTrackEnergy[y in RET_CAP_ENERGY], sum(EP[:vCAPTRACKENERGY][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackEnergy[y in RET_CAP_ENERGY], cum_min_retired_energy_cap_mw(gen[y])) @@ -390,10 +390,10 @@ function endogenous_retirement_energy!(EP::Model, end function endogenous_retirement_vre_stor_dc!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (VRE-Storage DC) Module") gen = inputs["RESOURCES"] @@ -425,7 +425,7 @@ function endogenous_retirement_vre_stor_dc!(EP::Model, @expression(EP, eNewCapTrackDC[y in RET_CAP_DC], sum(EP[:vCAPTRACKDC][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackDC[y in RET_CAP_DC], cum_min_retired_cap_inverter_mw(gen[y])) @@ -456,10 +456,10 @@ function endogenous_retirement_vre_stor_dc!(EP::Model, end function endogenous_retirement_vre_stor_solar!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (VRE-Storage Solar) Module") gen = inputs["RESOURCES"] @@ -491,7 +491,7 @@ function endogenous_retirement_vre_stor_solar!(EP::Model, @expression(EP, eNewCapTrackSolar[y in RET_CAP_SOLAR], sum(EP[:vCAPTRACKSOLAR][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackSolar[y in RET_CAP_SOLAR], cum_min_retired_cap_solar_mw(gen[y])) @@ -522,10 +522,10 @@ function endogenous_retirement_vre_stor_solar!(EP::Model, end function endogenous_retirement_vre_stor_wind!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (VRE-Storage Wind) Module") gen = inputs["RESOURCES"] @@ -557,7 +557,7 @@ function endogenous_retirement_vre_stor_wind!(EP::Model, @expression(EP, eNewCapTrackWind[y in RET_CAP_WIND], sum(EP[:vCAPTRACKWIND][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackWind[y in RET_CAP_WIND], cum_min_retired_cap_wind_mw(gen[y])) @@ -588,10 +588,10 @@ function endogenous_retirement_vre_stor_wind!(EP::Model, end function endogenous_retirement_vre_stor_stor!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (VRE-Storage Storage) Module") gen = inputs["RESOURCES"] @@ -623,7 +623,7 @@ function endogenous_retirement_vre_stor_stor!(EP::Model, @expression(EP, eNewCapTrackEnergy_VS[y in RET_CAP_STOR], sum(EP[:vCAPTRACKENERGY_VS][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackEnergy_VS[y in RET_CAP_STOR], cum_min_retired_energy_cap_mw(gen[y])) @@ -654,10 +654,10 @@ function endogenous_retirement_vre_stor_stor!(EP::Model, end function endogenous_retirement_vre_stor_discharge_dc!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (VRE-Storage Discharge DC) Module") gen = inputs["RESOURCES"] @@ -691,7 +691,7 @@ function endogenous_retirement_vre_stor_discharge_dc!(EP::Model, @expression(EP, eNewCapTrackDischargeDC[y in RET_CAP_DISCHARGE_DC], sum(EP[:vCAPTRACKDISCHARGEDC][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackDischargeDC[y in RET_CAP_DISCHARGE_DC], cum_min_retired_cap_discharge_dc_mw(gen[y])) @@ -723,10 +723,10 @@ function endogenous_retirement_vre_stor_discharge_dc!(EP::Model, end function endogenous_retirement_vre_stor_charge_dc!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (VRE-Storage Charge DC) Module") gen = inputs["RESOURCES"] @@ -757,7 +757,7 @@ function endogenous_retirement_vre_stor_charge_dc!(EP::Model, @expression(EP, eNewCapTrackChargeDC[y in RET_CAP_CHARGE_DC], sum(EP[:vCAPTRACKCHARGEDC][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackChargeDC[y in RET_CAP_CHARGE_DC], cum_min_retired_cap_charge_dc_mw(gen[y])) @@ -788,10 +788,10 @@ function endogenous_retirement_vre_stor_charge_dc!(EP::Model, end function endogenous_retirement_vre_stor_discharge_ac!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (VRE-Storage Discharge AC) Module") gen = inputs["RESOURCES"] @@ -824,7 +824,7 @@ function endogenous_retirement_vre_stor_discharge_ac!(EP::Model, @expression(EP, eNewCapTrackDischargeAC[y in RET_CAP_DISCHARGE_AC], sum(EP[:vCAPTRACKDISCHARGEAC][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackDischargeAC[y in RET_CAP_DISCHARGE_AC], cum_min_retired_cap_discharge_ac_mw(gen[y])) @@ -856,10 +856,10 @@ function endogenous_retirement_vre_stor_discharge_ac!(EP::Model, end function endogenous_retirement_vre_stor_charge_ac!(EP::Model, - inputs::Dict, - num_stages::Int, - cur_stage::Int, - stage_lens::Array{Int, 1}) + inputs::Dict, + num_stages::Int, + cur_stage::Int, + stage_lens::Array{Int, 1}) println("Endogenous Retirement (VRE-Storage Charge AC) Module") gen = inputs["RESOURCES"] @@ -890,7 +890,7 @@ function endogenous_retirement_vre_stor_charge_ac!(EP::Model, @expression(EP, eNewCapTrackChargeAC[y in RET_CAP_CHARGE_AC], sum(EP[:vCAPTRACKCHARGEAC][y, p] - for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) + for p in 1:get_retirement_stage(cur_stage, lifetime(gen[y]), stage_lens))) @expression(EP, eMinRetCapTrackChargeAC[y in RET_CAP_CHARGE_AC], cum_min_retired_cap_charge_ac_mw(gen[y])) diff --git a/src/multi_stage/write_multi_stage_costs.jl b/src/multi_stage/write_multi_stage_costs.jl index 20bed3a9f7..92a8c9a71c 100644 --- a/src/multi_stage/write_multi_stage_costs.jl +++ b/src/multi_stage/write_multi_stage_costs.jl @@ -30,7 +30,13 @@ function write_multi_stage_costs(outpath::String, settings_d::Dict, inputs_dict: if myopic DF = 1 # DF=1 because we do not apply discount factor in myopic case else - DF = 1 / (1 + wacc)^(stage_lens[p] * (p - 1)) # Discount factor applied to ALL costs in each stage + cum_stage_length = 0 + if p > 1 + for stage_counter in 1:(p - 1) + cum_stage_length += stage_lens[stage_counter] + end + end + DF = 1 / (1 + wacc)^(cum_stage_length) # Discount factor applied to ALL costs in each stage end df_costs[!, Symbol("TotalCosts_p$p")] = DF .* costs_d[p][!, Symbol("Total")] end @@ -39,13 +45,13 @@ function write_multi_stage_costs(outpath::String, settings_d::Dict, inputs_dict: for cost in ["cVar", "cNSE", "cStart", "cUnmetRsv", "cUnmetPolicyPenalty"] if cost in df_costs[!, :Costs] df_costs[df_costs[!, :Costs] .== cost, 2:end] = transpose(OPEXMULTS) .* - df_costs[df_costs[!, :Costs] .== cost, - 2:end] + df_costs[df_costs[!, :Costs] .== cost, 2:end] end end # Remove "cTotal" from results (as this includes Cost-to-Go) df_costs = df_costs[df_costs[!, :Costs] .!= "cTotal", :] + @warn("The cost calculation of the multi-stage GenX is approximate currently, and we will be refining it more in one of the future releases.") CSV.write(joinpath(outpath, "costs_multi_stage.csv"), df_costs) end diff --git a/src/multi_stage/write_multi_stage_outputs.jl b/src/multi_stage/write_multi_stage_outputs.jl new file mode 100644 index 0000000000..4e6d5612d1 --- /dev/null +++ b/src/multi_stage/write_multi_stage_outputs.jl @@ -0,0 +1,30 @@ +@doc raw""" + write_multi_stage_outputs(stats_d::Dict, + outpath::String, + settings_d::Dict, + inputs_dict::Dict) + +This function calls various methods which write multi-stage modeling outputs as .csv files. + +# Arguments: + * stats\_d: Dictionary which contains the run time, upper bound, and lower bound of each DDP iteration. + * outpath: String which represents the path to the Results directory. + * settings\_d: Dictionary containing settings configured in the GenX settings `genx_settings.yml` file as well as the multi-stage settings file `multi_stage_settings.yml`. + * inputs\_dict: Dictionary containing the input data for the multi-stage model. +""" +function write_multi_stage_outputs(stats_d::Dict, + outpath::String, + settings_d::Dict, + inputs_dict::Dict) + multi_stage_settings_d = settings_d["MultiStageSettingsDict"] + + write_multi_stage_capacities_discharge(outpath, multi_stage_settings_d) + write_multi_stage_capacities_charge(outpath, multi_stage_settings_d) + write_multi_stage_capacities_energy(outpath, multi_stage_settings_d) + if settings_d["NetworkExpansion"] == 1 + write_multi_stage_network_expansion(outpath, multi_stage_settings_d) + end + write_multi_stage_costs(outpath, multi_stage_settings_d, inputs_dict) + multi_stage_settings_d["Myopic"] == 0 && write_multi_stage_stats(outpath, stats_d) + write_multi_stage_settings(outpath, settings_d) +end diff --git a/src/multi_stage/write_multi_stage_stats.jl b/src/multi_stage/write_multi_stage_stats.jl index b0c089a9d0..7f37b04043 100644 --- a/src/multi_stage/write_multi_stage_stats.jl +++ b/src/multi_stage/write_multi_stage_stats.jl @@ -1,3 +1,8 @@ +_get_multi_stage_stats_filename() = "stats_multi_stage.csv" +function _get_multi_stage_stats_header() + ["Iteration_Number", "Seconds", "Upper_Bound", "Lower_Bound", "Relative_Gap"] +end + @doc raw""" write_multi_stage_stats(outpath::String, stats_d::Dict) @@ -9,6 +14,11 @@ inputs: * stats\_d – Dictionary which contains the run time, upper bound, and lower bound of each DDP iteration. """ function write_multi_stage_stats(outpath::String, stats_d::Dict) + filename = _get_multi_stage_stats_filename() + + # don't overwrite existing file + isfile(joinpath(outpath, filename)) && return nothing + times_a = stats_d["TIMES"] # Time (seconds) of each iteration upper_bounds_a = stats_d["UPPER_BOUNDS"] # Upper bound of each iteration lower_bounds_a = stats_d["LOWER_BOUNDS"] # Lower bound of each iteration @@ -19,11 +29,76 @@ function write_multi_stage_stats(outpath::String, stats_d::Dict) realtive_gap_a = (upper_bounds_a .- lower_bounds_a) ./ lower_bounds_a # Construct dataframe where first column is iteration number, second is iteration time - df_stats = DataFrame(Iteration_Number = iteration_count_a, - Seconds = times_a, - Upper_Bound = upper_bounds_a, - Lower_Bound = lower_bounds_a, - Relative_Gap = realtive_gap_a) + header = _get_multi_stage_stats_header() + df_stats = DataFrame(header .=> + [iteration_count_a, times_a, upper_bounds_a, lower_bounds_a, realtive_gap_a]) + + CSV.write(joinpath(outpath, filename), df_stats) + return nothing +end + +@doc raw""" + create_multi_stage_stats_file(outpath::String) + +Create an empty CSV file in the specified output directory with the filename `stats_multi_stage.csv`. +The file contains the columns defined in `_get_multi_stage_stats_header()`. +The function first generates the filename and header using `_get_multi_stage_stats_filename()` and +`_get_multi_stage_stats_header()` respectively. It then creates a DataFrame with column names as headers and +writes it into a CSV file in the specified output directory. + +# Arguments +- `outpath::String`: The output directory where the statistics file will be written. + +# Returns +- Nothing. A CSV file is written to the `outpath`. +""" +function create_multi_stage_stats_file(outpath::String) + filename = _get_multi_stage_stats_filename() + header = _get_multi_stage_stats_header() + df_stats = DataFrame([col_name => Float64[] for col_name in header]) + CSV.write(joinpath(outpath, filename), df_stats) +end + +@doc raw""" + update_multi_stage_stats_file(outpath::String, ic::Int64, upper_bound::Float64, lower_bound::Float64, iteration_time::Float64; new_row::Bool=false) + +Update a multi-stage statistics file. + +# Arguments +- `outpath::String`: The output directory where the statistics file will be written. +- `ic::Int64`: The iteration count. +- `upper_bound::Float64`: The upper bound value. +- `lower_bound::Float64`: The lower bound value. +- `iteration_time::Float64`: The iteration time value. +- `new_row::Bool=false`: Optional argument to determine whether to append a new row (if true) or update the current row (if false). + +The function first checks if the file exists. If it does not, it creates a new one. +Then, it reads the statistics from the existing file into a DataFrame. +It calculates the relative gap based on the upper and lower bounds, and either appends a new row or updates the current row based on the `new_row` argument. +Finally, it writes the updated DataFrame back to the file. + +# Returns +- Nothing. A CSV file is updated or created at the `outpath`. +""" +function update_multi_stage_stats_file(outpath::String, ic::Int64, upper_bound::Float64, + lower_bound::Float64, iteration_time::Float64; new_row::Bool = false) + filename = _get_multi_stage_stats_filename() + + # If the file does not exist, create it + if !isfile(joinpath(outpath, filename)) + create_multi_stage_stats_file(outpath) + end + + df_stats = CSV.read(joinpath(outpath, filename), DataFrame, types = Float64) + + relative_gap = (upper_bound - lower_bound) / lower_bound + + new_values = [ic, iteration_time, upper_bound, lower_bound, relative_gap] + + # If new_row is true, append the new values to the end of the dataframe + # otherwise, update the row at index ic + new_row ? push!(df_stats, new_values) : (df_stats[ic, :] = new_values) - CSV.write(joinpath(outpath, "stats_multi_stage.csv"), df_stats) + CSV.write(joinpath(outpath, filename), df_stats) + return nothing end diff --git a/src/time_domain_reduction/time_domain_reduction.jl b/src/time_domain_reduction/time_domain_reduction.jl index 71b39419b3..8f7f1e7a61 100644 --- a/src/time_domain_reduction/time_domain_reduction.jl +++ b/src/time_domain_reduction/time_domain_reduction.jl @@ -92,7 +92,8 @@ function parse_data(myinputs) end all_col_names = [demand_col_names; var_col_names; fuel_col_names] all_profiles = [demand_profiles..., var_profiles..., fuel_profiles...] - return demand_col_names, var_col_names, solar_col_names, wind_col_names, fuel_col_names, + return demand_col_names, + var_col_names, solar_col_names, wind_col_names, fuel_col_names, all_col_names, demand_profiles, var_profiles, solar_profiles, wind_profiles, fuel_profiles, all_profiles, @@ -183,7 +184,8 @@ function parse_multi_stage_data(inputs_dict) all_col_names = [demand_col_names; var_col_names; fuel_col_names] all_profiles = [demand_profiles..., var_profiles..., fuel_profiles...] - return demand_col_names, var_col_names, solar_col_names, wind_col_names, fuel_col_names, + return demand_col_names, + var_col_names, solar_col_names, wind_col_names, fuel_col_names, all_col_names, demand_profiles, var_profiles, solar_profiles, wind_profiles, fuel_profiles, all_profiles, @@ -234,11 +236,11 @@ K-Means: [https://juliastats.org/Clustering.jl/dev/kmeans.html](https://juliasta K-Medoids: [https://juliastats.org/Clustering.jl/stable/kmedoids.html](https://juliastats.org/Clustering.jl/stable/kmedoids.html) """ function cluster(ClusterMethod, - ClusteringInputDF, - NClusters, - nIters, - v = false, - random = true) + ClusteringInputDF, + NClusters, + nIters, + v = false, + random = true) if ClusterMethod == "kmeans" DistMatrix = pairwise(Euclidean(), Matrix(ClusteringInputDF), dims = 2) R = kmeans(Matrix(ClusteringInputDF), NClusters, init = :kmcen) @@ -338,7 +340,7 @@ system to be included among the extreme periods. They would select """ function get_extreme_period(DF, GDF, profKey, typeKey, statKey, - ConstCols, demand_col_names, solar_col_names, wind_col_names, v = false) + ConstCols, demand_col_names, solar_col_names, wind_col_names, v = false) if v println(profKey, " ", typeKey, " ", statKey) end @@ -477,15 +479,15 @@ hours that one hour in representative period $m$ represents in the original prof """ function get_demand_multipliers(ClusterOutputData, - InputData, - M, - W, - DemandCols, - TimestepsPerRepPeriod, - NewColNames, - NClusters, - Ncols, - v = false) + InputData, + M, + W, + DemandCols, + TimestepsPerRepPeriod, + NewColNames, + NClusters, + Ncols, + v = false) # Compute original zonal total demands zone_sums = Dict() for demandcol in DemandCols @@ -496,8 +498,8 @@ function get_demand_multipliers(ClusterOutputData, cluster_zone_sums = Dict() for m in 1:NClusters clustered_lp_DF = DataFrame(Dict(NewColNames[i] => ClusterOutputData[!, m][(TimestepsPerRepPeriod * (i - 1) + 1):(TimestepsPerRepPeriod * i)] - for i in 1:Ncols - if (Symbol(NewColNames[i]) in DemandCols))) + for i in 1:Ncols + if (Symbol(NewColNames[i]) in DemandCols))) cluster_zone_sums[m] = Dict() for demandcol in DemandCols cluster_zone_sums[m][demandcol] = sum(clustered_lp_DF[:, demandcol]) @@ -636,11 +638,11 @@ and wind profiles for co-located resources will be separated into different CSV after the clustering of the inputs has occurred. """ function cluster_inputs(inpath, - settings_path, - mysetup, - stage_id = -99, - v = false; - random = true) + settings_path, + mysetup, + stage_id = -99, + v = false; + random = true) if v println(now()) end @@ -762,7 +764,8 @@ function cluster_inputs(inpath, end # Remove Constant Columns - Add back later in final output - all_profiles, all_col_names, ConstData, ConstCols, ConstIdx = RemoveConstCols(all_profiles, + all_profiles, all_col_names, ConstData, ConstCols, ConstIdx = RemoveConstCols( + all_profiles, all_col_names, v) @@ -774,7 +777,7 @@ function cluster_inputs(inpath, # Put it together! InputData = DataFrame(Dict(all_col_names[c] => all_profiles[c] - for c in 1:length(all_col_names))) + for c in 1:length(all_col_names))) InputData = convert.(Float64, InputData) if v println("Demand (MW) and Capacity Factor Profiles: ") @@ -790,27 +793,30 @@ function cluster_inputs(inpath, # Normalize/standardize data based on user-provided method if ScalingMethod == "N" - normProfiles = [StatsBase.transform(fit(UnitRangeTransform, - InputData[:, c]; - dims = 1, - unit = true), - InputData[:, c]) for c in 1:length(OldColNames)] + normProfiles = [StatsBase.transform( + fit(UnitRangeTransform, + InputData[:, c]; + dims = 1, + unit = true), + InputData[:, c]) for c in 1:length(OldColNames)] elseif ScalingMethod == "S" - normProfiles = [StatsBase.transform(fit(ZScoreTransform, InputData[:, c]; dims = 1), - InputData[:, c]) for c in 1:length(OldColNames)] + normProfiles = [StatsBase.transform( + fit(ZScoreTransform, InputData[:, c]; dims = 1), + InputData[:, c]) for c in 1:length(OldColNames)] else println("ERROR InvalidScalingMethod: Use N for Normalization or S for Standardization.") println("CONTINUING using 0->1 normalization...") - normProfiles = [StatsBase.transform(fit(UnitRangeTransform, - InputData[:, c]; - dims = 1, - unit = true), - InputData[:, c]) for c in 1:length(OldColNames)] + normProfiles = [StatsBase.transform( + fit(UnitRangeTransform, + InputData[:, c]; + dims = 1, + unit = true), + InputData[:, c]) for c in 1:length(OldColNames)] end # Compile newly normalized/standardized profiles AnnualTSeriesNormalized = DataFrame(Dict(OldColNames[c] => normProfiles[c] - for c in 1:length(OldColNames))) + for c in 1:length(OldColNames))) # Optional pre-scaling of demand in order to give it more preference in clutering algorithm if DemandWeight != 1 # If we want to value demand more/less than capacity factors. Assume nonnegative. LW=1 means no scaling. @@ -887,7 +893,8 @@ function cluster_inputs(inpath, if v print(geoKey, " ") end - (stat, group_idx) = get_extreme_period(select(InputData, + (stat, group_idx) = get_extreme_period( + select(InputData, [:Group; Symbol.(z_cols_type)]), select(cgdf, [:Group; Symbol.(z_cols_type)]), profKey, @@ -934,18 +941,20 @@ function cluster_inputs(inpath, # from 8760 (# hours) by n (# profiles) DF to # 168*n (n period-stacked profiles) by 52 (# periods) DF DFsToConcat = [stack(InputData[isequal.(InputData.Group, w), :], OldColNames)[!, - :value] for w in 1:NumDataPoints if w <= NumDataPoints] + :value] for w in 1:NumDataPoints if w <= NumDataPoints] ModifiedData = DataFrame(Dict(Symbol(i) => DFsToConcat[i] for i in 1:NumDataPoints)) AnnualTSeriesNormalized[:, :Group] .= (1:Nhours) .÷ (TimestepsPerRepPeriod + 0.0001) .+ 1 - DFsToConcatNorm = [stack(AnnualTSeriesNormalized[isequal.(AnnualTSeriesNormalized.Group, - w), - :], - OldColNames)[!, - :value] for w in 1:NumDataPoints if w <= NumDataPoints] + DFsToConcatNorm = [stack( + AnnualTSeriesNormalized[ + isequal.(AnnualTSeriesNormalized.Group, + w), + :], + OldColNames)[!, + :value] for w in 1:NumDataPoints if w <= NumDataPoints] ModifiedDataNormalized = DataFrame(Dict(Symbol(i) => DFsToConcatNorm[i] - for i in 1:NumDataPoints)) + for i in 1:NumDataPoints)) # Remove extreme periods from normalized data before clustering NClusters = MinPeriods @@ -1129,14 +1138,14 @@ function cluster_inputs(inpath, for m in 1:NClusters rpDF = DataFrame(Dict(NewColNames[i] => ClusterOutputData[!, m][(TimestepsPerRepPeriod * (i - 1) + 1):(TimestepsPerRepPeriod * i)] - for i in 1:Ncols)) + for i in 1:Ncols)) gvDF = DataFrame(Dict(NewColNames[i] => ClusterOutputData[!, m][(TimestepsPerRepPeriod * (i - 1) + 1):(TimestepsPerRepPeriod * i)] - for i in 1:Ncols if (Symbol(NewColNames[i]) in VarCols))) + for i in 1:Ncols if (Symbol(NewColNames[i]) in VarCols))) dmDF = DataFrame(Dict(NewColNames[i] => ClusterOutputData[!, m][(TimestepsPerRepPeriod * (i - 1) + 1):(TimestepsPerRepPeriod * i)] - for i in 1:Ncols if (Symbol(NewColNames[i]) in DemandCols))) + for i in 1:Ncols if (Symbol(NewColNames[i]) in DemandCols))) if IncludeFuel fpDF = DataFrame(Dict(NewColNames[i] => ClusterOutputData[!, m][(TimestepsPerRepPeriod * (i - 1) + 1):(TimestepsPerRepPeriod * i)] - for i in 1:Ncols if (Symbol(NewColNames[i]) in FuelCols))) + for i in 1:Ncols if (Symbol(NewColNames[i]) in FuelCols))) end if !IncludeFuel fpDF = DataFrame(Placeholder = 1:TimestepsPerRepPeriod) @@ -1184,7 +1193,7 @@ function cluster_inputs(inpath, InputDataTest = InputData[(InputData.Group .<= NumDataPoints * 1.0), :] ClusterDataTest = vcat([rpDFs[a] for a in A]...) # To compare fairly, demand is not scaled here RMSE = Dict(c => rmse_score(InputDataTest[:, c], ClusterDataTest[:, c]) - for c in OldColNames) + for c in OldColNames) ##### Step 6: Print to File @@ -1198,10 +1207,10 @@ function cluster_inputs(inpath, end groups_per_stage = round.(Int, size(A, 1) * relative_lengths) group_ranges = [if i == 1 - 1:groups_per_stage[1] - else - (sum(groups_per_stage[1:(i - 1)]) + 1):sum(groups_per_stage[1:i]) - end + 1:groups_per_stage[1] + else + (sum(groups_per_stage[1:(i - 1)]) + 1):sum(groups_per_stage[1:i]) + end for i in 1:size(relative_lengths, 1)] Stage_Weights = Dict() @@ -1219,10 +1228,10 @@ function cluster_inputs(inpath, # Stage-specific weights and mappings cmap = countmap(A[group_ranges[per]]) # Count number of each rep. period in the planning stage weight_props = [if i in keys(cmap) - cmap[i] / N[i] - else - 0 - end + cmap[i] / N[i] + else + 0 + end for i in 1:size(M, 1)] # Proportions of each rep. period associated with each planning stage Stage_Weights[per] = weight_props .* W # Total hours that each rep. period represents within the planning stage Stage_PeriodMaps[per] = PeriodMap[group_ranges[per], :] @@ -1242,7 +1251,8 @@ function cluster_inputs(inpath, # Save output data to stage-specific locations ### TDR_Results/Demand_data_clustered.csv - demand_in = get_demand_dataframe(joinpath(inpath, "inputs", "inputs_p$per"), + demand_in = get_demand_dataframe( + joinpath(inpath, "inputs", "inputs_p$per"), mysetup["SystemFolder"]) demand_in[!, :Sub_Weights] = demand_in[!, :Sub_Weights] * 1.0 demand_in[1:length(Stage_Weights[per]), :Sub_Weights] .= Stage_Weights[per] @@ -1272,7 +1282,7 @@ function cluster_inputs(inpath, ### TDR_Results/Generators_variability.csv # Reset column ordering, add time index, and solve duplicate column name trouble with CSV.write's header kwarg GVColMap = Dict(RESOURCE_ZONES[i] => RESOURCES[i] - for i in 1:length(inputs_dict[1]["RESOURCE_NAMES"])) + for i in 1:length(inputs_dict[1]["RESOURCE_NAMES"])) GVColMap["Time_Index"] = "Time_Index" GVOutputData = GVOutputData[!, Symbol.(RESOURCE_ZONES)] insertcols!(GVOutputData, 1, :Time_Index => 1:size(GVOutputData, 1)) @@ -1402,7 +1412,7 @@ function cluster_inputs(inpath, # Reset column ordering, add time index, and solve duplicate column name trouble with CSV.write's header kwarg GVColMap = Dict(RESOURCE_ZONES[i] => RESOURCES[i] - for i in 1:length(myinputs["RESOURCE_NAMES"])) + for i in 1:length(myinputs["RESOURCE_NAMES"])) GVColMap["Time_Index"] = "Time_Index" GVOutputData = GVOutputData[!, Symbol.(RESOURCE_ZONES)] insertcols!(GVOutputData, 1, :Time_Index => 1:size(GVOutputData, 1)) @@ -1453,12 +1463,14 @@ function cluster_inputs(inpath, "Vre_and_stor_solar_variability.csv") WindVar_Outfile = joinpath(TimeDomainReductionFolder, "Vre_and_stor_wind_variability.csv") - CSV.write(joinpath(inpath, + CSV.write( + joinpath(inpath, "inputs", input_stage_directory, SolarVar_Outfile), solar_var) - CSV.write(joinpath(inpath, + CSV.write( + joinpath(inpath, "inputs", input_stage_directory, WindVar_Outfile), @@ -1494,7 +1506,8 @@ function cluster_inputs(inpath, if v println("Writing .yml settings...") end - YAML.write_file(joinpath(inpath, "inputs", input_stage_directory, YAML_Outfile), + YAML.write_file( + joinpath(inpath, "inputs", input_stage_directory, YAML_Outfile), myTDRsetup) end else @@ -1517,7 +1530,8 @@ function cluster_inputs(inpath, demand_in[!, :Time_Index] .= Time_Index_M for c in DemandCols - new_col = Union{Float64, Missings.Missing}[missing for i in 1:size(demand_in, 1)] + new_col = Union{Float64, Missings.Missing}[missing + for i in 1:size(demand_in, 1)] new_col[1:size(DMOutputData, 1)] = DMOutputData[!, c] demand_in[!, c] .= new_col end @@ -1532,7 +1546,7 @@ function cluster_inputs(inpath, # Reset column ordering, add time index, and solve duplicate column name trouble with CSV.write's header kwarg GVColMap = Dict(RESOURCE_ZONES[i] => RESOURCES[i] - for i in 1:length(myinputs["RESOURCE_NAMES"])) + for i in 1:length(myinputs["RESOURCE_NAMES"])) GVColMap["Time_Index"] = "Time_Index" GVOutputData = GVOutputData[!, Symbol.(RESOURCE_ZONES)] insertcols!(GVOutputData, 1, :Time_Index => 1:size(GVOutputData, 1)) diff --git a/src/write_outputs/capacity_reserve_margin/effective_capacity.jl b/src/write_outputs/capacity_reserve_margin/effective_capacity.jl index 92980eabf3..88a5e9f2b0 100644 --- a/src/write_outputs/capacity_reserve_margin/effective_capacity.jl +++ b/src/write_outputs/capacity_reserve_margin/effective_capacity.jl @@ -8,10 +8,10 @@ Effective capacity in a capacity reserve margin zone for certain resources in the given timesteps. """ function thermal_plant_effective_capacity(EP, - inputs, - resources::Vector{Int}, - capres_zone::Int, - timesteps::Vector{Int})::Matrix{Float64} + inputs, + resources::Vector{Int}, + capres_zone::Int, + timesteps::Vector{Int})::Matrix{Float64} eff_cap = thermal_plant_effective_capacity.(Ref(EP), Ref(inputs), resources, @@ -27,10 +27,10 @@ function thermal_plant_effective_capacity(EP::Model, inputs::Dict, y, capres_zon end function thermal_plant_effective_capacity(EP::Model, - inputs::Dict, - r_id::Int, - capres_zone::Int, - timesteps::Vector{Int})::Vector{Float64} + inputs::Dict, + r_id::Int, + capres_zone::Int, + timesteps::Vector{Int})::Vector{Float64} y = r_id gen = inputs["RESOURCES"] capresfactor = derating_factor(gen[y], tag = capres_zone) diff --git a/src/write_outputs/capacity_reserve_margin/write_capacity_value.jl b/src/write_outputs/capacity_reserve_margin/write_capacity_value.jl index 6301fa177c..b1aa431bb6 100644 --- a/src/write_outputs/capacity_reserve_margin/write_capacity_value.jl +++ b/src/write_outputs/capacity_reserve_margin/write_capacity_value.jl @@ -141,9 +141,9 @@ function write_capacity_value(path::AbstractString, inputs::Dict, setup::Dict, E Reserve = fill(Symbol("CapRes_$i"), G)) temp_dfCapValue = hcat(temp_dfCapValue, DataFrame(capvalue, :auto)) auxNew_Names = [Symbol("Resource"); - Symbol("Zone"); - Symbol("Reserve"); - [Symbol("t$t") for t in 1:T]] + Symbol("Zone"); + Symbol("Reserve"); + [Symbol("t$t") for t in 1:T]] rename!(temp_dfCapValue, auxNew_Names) append!(dfCapValue, temp_dfCapValue) end @@ -161,9 +161,9 @@ This is equal to the dual variable of the capacity constraint. Returns a vector, with units of $/MW """ function capacity_reserve_margin_price(EP::Model, - inputs::Dict, - setup::Dict, - capres_zone::Int)::Vector{Float64} + inputs::Dict, + setup::Dict, + capres_zone::Int)::Vector{Float64} ω = inputs["omega"] scale_factor = setup["ParameterScale"] == 1 ? ModelScalingFactor : 1 return dual.(EP[:cCapacityResMargin][capres_zone, :]) ./ ω * scale_factor diff --git a/src/write_outputs/capacity_reserve_margin/write_reserve_margin_revenue.jl b/src/write_outputs/capacity_reserve_margin/write_reserve_margin_revenue.jl index 5036b0759c..707147556b 100644 --- a/src/write_outputs/capacity_reserve_margin/write_reserve_margin_revenue.jl +++ b/src/write_outputs/capacity_reserve_margin/write_reserve_margin_revenue.jl @@ -9,9 +9,9 @@ Function for reporting the capacity revenue earned by each generator listed in t As a reminder, GenX models the capacity reserve margin (aka capacity market) at the time-dependent level, and each constraint either stands for an overall market or a locality constraint. """ function write_reserve_margin_revenue(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) scale_factor = setup["ParameterScale"] == 1 ? ModelScalingFactor : 1 gen = inputs["RESOURCES"] @@ -73,21 +73,27 @@ function write_reserve_margin_revenue(path::AbstractString, gen_VRE_STOR = gen.VreStorage tempresrev[VRE_STOR] = derating_factor.(gen_VRE_STOR, tag = i) .* ((value.(EP[:vP][VRE_STOR, :])) * weighted_price) - tempresrev[VRE_STOR_STOR] .-= derating_factor.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_discharge .!= 0) .| (gen_VRE_STOR.stor_dc_charge .!= 0) .| (gen_VRE_STOR.stor_ac_discharge .!= 0) .| (gen_VRE_STOR.stor_ac_charge .!= 0)], + tempresrev[VRE_STOR_STOR] .-= derating_factor.( + gen_VRE_STOR[(gen_VRE_STOR.stor_dc_discharge .!= 0) .| (gen_VRE_STOR.stor_dc_charge .!= 0) .| (gen_VRE_STOR.stor_ac_discharge .!= 0) .| (gen_VRE_STOR.stor_ac_charge .!= 0)], tag = i) .* (value.(EP[:vCHARGE_VRE_STOR][VRE_STOR_STOR, :]).data * weighted_price) - tempresrev[DC_DISCHARGE] .+= derating_factor.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_discharge .!= 0)], + tempresrev[DC_DISCHARGE] .+= derating_factor.( + gen_VRE_STOR[(gen_VRE_STOR.stor_dc_discharge .!= 0)], tag = i) .* ((value.(EP[:vCAPRES_DC_DISCHARGE][DC_DISCHARGE, - :]).data .* etainverter.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_discharge .!= 0)])) * + :]).data .* + etainverter.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_discharge .!= 0)])) * weighted_price) - tempresrev[AC_DISCHARGE] .+= derating_factor.(gen_VRE_STOR[(gen_VRE_STOR.stor_ac_discharge .!= 0)], + tempresrev[AC_DISCHARGE] .+= derating_factor.( + gen_VRE_STOR[(gen_VRE_STOR.stor_ac_discharge .!= 0)], tag = i) .* ((value.(EP[:vCAPRES_AC_DISCHARGE][AC_DISCHARGE, :]).data) * weighted_price) - tempresrev[DC_CHARGE] .-= derating_factor.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_charge .!= 0)], + tempresrev[DC_CHARGE] .-= derating_factor.( + gen_VRE_STOR[(gen_VRE_STOR.stor_dc_charge .!= 0)], tag = i) .* ((value.(EP[:vCAPRES_DC_CHARGE][DC_CHARGE, :]).data ./ etainverter.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_charge .!= 0)])) * weighted_price) - tempresrev[AC_CHARGE] .-= derating_factor.(gen_VRE_STOR[(gen_VRE_STOR.stor_ac_charge .!= 0)], + tempresrev[AC_CHARGE] .-= derating_factor.( + gen_VRE_STOR[(gen_VRE_STOR.stor_ac_charge .!= 0)], tag = i) .* ((value.(EP[:vCAPRES_AC_CHARGE][AC_CHARGE, :]).data) * weighted_price) end diff --git a/src/write_outputs/capacity_reserve_margin/write_reserve_margin_slack.jl b/src/write_outputs/capacity_reserve_margin/write_reserve_margin_slack.jl index c3f2ebf2c4..f6d71e4fb4 100644 --- a/src/write_outputs/capacity_reserve_margin/write_reserve_margin_slack.jl +++ b/src/write_outputs/capacity_reserve_margin/write_reserve_margin_slack.jl @@ -1,7 +1,7 @@ function write_reserve_margin_slack(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) NCRM = inputs["NCapacityReserveMargin"] T = inputs["T"] # Number of time steps (hours) dfResMar_slack = DataFrame(CRM_Constraint = [Symbol("CapRes_$res") for res in 1:NCRM], diff --git a/src/write_outputs/capacity_reserve_margin/write_reserve_margin_w.jl b/src/write_outputs/capacity_reserve_margin/write_reserve_margin_w.jl index 025f5cd4be..74b4efa7fe 100644 --- a/src/write_outputs/capacity_reserve_margin/write_reserve_margin_w.jl +++ b/src/write_outputs/capacity_reserve_margin/write_reserve_margin_w.jl @@ -8,7 +8,7 @@ function write_reserve_margin_w(path::AbstractString, inputs::Dict, setup::Dict, end dfResMar_w = hcat(dfResMar_w, DataFrame(temp_ResMar_w, :auto)) auxNew_Names_res = [Symbol("Constraint"); - [Symbol("CapRes_$i") for i in 1:inputs["NCapacityReserveMargin"]]] + [Symbol("CapRes_$i") for i in 1:inputs["NCapacityReserveMargin"]]] rename!(dfResMar_w, auxNew_Names_res) CSV.write(joinpath(path, "ReserveMargin_w.csv"), dfResMar_w) end diff --git a/src/write_outputs/co2_cap/write_co2_cap.jl b/src/write_outputs/co2_cap/write_co2_cap.jl index a78ddf2f23..19cba87d71 100644 --- a/src/write_outputs/co2_cap/write_co2_cap.jl +++ b/src/write_outputs/co2_cap/write_co2_cap.jl @@ -5,7 +5,8 @@ Function for reporting carbon price associated with carbon cap constraints. """ function write_co2_cap(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) - dfCO2Price = DataFrame(CO2_Cap = [Symbol("CO2_Cap_$cap") for cap in 1:inputs["NCO2Cap"]], + dfCO2Price = DataFrame( + CO2_Cap = [Symbol("CO2_Cap_$cap") for cap in 1:inputs["NCO2Cap"]], CO2_Price = (-1) * (dual.(EP[:cCO2Emissions_systemwide]))) if setup["ParameterScale"] == 1 dfCO2Price.CO2_Price .*= ModelScalingFactor # Convert Million$/kton to $/ton diff --git a/src/write_outputs/energy_share_requirement/write_esr_revenue.jl b/src/write_outputs/energy_share_requirement/write_esr_revenue.jl index 246853eaec..3bb98c6e1f 100644 --- a/src/write_outputs/energy_share_requirement/write_esr_revenue.jl +++ b/src/write_outputs/energy_share_requirement/write_esr_revenue.jl @@ -4,11 +4,11 @@ Function for reporting the renewable/clean credit revenue earned by each generator listed in the input file. GenX will print this file only when RPS/CES is modeled and the shadow price can be obtained form the solver. Each row corresponds to a generator, and each column starting from the 6th to the second last is the total revenue earned from each RPS constraint. The revenue is calculated as the total annual generation (if elgible for the corresponding constraint) multiplied by the RPS/CES price. The last column is the total revenue received from all constraint. The unit is \$. """ function write_esr_revenue(path::AbstractString, - inputs::Dict, - setup::Dict, - dfPower::DataFrame, - dfESR::DataFrame, - EP::Model) + inputs::Dict, + setup::Dict, + dfPower::DataFrame, + dfESR::DataFrame, + EP::Model) gen = inputs["RESOURCES"] regions = region.(gen) clusters = cluster.(gen) @@ -66,14 +66,16 @@ function write_esr_revenue(path::AbstractString, dfESRRev[SOLAR_WIND, esr_col] = (((value.(EP[:vP_WIND][SOLAR_WIND, :]).data * weight) .* - esr_vrestor.(gen_VRE_STOR[solar_and_wind_resources], + esr_vrestor.( + gen_VRE_STOR[solar_and_wind_resources], tag = i) * price) + (value.(EP[:vP_SOLAR][SOLAR_WIND, :]).data .* etainverter.(gen_VRE_STOR[solar_and_wind_resources]) * weight) .* - esr_vrestor.(gen_VRE_STOR[solar_and_wind_resources], + esr_vrestor.( + gen_VRE_STOR[solar_and_wind_resources], tag = i) * price) end end diff --git a/src/write_outputs/hydrogen/write_hourly_matching_prices.jl b/src/write_outputs/hydrogen/write_hourly_matching_prices.jl index 00d0ce3220..794eb9111b 100644 --- a/src/write_outputs/hydrogen/write_hourly_matching_prices.jl +++ b/src/write_outputs/hydrogen/write_hourly_matching_prices.jl @@ -1,7 +1,7 @@ function write_hourly_matching_prices(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) T = inputs["T"] # Number of time steps (hours) Z = inputs["Z"] # Number of zones scale_factor = setup["ParameterScale"] == 1 ? ModelScalingFactor : 1 @@ -10,8 +10,9 @@ function write_hourly_matching_prices(path::AbstractString, dfHourlyMatchPrices = DataFrame(Zone = 1:Z) # The unit is $/MWh # Dividing dual variable for each hour with corresponding hourly weight to retrieve marginal cost of the constraint dfHourlyMatchPrices = hcat(dfHourlyMatchPrices, - DataFrame(dual.(EP[:cHourlyMatching]).data ./ transpose(inputs["omega"]) * - scale_factor, + DataFrame( + dual.(EP[:cHourlyMatching]).data ./ transpose(inputs["omega"]) * + scale_factor, :auto)) auxNew_Names = [Symbol("Zone"); [Symbol("t$t") for t in 1:T]] diff --git a/src/write_outputs/long_duration_storage/write_opwrap_lds_stor_init.jl b/src/write_outputs/long_duration_storage/write_opwrap_lds_stor_init.jl index 81587b655a..66c8ed7efb 100644 --- a/src/write_outputs/long_duration_storage/write_opwrap_lds_stor_init.jl +++ b/src/write_outputs/long_duration_storage/write_opwrap_lds_stor_init.jl @@ -1,7 +1,7 @@ function write_opwrap_lds_stor_init(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) ## Extract data frames from input dictionary gen = inputs["RESOURCES"] zones = zone_id.(gen) diff --git a/src/write_outputs/min_max_capacity_requirement/write_maximum_capacity_requirement.jl b/src/write_outputs/min_max_capacity_requirement/write_maximum_capacity_requirement.jl index 997057955e..8d1b3450ee 100644 --- a/src/write_outputs/min_max_capacity_requirement/write_maximum_capacity_requirement.jl +++ b/src/write_outputs/min_max_capacity_requirement/write_maximum_capacity_requirement.jl @@ -1,10 +1,11 @@ function write_maximum_capacity_requirement(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) NumberOfMaxCapReqs = inputs["NumberOfMaxCapReqs"] - dfMaxCapPrice = DataFrame(Constraint = [Symbol("MaxCapReq_$maxcap") - for maxcap in 1:NumberOfMaxCapReqs], + dfMaxCapPrice = DataFrame( + Constraint = [Symbol("MaxCapReq_$maxcap") + for maxcap in 1:NumberOfMaxCapReqs], Price = -dual.(EP[:cZoneMaxCapReq])) scale_factor = setup["ParameterScale"] == 1 ? ModelScalingFactor : 1 diff --git a/src/write_outputs/min_max_capacity_requirement/write_minimum_capacity_requirement.jl b/src/write_outputs/min_max_capacity_requirement/write_minimum_capacity_requirement.jl index 346213ae61..bae7d17ee9 100644 --- a/src/write_outputs/min_max_capacity_requirement/write_minimum_capacity_requirement.jl +++ b/src/write_outputs/min_max_capacity_requirement/write_minimum_capacity_requirement.jl @@ -1,10 +1,11 @@ function write_minimum_capacity_requirement(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) NumberOfMinCapReqs = inputs["NumberOfMinCapReqs"] - dfMinCapPrice = DataFrame(Constraint = [Symbol("MinCapReq_$mincap") - for mincap in 1:NumberOfMinCapReqs], + dfMinCapPrice = DataFrame( + Constraint = [Symbol("MinCapReq_$mincap") + for mincap in 1:NumberOfMinCapReqs], Price = dual.(EP[:cZoneMinCapReq])) scale_factor = setup["ParameterScale"] == 1 ? ModelScalingFactor : 1 diff --git a/src/write_outputs/reserves/write_operating_reserve_price_revenue.jl b/src/write_outputs/reserves/write_operating_reserve_price_revenue.jl index c3d4f389a4..79ab9b8cbe 100644 --- a/src/write_outputs/reserves/write_operating_reserve_price_revenue.jl +++ b/src/write_outputs/reserves/write_operating_reserve_price_revenue.jl @@ -8,9 +8,9 @@ Function for reporting the operating reserve and regulation revenue earned by ge As a reminder, GenX models the operating reserve and regulation at the time-dependent level, and each constraint either stands for an overall market or a locality constraint. """ function write_operating_reserve_regulation_revenue(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) scale_factor = setup["ParameterScale"] == 1 ? ModelScalingFactor : 1 gen = inputs["RESOURCES"] diff --git a/src/write_outputs/reserves/write_rsv.jl b/src/write_outputs/reserves/write_rsv.jl index f48b40c034..ba38ccb727 100644 --- a/src/write_outputs/reserves/write_rsv.jl +++ b/src/write_outputs/reserves/write_rsv.jl @@ -18,9 +18,9 @@ function write_rsv(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) total_unmet = sum(unmet_vec) dfRsv = hcat(dfRsv, DataFrame(rsv, :auto)) auxNew_Names = [Symbol("Resource"); - Symbol("Zone"); - Symbol("AnnualSum"); - [Symbol("t$t") for t in 1:T]] + Symbol("Zone"); + Symbol("AnnualSum"); + [Symbol("t$t") for t in 1:T]] rename!(dfRsv, auxNew_Names) total = DataFrame(["Total" 0 sum(dfRsv.AnnualSum) zeros(1, T)], :auto) diff --git a/src/write_outputs/transmission/write_transmission_flows.jl b/src/write_outputs/transmission/write_transmission_flows.jl index 199f61e72f..b143496b69 100644 --- a/src/write_outputs/transmission/write_transmission_flows.jl +++ b/src/write_outputs/transmission/write_transmission_flows.jl @@ -1,7 +1,7 @@ function write_transmission_flows(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) # Transmission related values T = inputs["T"] # Number of time steps (hours) L = inputs["L"] # Number of transmission lines diff --git a/src/write_outputs/transmission/write_transmission_losses.jl b/src/write_outputs/transmission/write_transmission_losses.jl index 6587cac4ed..8dc15f7b28 100644 --- a/src/write_outputs/transmission/write_transmission_losses.jl +++ b/src/write_outputs/transmission/write_transmission_losses.jl @@ -1,7 +1,7 @@ function write_transmission_losses(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) T = inputs["T"] # Number of time steps (hours) L = inputs["L"] # Number of transmission lines LOSS_LINES = inputs["LOSS_LINES"] diff --git a/src/write_outputs/write_capacityfactor.jl b/src/write_outputs/write_capacityfactor.jl index d2b8e94f20..17a019c2ff 100644 --- a/src/write_outputs/write_capacityfactor.jl +++ b/src/write_outputs/write_capacityfactor.jl @@ -67,10 +67,12 @@ function write_capacityfactor(path::AbstractString, inputs::Dict, setup::Dict, E # Capacity factor for electrolyzers is based on vUSE variable not vP if (!isempty(ELECTROLYZER)) dfCapacityfactor.AnnualSum[ELECTROLYZER] .= value.(EP[:vUSE][ELECTROLYZER, - :]).data * inputs["omega"] * scale_factor + :]).data * inputs["omega"] * + scale_factor dfCapacityfactor.CapacityFactor[ELECTROLYZER] .= (dfCapacityfactor.AnnualSum[ELECTROLYZER] ./ dfCapacityfactor.Capacity[ELECTROLYZER]) / - sum(inputs["omega"][t] for t in 1:T) + sum(inputs["omega"][t] + for t in 1:T) end CSV.write(joinpath(path, "capacityfactor.csv"), dfCapacityfactor) diff --git a/src/write_outputs/write_co2.jl b/src/write_outputs/write_co2.jl index 1a6d0eab6b..729562425e 100644 --- a/src/write_outputs/write_co2.jl +++ b/src/write_outputs/write_co2.jl @@ -10,9 +10,9 @@ function write_co2(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) end function write_co2_emissions_plant(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) gen = inputs["RESOURCES"] G = inputs["G"] # Number of resources (generators, storage, DR, and DERs) diff --git a/src/write_outputs/write_costs.jl b/src/write_outputs/write_costs.jl index 055979a126..0124c67828 100644 --- a/src/write_outputs/write_costs.jl +++ b/src/write_outputs/write_costs.jl @@ -22,7 +22,7 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) "cUnmetRsv", "cNetworkExp", "cUnmetPolicyPenalty", - "cCO2", + "cCO2" ] if !isempty(VRE_STOR) push!(cost_list, "cGridConnection") @@ -70,7 +70,7 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) 0.0, 0.0, 0.0, - 0.0, + 0.0 ] else total_cost = [ @@ -83,7 +83,7 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) 0.0, 0.0, 0.0, - 0.0, + 0.0 ] end @@ -298,7 +298,7 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) "-", "-", "-", - tempCCO2, + tempCCO2 ] if !isempty(VRE_STOR) push!(temp_cost_list, "-") diff --git a/src/write_outputs/write_emissions.jl b/src/write_outputs/write_emissions.jl index e244df74c9..b53641bc6b 100644 --- a/src/write_outputs/write_emissions.jl +++ b/src/write_outputs/write_emissions.jl @@ -32,8 +32,8 @@ function write_emissions(path::AbstractString, inputs::Dict, setup::Dict, EP::Mo DataFrame(tempCO2Price, :auto), DataFrame(AnnualSum = Array{Float64}(undef, Z))) auxNew_Names = [Symbol("Zone"); - [Symbol("CO2_Price_$cap") for cap in 1:inputs["NCO2Cap"]]; - Symbol("AnnualSum")] + [Symbol("CO2_Price_$cap") for cap in 1:inputs["NCO2Cap"]]; + Symbol("AnnualSum")] rename!(dfEmissions, auxNew_Names) else dfEmissions = DataFrame(Zone = 1:Z, AnnualSum = Array{Float64}(undef, Z)) @@ -48,10 +48,11 @@ function write_emissions(path::AbstractString, inputs::Dict, setup::Dict, EP::Mo if setup["WriteOutputs"] == "annual" total = DataFrame(["Total" sum(dfEmissions.AnnualSum)], [:Zone; :AnnualSum]) if setup["CO2Cap"] >= 1 - total = DataFrame(["Total" zeros(1, inputs["NCO2Cap"]) sum(dfEmissions.AnnualSum)], + total = DataFrame( + ["Total" zeros(1, inputs["NCO2Cap"]) sum(dfEmissions.AnnualSum)], [:Zone; - [Symbol("CO2_Price_$cap") for cap in 1:inputs["NCO2Cap"]]; - :AnnualSum]) + [Symbol("CO2_Price_$cap") for cap in 1:inputs["NCO2Cap"]]; + :AnnualSum]) end dfEmissions = vcat(dfEmissions, total) CSV.write(joinpath(path, "emissions.csv"), dfEmissions) @@ -60,11 +61,12 @@ function write_emissions(path::AbstractString, inputs::Dict, setup::Dict, EP::Mo DataFrame(emissions_by_zone * scale_factor, :auto)) if setup["CO2Cap"] >= 1 auxNew_Names = [Symbol("Zone"); - [Symbol("CO2_Price_$cap") for cap in 1:inputs["NCO2Cap"]]; - Symbol("AnnualSum"); - [Symbol("t$t") for t in 1:T]] + [Symbol("CO2_Price_$cap") for cap in 1:inputs["NCO2Cap"]]; + Symbol("AnnualSum"); + [Symbol("t$t") for t in 1:T]] rename!(dfEmissions, auxNew_Names) - total = DataFrame(["Total" zeros(1, inputs["NCO2Cap"]) sum(dfEmissions[!, + total = DataFrame( + ["Total" zeros(1, inputs["NCO2Cap"]) sum(dfEmissions[!, :AnnualSum]) fill(0.0, (1, T))], :auto) for t in 1:T @@ -73,10 +75,11 @@ function write_emissions(path::AbstractString, inputs::Dict, setup::Dict, EP::Mo end else auxNew_Names = [Symbol("Zone"); - Symbol("AnnualSum"); - [Symbol("t$t") for t in 1:T]] + Symbol("AnnualSum"); + [Symbol("t$t") for t in 1:T]] rename!(dfEmissions, auxNew_Names) - total = DataFrame(["Total" sum(dfEmissions[!, :AnnualSum]) fill(0.0, + total = DataFrame( + ["Total" sum(dfEmissions[!, :AnnualSum]) fill(0.0, (1, T))], :auto) for t in 1:T @@ -108,8 +111,8 @@ function write_emissions(path::AbstractString, inputs::Dict, setup::Dict, EP::Mo dfEmissions = hcat(dfEmissions, DataFrame(emissions_by_zone * scale_factor, :auto)) auxNew_Names = [Symbol("Zone"); - Symbol("AnnualSum"); - [Symbol("t$t") for t in 1:T]] + Symbol("AnnualSum"); + [Symbol("t$t") for t in 1:T]] rename!(dfEmissions, auxNew_Names) total = DataFrame(["Total" sum(dfEmissions[!, :AnnualSum]) fill(0.0, (1, T))], :auto) diff --git a/src/write_outputs/write_fuel_consumption.jl b/src/write_outputs/write_fuel_consumption.jl index 47d70301cf..4ff094289d 100644 --- a/src/write_outputs/write_fuel_consumption.jl +++ b/src/write_outputs/write_fuel_consumption.jl @@ -15,9 +15,9 @@ function write_fuel_consumption(path::AbstractString, inputs::Dict, setup::Dict, end function write_fuel_consumption_plant(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) gen = inputs["RESOURCES"] HAS_FUEL = inputs["HAS_FUEL"] @@ -41,15 +41,10 @@ function write_fuel_consumption_plant(path::AbstractString, tempannualsum_fuel_heat_multi_total = zeros(length(HAS_FUEL)) tempannualsum_fuel_cost_multi = zeros(length(HAS_FUEL)) for g in MULTI_FUELS - tempannualsum_fuel_heat_multi_generation[findfirst(x -> x == g, HAS_FUEL)] = value.(EP[:ePlantFuelConsumptionYear_multi_generation][g, - i]) - tempannualsum_fuel_heat_multi_start[findfirst(x -> x == g, HAS_FUEL)] = value.(EP[:ePlantFuelConsumptionYear_multi_start][g, - i]) - tempannualsum_fuel_heat_multi_total[findfirst(x -> x == g, HAS_FUEL)] = value.(EP[:ePlantFuelConsumptionYear_multi][g, - i]) - tempannualsum_fuel_cost_multi[findfirst(x -> x == g, HAS_FUEL)] = value.(EP[:ePlantCFuelOut_multi][g, - i]) + value.(EP[:ePlantCFuelOut_multi_start][g, - i]) + tempannualsum_fuel_heat_multi_generation[findfirst(x -> x == g, HAS_FUEL)] = value.(EP[:ePlantFuelConsumptionYear_multi_generation][g,i]) + tempannualsum_fuel_heat_multi_start[findfirst(x -> x == g, HAS_FUEL)] = value.(EP[:ePlantFuelConsumptionYear_multi_start][g,i]) + tempannualsum_fuel_heat_multi_total[findfirst(x -> x == g, HAS_FUEL)] = value.(EP[:ePlantFuelConsumptionYear_multi][g,i]) + tempannualsum_fuel_cost_multi[findfirst(x -> x == g, HAS_FUEL)] = value.(EP[:ePlantCFuelOut_multi][g,i]) + value.(EP[:ePlantCFuelOut_multi_start][g,i]) end if setup["ParameterScale"] == 1 tempannualsum_fuel_heat_multi_generation *= ModelScalingFactor @@ -74,9 +69,9 @@ function write_fuel_consumption_plant(path::AbstractString, end function write_fuel_consumption_ts(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) T = inputs["T"] # Number of time steps (hours) HAS_FUEL = inputs["HAS_FUEL"] @@ -103,9 +98,9 @@ function write_fuel_consumption_ts(path::AbstractString, end function write_fuel_consumption_tot(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) # types of fuel fuel_types = inputs["fuels"] fuel_number = length(fuel_types) diff --git a/src/write_outputs/write_net_revenue.jl b/src/write_outputs/write_net_revenue.jl index a647a81dbd..ef42e9ca87 100644 --- a/src/write_outputs/write_net_revenue.jl +++ b/src/write_outputs/write_net_revenue.jl @@ -4,20 +4,20 @@ Function for writing net revenue of different generation technologies. """ function write_net_revenue(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model, - dfCap::DataFrame, - dfESRRev::DataFrame, - dfResRevenue::DataFrame, - dfChargingcost::DataFrame, - dfPower::DataFrame, - dfEnergyRevenue::DataFrame, - dfSubRevenue::DataFrame, - dfRegSubRevenue::DataFrame, - dfVreStor::DataFrame, - dfOpRegRevenue::DataFrame, - dfOpRsvRevenue::DataFrame) + inputs::Dict, + setup::Dict, + EP::Model, + dfCap::DataFrame, + dfESRRev::DataFrame, + dfResRevenue::DataFrame, + dfChargingcost::DataFrame, + dfPower::DataFrame, + dfEnergyRevenue::DataFrame, + dfSubRevenue::DataFrame, + dfRegSubRevenue::DataFrame, + dfVreStor::DataFrame, + dfOpRegRevenue::DataFrame, + dfOpRsvRevenue::DataFrame) gen = inputs["RESOURCES"] zones = zone_id.(gen) regions = region.(gen) @@ -114,14 +114,13 @@ function write_net_revenue(path::AbstractString, end if !isempty(DC_DISCHARGE) dfNetRevenue.Var_OM_cost_out[DC_DISCHARGE] += var_om_cost_per_mwh_discharge_dc.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_discharge .!= 0)]) .* - (value.(EP[:vP_DC_DISCHARGE][DC_DISCHARGE, - :]).data .* etainverter.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_discharge .!= 0)]) * + (value.(EP[:vP_DC_DISCHARGE][DC_DISCHARGE,:]).data .* + etainverter.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_discharge .!= 0)]) * inputs["omega"]) end if !isempty(AC_DISCHARGE) dfNetRevenue.Var_OM_cost_out[AC_DISCHARGE] += var_om_cost_per_mwh_discharge_ac.(gen_VRE_STOR[(gen_VRE_STOR.stor_ac_discharge .!= 0)]) .* - (value.(EP[:vP_AC_DISCHARGE][AC_DISCHARGE, - :]).data * inputs["omega"]) + (value.(EP[:vP_AC_DISCHARGE][AC_DISCHARGE,:]).data * inputs["omega"]) end end if setup["ParameterScale"] == 1 @@ -147,8 +146,8 @@ function write_net_revenue(path::AbstractString, if !isempty(VRE_STOR) if !isempty(DC_CHARGE) dfNetRevenue.Var_OM_cost_in[DC_CHARGE] += var_om_cost_per_mwh_charge_dc.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_charge .!= 0)]) .* - (value.(EP[:vP_DC_CHARGE][DC_CHARGE, - :]).data ./ etainverter.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_charge .!= 0)]) * + (value.(EP[:vP_DC_CHARGE][DC_CHARGE,:]).data ./ + etainverter.(gen_VRE_STOR[(gen_VRE_STOR.stor_dc_charge .!= 0)]) * inputs["omega"]) end if !isempty(AC_CHARGE) diff --git a/src/write_outputs/write_nse.jl b/src/write_outputs/write_nse.jl index 03605c1dd4..cf34845888 100644 --- a/src/write_outputs/write_nse.jl +++ b/src/write_outputs/write_nse.jl @@ -26,9 +26,9 @@ function write_nse(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) else # setup["WriteOutputs"] == "full" dfNse = hcat(dfNse, DataFrame(nse, :auto)) auxNew_Names = [Symbol("Segment"); - Symbol("Zone"); - Symbol("AnnualSum"); - [Symbol("t$t") for t in 1:T]] + Symbol("Zone"); + Symbol("AnnualSum"); + [Symbol("t$t") for t in 1:T]] rename!(dfNse, auxNew_Names) total = DataFrame(["Total" 0 sum(dfNse[!, :AnnualSum]) fill(0.0, (1, T))], :auto) diff --git a/src/write_outputs/write_outputs.jl b/src/write_outputs/write_outputs.jl index 02143a86bd..892e672d17 100644 --- a/src/write_outputs/write_outputs.jl +++ b/src/write_outputs/write_outputs.jl @@ -273,7 +273,8 @@ function write_outputs(EP::Model, path::AbstractString, setup::Dict, inputs::Dic if output_settings_d["WriteEnergyRevenue"] || output_settings_d["WriteNetRevenue"] - elapsed_time_energy_rev = @elapsed dfEnergyRevenue = write_energy_revenue(path, + elapsed_time_energy_rev = @elapsed dfEnergyRevenue = write_energy_revenue( + path, inputs, setup, EP) @@ -283,7 +284,8 @@ function write_outputs(EP::Model, path::AbstractString, setup::Dict, inputs::Dic if output_settings_d["WriteChargingCost"] || output_settings_d["WriteNetRevenue"] - elapsed_time_charging_cost = @elapsed dfChargingcost = write_charging_cost(path, + elapsed_time_charging_cost = @elapsed dfChargingcost = write_charging_cost( + path, inputs, setup, EP) @@ -293,7 +295,8 @@ function write_outputs(EP::Model, path::AbstractString, setup::Dict, inputs::Dic if output_settings_d["WriteSubsidyRevenue"] || output_settings_d["WriteNetRevenue"] - elapsed_time_subsidy = @elapsed dfSubRevenue, dfRegSubRevenue = write_subsidy_revenue(path, + elapsed_time_subsidy = @elapsed dfSubRevenue, dfRegSubRevenue = write_subsidy_revenue( + path, inputs, setup, EP) @@ -361,7 +364,8 @@ function write_outputs(EP::Model, path::AbstractString, setup::Dict, inputs::Dic if output_settings_d["WriteReserveMarginRevenue"] || output_settings_d["WriteNetRevenue"] - elapsed_time_res_rev = @elapsed dfResRevenue = write_reserve_margin_revenue(path, + elapsed_time_res_rev = @elapsed dfResRevenue = write_reserve_margin_revenue( + path, inputs, setup, EP) @@ -392,7 +396,8 @@ function write_outputs(EP::Model, path::AbstractString, setup::Dict, inputs::Dic dfOpRegRevenue = DataFrame() dfOpRsvRevenue = DataFrame() if setup["OperationalReserves"] == 1 && has_duals(EP) - elapsed_time_op_res_rev = @elapsed dfOpRegRevenue, dfOpRsvRevenue = write_operating_reserve_regulation_revenue(path, + elapsed_time_op_res_rev = @elapsed dfOpRegRevenue, dfOpRsvRevenue = write_operating_reserve_regulation_revenue( + path, inputs, setup, EP) @@ -436,7 +441,8 @@ function write_outputs(EP::Model, path::AbstractString, setup::Dict, inputs::Dic end if setup["HydrogenHourlyMatching"] == 1 && output_settings_d["WriteHourlyMatchingPrices"] - elapsed_time_hourly_matching_prices = @elapsed write_hourly_matching_prices(path, + elapsed_time_hourly_matching_prices = @elapsed write_hourly_matching_prices( + path, inputs, setup, EP) @@ -488,14 +494,14 @@ end Internal function for writing full time series outputs. This function wraps the instructions for creating the full time series output files. """ function write_fulltimeseries(fullpath::AbstractString, - dataOut::Matrix{Float64}, - dfOut::DataFrame) + dataOut::Matrix{Float64}, + dfOut::DataFrame) T = size(dataOut, 2) dfOut = hcat(dfOut, DataFrame(dataOut, :auto)) auxNew_Names = [Symbol("Resource"); - Symbol("Zone"); - Symbol("AnnualSum"); - [Symbol("t$t") for t in 1:T]] + Symbol("Zone"); + Symbol("AnnualSum"); + [Symbol("t$t") for t in 1:T]] rename!(dfOut, auxNew_Names) total = DataFrame(["Total" 0 sum(dfOut[!, :AnnualSum]) fill(0.0, (1, T))], auxNew_Names) total[!, 4:(T + 3)] .= sum(dataOut, dims = 1) diff --git a/src/write_outputs/write_power_balance.jl b/src/write_outputs/write_power_balance.jl index b945755c94..c574ae67aa 100644 --- a/src/write_outputs/write_power_balance.jl +++ b/src/write_outputs/write_power_balance.jl @@ -36,15 +36,14 @@ function write_power_balance(path::AbstractString, inputs::Dict, setup::Dict, EP STOR_ALL_ZONE = intersect(resources_in_zone_by_rid(gen, z), STOR_ALL) powerbalance[(z - 1) * L + 2, :] = sum(value.(EP[:vP][STOR_ALL_ZONE, :]), dims = 1) - powerbalance[(z - 1) * L + 3, :] = (-1) * - sum((value.(EP[:vCHARGE][STOR_ALL_ZONE, - :]).data), + powerbalance[(z - 1) * L + 3, :] = (-1) * sum( + (value.(EP[:vCHARGE][STOR_ALL_ZONE,:]).data), dims = 1) end if !isempty(intersect(resources_in_zone_by_rid(gen, z), FLEX)) FLEX_ZONE = intersect(resources_in_zone_by_rid(gen, z), FLEX) - powerbalance[(z - 1) * L + 4, :] = sum((value.(EP[:vCHARGE_FLEX][FLEX_ZONE, - :]).data), + powerbalance[(z - 1) * L + 4, :] = sum( + (value.(EP[:vCHARGE_FLEX][FLEX_ZONE,:]).data), dims = 1) powerbalance[(z - 1) * L + 5, :] = (-1) * sum(value.(EP[:vP][FLEX_ZONE, :]), dims = 1) @@ -61,9 +60,8 @@ function write_power_balance(path::AbstractString, inputs::Dict, setup::Dict, EP powerbalance[(z - 1) * L + 10, :] = (((-1) * inputs["pD"][:, z]))' # Transpose if !isempty(ELECTROLYZER) ELECTROLYZER_ZONE = intersect(resources_in_zone_by_rid(gen, z), ELECTROLYZER) - powerbalance[(z - 1) * L + 11, :] = (-1) * - sum(value.(EP[:vUSE][ELECTROLYZER_ZONE, - :].data), + powerbalance[(z - 1) * L + 11, :] = (-1) * sum( + value.(EP[:vUSE][ELECTROLYZER_ZONE,:].data), dims = 1) end # VRE storage discharge and charge @@ -91,9 +89,9 @@ function write_power_balance(path::AbstractString, inputs::Dict, setup::Dict, EP else # setup["WriteOutputs"] == "full" dfPowerBalance = hcat(dfPowerBalance, DataFrame(powerbalance, :auto)) auxNew_Names = [Symbol("BalanceComponent"); - Symbol("Zone"); - Symbol("AnnualSum"); - [Symbol("t$t") for t in 1:T]] + Symbol("Zone"); + Symbol("AnnualSum"); + [Symbol("t$t") for t in 1:T]] rename!(dfPowerBalance, auxNew_Names) CSV.write(joinpath(path, "power_balance.csv"), dftranspose(dfPowerBalance, false), diff --git a/src/write_outputs/write_storagedual.jl b/src/write_outputs/write_storagedual.jl index 17df67d47b..9f374b80e8 100644 --- a/src/write_outputs/write_storagedual.jl +++ b/src/write_outputs/write_storagedual.jl @@ -29,32 +29,40 @@ function write_storagedual(path::AbstractString, inputs::Dict, setup::Dict, EP:: if !isempty(STOR_ALL) STOR_ALL_NONLDS = setdiff(STOR_ALL, inputs["STOR_LONG_DURATION"]) STOR_ALL_LDS = intersect(STOR_ALL, inputs["STOR_LONG_DURATION"]) - dual_values[STOR_ALL, INTERIOR_SUBPERIODS] = (dual.(EP[:cSoCBalInterior][INTERIOR_SUBPERIODS, + dual_values[STOR_ALL, INTERIOR_SUBPERIODS] = (dual.(EP[:cSoCBalInterior][ + INTERIOR_SUBPERIODS, STOR_ALL]).data ./ inputs["omega"][INTERIOR_SUBPERIODS])' - dual_values[STOR_ALL_NONLDS, START_SUBPERIODS] = (dual.(EP[:cSoCBalStart][START_SUBPERIODS, + dual_values[STOR_ALL_NONLDS, START_SUBPERIODS] = (dual.(EP[:cSoCBalStart][ + START_SUBPERIODS, STOR_ALL_NONLDS]).data ./ inputs["omega"][START_SUBPERIODS])' if !isempty(STOR_ALL_LDS) if inputs["REP_PERIOD"] > 1 - dual_values[STOR_ALL_LDS, START_SUBPERIODS] = (dual.(EP[:cSoCBalLongDurationStorageStart][1:REP_PERIOD, + dual_values[STOR_ALL_LDS, START_SUBPERIODS] = (dual.(EP[:cSoCBalLongDurationStorageStart][ + 1:REP_PERIOD, STOR_ALL_LDS]).data ./ inputs["omega"][START_SUBPERIODS])' else - dual_values[STOR_ALL_LDS, START_SUBPERIODS] = (dual.(EP[:cSoCBalStart][START_SUBPERIODS, + dual_values[STOR_ALL_LDS, START_SUBPERIODS] = (dual.(EP[:cSoCBalStart][ + START_SUBPERIODS, STOR_ALL_LDS]).data ./ inputs["omega"][START_SUBPERIODS])' end end end if !isempty(VRE_STOR) - dual_values[VS_STOR, INTERIOR_SUBPERIODS] = ((dual.(EP[:cSoCBalInterior_VRE_STOR][VS_STOR, + dual_values[VS_STOR, INTERIOR_SUBPERIODS] = ((dual.(EP[:cSoCBalInterior_VRE_STOR][ + VS_STOR, INTERIOR_SUBPERIODS]).data)' ./ inputs["omega"][INTERIOR_SUBPERIODS])' - dual_values[VS_NONLDS, START_SUBPERIODS] = ((dual.(EP[:cSoCBalStart_VRE_STOR][VS_NONLDS, + dual_values[VS_NONLDS, START_SUBPERIODS] = ((dual.(EP[:cSoCBalStart_VRE_STOR][ + VS_NONLDS, START_SUBPERIODS]).data)' ./ inputs["omega"][START_SUBPERIODS])' if !isempty(VS_LDS) if inputs["REP_PERIOD"] > 1 - dual_values[VS_LDS, START_SUBPERIODS] = ((dual.(EP[:cVreStorSoCBalLongDurationStorageStart][VS_LDS, + dual_values[VS_LDS, START_SUBPERIODS] = ((dual.(EP[:cVreStorSoCBalLongDurationStorageStart][ + VS_LDS, 1:REP_PERIOD]).data)' ./ inputs["omega"][START_SUBPERIODS])' else - dual_values[VS_LDS, START_SUBPERIODS] = ((dual.(EP[:cSoCBalStart_VRE_STOR][VS_LDS, + dual_values[VS_LDS, START_SUBPERIODS] = ((dual.(EP[:cSoCBalStart_VRE_STOR][ + VS_LDS, START_SUBPERIODS]).data)' ./ inputs["omega"][START_SUBPERIODS])' end end diff --git a/src/write_outputs/write_subsidy_revenue.jl b/src/write_outputs/write_subsidy_revenue.jl index 5be5d66b48..f74bede14c 100644 --- a/src/write_outputs/write_subsidy_revenue.jl +++ b/src/write_outputs/write_subsidy_revenue.jl @@ -67,7 +67,8 @@ function write_subsidy_revenue(path::AbstractString, inputs::Dict, setup::Dict, if !isempty(MIN_CAP_GEN_SOLAR) dfRegSubRevenue.SubsidyRevenue[MIN_CAP_GEN_SOLAR] .+= ((value.(EP[:eTotalCap_SOLAR][MIN_CAP_GEN_SOLAR]).data) .* - etainverter.(gen[ids_with_policy(gen, + etainverter.(gen[ids_with_policy( + gen, min_cap_solar, tag = mincap)]) * diff --git a/src/write_outputs/write_vre_stor.jl b/src/write_outputs/write_vre_stor.jl index ef261d5196..5e54303ebf 100644 --- a/src/write_outputs/write_vre_stor.jl +++ b/src/write_outputs/write_vre_stor.jl @@ -263,7 +263,7 @@ function write_vre_stor_capacity(path::AbstractString, inputs::Dict, setup::Dict :StartDischargeACCap, :RetDischargeACCap, :NewDischargeACCap, - :EndDischargeACCap, + :EndDischargeACCap ] dfCap[!, columns_to_scale] .*= ModelScalingFactor end @@ -365,9 +365,9 @@ end Function for writing the vre-storage discharging decision variables/expressions. """ function write_vre_stor_discharge(path::AbstractString, - inputs::Dict, - setup::Dict, - EP::Model) + inputs::Dict, + setup::Dict, + EP::Model) gen = inputs["RESOURCES"] gen_VRE_STOR = gen.VreStorage T = inputs["T"] diff --git a/test/runtests.jl b/test/runtests.jl index ff624f632e..f2fc3896ac 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -56,3 +56,13 @@ end include("test_retrofit.jl") end end + +# Test writing outputs +@testset "Writing outputs " begin + for test_file in filter!(x -> endswith(x, ".jl"), readdir("writing_outputs")) + include("writing_outputs/$test_file") + end +end + +# Remove temporary files +isdir(results_path) && rm(results_path, force = true, recursive = true) diff --git a/test/test_load_resource_data.jl b/test/test_load_resource_data.jl index d03674fffe..00d979b197 100644 --- a/test/test_load_resource_data.jl +++ b/test/test_load_resource_data.jl @@ -29,8 +29,8 @@ function test_ids_with_positive(attr::Symbol, gen, dfGen) end function prepare_inputs_true(test_path::AbstractString, - in_filenames::InputsTrue, - setup::Dict) + in_filenames::InputsTrue, + setup::Dict) gen_filename = in_filenames.gen_filename inputs_filename = in_filenames.inputs_filename dfGen = GenX.load_dataframe(joinpath(test_path, gen_filename)) @@ -103,17 +103,21 @@ function test_load_scaled_resources_data(gen, dfGen) max_fuels = maximum(GenX.num_fuels.(gen)) for i in 1:max_fuels @test findall(g -> GenX.max_cofire_cols(g, tag = i) < 1, gen[MULTI_FUELS]) == - dfGen[dfGen[!, Symbol(string("fuel", i, "_max_cofire_level"))] .< 1, :][!, + dfGen[dfGen[!, Symbol(string("fuel", i, "_max_cofire_level"))] .< 1, :][ + !, :r_id] @test findall(g -> GenX.max_cofire_start_cols(g, tag = i) < 1, - gen[MULTI_FUELS]) == dfGen[dfGen[!, Symbol(string("fuel", i, "_max_cofire_level_start"))] .< 1, + gen[MULTI_FUELS]) == dfGen[ + dfGen[!, Symbol(string("fuel", i, "_max_cofire_level_start"))] .< 1, :][!, :r_id] @test findall(g -> GenX.min_cofire_cols(g, tag = i) > 0, gen[MULTI_FUELS]) == - dfGen[dfGen[!, Symbol(string("fuel", i, "_min_cofire_level"))] .> 0, :][!, + dfGen[dfGen[!, Symbol(string("fuel", i, "_min_cofire_level"))] .> 0, :][ + !, :r_id] @test findall(g -> GenX.min_cofire_start_cols(g, tag = i) > 0, - gen[MULTI_FUELS]) == dfGen[dfGen[!, Symbol(string("fuel", i, "_min_cofire_level_start"))] .> 0, + gen[MULTI_FUELS]) == dfGen[ + dfGen[!, Symbol(string("fuel", i, "_min_cofire_level_start"))] .> 0, :][!, :r_id] @test GenX.fuel_cols.(gen, tag = i) == dfGen[!, Symbol(string("fuel", i))] diff --git a/test/test_multistage.jl b/test/test_multistage.jl index 9798f3789b..96b135ad08 100644 --- a/test/test_multistage.jl +++ b/test/test_multistage.jl @@ -115,7 +115,8 @@ function test_update_cumulative_min_ret!() for t in 1:3 inpath_sub = joinpath(test_path, "cum_min_ret", string("inputs_p", t)) - true_min_retirements[t] = CSV.read(joinpath(inpath_sub, + true_min_retirements[t] = CSV.read( + joinpath(inpath_sub, "resources", "Resource_multistage_data.csv"), DataFrame) diff --git a/test/utilities.jl b/test/utilities.jl index 300be0f613..a98440ee7b 100644 --- a/test/utilities.jl +++ b/test/utilities.jl @@ -12,9 +12,12 @@ struct CSVFileNotFound <: Exception end Base.showerror(io::IO, e::CSVFileNotFound) = print(io, e.filefullpath, " not found") +const results_path = "results" +!isdir(results_path) && mkdir(results_path) + function run_genx_case_testing(test_path::AbstractString, - test_setup::Dict, - optimizer::Any = HiGHS.Optimizer) + test_setup::Dict, + optimizer::Any = HiGHS.Optimizer) # Merge the genx_setup with the default settings settings = GenX.default_settings() merge!(settings, test_setup) @@ -34,8 +37,8 @@ function run_genx_case_testing(test_path::AbstractString, end function run_genx_case_conflict_testing(test_path::AbstractString, - test_setup::Dict, - optimizer::Any = HiGHS.Optimizer) + test_setup::Dict, + optimizer::Any = HiGHS.Optimizer) # Merge the genx_setup with the default settings settings = GenX.default_settings() @@ -55,8 +58,8 @@ function run_genx_case_conflict_testing(test_path::AbstractString, end function run_genx_case_simple_testing(test_path::AbstractString, - genx_setup::Dict, - optimizer::Any) + genx_setup::Dict, + optimizer::Any) # Run the case OPTIMIZER = configure_solver(test_path, optimizer) inputs = load_inputs(genx_setup, test_path) @@ -66,8 +69,8 @@ function run_genx_case_simple_testing(test_path::AbstractString, end function run_genx_case_multistage_testing(test_path::AbstractString, - genx_setup::Dict, - optimizer::Any) + genx_setup::Dict, + optimizer::Any) # Run the case OPTIMIZER = configure_solver(test_path, optimizer) @@ -90,13 +93,13 @@ function run_genx_case_multistage_testing(test_path::AbstractString, # Step 2) Generate model model_dict[t] = generate_model(genx_setup, inputs_dict[t], OPTIMIZER) end - model_dict, _, inputs_dict = run_ddp(model_dict, genx_setup, inputs_dict) + model_dict, _, inputs_dict = run_ddp(results_path, model_dict, genx_setup, inputs_dict) return model_dict, inputs_dict, OPTIMIZER end function write_testlog(test_path::AbstractString, - message::AbstractString, - test_result::TestResult) + message::AbstractString, + test_result::TestResult) # Save the results to a log file # Format: datetime, message, test result @@ -119,9 +122,9 @@ function write_testlog(test_path::AbstractString, end function write_testlog(test_path::AbstractString, - obj_test::Real, - optimal_tol::Real, - test_result::TestResult) + obj_test::Real, + optimal_tol::Real, + test_result::TestResult) # Save the results to a log file # Format: datetime, objective value ± tolerance, test result message = "$obj_test ± $optimal_tol" @@ -129,9 +132,9 @@ function write_testlog(test_path::AbstractString, end function write_testlog(test_path::AbstractString, - obj_test::Vector{<:Real}, - optimal_tol::Vector{<:Real}, - test_result::TestResult) + obj_test::Vector{<:Real}, + optimal_tol::Vector{<:Real}, + test_result::TestResult) # Save the results to a log file # Format: datetime, [objective value ± tolerance], test result @assert length(obj_test) == length(optimal_tol) diff --git a/test/writing_outputs/test_writing_stats_ms.jl b/test/writing_outputs/test_writing_stats_ms.jl new file mode 100644 index 0000000000..305fa04a77 --- /dev/null +++ b/test/writing_outputs/test_writing_stats_ms.jl @@ -0,0 +1,107 @@ +module TestWritingStatsMs + +using Test +using CSV, DataFrames +using GenX + +# create temporary directory for testing +mkpath("writing_outputs/multi_stage_stats_tmp") +outpath = "writing_outputs/multi_stage_stats_tmp" +filename = GenX._get_multi_stage_stats_filename() + +function test_header() + # Note: if this test fails, it means that the header in the function _get_multi_stage_stats_header() has been changed. + # Make sure to check that the code is consistent with the new header, and update the test accordingly. + header = GenX._get_multi_stage_stats_header() + @test header == + ["Iteration_Number", "Seconds", "Upper_Bound", "Lower_Bound", "Relative_Gap"] +end + +function test_skip_existing_file() + touch(joinpath(outpath, filename)) + # If the file already exists, don't overwrite it + write_multi_stage_stats = GenX.write_multi_stage_stats(outpath, Dict()) + @test isnothing(write_multi_stage_stats) + rm(joinpath(outpath, filename)) +end + +function test_write_multi_stage_stats(iter::Int64 = 10) + # test writing stats to file for `iter` number of iterations + times_a, upper_bounds_a, lower_bounds_a = rand(iter), rand(iter), rand(iter) + stats_d = Dict("TIMES" => times_a, "UPPER_BOUNDS" => upper_bounds_a, + "LOWER_BOUNDS" => lower_bounds_a) + + @test isnothing(GenX.write_multi_stage_stats(outpath, stats_d)) + df_stats = CSV.read(joinpath(outpath, filename), DataFrame) + header = GenX._get_multi_stage_stats_header() + @test size(df_stats) == (iter, length(header)) + for i in 1:iter + test_stats_d(df_stats, i, times_a[i], upper_bounds_a[i], lower_bounds_a[i], + (upper_bounds_a[i] - lower_bounds_a[i]) / lower_bounds_a[i]) + end + rm(joinpath(outpath, filename)) +end + +function test_create_multi_stage_stats_file() + GenX.create_multi_stage_stats_file(outpath) + df_stats = CSV.read(joinpath(outpath, filename), DataFrame) + @test size(df_stats, 1) == 0 + @test size(df_stats, 2) == 5 + @test names(df_stats) == GenX._get_multi_stage_stats_header() + rm(joinpath(outpath, filename)) +end + +function test_update_multi_stage_stats_file(iter::Int64 = 10) + # test updating the stats file with new values + header = GenX._get_multi_stage_stats_header() + GenX.create_multi_stage_stats_file(outpath) + lower_bound = rand() + iteration_time = rand() + for i in 1:iter + # upper bound is updated + upper_bound = rand() + GenX.update_multi_stage_stats_file( + outpath, i, upper_bound, lower_bound, iteration_time, new_row = true) + df_stats = CSV.read(joinpath(outpath, filename), DataFrame) + test_stats_d(df_stats, i, iteration_time, upper_bound, lower_bound, + (upper_bound - lower_bound) / lower_bound) + # lower bound is updated + lower_bound = rand() + GenX.update_multi_stage_stats_file( + outpath, i, upper_bound, lower_bound, iteration_time) + df_stats = CSV.read(joinpath(outpath, filename), DataFrame) + test_stats_d(df_stats, i, iteration_time, upper_bound, lower_bound, + (upper_bound - lower_bound) / lower_bound) + # iteration time is updated + iteration_time = rand() + GenX.update_multi_stage_stats_file( + outpath, i, upper_bound, lower_bound, iteration_time) + df_stats = CSV.read(joinpath(outpath, filename), DataFrame) + test_stats_d(df_stats, i, iteration_time, upper_bound, lower_bound, + (upper_bound - lower_bound) / lower_bound) + # test size + @test size(df_stats) == (i, length(header)) + end + rm(joinpath(outpath, filename)) +end + +function test_stats_d(df_stats, i, iteration_time, upper_bound, lower_bound, relative_gap) + header = GenX._get_multi_stage_stats_header() + @test df_stats[i, header[1]] == i + @test df_stats[i, header[2]] == iteration_time + @test df_stats[i, header[3]] == upper_bound + @test df_stats[i, header[4]] == lower_bound + @test df_stats[i, header[5]] == relative_gap +end + +@testset "Test writing multi-stage stats" begin + test_header() + test_skip_existing_file() + test_write_multi_stage_stats() + test_create_multi_stage_stats_file() + test_update_multi_stage_stats_file() +end + +rm(outpath) + +end # module TestWritingStatsMs