import Pkg
Pkg.activate(@__DIR__)
Pkg.instantiate()
BEE 4750 Homework 4: Linear Programming and Capacity Expansion
If you are enrolled in the course, make sure that you use the GitHub Classroom link provided in Ed Discussion, or you may not be able to get help if you run into problems.
Otherwise, you can find the Github repository here.
Overview
Instructions
- Problem 1 asks you to formulate and solve a resource allocation problem using linear programming.
- Problem 2 asks you to formulate, solve, and analyze a standard generating capacity expansion problem.
- Problem 3 (5750 only) asks you to add a CO2 constraint to the capacity expansion problem and identify changes in the resulting solution.
Load Environment
The following code loads the environment and makes sure all needed packages are installed. This should be at the start of most Julia scripts.
using JuMP
using HiGHS
using DataFrames
using Plots
using Measures
using CSV
using MarkdownTables
Problems (Total: 50/60 Points)
Problem 1 (20 points)
A farmer has access to a pesticide which can be used on corn, soybeans, and wheat fields, and costs $70/kg-yr to apply. The crop yields the farmer can obtain following crop yields by applying varying rates of pesticides to the field are shown in Table 1.
Application Rate (kg/ha) | Soybean (kg/ha) | Wheat (kg/ha) | Corn (kg/ha) |
---|---|---|---|
0 | 2900 | 3500 | 5900 |
1 | 3800 | 4100 | 6700 |
2 | 4400 | 4200 | 7900 |
The costs of production, excluding pesticides, for each crop, and selling prices, are shown in Table 2.
Crop | Production Cost ($/ha-yr) | Selling Price ($/kg) |
---|---|---|
Soybeans | 350 | 0.36 |
Wheat | 280 | 0.27 |
Corn | 390 | 0.22 |
Recently, environmental authorities have declared that farms cannot have an average application rate on soybeans, wheat, and corn which exceeds 0.8, 0.7, and 0.6 kg/ha, respectively. The farmer has asked you for advice on how they should plant crops and apply pesticides to maximize profits over 130 total ha while remaining in regulatory compliance if demand for each crop (which is the maximum the market would buy) this year is 250,000 kg?
In this problem:
- Formulate a linear program for this resource allocation problem, including clear definitions of decision variable(s) (including units), objective function(s), and constraint(s) (make sure to explain functions and constraints with any needed derivations and explanations). Tip: Make sure that all of your constraints are linear.
- Implement the program in
JuMP.jl
and find the solution. How many ha should the farmer dedicate to each crop and with what pesticide application rate(s)? How much profit will the farmer expect to make? - The farmer has an opportunity to buy an extra 10 ha of land. How much extra profit would this land be worth to the farmer? Discuss why this value makes sense and whether you would recommend the farmer should make the purchase.
Problem 2 (30 points)
For this problem, we will use hourly load (demand) data from 2013 in New York’s Zone C (which includes Ithaca). The load data is loaded and plotted below in Figure 1.
# load the data, pull Zone C, and reformat the DataFrame
= DataFrame(CSV.File("data/2013_hourly_load_NY.csv"))
NY_demand rename!(NY_demand, :"Time Stamp" => :Date)
= NY_demand[:, [:Date, :C]]
demand rename!(demand, :C => :Demand)
:, :Hour] = 1:nrow(demand)
demand[
# plot demand
plot(demand.Hour, demand.Demand, xlabel="Hour of Year", ylabel="Demand (MWh)", label=:false)
Next, we load the generator data, shown in Table 3. This data includes fixed costs ($/MW installed), variable costs ($/MWh generated), and CO2 emissions intensity (tCO2/MWh generated).
= DataFrame(CSV.File("data/generators.csv")) gens
Row | Plant | FixedCost | VarCost | Emissions |
---|---|---|---|---|
String15 | Int64 | Int64 | Float64 | |
1 | Geothermal | 450000 | 0 | 0.0 |
2 | Coal | 220000 | 24 | 1.0 |
3 | NG CCGT | 82000 | 30 | 0.43 |
4 | NG CT | 65000 | 40 | 0.55 |
5 | Wind | 91000 | 0 | 0.0 |
6 | Solar | 70000 | 0 | 0.0 |
Finally, we load the hourly solar and wind capacity factors, which are plotted in Figure 2. These tell us the fraction of installed capacity which is expected to be available in a given hour for generation (typically based on the average meteorology).
# load capacify factors into a DataFrame
= DataFrame(CSV.File("data/wind_solar_capacity_factors.csv"))
cap_factor
# plot January capacity factors
= plot(cap_factor.Wind[1:(24*31)], label="Wind")
p1 plot!(cap_factor.Solar[1:(24*31)], label="Solar")
xaxis!("Hour of the Month")
yaxis!("Capacity Factor")
= plot(cap_factor.Wind[4344:4344+(24*31)], label="Wind")
p2 plot!(cap_factor.Solar[4344:4344+(24*31)], label="Solar")
xaxis!("Hour of the Month")
yaxis!("Capacity Factor")
display(p1)
display(p2)
You have been asked to develop a generating capacity expansion plan for the utility in Riley County, NY, which currently has no existing electrical generation infrastructure. The utility can build any of the following plant types: geothermal, coal, natural gas combined cycle gas turbine (CCGT), natural gas combustion turbine (CT), solar, and wind.
While coal, CCGT, and CT plants can generate at their full installed capacity, geothermal plants operate at maximum 85% capacity, and solar and wind available capacities vary by the hour depend on the expected meteorology. The utility will also penalize any non-served demand at a rate of $10,000/MWh.
In this problem:
- Formulate a linear program for this capacity expansion problem, including clear definitions of decision variable(s) (including units), objective function(s), and constraint(s) (make sure to explain functions and constraints with any needed derivations and explanations).
- Implement your linear program in
JuMP.jl
. Find the optimal solution. How much should the utility build of each type of generating plant? What will the total cost be? How much energy will be non-served? - What fraction of annual generation does each plant type produce? How does this compare to the breakdown of built capacity that you found in Problem 1.5? Do these results make sense given the generator data?
- Make a plot of the electricity price in each hour. Discuss any trends that you see.
Use round(x; digits=n)
to report values to the appropriate precision! If your number is on a different order of magnitude and you want to round to a certain number of significant digits, you can use round(x; sigdigits=n)
.
value.(x)
will report the values of a JuMP
variable x
, but it will return a special container which holds other information about x
that is useful for JuMP
. This means that you can’t use this output directly for further calculations. To just extract the values, use value.(x).data
.
The output of specifying model components (variable or constraints) can be quite large for this problem because of the number of time periods. If you end a cell with an @variable
or @constraint
command, I highly recommend suppressing output by adding a semi-colon after the last command, or you might find that your notebook crashes.
Problem 3 (10 points)
This problem is only required for students in BEE 5750.
The NY state legislature is considering enacting an annual CO2 limit, which for the utility would limit the emissions in its footprint to 1.5 MtCO2/yr.
In this problem:
- Reformulate your linear program from Problem 2 with any necessary changes to capture the CO2 limit.
- Implement the new optimization problem and find the optimal solution. How much should the utility build of each type of generating plant? What is different from your plan from Problem 1? Do these changes make sense?
- What would the value to the utility be of allowing it to emit an additional 1000 tCO2/yr? An additional 5000?
References
List any external references consulted, including classmates.