Design and Scheduling (MILP)#

The primary motivation behind developing gana was to be able to write multiscale models seamlessly.

Here is how you write one

# !pip install gana # uncomment and run to install Gana, if not in environment
from gana import Prg, I, V, P, inf

Problem Definition#

We have a solar photovoltaic array and battery.

The sun shines (intermittently) and power flows. We need to meet a certain demand for power.

The problem is modeled using the resource task network (RTN) methodology as proposed in Pantelides, C.C., 1994, July. Unified frameworks for optimal process planning and scheduling. In Proceedings on the second conference on foundations of computer aided operations (pp. 253-274).

p = Prg()
p.y = I(size=1, tag='year')
p.q = I(size=3, tag='three quarters')

Resource indices#

p.res_cons = I('solar')
p.res_dem = I('power')
p.res_stg = I('charge')
p.res = p.res_cons | p.res_dem | p.res_stg

Process indices#

p.pro_var = I('pv')
p.pro_cer = I('li', 'li_d')
p.pro = p.pro_var | p.pro_cer
p.pro.show(True)
\[\displaystyle {pro} = \{ {{pv}, {li}, {li\_d}} \}\]

Parameters#

p.dm_fac = P(p.power, p.q, _=[0.5, 1, 0.5], tag='demand factor (variability)')
p.pv_fac = P(p.pv, p.q, _=[1, 0, 0.5], tag='intermittency factor for solar')
p.demand = P(p.res_dem, p.q, _=[100] * 3, tag='nominal demand')
p.capex = P(
    p.pro, p.y, _=[5000, 1000, 0], tag='capital expenditure for solar processes'
)
p.fopex = P(p.pro, p.y, _=[500, 100, 0], tag='fixed operating expenditure')
p.vopex = P(p.pro, p.y, _=[10, 50, 0], tag='variable operating expenditure')

Variables#

p.cap_p = V(p.pro, p.y, tag='nameplate production capacity')
p.cap_s = V(p.res_stg, p.y, tag='nameplate storage capacity')
p.sell = V(p.res_dem, p.q, tag='amount of power sold')
p.con = V(p.res_cons, p.q, tag='amount of solar consumed')
p.inv = V(p.res_stg, p.q, tag='charge inventory')
p.prod = V(p.pro, p.q, tag='production capacity utilization')
p.ex_cap = V(p.pro, p.y, tag='capital expenditure')
p.ex_fop = V(p.pro, p.y, tag='fixed operating expenditure')
p.ex_vop = V(p.pro, p.y, tag='variable operating expenditure')

Decision Constraints#

These lend degrees of freedom to the model

p.con_capmax = p.cap_p(p.pro, p.y) <= 200
p.con_capstg = p.cap_s(p.charge, p.y) <= 200
p.con_consmax = p.con(p.res_cons, p.q) <= 200
p.con_sell = p.sell(p.power, p.q) >= p.dm_fac(p.power, p.q) * p.demand(p.power, p.q)

Multiscale Decision Constraints#

Here the constraint LHS is index by quarter and the RHS by year.

Thus, temporally disparate phenomena (scheduling and capacity sizing) are modeled

p.con_pv = p.prod(p.pv, p.q) <= p.pv_fac(p.pv, p.q) * p.cap_p(p.pv, p.y)
p.con_pv.show(True)
\[\displaystyle [10]\text{ }{\mathbf{prod}}_{{pv},{{{q}_{0}}}} - 1.0 \cdot {\mathbf{cap\_p}}_{{pv},{{{y}_{0}}}} \leq 0\]
\[\displaystyle [11]\text{ }{\mathbf{prod}}_{{pv},{{{q}_{1}}}} - 0 \cdot {\mathbf{cap\_p}}_{{pv},{{{y}_{0}}}} \leq 0\]
\[\displaystyle [12]\text{ }{\mathbf{prod}}_{{pv},{{{q}_{2}}}} - 0.5 \cdot {\mathbf{cap\_p}}_{{pv},{{{y}_{0}}}} \leq 0\]
p.con_prod = p.prod(p.pro_cer, p.q) <= p.cap_p(p.pro_cer, p.y)
p.con_inv = p.inv(p.charge, p.q) <= p.cap_s(p.charge, p.y)

Calculations#

p.con_vopex = p.ex_vop(p.pro, p.y) == p.vopex(p.pro, p.y) * sum(
    p.prod(p.pro, q) for q in p.q
)
p.con_capex = p.ex_cap(p.pro, p.y) == p.capex(p.pro, p.y) * p.cap_p(p.pro, p.y)
p.con_fopex = p.ex_fop(p.pro, p.y) == p.fopex(p.pro, p.y) * p.cap_p(p.pro, p.y)

Balances#

p.con_solar = p.prod(p.pv, p.q) == p.con(p.solar, p.q)
p.con_power = (
    sum(p.prod(i, p.q) for i in p.pro_var)
    - p.prod(p.li, p.q)
    + p.prod(p.li_d, p.q)
    - p.sell(p.power, p.q)
    == 0
)
p.con_charge = (
    p.prod(p.li, p.q)
    - p.prod(p.li_d, p.q)
    + p.inv(p.charge, p.q - 1)
    - p.inv(p.charge, p.q)
    == 0
)

Objective#

p.o = inf(sum(p.ex_cap) + sum(p.ex_vop) + sum(p.ex_fop))

Display Program#

p.show()

Mathematical Program for prog



Index Sets

\[\displaystyle {y} = \{ {{{{y}_{0}}}} \}\]
\[\displaystyle {q} = \{ {{{{q}_{0}}}, {{{q}_{1}}}, {{{q}_{2}}}} \}\]
\[\displaystyle {res\_cons} = \{ {{solar}} \}\]
\[\displaystyle {res\_dem} = \{ {{power}} \}\]
\[\displaystyle {res\_stg} = \{ {{charge}} \}\]
\[\displaystyle {res} = \{ {{solar}, {power}, {charge}} \}\]
\[\displaystyle {pro\_var} = \{ {{pv}} \}\]
\[\displaystyle {pro\_cer} = \{ {{li}, {li\_d}} \}\]
\[\displaystyle {pro} = \{ {{pv}, {li}, {li\_d}} \}\]



Objective

\[\displaystyle min \hspace{0.2cm} {\mathbf{ex\_cap}}_{{pv},{{{y}_{0}}}} + {\mathbf{ex\_cap}}_{{li},{{{y}_{0}}}} + {\mathbf{ex\_cap}}_{{li\_d},{{{y}_{0}}}} + {\mathbf{ex\_vop}}_{{pv},{{{y}_{0}}}} + {\mathbf{ex\_vop}}_{{li},{{{y}_{0}}}} + {\mathbf{ex\_vop}}_{{li\_d},{{{y}_{0}}}} + {\mathbf{ex\_fop}}_{{pv},{{{y}_{0}}}} + {\mathbf{ex\_fop}}_{{li},{{{y}_{0}}}} + {\mathbf{ex\_fop}}_{{li\_d},{{{y}_{0}}}}\]



s.t.



Inequality Constraint Sets

\[\displaystyle [0]\text{ }{\mathbf{{\mathbf{cap\_p}}}}_{{pro},{y}} - 200 \leq 0\]
\[\displaystyle [1]\text{ }{\mathbf{{\mathbf{{\mathbf{cap\_s}}}}}}_{{charge},{y}} - 200 \leq 0\]
\[\displaystyle [2]\text{ }{\mathbf{{\mathbf{con}}}}_{{res\_cons},{q}} - 200 \leq 0\]
\[\displaystyle [3]\text{ }-{\mathbf{{\mathbf{{\mathbf{sell}}}}}}_{{power},{q}} + {\mathrm{{\mathrm{dm\_fac}}}}_{{power},{q}} \leq 0\]
\[\displaystyle [4]\text{ }{\mathbf{{\mathbf{{\mathbf{prod}}}}}}_{{pv},{q}} - {\mathrm{pv\_fac}}_{{pv},{q}} \cdot {\mathbf{{\mathbf{{\mathbf{cap\_p}}}}}}_{{pv},{y}} \leq 0\]
\[\displaystyle [5]\text{ }{\mathbf{{\mathbf{{\mathbf{prod}}}}}}_{{pro\_cer},{q}} - {\mathbf{{\mathbf{{\mathbf{cap\_p}}}}}}_{{pro\_cer},{y}} \leq 0\]
\[\displaystyle [6]\text{ }{\mathbf{{\mathbf{{\mathbf{inv}}}}}}_{{charge},{q}} - {\mathbf{{\mathbf{{\mathbf{cap\_s}}}}}}_{{charge},{y}} \leq 0\]



Equality Constraint Sets

\[\displaystyle [7]\text{ }{\mathbf{{\mathbf{ex\_vop}}}}_{{pro},{y}} - ({\mathrm{vopex}}_{{pro},{y}} \cdot {\mathbf{{\mathbf{{\mathbf{{\mathbf{prod}}}}}}}}_{{pro},{{{q}_{0}}}} + {\mathrm{vopex}}_{{pro},{y}} \cdot {\mathbf{{\mathbf{{\mathbf{{\mathbf{prod}}}}}}}}_{{pro},{{{q}_{1}}}} + {\mathrm{vopex}}_{{pro},{y}} \cdot {\mathbf{{\mathbf{{\mathbf{{\mathbf{prod}}}}}}}}_{{pro},{{{q}_{2}}}}) = 0\]
\[\displaystyle [8]\text{ }{\mathbf{{\mathbf{ex\_cap}}}}_{{pro},{y}} - {\mathrm{capex}}_{{pro},{y}} \cdot {\mathbf{{\mathbf{cap\_p}}}}_{{pro},{y}} = 0\]
\[\displaystyle [9]\text{ }{\mathbf{{\mathbf{ex\_fop}}}}_{{pro},{y}} - {\mathrm{fopex}}_{{pro},{y}} \cdot {\mathbf{{\mathbf{cap\_p}}}}_{{pro},{y}} = 0\]
\[\displaystyle [10]\text{ }{\mathbf{{\mathbf{{\mathbf{prod}}}}}}_{{pv},{q}} - {\mathbf{{\mathbf{{\mathbf{con}}}}}}_{{solar},{q}} = 0\]
\[\displaystyle [11]\text{ }{\mathbf{{\mathbf{{\mathbf{prod}}}}}}_{{pv},{q}} - {\mathbf{{\mathbf{{\mathbf{prod}}}}}}_{{li},{q}} + {\mathbf{{\mathbf{{\mathbf{prod}}}}}}_{{li\_d},{q}} - {\mathbf{{\mathbf{{\mathbf{sell}}}}}}_{{power},{q}} = 0\]
\[\displaystyle [12]\text{ }{\mathbf{{\mathbf{{\mathbf{prod}}}}}}_{{li},{q}} - {\mathbf{{\mathbf{{\mathbf{prod}}}}}}_{{li\_d},{q}} + {\mathbf{{\mathbf{{\mathbf{inv}}}}}}_{{charge},{{q}-1}} - {\mathbf{{\mathbf{{\mathbf{inv}}}}}}_{{charge},{q}} = 0\]



Functions

\[\displaystyle [0]\text{ }{\mathbf{ex\_cap}}_{{pv},{{{y}_{0}}}} + {\mathbf{ex\_cap}}_{{li},{{{y}_{0}}}} + {\mathbf{ex\_cap}}_{{li\_d},{{{y}_{0}}}} + {\mathbf{ex\_vop}}_{{pv},{{{y}_{0}}}} + {\mathbf{ex\_vop}}_{{li},{{{y}_{0}}}} + {\mathbf{ex\_vop}}_{{li\_d},{{{y}_{0}}}} + {\mathbf{ex\_fop}}_{{pv},{{{y}_{0}}}} + {\mathbf{ex\_fop}}_{{li},{{{y}_{0}}}} + {\mathbf{ex\_fop}}_{{li\_d},{{{y}_{0}}}}\]

Solution#

p.opt()
📝  Generated prog.mps                                                       ⏱ 0.0014 s
Set parameter Username
Academic license - for non-commercial use only - expires 2026-08-01
Read MPS format model from file prog.mps
Reading time = 0.00 seconds
PROG: 40 rows, 31 columns, 81 nonzeros
📝  Generated gurobipy model. See .formulation                               ⏱ 0.0142 s
Gurobi Optimizer version 12.0.3 build v12.0.3rc0 (win64 - Windows 11.0 (26100.2))

CPU model: 13th Gen Intel(R) Core(TM) i7-13700, instruction set [SSE2|AVX|AVX2]
Thread count: 16 physical cores, 24 logical processors, using up to 24 threads

Optimize a model with 40 rows, 31 columns and 81 nonzeros
Model fingerprint: 0x1ce81a32
Coefficient statistics:
  Matrix range     [5e-01, 5e+03]
  Objective range  [1e+00, 1e+00]
  Bounds range     [0e+00, 0e+00]
  RHS range        [5e+01, 2e+02]
Presolve removed 30 rows and 20 columns
Presolve time: 0.00s
Presolved: 10 rows, 11 columns, 25 nonzeros

Iteration    Objective       Primal Inf.    Dual Inf.      Time
       0    4.9949000e+02   3.748738e+01   0.000000e+00      0s
       8    9.4200000e+05   0.000000e+00   0.000000e+00      0s

Solved in 8 iterations and 0.01 seconds (0.00 work units)
Optimal objective  9.420000000e+05
📝  Generated Solution object for prog. See .solution                        ⏱ 0.0001 s
✅  prog optimized using gurobi. Display using .output()                     ⏱ 0.0220 s