L3-Cashflow
Optimization modeling in Julia¶
In this notebook we describe how to formulate and solve an optimization problem using mathematical programming tools in Julia.
Copyright By PowCoder代写 加微信 powcoder
Step 0: Setup¶
Step 0.1: load packages¶
We can load packages using the using keyword:
using JuMP, Gurobi
Step 0.2: set up problem data¶
Recall that we are trying to solve the following optimization problem:
\begin{aligned}
\max z_6\\
\text{ s.t.}\\
-x_1 &- y_1 & & & & & + & z_1 & =&& -150\\
-x_2 &- y_2 & + & 1.01x_1 & – & 1.003z_1 & + & z_2 &=&& -100\\
-x_3 &- y_3 & + & 1.01x_2 & – & 1.003z_2 & + & z_3 &=&& 200\\
-x_4 &+ 1.02y_1 & + & 1.01x_3 & – & 1.003z_3 & + & z_4 &=&& -200\\
-x_5 &+ 1.02y_2 & + & 1.01x_4 & – & 1.003z_4 & + & z_5 &=&& 50\\
&+ 1.02y_3 & + & 1.01x_5 & – & 1.003z_5 & + & z_6 &=&& 300\\
& & & & & & & x_i& \le&& 100 &\quad\forall i=1,\ldots,5\\
& & & & & & & x_i& \ge&& 0 &\quad\forall i=1,\ldots,5\\
& & & & & & & y_i& \ge&& 0 &\quad\forall i=1,\ldots,3\\
& & & & & & & z_i& \ge&& 0 &\quad\forall i=1,\ldots,6\\
\end{aligned}
Accordingly, we will define the cash flow requirements for each month as a vector:
cash_flow_requirements = [-150, -100, 200, -200, 50, 300]
Step 1: Formulate the optimization problem in JuMP¶
Step 1.1: create a JuMP model object¶
A JuMP model is a container that holds all the ingredients of our optimization model
model = JuMP.Model()
At this point our model is quite sad and empty, so let’s build it up!
Step 1.2: attach a solver to our JuMP model¶
For clarity, in JuMP a “solver” refers to an actual third-party solver software (e.g. Gurobi) that we use to solve optimization problems, while an “optimizer” refers to the Julia interface for the solver software.
set_optimizer(model, Gurobi.Optimizer)
Step 1.3: decision variables¶
We’re now ready to define decision variables and their associated bounds.
@variable(model, 0 <= x[1:5] <= 100)
@variable(model, y[1:3] >= 0)
@variable(model, z[1:6] >= 0)
The above output is not very legible, so we can actually print out our model in math notation as follows:
print(model)
Note that the above print function is a bit dangerous when the optimization problem is very large, as writing down each constraint in this way may create huge amounts of input.
Step 1.4: Constraints¶
Now we’re ready to add the constraints of the optimization problem to our model container!
# Note that Julia indexes arrays from 1
# Remember to use == for equality constraints and not =
@constraint(model, January, -x[1] – y[1] + z[1] == cash_flow_requirements[1])
@constraint(model, February, -x[2] – y[2] + 1.01x[1] – 1.003z[1] + z[2] == cash_flow_requirements[2])
@constraint(model, March, -x[3] – y[3] + 1.01x[2] – 1.003z[2] + z[3] == cash_flow_requirements[3])
@constraint(model, April, -x[4] + 1.02y[1] + 1.01x[3] – 1.003z[3] + z[4] == cash_flow_requirements[4])
@constraint(model, May, -x[5] + 1.02y[2] + 1.01x[4] – 1.003z[4] + z[5] == cash_flow_requirements[5])
@constraint(model, June, 1.02y[3] + 1.01x[5] – 1.003z[5] + z[6] == cash_flow_requirements[6])
Let’s take a look at our model!
print(model)
Step 1.5: Objective¶
Notice that the first line of the problem description is “feasibility”. That’s because we haven’t yet defined an objective! Let’s do it now:
@objective(model, Max, z[6])
Step 2: solving the problem and inspecting the solution¶
We’re now ready to solve the problem! We can call on Gurobi to do its magic using the optimize! function
# In Julia, an exclamation mark at the end of a function name means it is changing/mutating its argument
optimize!(model)
We can take a look at the objective value by reading the log, or by extracting it directly:
objective_value(model)
And similarly we can obtain the optimal variable values:
value(y[1])
Oops! Calling value on a vector of variables doesn’t work directly, but we can get around that using a list comprehension:
optimal_x_values = [value(x[i]) for i=eachindex(x)]
Alternatively, we can use unique Julia syntax to “vectorize” the value function:
optimal_x_values = value.(x)
Step 3: sensitivity report¶
Much like in Excel, we can query the sensitivity information we need:
# @show is a nifty little tool that executes a line of code and also displays the value – great for debugging!
@show shadow_price(January)
@show shadow_price(March); # the semicolon suppresses output to the console
@show reduced_cost(x[1]);
If we also want the allowable ranges, we have to do a little bit more work:
report = lp_sensitivity_report(model)
@show report[January]; # we see the maximum allowable decrease and increase
Appendix: some JuMP syntax that may be useful¶
model = Model(Gurobi.Optimizer)
# creating “ragged” arrays of variables
@variable(model, x[i=1:5, j=1:i^2] >= i^2);
# creating multiple constraints and summing over variables
@constraint(model, sumconstraint[i=1:5], sum(x[i, j] for j = 1:i^2 if j % 2 == 0) >= 0);
A = [1 0 ; 0 1]
model = Model(Gurobi.Optimizer)
@variable(model, x[1:2] >= 0)
b = [3, 2]
@constraint(model, A * x .<= b) 程序代写 CS代考 加微信: powcoder QQ: 1823890830 Email: powcoder@163.com