diff --git a/CHANGELOG.md b/CHANGELOG.md index 90cdf5dbef..4159ae7dfb 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,6 +8,11 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## Unreleased ### Added +- Feature CO2 and fuel module (#536) + Adds a fuel module which enables modeling of fuel usage via (1) a constant heat rate and (2) + piecewise-linear approximation of heat rate curves. + Adds a CO2 module that determines the CO2 emissions based on fuel consumption, CO2 capture + fraction, and whether the feedstock is biomass. - Feature electrolysis basic (#525) Adds hydrogen electrolyzer model which enables the addition of hydrogen electrolyzer demands along with optional clean supply constraints. diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/CO2_cap.csv b/Example_Systems/PiecewiseFuel_CO2_Example/CO2_cap.csv new file mode 100644 index 0000000000..4971b2515f --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/CO2_cap.csv @@ -0,0 +1,2 @@ +,Network_zones,CO_2_Cap_Zone_1,CO_2_Max_tons_MWh_1,CO_2_Max_Mtons_1 +NE,z1,1,0,0 diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/Demand_data.csv b/Example_Systems/PiecewiseFuel_CO2_Example/Demand_data.csv new file mode 100644 index 0000000000..01f65be662 --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/Demand_data.csv @@ -0,0 +1,25 @@ +Voll,Demand_Segment,Cost_of_Demand_Curtailment_per_MW,Max_Demand_Curtailment,$/MWh,Rep_Periods,Timesteps_per_Rep_Period,Sub_Weights,Time_Index,Demand_MW_z1 +2000,1,1,1,2000,1,24,8760,1,10000 +,2,0.9,0.04,1800,,,,2,9500 +,3,0.55,0.024,1100,,,,3,9000 +,4,0.2,0.003,400,,,,4,8500 +,,,,,,,,5,8000 +,,,,,,,,6,7500 +,,,,,,,,7,8000 +,,,,,,,,8,8500 +,,,,,,,,9,9000 +,,,,,,,,10,9500 +,,,,,,,,11,10000 +,,,,,,,,12,10500 +,,,,,,,,13,11000 +,,,,,,,,14,11500 +,,,,,,,,15,12000 +,,,,,,,,16,12500 +,,,,,,,,17,13000 +,,,,,,,,18,14000 +,,,,,,,,19,15000 +,,,,,,,,20,16000 +,,,,,,,,21,17000 +,,,,,,,,22,16000 +,,,,,,,,23,15000 +,,,,,,,,24,14000 diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/Fuels_data.csv b/Example_Systems/PiecewiseFuel_CO2_Example/Fuels_data.csv new file mode 100644 index 0000000000..d5b36b7fbf --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/Fuels_data.csv @@ -0,0 +1,26 @@ +Time_Index,NG,Biomass +0,0.05306,0.095 +1,4,8 +2,4,8 +3,4,8 +4,4,8 +5,4,8 +6,4,8 +7,4,8 +8,4,8 +9,4,8 +10,4,8 +11,4,8 +12,4,8 +13,4,8 +14,4,8 +15,4,8 +16,4,8 +17,4,8 +18,4,8 +19,4,8 +20,4,8 +21,4,8 +22,4,8 +23,4,8 +24,4,8 diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/Generators_data.csv b/Example_Systems/PiecewiseFuel_CO2_Example/Generators_data.csv new file mode 100644 index 0000000000..b0c67824f8 --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/Generators_data.csv @@ -0,0 +1,4 @@ +Resource,Zone,THERM,MUST_RUN,STOR,FLEX,HYDRO,VRE,LDS,Num_VRE_Bins,New_Build,Existing_Cap_MW,Existing_Cap_MWh,Existing_Charge_Cap_MW,Max_Cap_MW,Max_Cap_MWh,Max_Charge_Cap_MW,Min_Cap_MW,Inv_Cost_per_MWyr,Inv_Cost_per_MWhyr,Inv_Cost_Charge_per_MWyr,Fixed_OM_Cost_per_MWyr,Fixed_OM_Cost_per_MWhyr,Fixed_OM_Cost_Charge_per_MWyr,Var_OM_Cost_per_MWh,Var_OM_Cost_per_MWh_In,Heat_Rate_MMBTU_per_MWh,Fuel,Cap_Size,Start_Cost_per_MW,Start_Fuel_MMBTU_per_MW,Up_Time,Down_Time,Ramp_Up_Percentage,Ramp_Dn_Percentage,Min_Power,Eff_Up,Eff_Down,Resource_Type,region,cluster,PWFU_Fuel_Usage_Zero_Load_MMBTU_per_h,PWFU_Heat_Rate_MMBTU_per_MWh_1,PWFU_Heat_Rate_MMBTU_per_MWh_2,PWFU_Load_Point_MW_1,PWFU_Load_Point_MW_2,CO2_Capture_Fraction,CO2_Capture_Fraction_Startup,CCS_Disposal_Cost_per_Metric_Ton,Biomass +onshore_wind,1,0,0,0,0,0,1,0,1,0,15000,0,0,-1,-1,-1,0,0,0,0,0,0,0,0,0,0,None,100,0,0,0,0,1,1,0,1,1,onshore_wind_turbine,NE,1,0,0,0,0,0,0,0,0,0 +natural_gas_combined_cycle_ccs,1,1,0,0,0,0,0,0,0,0,15000,0,0,-1,-1,-1,0,0,0,0,0,0,0,5,0,0,NG,250,91,2,6,6,0.64,0.64,0.468,1,1,natural_gas_fired_combined_cycle_ccs,NE,1,400,6,7.2,160,250,0.9,0.6,20,0 +biomass_ccs,1,1,0,0,0,0,0,0,0,0,300,0,0,-1,-1,-1,0,0,0,0,0,0,0,10,0,10,Biomass,300,91,2,6,6,0.64,0.64,0.468,1,1,biomass_ccs,NE,1,0,0,0,0,0,0.9,0.6,20,1 diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/Generators_variability.csv b/Example_Systems/PiecewiseFuel_CO2_Example/Generators_variability.csv new file mode 100644 index 0000000000..9d2a2871bd --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/Generators_variability.csv @@ -0,0 +1,25 @@ +Time_Index,onshore_wind +1,0.889717042 +2,0.877715468 +3,0.903424203 +4,0.895153165 +5,0.757258117 +6,0.630928695 +7,0.557177782 +8,0.6072492 +9,0.423417866 +10,0.007470775 +11,0.002535942 +12,0.002153709 +13,0.00445132 +14,0.007711587 +15,0.100848213 +16,0.201802149 +17,0.141933054 +18,0.567022562 +19,0.946024895 +20,0.923394203 +21,0.953386247 +22,0.929205418 +23,0.849528909 +24,0.665570974 diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/PiecewiseFuelUsage_data_description.pptx b/Example_Systems/PiecewiseFuel_CO2_Example/PiecewiseFuelUsage_data_description.pptx new file mode 100644 index 0000000000..76b8a24af7 Binary files /dev/null and b/Example_Systems/PiecewiseFuel_CO2_Example/PiecewiseFuelUsage_data_description.pptx differ diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/README.md b/Example_Systems/PiecewiseFuel_CO2_Example/README.md new file mode 100644 index 0000000000..0eeec17c6e --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/README.md @@ -0,0 +1,17 @@ +# PiecewiseFuel_CO2 + +**PiecewiseFuel_CO2** is a 24 hr example and contains only one zone. It is condensed for easy and quick testing of CO2, biomass, and piecewise fuel usage related functions of the GenX. This testing +system only includes natural gas ccs, wind, and biomass with ccs, all set at a fixed initial +capacity, and does not allow for building additional capacity. For natural gas ccs generator, we provide picewise fuel usage (PWFU) parameters to represent the fuel consumption at differernt load point. Please refer to "PiecewiseFuelUsage_data_description.pptx" for a sylized visual representation of PWFU segments and corresponding data requirements. When UC >= 1 and PWFU parameters are provided in "Generator_data.csv", the standard heat rate (i.e., Heat_Rate_MMBTU_per_MWh) will not be used unless UC = 0. + +To run the model, first navigate to the example directory at `GenX/Example_Systems/PiecewiseFuel_CO2`: + +`cd("Example_Systems/PiecewiseFuel_CO2")` + +Next, ensure that your settings in `GenX_settings.yml` are correct. The linear clustering unit commitment method (UC = 2) is used. The default settings use the solver HiGHS (`Solver: HiGHS`). A mass-based carbon cap of 0 t CO2 (net-zero) is specified in the `CO2_cap.csv` input file. + +Once the settings are confirmed, run the model with the `Run.jl` script in the example directory: + +`include("Run.jl")` + +Once the model has completed, results will write to the `Results` directory. diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/Run.jl b/Example_Systems/PiecewiseFuel_CO2_Example/Run.jl new file mode 100644 index 0000000000..b44ca23ec1 --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/Run.jl @@ -0,0 +1,3 @@ +using GenX + +run_genx_case!(dirname(@__FILE__)) diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/Settings/genx_settings.yml b/Example_Systems/PiecewiseFuel_CO2_Example/Settings/genx_settings.yml new file mode 100644 index 0000000000..53646d1cd3 --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/Settings/genx_settings.yml @@ -0,0 +1,101 @@ +## Model settings +# -------------- + +# Number of segments used in piecewise linear approximation of transmission losses; +# 1 = linear +# >2 = piecewise quadratic +Trans_Loss_Segments: 1 + +# Unit committment of thermal power plants; +# 0 = not active +# 1 = active using integer clustering +# 2 = active using linearized clustering +UCommit: 2 + +## Policy settings +# --------------- + +# Regulation (primary) and operating (secondary) reserves +# 0 = not active +# 1 = active systemwide +# Reserves: 0 + +# Minimum qualifying renewables penetration +# 0 = not active +# 1 = active systemwide +# EnergyShareRequirement: 0 + +# Number of capacity reserve margin constraints +# 0 = not active +# 1 = active systemwide +# CapacityReserveMargin: 0 + +# CO2 emissions cap +# 0 = not active (no CO2 emission limit) +# 1 = mass-based emission limit constraint +# 2 = demand + rate-based emission limit constraint +# 3 = generation + rate-based emission limit constraint +CO2Cap: 1 + +# Let lost energy be accounted for in Energy Share Requirement and CO2 constraints +# 0 = not active (DO NOT account for energy lost); +# 1 = active systemwide (DO account for energy lost) +# StorageLosses: 0 + +# Activate minimum technology carveout constraints +# 0 = not active +# 1 = active +# MinCapReq: 0 + +# Activate maximum technology limits +# 0 = not active +# 1 = active +# MaxCapReq: 0 +# +# Transmission network expansion +# 0 = not active +# 1 = active systemwide +# NetworkExpansion: 1 + +# Hydrogen electrolyzer hourly supply matching required? +# 0 = constraint ot active +# 1 = constraint active for all electrolyzer resources +#HydrogenHourlyMatching: 0 + +## Solver settings +# ---------------- + +# Available solvers: Gurobi, CPLEX, Clps +#Solver: HiGHS + +# Internally scale demand, capacity and power variables in GW rather than MW. +# 0 = not active +# 1 = active systemwide +ParameterScale: 1 + +# Time domain reduce (i.e. cluster) inputs based on demand, resource availability profiles, and fuel prices; +# 0 = not active (use input data as provided) +# 1 = active (cluster input data, or use data that has already been clustered) +#TimeDomainReduction: 0 + +# Directory where results from time domain reduction will be saved. +# If results already exist, these will be used without running time domain reduction script again. +#TimeDomainReductionFolder: "TDR_Results" + +## Output settings +# ---------------- + +# Write shadow prices of LP or relaxed MILP +# 0 = not active +# 1 = active +WriteShadowPrices: 1 + +# Overwrite existing results in output folder or create a new one; +# 0 = create new folder +# 1 = overwrite existing results +# OverwriteResults: 0 + +# Write the model formulation as an output; +# 0 = active +# 1 = not active +# PrintModel: 0 diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/Settings/gurobi_settings.yml b/Example_Systems/PiecewiseFuel_CO2_Example/Settings/gurobi_settings.yml new file mode 100644 index 0000000000..0aa53a56d3 --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/Settings/gurobi_settings.yml @@ -0,0 +1,15 @@ +# Gurobi Solver Parameters +# Common solver settings +Feasib_Tol: 1.0e-05 # Constraint (primal) feasibility tolerances. +Optimal_Tol: 1e-5 # Dual feasibility tolerances. +TimeLimit: 110000 # Limits total time solver. +Pre_Solve: 1 # Controls presolve level. +Method: 2 # Algorithm used to solve continuous models (including MIP root relaxation). + +#Gurobi-specific solver settings +MIPGap: 1e-3 # Relative (p.u. of optimal) mixed integer optimality tolerance for MIP problems (ignored otherwise). +BarConvTol: 1.0e-08 # Barrier convergence tolerance (determines when barrier terminates). +NumericFocus: 0 # Numerical precision emphasis. +Crossover: 0 # Barrier crossver strategy. +PreDual: 0 # Decides whether presolve should pass the primal or dual linear programming problem to the LP optimization algorithm. +AggFill: 10 # Allowed fill during presolve aggregation. diff --git a/Example_Systems/PiecewiseFuel_CO2_Example/Settings/highs_settings.yml b/Example_Systems/PiecewiseFuel_CO2_Example/Settings/highs_settings.yml new file mode 100644 index 0000000000..bb4b286074 --- /dev/null +++ b/Example_Systems/PiecewiseFuel_CO2_Example/Settings/highs_settings.yml @@ -0,0 +1,13 @@ +# HiGHS Solver Parameters +# Common solver settings +Feasib_Tol: 1.0e-05 # Primal feasibility tolerance # [type: double, advanced: false, range: [1e-10, inf], default: 1e-07] +Optimal_Tol: 1.0e-05 # Dual feasibility tolerance # [type: double, advanced: false, range: [1e-10, inf], default: 1e-07] +TimeLimit: 1.0e23 # Time limit # [type: double, advanced: false, range: [0, inf], default: inf] +Pre_Solve: choose # Presolve option: "off", "choose" or "on" # [type: string, advanced: false, default: "choose"] +Method: choose #HiGHS-specific solver settings # Solver option: "simplex", "choose" or "ipm" # [type: string, advanced: false, default: "choose"] + +#highs-specific solver settings + +# run the crossover routine for ipx +# [type: string, advanced: "on", range: {"off", "on"}, default: "off"] +run_crossover: "off" diff --git a/docs/src/core.md b/docs/src/core.md index af2b19cca2..00bca332e2 100644 --- a/docs/src/core.md +++ b/docs/src/core.md @@ -29,8 +29,14 @@ Pages = ["transmission.jl"] Modules = [GenX] Pages = ["ucommit.jl"] ``` -## Emissions +## CO2 ```@autodocs Modules = [GenX] -Pages = ["emissions.jl"] +Pages = ["co2.jl"] +``` + +## Fuel +```@autodocs +Modules = [GenX] +Pages = ["fuel.jl"] ``` diff --git a/docs/src/data_documentation.md b/docs/src/data_documentation.md index e50050351a..068b584302 100644 --- a/docs/src/data_documentation.md +++ b/docs/src/data_documentation.md @@ -384,11 +384,24 @@ This file contains cost and performance parameters for various generators and ot |MinCapTag\_*| Eligibility of resources to participate in Minimum Technology Carveout constraint. \* corresponds to the ith row of the file `Minimum_capacity_requirement.csv`. Note that this eligibility must be 0 for co-located VRE-STOR resources (policy inputs are read from the specific VRE-STOR dataframe).| |**MaxCapReq = 1**| |MaxCapTag\_*| Eligibility of resources to participate in Maximum Technology Carveout constraint. \* corresponds to the ith row of the file `Maximum_capacity_requirement.csv`. Note that this eligibility must be 0 for co-located VRE-STOR resources (policy inputs are read from the specific VRE-STOR dataframe).| +|**PiecewiseFuelUsage-related parameters required if any resources have nonzero PWFU fuel usage, heat rates, and load points**| +|PWFU\_Fuel\_Usage\_Zero\_Load\_MMBTU\_per\_h|The fuel usage (MMBTU/h) for the first PWFU segemnt (y-intercept) at zero load.| +|PWFU\_Heat\_Rate\_MMBTU\_per\_MWh\_*i| The slope of fuel usage function of the segment i.| +|PWFU\_Load\_Point\_MW\_*i| The end of segment i (MW).| |**Electrolyzer related parameters required if the set ELECTROLYZER is not empty**| |Hydrogen_MWh_Per_Tonne| Electrolyzer efficiency in megawatt-hours (MWh) of electricity per metric tonne of hydrogen produced (MWh/t)| |Electrolyzer_Min_kt| Minimum annual quantity of hydrogen that must be produced by electrolyzer in kilotonnes (kt)| |Hydrogen_Price_Per_Tonne| Price (or value) of hydrogen per metric tonne ($/t)| |Qualified_Hydrogen_Supply| {0,1}, Indicates that generator or storage resources is eligible to supply electrolyzers in the same zone (used for hourly clean supply constraint)| +|**CO2-related parameters required if any resources have nonzero CO2_Capture_Fraction**| +|CO2\_Capture\_Fraction |[0,1], The CO2 capture fraction of CCS-equipped power plants during steady state operation. This value should be 0 for generators without CCS. | +|CO2\_Capture\_Fraction\_Startup |[0,1], The CO2 capture fraction of CCS-equipped power plants during the startup events. This value should be 0 for generators without CCS | +|Biomass | {0, 1}, Flag to indicate if generator uses biomass as feedstock (optional input column).| +||Biomass = 0: Not part of set (default). | +||Biomass = 1: Uses biomass as fuel.| +|CCS\_Disposal\_Cost\_per\_Metric_Ton | Cost associated with CCS disposal ($/tCO2), including pipeline, injection and storage costs of CCS-equipped generators.| + + ### 2.2 Optional inputs files diff --git a/src/configure_settings/configure_settings.jl b/src/configure_settings/configure_settings.jl index 352facdae4..d95e3f390f 100644 --- a/src/configure_settings/configure_settings.jl +++ b/src/configure_settings/configure_settings.jl @@ -23,7 +23,7 @@ function default_settings() "MethodofMorris" => 0, "IncludeLossesInESR" => 0, "EnableJuMPStringNames" => false, - "HydrogenHourlyMatching" => 0, + "HydrogenHourlyMatching" => 0 ) end diff --git a/src/load_inputs/load_fuels_data.jl b/src/load_inputs/load_fuels_data.jl index 4dee038a4a..fe15ea824c 100644 --- a/src/load_inputs/load_fuels_data.jl +++ b/src/load_inputs/load_fuels_data.jl @@ -1,6 +1,5 @@ @doc raw""" - load_fuels_data!(setup::Dict, path::AbstractString, inputs::Dict) - +load_fuels_data!(setup::Dict, path::AbstractString, inputs::Dict) Read input parameters related to fuel costs and CO$_2$ content of fuels """ function load_fuels_data!(setup::Dict, path::AbstractString, inputs::Dict) @@ -32,9 +31,10 @@ function load_fuels_data!(setup::Dict, path::AbstractString, inputs::Dict) scale_factor = setup["ParameterScale"] == 1 ? ModelScalingFactor : 1 for i = 1:length(fuels) + # fuel cost is in $/MMBTU w/o scaling, $/Billon BTU w/ scaling fuel_costs[fuels[i]] = costs[:,i] / scale_factor - # fuel_CO2 is kton/MMBTU with scaling, or ton/MMBTU without scaling. - fuel_CO2[fuels[i]] = CO2_content[i] / scale_factor + # No need to scale fuel_CO2, fuel_CO2 is ton/MMBTU or kton/Billion BTU + fuel_CO2[fuels[i]] = CO2_content[i] end inputs["fuels"] = fuels diff --git a/src/load_inputs/load_generators_data.jl b/src/load_inputs/load_generators_data.jl index 106feaf964..0ee509501c 100644 --- a/src/load_inputs/load_generators_data.jl +++ b/src/load_inputs/load_generators_data.jl @@ -185,44 +185,39 @@ function load_generators_data!(setup::Dict, path::AbstractString, inputs_gen::Di end end - if setup["UCommit"]>=1 - # Fuel consumed on start-up (million BTUs per MW per start) if unit commitment is modelled - start_fuel = convert(Array{Float64}, gen_in[!,:Start_Fuel_MMBTU_per_MW]) - # Fixed cost per start-up ($ per MW per start) if unit commitment is modelled + if setup["UCommit"] >= 1 start_cost = convert(Array{Float64}, gen_in[!,:Start_Cost_per_MW]) inputs_gen["C_Start"] = zeros(Float64, G, inputs_gen["T"]) - gen_in[!,:CO2_per_Start] = zeros(Float64, G) end - # Heat rate of all resources (million BTUs/MWh) - heat_rate = convert(Array{Float64}, gen_in[!,:Heat_Rate_MMBTU_per_MWh]) - # Fuel used by each resource - fuel_type = gen_in[!,:Fuel] - # Maximum fuel cost in $ per MWh and CO2 emissions in tons per MWh - inputs_gen["C_Fuel_per_MWh"] = zeros(Float64, G, inputs_gen["T"]) - gen_in[!,:CO2_per_MWh] = zeros(Float64, G) + # scale the start costs for g in 1:G - # NOTE: When Setup[ParameterScale] =1, fuel costs are scaled in fuels_data.csv, so no if condition needed to scale C_Fuel_per_MWh - inputs_gen["C_Fuel_per_MWh"][g,:] = fuel_costs[fuel_type[g]].*heat_rate[g] - gen_in[g,:CO2_per_MWh] = fuel_CO2[fuel_type[g]]*heat_rate[g] - gen_in[g,:CO2_per_MWh] *= scale_factor - # kton/MMBTU * MMBTU/MWh = kton/MWh, to get kton/GWh, we need to mutiply 1000 if g in inputs_gen["COMMIT"] - # Start-up cost is sum of fixed cost per start plus cost of fuel consumed on startup. - # CO2 from fuel consumption during startup also calculated - - inputs_gen["C_Start"][g,:] = gen_in[g,:Cap_Size] * (fuel_costs[fuel_type[g]] .* start_fuel[g] .+ start_cost[g]) - # No need to re-scale C_Start since Cap_size, fuel_costs and start_cost are scaled When Setup[ParameterScale] =1 - Dharik - gen_in[g,:CO2_per_Start] = gen_in[g,:Cap_Size]*(fuel_CO2[fuel_type[g]]*start_fuel[g]) - gen_in[g,:CO2_per_Start] *= scale_factor - # Setup[ParameterScale] =1, gen_in[g,:Cap_Size] is GW, fuel_CO2[fuel_type[g]] is ktons/MMBTU, start_fuel is MMBTU/MW, - # thus the overall is MTons/GW, and thus gen_in[g,:CO2_per_Start] is Mton, to get kton, change we need to multiply 1000 - # Setup[ParameterScale] =0, gen_in[g,:Cap_Size] is MW, fuel_CO2[fuel_type[g]] is tons/MMBTU, start_fuel is MMBTU/MW, - # thus the overall is MTons/GW, and thus gen_in[g,:CO2_per_Start] is ton + # Start-up cost is sum of fixed cost per start startup. + inputs_gen["C_Start"][g,:] .= gen_in[g,:Cap_Size] * ( start_cost[g]) end end load_vre_stor_data!(inputs_gen, setup, path) + + + # write zeros if col names are not in the gen_in dataframe + required_cols_for_co2 = ["Biomass", "CO2_Capture_Fraction", "CO2_Capture_Fraction_Startup", "CCS_Disposal_Cost_per_Metric_Ton"] + for col in required_cols_for_co2 + ensure_column!(gen_in, col, 0) + end + + # Scale CCS_Disposal_Cost_per_Metric_Ton for CCS units + gen_in.CCS_Disposal_Cost_per_Metric_Ton /= scale_factor + + # get R_ID when fuel is not None + inputs_gen["HAS_FUEL"] = gen_in[(gen_in[!,:Fuel] .!= "None"),:R_ID] + + # Piecewise fuel usage option + if setup["UCommit"] > 0 + process_piecewisefuelusage!(inputs_gen, scale_factor) + end + println(filename * " Successfully Read!") end @@ -501,4 +496,107 @@ function load_vre_stor_data!(inputs_gen::Dict, setup::Dict, path::AbstractString inputs_gen["dfVRE_STOR"] = DataFrame() end summarize_errors(error_strings) +end + + +function process_piecewisefuelusage!(inputs::Dict, scale_factor) + gen_in = inputs["dfGen"] + inputs["PWFU_Num_Segments"] = 0 + inputs["THERM_COMMIT_PWFU"] = Int64[] + + if any(occursin.(Ref("PWFU_"), names(gen_in))) + heat_rate_mat = extract_matrix_from_dataframe(gen_in, "PWFU_Heat_Rate_MMBTU_per_MWh") + load_point_mat = extract_matrix_from_dataframe(gen_in, "PWFU_Load_Point_MW") + + # check data input + validate_piecewisefuelusage(heat_rate_mat, load_point_mat) + + # determine if a generator contains piecewise fuel usage segment based on non-zero heatrate + gen_in.HAS_PWFU = any(heat_rate_mat .!= 0 , dims = 2)[:] + num_segments = size(heat_rate_mat)[2] + + # translate the inital fuel usage, heat rate, and load points into intercept for each segment + fuel_usage_zero_load = gen_in[!,"PWFU_Fuel_Usage_Zero_Load_MMBTU_per_h"] + # construct a matrix for intercept + intercept_mat = zeros(size(heat_rate_mat)) + # PWFU_Fuel_Usage_MMBTU_per_h is always the intercept of the first segment + intercept_mat[:,1] = fuel_usage_zero_load + + # create a function to compute intercept if we have more than one segment + function calculate_intercepts(slope, intercept_1, load_point) + m, n = size(slope) + # Initialize the intercepts matrix with zeros + intercepts = zeros(m, n) + # The first segment's intercepts should be intercept_1 vector + intercepts[:, 1] = intercept_1 + # Calculate intercepts for the other segments using the load points (i.e., intersection points) + for j in 1:n-1 + for i in 1:m + current_slope = slope[i, j+1] + previous_slope = slope[i, j] + # If the current slope is 0, then skip the calculation and return 0 + if current_slope == 0 + intercepts[i, j+1] = 0.0 + else + # y = a*x + b; => b = y - ax + # Calculate y-coordinate of the intersection + y = previous_slope * load_point[i, j] + intercepts[i, j] + # determine the new intercept + b = y - current_slope * load_point[i, j] + intercepts[i, j+1] = b + end + end + end + return intercepts + end + + if num_segments > 1 + # determine the intercept for the rest of segment if num_segments > 1 + intercept_mat = calculate_intercepts(heat_rate_mat, fuel_usage_zero_load, load_point_mat) + end + + # 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),"PWFU_Heat_Rate_MMBTU_per_MWh"),names(gen_in))) + slope_df = DataFrame(heat_rate_mat, Symbol.(slope_cols)) + PWFU_data = hcat(slope_df, intercept_df) + # no need to scale sclope, but intercept should be scaled when parameterscale is on (MMBTU -> billion BTU) + PWFU_data[!, intercept_cols] ./= scale_factor + + inputs["slope_cols"] = slope_cols + inputs["intercept_cols"] = intercept_cols + inputs["PWFU_data"] = PWFU_data + inputs["PWFU_Num_Segments"] =num_segments + inputs["THERM_COMMIT_PWFU"] = intersect(gen_in[gen_in.THERM.==1,:R_ID], gen_in[gen_in.HAS_PWFU,:R_ID]) + end +end + +function validate_piecewisefuelusage(heat_rate_mat, load_point_mat) + # it's possible to construct piecewise fuel consumption with n of heat rate and n-1 of load point. + # if a user feed n of heat rate and more than n of load point, throw a error message, and then use + # n of heat rate and n-1 load point to construct the piecewise fuel usage fuction + if size(heat_rate_mat)[2] < size(load_point_mat)[2] + @error """ The numbers of heatrate data are less than load points, we found $(size(heat_rate_mat)[2]) of heat rate, + and $(size(load_point_mat)[2]) of load points. We will just use $(size(heat_rate_mat)[2]) of heat rate, and $(size(heat_rate_mat)[2]-1) + load point to create piecewise fuel usage + """ + end + + # check if values for piecewise fuel consumption make sense. Negative heat rate or load point are not allowed + if any(heat_rate_mat .< 0) | any(load_point_mat .< 0) + @error """ Neither heat rate nor load point can be negative + """ + error("Invalid inputs detected for piecewise fuel usage") + end + # for non-zero values, heat rates and load points should follow an increasing trend + if any([any(diff(filter(!=(0), row)) .< 0) for row in eachrow(heat_rate_mat)]) + @error """ Heat rates should follow an increasing trend + """ + error("Invalid inputs detected for piecewise fuel usage") + elseif any([any(diff(filter(!=(0), row)) .< 0) for row in eachrow(load_point_mat)]) + @error """load points should follow an increasing trend + """ + error("Invalid inputs detected for piecewise fuel usage") + end end \ No newline at end of file diff --git a/src/model/core/co2.jl b/src/model/core/co2.jl new file mode 100644 index 0000000000..3d4e14edf0 --- /dev/null +++ b/src/model/core/co2.jl @@ -0,0 +1,105 @@ +@doc raw""" + co2!(EP::Model, inputs::Dict) + +This function creates expressions to account for CO2 emissions as well as captured and sequestrated +CO2 from thermal generators. It also has the capability to model the negative CO2 emissions +from bioenergy with carbon capture and storage. + +***** Expressions ***** + +For thermal generators which combust fuels (e.g., coal, natural gas, and biomass), the net CO2 +emission to the environment is a function of fuel consumption, CO2 emission factor, CO2 capture +fraction, and whether the feedstock is biomass. Biomass is a factor in this equation because +biomass generators are assumed to generate zero net CO2 emissions, or negative net CO2 emissions +in the case that the CO2 they emit is captured and sequestered underground. + +If a user wishes to represent a generator that combusts biomass, then in Generators_data.csv, +the "Biomass" column (boolean, 1 or 0), which represents if a generator $y$ uses biomass or not, should be set to 1. +The CO2 emissions from such a generator will be assumed to be zero without CCS and negative with CCS. + +The CO2 emissions from generator $y$ at time $t$ are determined by total fuel +consumption (MMBTU) multiplied by the CO2 content of the fuel (tCO2/MMBTU), and by +(1 - Biomass [0 or 1] - CO2 capture fraction [a fraction, between 0 - 1]). +The CO2 capture fraction could be differernt during the steady-state and startup events +(generally startup events have a lower CO2 capture fraction), so we use distinct CO2 capture fractions +to determine the emissions. +In short, the CO2 emissions for a generator depend on the CO2 emission factor from fuel combustion, +the CO2 capture fraction, and whether the generator uses biomass. + +```math +\begin{aligned} +eEmissionsByPlant_{g,t} = (1-Biomass_y- CO2\_Capture\_Fraction_y) * vFuel_{y,t} * CO2_{content} + (1-Biomass_y- CO2\_Capture\_Fraction\_Startup_y) * eStartFuel_{y,t} * CO2_{content} +\hspace{1cm} \forall y \in G, \forall t \in T, Biomass_y \in {{0,1}} +\end{aligned} +``` + +Where $Biomass_y$ represents a binary variable (1 or 0) that determines if the generator $y$ +uses biomass, and $CO2\_Capture\_Fraction_y$ represents a fraction for CO2 capture rate. + +In addition to CO2 emissions, for generators with a non-zero CO2 capture rate, we also +determine the amount of CO2 being captured and sequestered. The CO2 emissions from +generator $y$ at time $t$, denoted by $eEmissionsCaptureByPlant_{g,t}$, are determined by +total fuel consumption (MMBTU) multiplied by the $CO_2$ content of the fuel (tCO2/MMBTU), +times CO2 capture rate. + +```math +\begin{aligned} +eEmissionsCaptureByPlant_{g,t} = CO2\_Capture\_Fraction_y * vFuel_{y,t} * CO2_{content} + CO2\_Capture\_Fraction\_Startup_y * eStartFuel_{y,t} * CO2_{content} +\hspace{1cm} \forall y \in G, \forall t \in T +\end{aligned} +``` + +""" +function co2!(EP::Model, inputs::Dict) + + println("CO2 Module") + + dfGen = inputs["dfGen"] + G = inputs["G"] # Number of resources (generators, storage, DR, and DERs) + T = inputs["T"] # Number of time steps (hours) + Z = inputs["Z"] # Number of zones + fuel_CO2 = inputs["fuel_CO2"] # CO2 content of fuel (t CO2/MMBTU or ktCO2/Billion BTU) + ### Expressions ### + # CO2 emissions from power plants in "Generators_data.csv" + # If all the CO2 capture fractions from Generators_data are zeros, the CO2 emissions from thermal generators are determined by fuel consumption times CO2 content per MMBTU + if all(dfGen.CO2_Capture_Fraction .==0) + @expression(EP, eEmissionsByPlant[y=1:G, t=1:T], + ((1-dfGen[y, :Biomass]) *(EP[:vFuel][y, t] + EP[:eStartFuel][y, t]) * fuel_CO2[dfGen[y,:Fuel]])) + else + @info "Using the CO2 module to determine the CO2 emissions of CCS-equipped plants" + # CO2_Capture_Fraction refers to the CO2 capture rate of CCS equiped power plants at a steady state + # CO2_Capture_Fraction_Startup refers to the CO2 capture rate of CCS equiped power plants during startup events + @expression(EP, eEmissionsByPlant[y=1:G, t=1:T], + (1-dfGen[y, :Biomass] - dfGen[y, :CO2_Capture_Fraction]) * EP[:vFuel][y, t] * fuel_CO2[dfGen[y,:Fuel]]+ + (1-dfGen[y, :Biomass] - dfGen[y, :CO2_Capture_Fraction_Startup]) * EP[:eStartFuel][y, t] * fuel_CO2[dfGen[y,:Fuel]]) + + # CO2 captured from power plants in "Generators_data.csv" + @expression(EP, eEmissionsCaptureByPlant[y=1:G, t=1:T], + dfGen[y, :CO2_Capture_Fraction] * EP[:vFuel][y, t] * fuel_CO2[dfGen[y,:Fuel]]+ + dfGen[y, :CO2_Capture_Fraction_Startup] * EP[:eStartFuel][y, t] * fuel_CO2[dfGen[y,:Fuel]]) + + @expression(EP, eEmissionsCaptureByPlantYear[y=1:G], + sum(inputs["omega"][t] * eEmissionsCaptureByPlant[y, 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=1:G], + sum(inputs["omega"][t] * eEmissionsCaptureByPlant[y, t] * + dfGen[y, :CCS_Disposal_Cost_per_Metric_Ton] for t in 1:T)) + + @expression(EP, eZonalCCO2Sequestration[z=1:Z], + sum(ePlantCCO2Sequestration[y] + for y in dfGen[(dfGen[!, :Zone].==z), :R_ID])) + + @expression(EP, eTotaleCCO2Sequestration, + sum(eZonalCCO2Sequestration[z] for z in 1:Z)) + + add_to_expression!(EP[:eObj], EP[:eTotaleCCO2Sequestration]) + end + + # emissions by zone + @expression(EP, eEmissionsByZone[z = 1:Z, t = 1:T], + sum(eEmissionsByPlant[y, t] for y in dfGen[(dfGen[!, :Zone].==z), :R_ID])) + return EP + +end diff --git a/src/model/core/discharge/discharge.jl b/src/model/core/discharge/discharge.jl index e226e0302f..977ea4f04c 100644 --- a/src/model/core/discharge/discharge.jl +++ b/src/model/core/discharge/discharge.jl @@ -1,11 +1,11 @@ @doc raw""" discharge(EP::Model, inputs::Dict, setup::Dict) This module defines the power decision variable $\Theta_{y,t} \forall y \in \mathcal{G}, t \in \mathcal{T}$, representing energy injected into the grid by resource $y$ by at time period $t$. -This module additionally defines contributions to the objective function from variable costs of generation (variable O&M plus fuel cost) from all resources $y \in \mathcal{G}$ over all time periods $t \in \mathcal{T}$: +This module additionally defines contributions to the objective function from variable costs of generation (variable O&M) from all resources $y \in \mathcal{G}$ over all time periods $t \in \mathcal{T}$: ```math \begin{aligned} Obj_{Var\_gen} = - \sum_{y \in \mathcal{G} } \sum_{t \in \mathcal{T}}\omega_{t}\times(\pi^{VOM}_{y} + \pi^{FUEL}_{y})\times \Theta_{y,t} + \sum_{y \in \mathcal{G} } \sum_{t \in \mathcal{T}}\omega_{t}\times(\pi^{VOM}_{y})\times \Theta_{y,t} \end{aligned} ``` """ @@ -27,9 +27,8 @@ function discharge!(EP::Model, inputs::Dict, setup::Dict) ## Objective Function Expressions ## - # Variable costs of "generation" for resource "y" during hour "t" = variable O&M plus fuel cost - @expression(EP, eCVar_out[y=1:G,t=1:T], (inputs["omega"][t]*(dfGen[y,:Var_OM_Cost_per_MWh]+inputs["C_Fuel_per_MWh"][y,t])*vP[y,t])) - #@expression(EP, eCVar_out[y=1:G,t=1:T], (round(inputs["omega"][t]*(dfGen[y,:Var_OM_Cost_per_MWh]+inputs["C_Fuel_per_MWh"][y,t]), digits=RD)*vP[y,t])) + # Variable costs of "generation" for resource "y" during hour "t" = variable O&M + @expression(EP, eCVar_out[y=1:G,t=1:T], (inputs["omega"][t]*(dfGen[y,:Var_OM_Cost_per_MWh]*vP[y,t]))) # Sum individual resource contributions to variable discharging costs to get total variable discharging costs @expression(EP, eTotalCVarOutT[t=1:T], sum(eCVar_out[y,t] for y in 1:G)) @expression(EP, eTotalCVarOut, sum(eTotalCVarOutT[t] for t in 1:T)) diff --git a/src/model/core/emissions.jl b/src/model/core/emissions.jl deleted file mode 100644 index 0bb007463f..0000000000 --- a/src/model/core/emissions.jl +++ /dev/null @@ -1,26 +0,0 @@ -@doc raw""" - emissions(EP::Model, inputs::Dict) - -This function creates expression to add the CO2 emissions by plants in each zone, which is subsequently added to the total emissions -""" -function emissions!(EP::Model, inputs::Dict) - - println("Emissions Module (for CO2 Policy modularization)") - - dfGen = inputs["dfGen"] - - G = inputs["G"] # Number of resources (generators, storage, DR, and DERs) - T = inputs["T"] # Number of time steps (hours) - Z = inputs["Z"] # Number of zones - - @expression(EP, eEmissionsByPlant[y=1:G,t=1:T], - - if y in inputs["COMMIT"] - dfGen[y,:CO2_per_MWh]*EP[:vP][y,t]+dfGen[y,:CO2_per_Start]*EP[:vSTART][y,t] - else - dfGen[y,:CO2_per_MWh]*EP[:vP][y,t] - end - ) - @expression(EP, eEmissionsByZone[z=1:Z, t=1:T], sum(eEmissionsByPlant[y,t] for y in dfGen[(dfGen[!,:Zone].==z),:R_ID])) - -end diff --git a/src/model/core/fuel.jl b/src/model/core/fuel.jl new file mode 100644 index 0000000000..edb424c260 --- /dev/null +++ b/src/model/core/fuel.jl @@ -0,0 +1,154 @@ + +@doc raw""" + fuel!(EP::Model, inputs::Dict, setup::Dict) + +This function creates expressions to account for total fuel consumption (e.g., coal, +natural gas, hydrogen, etc). It also has the capability to model heat rates that are +a function of load via a piecewise-linear approximation. + +***** Expressions ****** +Users have two options to model the fuel consumption as a function of power generation: +(1). Use a constant heat rate, regardless of the minimum load or maximum load; and +(2). Use the PiecewiseFuelUsage-related parameters to model the fuel consumption via a +piecewise-linear approximation of the heat rate curves. By using this option, users can represent +the fact that most generators have a decreasing heat rate as a function of load. + +(1). Constant heat rate. +The fuel consumption for power generation $vFuel_{y,t}$ is determined by power generation +($vP_{y,t}$) mutiplied by the corresponding heat rate ($Heat\_Rate_y$). +The fuel costs for power generation and start fuel for a plant $y$ at time $t$, +denoted by $eCFuelOut_{y,t}$ and $eFuelStart$, are determined by fuel consumption ($vFuel_{y,t}$ +and $eStartFuel$) multiplied by the fuel costs (\$/MMBTU) +(2). Piecewise-linear approximation +With this formulation, the heat rate of generators becomes a function of load. +In reality this relationship takes a nonlinear form, but we model it +through a piecewise-linear approximation: + +```math +\begin{aligned} +vFuel_{y,t} >= vP_{y,t} * h_{y,x} + U_{g,t}* f_{y,x} +\hspace{1cm} \forall y \in G, \forall t \in T, \forall x \in X +\end{aligned} +``` +Where $h_{y,x}$ represents the heat rate slope for generator $y$ in segment $x$ [MMBTU/MWh], + $f_{y,x}$ represents the heat rate intercept (MMBTU) for a generator $y$ in segment $x$ [MMBTU], +and $U_{y,t}$ represents the commitment status of a generator $y$ at time $t$. These parameters +are optional inputs to "Generators_data.csv". +When Unit commitment is on, if a user provides slope and intercept, the standard heat rate +(i.e., Heat_Rate_MMBTU_per_MWh) will not be used. When unit commitment is off, the model will +always use the standard heat rate. +The user should determine the slope and intercept parameters based on the Cap_Size of the plant. +For example, when a plant is operating at the full load (i.e., power output equal to the Cap_Size), +the fuel usage determined by the effective segment divided by Cap_Size should be equal to the +heat rate at full-load. + +Since fuel consumption and fuel costs are postive, the optimization will force the fuel usage +to be equal to the highest fuel usage segment for any given value of vP. +When the power output is zero, the commitment variable $U_{g,t}$ will bring the intercept +to be zero such that the fuel consumption is zero when thermal units are offline. + +In order to run piecewise fuel consumption module, +the unit commitment must be turned on (UC = 1 or 2), and users should provide PWFU_Slope_* and +PWFU_Intercept_* for at least one segment. +""" + +function fuel!(EP::Model, inputs::Dict, setup::Dict) + println("Fuel Module") + dfGen = inputs["dfGen"] + T = inputs["T"] # Number of time steps (hours) + Z = inputs["Z"] # Number of zones + G = inputs["G"] + HAS_FUEL = inputs["HAS_FUEL"] + THERM_COMMIT = inputs["THERM_COMMIT"] + fuels = inputs["fuels"] + NUM_FUEL = length(fuels) + + + # create variable for fuel consumption for output + @variable(EP, vFuel[y in 1:G, t = 1:T] >= 0) + + ### Expressions #### + # Fuel consumed on start-up (MMBTU or kMMBTU (scaled)) + # if unit commitment is modelled + @expression(EP, eStartFuel[y in 1:G, t = 1:T], + if y in THERM_COMMIT + (dfGen[y,:Cap_Size] * EP[:vSTART][y, t] * + dfGen[y,:Start_Fuel_MMBTU_per_MW]) + else + 0 + end) + + # fuel_cost is in $/MMBTU (M$/billion BTU if scaled) + # vFuel and eStartFuel is MMBTU (or billion BTU if scaled) + # eCFuel_start or eCFuel_out is $ or Million$ + # Start up fuel cost + @expression(EP, eCFuelStart[y = 1:G, t = 1:T], + (inputs["fuel_costs"][dfGen[y,:Fuel]][t] * EP[:eStartFuel][y, t])) + # plant level start-up fuel cost for output + @expression(EP, ePlantCFuelStart[y = 1:G], + sum(inputs["omega"][t] * EP[:eCFuelStart][y, t] for t in 1:T)) + # zonal level total fuel cost for output + @expression(EP, eZonalCFuelStart[z = 1:Z], + sum(EP[:ePlantCFuelStart][y] for y in dfGen[dfGen[!, :Zone].==z, :R_ID])) + + # Fuel cost for power generation + @expression(EP, eCFuelOut[y = 1:G, t = 1:T], + (inputs["fuel_costs"][dfGen[y,:Fuel]][t] * EP[:vFuel][y, t])) + # plant level start-up fuel cost for output + @expression(EP, ePlantCFuelOut[y = 1:G], + sum(inputs["omega"][t] * EP[:eCFuelOut][y, t] for t in 1:T)) + # zonal level total fuel cost for output + @expression(EP, eZonalCFuelOut[z = 1:Z], + sum(EP[:ePlantCFuelOut][y] for y in dfGen[dfGen[!, :Zone].==z, :R_ID])) + + + # system level total fuel cost for output + @expression(EP, eTotalCFuelOut, sum(eZonalCFuelOut[z] for z in 1:Z)) + @expression(EP, eTotalCFuelStart, sum(eZonalCFuelStart[z] for z in 1:Z)) + + + add_to_expression!(EP[:eObj], EP[:eTotalCFuelOut] + EP[:eTotalCFuelStart]) + + #fuel consumption (MMBTU or Billion BTU) + @expression(EP, eFuelConsumption[f in 1:NUM_FUEL, t in 1:T], + sum(EP[:vFuel][y, t] + EP[:eStartFuel][y,t] + for y in resources_with_fuel(dfGen, fuels[f]))) + + @expression(EP, eFuelConsumptionYear[f in 1:NUM_FUEL], + sum(inputs["omega"][t] * EP[:eFuelConsumption][f, t] for t in 1:T)) + + + ### Constraint ### + ### only apply constraint to generators with fuel type other than None + @constraint(EP, FuelCalculation[y in setdiff(HAS_FUEL, THERM_COMMIT), t = 1:T], + EP[:vFuel][y, t] - EP[:vP][y, t] * dfGen[y, :Heat_Rate_MMBTU_per_MWh] == 0) + + if !isempty(THERM_COMMIT) + # Only apply piecewise fuel consumption to thermal generators in THERM_COMMIT_PWFU set + THERM_COMMIT_PWFU = inputs["THERM_COMMIT_PWFU"] + # segemnt for piecewise fuel usage + if !isempty(THERM_COMMIT_PWFU) + segs = 1:inputs["PWFU_Num_Segments"] + PWFU_data = inputs["PWFU_data"] + slope_cols = inputs["slope_cols"] + intercept_cols = inputs["intercept_cols"] + segment_intercept(y, seg) = PWFU_data[y, intercept_cols[seg]] + segment_slope(y, seg) = PWFU_data[y, slope_cols[seg]] + # constraint for piecewise fuel consumption + @constraint(EP, 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))) + end + # constraint for fuel consumption at a constant heat rate + @constraint(EP, FuelCalculationCommit[y in setdiff(THERM_COMMIT,THERM_COMMIT_PWFU), t = 1:T], + EP[:vFuel][y, t] - EP[:vP][y, t] * dfGen[y, :Heat_Rate_MMBTU_per_MWh] == 0) + end + + return EP +end + + +function resources_with_fuel(df::DataFrame, fuel::AbstractString) + return df[df[!, :Fuel] .== fuel, :R_ID] +end + diff --git a/src/model/generate_model.jl b/src/model/generate_model.jl index 0ea10cc507..2577313df3 100644 --- a/src/model/generate_model.jl +++ b/src/model/generate_model.jl @@ -124,7 +124,10 @@ function generate_model(setup::Dict,inputs::Dict,OPTIMIZER::MOI.OptimizerWithAtt ucommit!(EP, inputs, setup) end - emissions!(EP, inputs) + fuel!(EP, inputs, setup) + + co2!(EP, inputs) + if setup["Reserves"] > 0 reserves!(EP, inputs, setup) @@ -185,6 +188,7 @@ function generate_model(setup::Dict,inputs::Dict,OPTIMIZER::MOI.OptimizerWithAtt end # Policies + # CO2 emissions limits if setup["CO2Cap"] > 0 co2_cap!(EP, inputs, setup) diff --git a/src/model/policies/co2_cap.jl b/src/model/policies/co2_cap.jl index 5238342e17..c8903ae028 100644 --- a/src/model/policies/co2_cap.jl +++ b/src/model/policies/co2_cap.jl @@ -67,7 +67,7 @@ Note that the generator-side rate-based constraint can be used to represent a fe """ function co2_cap!(EP::Model, inputs::Dict, setup::Dict) - println("C02 Policies Module") + println("CO2 Policies Module") dfGen = inputs["dfGen"] SEG = inputs["SEG"] # Number of lines diff --git a/src/write_outputs/write_co2.jl b/src/write_outputs/write_co2.jl new file mode 100644 index 0000000000..3cbf100fd0 --- /dev/null +++ b/src/write_outputs/write_co2.jl @@ -0,0 +1,63 @@ +@doc raw""" + write_co2(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) + +Function for reporting time-dependent CO2 emissions by zone. + +""" +function write_co2(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) + write_co2_emissions_plant(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) + write_co2_capture_plant(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) +end + + +function write_co2_emissions_plant(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) + dfGen = inputs["dfGen"] + G = inputs["G"] # Number of resources (generators, storage, DR, and DERs) + T = inputs["T"] # Number of time steps (hours) + Z = inputs["Z"] # Number of zones + + # CO2 emissions by plant + dfEmissions_plant = DataFrame(Resource=inputs["RESOURCES"], Zone=dfGen[!, :Zone], AnnualSum=zeros(G)) + emissions_plant = value.(EP[:eEmissionsByPlant]) + if setup["ParameterScale"] == 1 + emissions_plant *= ModelScalingFactor + end + dfEmissions_plant.AnnualSum .= emissions_plant * inputs["omega"] + dfEmissions_plant = hcat(dfEmissions_plant, DataFrame(emissions_plant, :auto)) + + auxNew_Names = [Symbol("Resource"); Symbol("Zone"); Symbol("AnnualSum"); [Symbol("t$t") for t = 1:T]] + rename!(dfEmissions_plant, auxNew_Names) + + total = DataFrame(["Total" 0 sum(dfEmissions_plant[!, :AnnualSum]) fill(0.0, (1, T))], auxNew_Names) + total[!, 4:T+3] .= sum(emissions_plant, dims=1) + dfEmissions_plant = vcat(dfEmissions_plant, total) + CSV.write(joinpath(path, "emissions_plant.csv"), dftranspose(dfEmissions_plant, false), writeheader=false) +end + +function write_co2_capture_plant(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) + dfGen = inputs["dfGen"] + G = inputs["G"] # Number of resources (generators, storage, DR, and DERs) + T = inputs["T"] # Number of time steps (hours) + Z = inputs["Z"] # Number of zones + + dfCapturedEmissions_plant = DataFrame(Resource=inputs["RESOURCES"], Zone=dfGen[!, :Zone], AnnualSum=zeros(G)) + if any(dfGen.CO2_Capture_Fraction .!= 0) + # Captured CO2 emissions by plant + emissions_captured_plant = zeros(G, T) + emissions_captured_plant = (value.(EP[:eEmissionsCaptureByPlant])) + if setup["ParameterScale"] == 1 + emissions_captured_plant *= ModelScalingFactor + end + dfCapturedEmissions_plant.AnnualSum .= emissions_captured_plant * inputs["omega"] + dfCapturedEmissions_plant = hcat(dfCapturedEmissions_plant, DataFrame(emissions_captured_plant, :auto)) + + auxNew_Names = [Symbol("Resource"); Symbol("Zone"); Symbol("AnnualSum"); [Symbol("t$t") for t = 1:T]] + rename!(dfCapturedEmissions_plant, auxNew_Names) + + total = DataFrame(["Total" 0 sum(dfCapturedEmissions_plant[!, :AnnualSum]) fill(0.0, (1, T))], auxNew_Names) + total[!, 4:T+3] .= sum(emissions_captured_plant, dims=1) + dfCapturedEmissions_plant = vcat(dfCapturedEmissions_plant, total) + + CSV.write(joinpath(path, "captured_emissions_plant.csv"), dftranspose(dfCapturedEmissions_plant, false), writeheader=false) + end +end diff --git a/src/write_outputs/write_costs.jl b/src/write_outputs/write_costs.jl index 39c4e02427..d13d7903e4 100644 --- a/src/write_outputs/write_costs.jl +++ b/src/write_outputs/write_costs.jl @@ -12,7 +12,7 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) VRE_STOR = inputs["VRE_STOR"] ELECTROLYZER = inputs["ELECTROLYZER"] - cost_list = ["cTotal", "cFix", "cVar", "cNSE", "cStart", "cUnmetRsv", "cNetworkExp", "cUnmetPolicyPenalty"] + cost_list = ["cTotal", "cFix", "cVar", "cFuel" ,"cNSE", "cStart", "cUnmetRsv", "cNetworkExp", "cUnmetPolicyPenalty", "cCO2"] if !isempty(VRE_STOR) push!(cost_list, "cGridConnection") end @@ -23,7 +23,9 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) cVar = value(EP[:eTotalCVarOut])+ (!isempty(inputs["STOR_ALL"]) ? value(EP[:eTotalCVarIn]) : 0.0) + (!isempty(inputs["FLEX"]) ? value(EP[:eTotalCVarFlexIn]) : 0.0) cFix = value(EP[:eTotalCFix]) + (!isempty(inputs["STOR_ALL"]) ? value(EP[:eTotalCFixEnergy]) : 0.0) + (!isempty(inputs["STOR_ASYMMETRIC"]) ? value(EP[:eTotalCFixCharge]) : 0.0) - + + cFuel = value.(EP[:eTotalCFuelOut]) + if !isempty(VRE_STOR) cFix += ((!isempty(inputs["VS_DC"]) ? value(EP[:eTotalCFixDC]) : 0.0) + (!isempty(inputs["VS_SOLAR"]) ? value(EP[:eTotalCFixSolar]) : 0.0) + (!isempty(inputs["VS_WIND"]) ? value(EP[:eTotalCFixWind]) : 0.0)) cVar += ((!isempty(inputs["VS_SOLAR"]) ? value(EP[:eTotalCVarOutSolar]) : 0.0) + (!isempty(inputs["VS_WIND"]) ? value(EP[:eTotalCVarOutWind]) : 0.0)) @@ -31,9 +33,9 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) cFix += ((!isempty(inputs["VS_STOR"]) ? value(EP[:eTotalCFixStor]) : 0.0) + (!isempty(inputs["VS_ASYM_DC_CHARGE"]) ? value(EP[:eTotalCFixCharge_DC]) : 0.0) + (!isempty(inputs["VS_ASYM_DC_DISCHARGE"]) ? value(EP[:eTotalCFixDischarge_DC]) : 0.0) + (!isempty(inputs["VS_ASYM_AC_CHARGE"]) ? value(EP[:eTotalCFixCharge_AC]) : 0.0) + (!isempty(inputs["VS_ASYM_AC_DISCHARGE"]) ? value(EP[:eTotalCFixDischarge_AC]) : 0.0)) cVar += (!isempty(inputs["VS_STOR"]) ? value(EP[:eTotalCVarStor]) : 0.0) end - total_cost = [objective_value(EP), cFix, cVar, value(EP[:eTotalCNSE]), 0.0, 0.0, 0.0, 0.0, 0.0] + total_cost =[objective_value(EP), cFix, cVar, cFuel, value(EP[:eTotalCNSE]), 0.0, 0.0, 0.0, 0.0, 0.0, 0.0] else - total_cost = [objective_value(EP), cFix, cVar, value(EP[:eTotalCNSE]), 0.0, 0.0, 0.0, 0.0] + total_cost = [objective_value(EP), cFix, cVar, cFuel, value(EP[:eTotalCNSE]), 0.0, 0.0, 0.0, 0.0, 0.0] end if !isempty(ELECTROLYZER) @@ -47,51 +49,58 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) end if setup["UCommit"]>=1 - dfCost[5,2] = value(EP[:eTotalCStart]) + dfCost[6,2] = value(EP[:eTotalCStart]) + value(EP[:eTotalCFuelStart]) end if setup["Reserves"]==1 - dfCost[6,2] = value(EP[:eTotalCRsvPen]) + dfCost[7,2] = value(EP[:eTotalCRsvPen]) end if setup["NetworkExpansion"] == 1 && Z > 1 - dfCost[7,2] = value(EP[:eTotalCNetworkExp]) + dfCost[8,2] = value(EP[:eTotalCNetworkExp]) end if haskey(inputs, "dfCapRes_slack") - dfCost[8,2] += value(EP[:eCTotalCapResSlack]) + dfCost[9,2] += value(EP[:eCTotalCapResSlack]) end if haskey(inputs, "dfESR_slack") - dfCost[8,2] += value(EP[:eCTotalESRSlack]) + dfCost[9,2] += value(EP[:eCTotalESRSlack]) end if haskey(inputs, "dfCO2Cap_slack") - dfCost[8,2] += value(EP[:eCTotalCO2CapSlack]) + dfCost[9,2] += value(EP[:eCTotalCO2CapSlack]) end if haskey(inputs, "MinCapPriceCap") - dfCost[8,2] += value(EP[:eTotalCMinCapSlack]) + dfCost[9,2] += value(EP[:eTotalCMinCapSlack]) end if !isempty(VRE_STOR) - dfCost[!,2][9] = value(EP[:eTotalCGrid]) * (setup["ParameterScale"] == 1 ? ModelScalingFactor^2 : 1) + dfCost[!,2][11] = value(EP[:eTotalCGrid]) * (setup["ParameterScale"] == 1 ? ModelScalingFactor^2 : 1) + end + + if any(dfGen.CO2_Capture_Fraction .!= 0) + dfCost[10,2] += value(EP[:eTotaleCCO2Sequestration]) end if setup["ParameterScale"] == 1 - dfCost[5,2] *= ModelScalingFactor^2 dfCost[6,2] *= ModelScalingFactor^2 dfCost[7,2] *= ModelScalingFactor^2 dfCost[8,2] *= ModelScalingFactor^2 + dfCost[9,2] *= ModelScalingFactor^2 + dfCost[10,2] *= ModelScalingFactor^2 end for z in 1:Z tempCTotal = 0.0 tempCFix = 0.0 tempCVar = 0.0 + tempCFuel = 0.0 tempCStart = 0.0 tempCNSE = 0.0 tempHydrogenValue = 0.0 + tempCCO2 = 0.0 Y_ZONE = dfGen[dfGen[!,:Zone].==z,:R_ID] STOR_ALL_ZONE = intersect(inputs["STOR_ALL"], Y_ZONE) @@ -107,6 +116,9 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) tempCVar = sum(value.(EP[:eCVar_out][Y_ZONE,:])) tempCTotal += tempCVar + tempCFuel = sum(value.(EP[:ePlantCFuelOut][Y_ZONE,:])) + tempCTotal += tempCFuel + if !isempty(STOR_ALL_ZONE) eCVar_in = sum(value.(EP[:eCVar_in][STOR_ALL_ZONE,:])) tempCVar += eCVar_in @@ -192,7 +204,7 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) end if setup["UCommit"] >= 1 && !isempty(COMMIT_ZONE) - eCStart = sum(value.(EP[:eCStart][COMMIT_ZONE,:])) + eCStart = sum(value.(EP[:eCStart][COMMIT_ZONE,:])) + sum(value.(EP[:ePlantCFuelStart][COMMIT_ZONE,:])) tempCStart += eCStart tempCTotal += eCStart end @@ -206,16 +218,22 @@ function write_costs(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) tempCNSE = sum(value.(EP[:eCNSE][:,:,z])) tempCTotal += tempCNSE + if any(dfGen.CO2_Capture_Fraction .!=0) + tempCCO2 = sum(value.(EP[:ePlantCCO2Sequestration][Y_ZONE,:])) + tempCTotal += tempCCO2 + end + if setup["ParameterScale"] == 1 tempCTotal *= ModelScalingFactor^2 tempCFix *= ModelScalingFactor^2 tempCVar *= ModelScalingFactor^2 + tempCFuel *= ModelScalingFactor^2 tempCNSE *= ModelScalingFactor^2 tempCStart *= ModelScalingFactor^2 tempHydrogenValue *= ModelScalingFactor^2 + tempCCO2 *= ModelScalingFactor^2 end - - temp_cost_list = [tempCTotal, tempCFix, tempCVar, tempCNSE, tempCStart, "-", "-", "-"] + temp_cost_list = [tempCTotal, tempCFix, tempCVar, tempCFuel,tempCNSE, tempCStart, "-", "-", "-", tempCCO2] if !isempty(VRE_STOR) push!(temp_cost_list, "-") end diff --git a/src/write_outputs/write_fuel_consumption.jl b/src/write_outputs/write_fuel_consumption.jl new file mode 100644 index 0000000000..ec88b3bae8 --- /dev/null +++ b/src/write_outputs/write_fuel_consumption.jl @@ -0,0 +1,60 @@ + +@doc raw""" + write_fuel_consumption(path::AbstractString, inputs::Dict, setup::Dict, EP::Model). +Write fuel consumption of each power plant. +""" +function write_fuel_consumption(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) + + write_fuel_consumption_plant(path::AbstractString,inputs::Dict, setup::Dict, EP::Model) + write_fuel_consumption_ts(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) + write_fuel_consumption_tot(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) +end + +function write_fuel_consumption_plant(path::AbstractString,inputs::Dict, setup::Dict, EP::Model) + dfGen = inputs["dfGen"] + G = inputs["G"] + HAS_FUEL = inputs["HAS_FUEL"] + # Fuel consumption cost by each resource, including start up fuel + dfPlantFuel = DataFrame(Resource = inputs["RESOURCES"][HAS_FUEL], + Fuel = dfGen[HAS_FUEL, :Fuel], + Zone = dfGen[HAS_FUEL,:Zone], + AnnualSum = zeros(length(HAS_FUEL))) + tempannualsum = value.(EP[:ePlantCFuelOut][HAS_FUEL]) + value.(EP[:ePlantCFuelStart][HAS_FUEL]) + + if setup["ParameterScale"] == 1 + tempannualsum *= ModelScalingFactor^2 # + end + dfPlantFuel.AnnualSum .+= tempannualsum + CSV.write(joinpath(path, "Fuel_cost_plant.csv"), dfPlantFuel) +end + + +function write_fuel_consumption_ts(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) + T = inputs["T"] # Number of time steps (hours) + HAS_FUEL = inputs["HAS_FUEL"] + # Fuel consumption by each resource per time step, unit is MMBTU + dfPlantFuel_TS = DataFrame(Resource = inputs["RESOURCES"][HAS_FUEL]) + tempts = value.(EP[:vFuel] + EP[:eStartFuel])[HAS_FUEL,:] + if setup["ParameterScale"] == 1 + tempts *= ModelScalingFactor # kMMBTU to MMBTU + end + dfPlantFuel_TS = hcat(dfPlantFuel_TS, + DataFrame(tempts, [Symbol("t$t") for t in 1:T])) + CSV.write(joinpath(path, "FuelConsumption_plant_MMBTU.csv"), + dftranspose(dfPlantFuel_TS, false), writeheader=false) +end + + +function write_fuel_consumption_tot(path::AbstractString, inputs::Dict, setup::Dict, EP::Model) + # types of fuel + fuel_types = inputs["fuels"] + fuel_number = length(fuel_types) + dfFuel = DataFrame(Fuel = fuel_types, + AnnualSum = zeros(fuel_number)) + tempannualsum = value.(EP[:eFuelConsumptionYear]) + if setup["ParameterScale"] == 1 + tempannualsum *= ModelScalingFactor # billion MMBTU to MMBTU + end + dfFuel.AnnualSum .+= tempannualsum + CSV.write(joinpath(path,"FuelConsumption_total_MMBTU.csv"), dfFuel) +end diff --git a/src/write_outputs/write_net_revenue.jl b/src/write_outputs/write_net_revenue.jl index b8d2f1f7d6..fd1d473c45 100644 --- a/src/write_outputs/write_net_revenue.jl +++ b/src/write_outputs/write_net_revenue.jl @@ -77,7 +77,7 @@ function write_net_revenue(path::AbstractString, inputs::Dict, setup::Dict, EP:: end # Add fuel cost to the dataframe - dfNetRevenue.Fuel_cost = (inputs["C_Fuel_per_MWh"] .* value.(EP[:vP])) * inputs["omega"] + dfNetRevenue.Fuel_cost = sum(value.(EP[:ePlantCFuelOut]), dims = 2) if setup["ParameterScale"] == 1 dfNetRevenue.Fuel_cost *= ModelScalingFactor^2 # converting Million US$ to US$ end @@ -102,8 +102,9 @@ function write_net_revenue(path::AbstractString, inputs::Dict, setup::Dict, EP:: # Add start-up cost to the dataframe dfNetRevenue.StartCost = zeros(nrow(dfNetRevenue)) if setup["UCommit"]>=1 && !isempty(COMMIT) - # if you don't use vec, dimension won't match - dfNetRevenue.StartCost[COMMIT] .= vec(sum(value.(EP[:eCStart][COMMIT, :]).data, dims = 2)) + start_costs = vec(sum(value.(EP[:eCStart][COMMIT, :]).data, dims = 2)) + start_fuel_costs = vec(value.(EP[:ePlantCFuelStart][COMMIT])) + dfNetRevenue.StartCost[COMMIT] .= start_costs + start_fuel_costs end if setup["ParameterScale"] == 1 dfNetRevenue.StartCost *= ModelScalingFactor^2 # converting Million US$ to US$ diff --git a/src/write_outputs/write_outputs.jl b/src/write_outputs/write_outputs.jl index 4a6bbb6f69..13e7df2f57 100644 --- a/src/write_outputs/write_outputs.jl +++ b/src/write_outputs/write_outputs.jl @@ -128,6 +128,14 @@ function write_outputs(EP::Model, path::AbstractString, setup::Dict, inputs::Dic println(elapsed_time_lds_dstor) end + elapsed_time_fuel_consumption = @elapsed write_fuel_consumption(path, inputs, setup, EP) + println("Time elapsed for writing fuel consumption is") + println(elapsed_time_fuel_consumption) + + elapsed_time_emissions = @elapsed write_co2(path, inputs, setup, EP) + println("Time elapsed for writing co2 is") + println(elapsed_time_emissions) + # Temporary! Suppress these outputs until we know that they are compatable with multi-stage modeling if setup["MultiStage"] == 0 dfPrice = DataFrame()