4. Networks 3: Production-Inventory Planning: Single Product¶

ISE 754, Fall 2024

Package Used: No new packages used.

Production-inventory planning models are one of the main uses of mathematical programming in industry. They provide a means to make complex decisions over a rolling planning horizon.

1. Single-Product with Variable Capacity: LP¶

The following production-inventory LP planning model is for a single product, produced in m stages with a rolling horizon of t periods (e.g., months). Determining if production capacity $u_{ij}$ should be available for each period, e.g., $u_{ij} = u_i$ or $u_{ij} = 0$ is used to plan production.

$ \begin{eqnarray*} \quad \mbox{Minimize} \quad \sum_{i=1}^{m} \sum_{j=1}^{t} c_i^p x_{ij} + \sum_{i=1}^{m} \sum_{j=2}^{t+1} c_i^i y_{ij} \\ \quad \mbox{subject to} \quad -x_{ij} + x_{(i+1)j} - y_{ij} + y_{i(j+1)} &=& 0 \mbox{,} &\quad i = 1, \dots, m-1; j = 1, \dots, t &\quad (a) \\ -x_{m,j} - y_{m,j} + y_{m(j+1)} &=& d_j \mbox{,} &\quad j = 1, \dots, t &\quad (b) \\ x_{ij} &\le& u_{ij} \mbox{,} &\quad i = 1, \dots, m; j = 1, \dots, t \\ y_{i,1} &=& y_i^0 \mbox{,} &\quad i = 1, \dots, m \\ y_{i(t+1)} &=& y_i^{t+1} \mbox{,} &\quad i = 1, \dots, m \\ x, y &\ge& 0 \mbox{,} \\ \end{eqnarray*} $

where,

$ \begin{eqnarray*} \quad m &=& \mbox{number of production stages} \\ t &=& \mbox{number of periods of production} \\ c_i^p &=& \mbox{production cost (dollar/ton) in stage } i \\ x_{ij} &=& \mbox{production (ton) at stage } i \mbox{ in period } j \\ c_i^i &=& \mbox{inventory cost (dollar/ton) in stage } i \\ y_{ij} &=& \mbox{stage-} i \mbox{ inventory (ton) from period } j-1 \mbox{ to } j \\ d_j &=& \mbox{demand (ton) in period } j \\ u_{ij} &=& \mbox{production capacity (ton) of stage } i \mbox{ in period } j \\ y_i^0 &=& \mbox{initial inventory (ton) of stage } i \\ y_i^t &=& \mbox{final inventory (ton) of stage } i \mbox{.} \\ \end{eqnarray*} $

In the objective function, inventory costs are summed from 2 to t + 1 because the period 1 initial inventory cost was accounted for in period 0 prior to period 1 (rolling horizon). Constraints (a) and (b) in the model are the flow-balance equations for each node, while the remaining constraints represent lower- and upper-bounds on the decision variables x and y.

Ex: Capacity Availability as Input to Model¶

A single product is produced in two stages. The planning period is three months, and the model is run each month using current demand forecasts and scheduled production capacity. In some periods, the capacity is input as zero, which indicates that it will be used for other purposes during that period. The following shows the production-inventory planning network:

image.png

In the following formulation, the flow-balance constraints are the only constraints in the model because the capacity constraints are handled as upper bounds on $x_{ij}$. Inventory costs $c_i^i$ can be determined by multiplying an annual inventory carrying rate, $h$, by the cummulative production cost:

$ \begin{eqnarray*} \quad c_i^i &=& h\sum_{k=1}^{i} c_k^p \\ h &=& \mbox{per-period inventory carrying rate (\$/\$-period) .} \\ \end{eqnarray*} $

Assuming an annual inventory carrying rate of 0.3 \$/\$-yr, the monthly rate is $ h = 0.3/12$. The inventory cost per ton at stage i is then the monthly rate times the cumulative cost to produce each ton from stage 1 through i:

$ \begin{eqnarray*} \quad c_1^i &=& \dfrac{0.3}{12} 200 = 5 \\ c_2^i &=& \dfrac{0.3}{12} \bigl( 200 + 800 \bigr) = 25 \end{eqnarray*} $.

In [1]:
u = [50  0 50;                     # capacity of stage m in period t (ton)
     60  0  0];
d = [20, 10, 15]                   # demand in period t (ton)
cᵖ = [200, 800]                    # production cost in stage m ($/ton)
h = 0.3/12                         # monthly inventory carrying rate (1/month)
@show cⁱ = cumsum(cᵖ)*h            # inventory cost for stage m ($/ton)
y⁰ = [0, 0]                        # initial inventory of stage m (ton)
yᵗ⁺¹ = [35, 0]                     # final inventory of stage m (ton)
m = length(cᵖ)                     # number of production stages = 2
t = length(d)                      # number of periods of production = 3
LBʸ = [y⁰ zeros(m, t-1) yᵗ⁺¹]      # same LB and UB to fix inventory level
UBʸ = [y⁰ fill(Inf, m, t-1) yᵗ⁺¹]

cⁱ = cumsum(cᵖ) * h = [5.0, 24.999999999999996]
Out[1]:
2×4 Matrix{Float64}:
 0.0  Inf  Inf  35.0
 0.0  Inf  Inf   0.0
In [2]:
using JuMP, Cbc

model = Model(Cbc.Optimizer)
@variable(model, 0 <= x[i=1:m, j=1:t] <= u[i,j] )
@variable(model, LBʸ[i,j] <= y[i=1:m, j=1:t+1] <= UBʸ[i,j] )
@objective(model, Min, 
    sum(cᵖ[i]*x[i,j] for i=1:m, j=1:t ) + 
    sum(cⁱ[i]*y[i,j] for i=1:m, j=2:t+1 ) )
# (a) flow-balance all but last stage
@constraint(model, [i=1:m-1, j=1:t], -x[i,j] + x[i+1,j] - y[i,j] + y[i,j+1] == 0 )
# (b) flow-balance for last stage
@constraint(model, [j=1:t], -x[m,j] - y[m,j] + y[m,j+1] == -d[j] )
print(model)
set_attribute(model, "logLevel", 0)
optimize!(model)
println(solution_summary(model).termination_status)
TCᵒ, xᵒ, yᵒ = objective_value(model), Array(value.(x)), Array(value.(y))
$$ \begin{aligned} \min\quad & 200 x_{1,1} + 200 x_{1,2} + 200 x_{1,3} + 800 x_{2,1} + 800 x_{2,2} + 800 x_{2,3} + 5 y_{1,2} + 5 y_{1,3} + 5 y_{1,4} + 24.999999999999996 y_{2,2} + 24.999999999999996 y_{2,3} + 24.999999999999996 y_{2,4}\\ \text{Subject to} \quad & -x_{1,1} + x_{2,1} - y_{1,1} + y_{1,2} = 0\\ & -x_{1,2} + x_{2,2} - y_{1,2} + y_{1,3} = 0\\ & -x_{1,3} + x_{2,3} - y_{1,3} + y_{1,4} = 0\\ & -x_{2,1} - y_{2,1} + y_{2,2} = -20\\ & -x_{2,2} - y_{2,2} + y_{2,3} = -10\\ & -x_{2,3} - y_{2,3} + y_{2,4} = -15\\ & x_{1,1} \geq 0\\ & x_{2,1} \geq 0\\ & x_{1,2} \geq 0\\ & x_{2,2} \geq 0\\ & x_{1,3} \geq 0\\ & x_{2,3} \geq 0\\ & y_{1,1} \geq 0\\ & y_{2,1} \geq 0\\ & y_{1,2} \geq 0\\ & y_{2,2} \geq 0\\ & y_{1,3} \geq 0\\ & y_{2,3} \geq 0\\ & y_{1,4} \geq 35\\ & y_{2,4} \geq 0\\ & x_{1,1} \leq 50\\ & x_{2,1} \leq 60\\ & x_{1,2} \leq 0\\ & x_{2,2} \leq 0\\ & x_{1,3} \leq 50\\ & x_{2,3} \leq 0\\ & y_{1,1} \leq 0\\ & y_{2,1} \leq 0\\ & y_{1,4} \leq 35\\ & y_{2,4} \leq 0\\ \end{aligned} $$ 
OPTIMAL
Presolve 0 (-6) rows, 0 (-14) columns and 0 (-21) elements
Optimal - objective value 53175
After Postsolve, objective 53175, infeasibilities - dual 0 (0), primal 0 (0)
Optimal objective 53175 - 0 iterations time 0.002, Presolve 0.00
Out[2]:
(53175.0, [45.0 0.0 35.0; 45.0 0.0 0.0], [0.0 -0.0 0.0 35.0; 0.0 25.0 15.0 0.0])

The following shows the optimal production and inventory flows (in blue) in the production-inventory planning network. Note that the flow is balanced at each node; that is, the flow into a node equals the combined production, inventory, and demand flow out of the node.

image.png

Optimal production flows (and required demands):

In [3]:
[xᵒ; d']
Out[3]:
3×3 Matrix{Float64}:
 45.0   0.0  35.0
 45.0   0.0   0.0
 20.0  10.0  15.0

Optimal inventory flows:

In [4]:
yᵒ
Out[4]:
2×4 Matrix{Float64}:
 0.0  -0.0   0.0  35.0
 0.0  25.0  15.0   0.0

Ex: Scheduled Maintenance¶

A single product is produced in a two-stage production process. Stage one has a capacity to produce 40 tons of the product per month, and there is a cost of \$200 for each ton of product produced. For stage two, capacity is 80 tons per month, and there is a cost of \$800 per ton. Due to scheduled maintenance, no capacity will be available during month five for stage 1 and during months three and five for stage 2. Demand for the next six months is estimated to be 25, 15, 10, 50, 25, and 15 tons per month. Currently, five tons of material is in storage at stage 1 and 35 tons at stage 2; at the end of the planning period, ten tons of material should be in storage at stage 1 and 15 at stage 2. Up to 30 tons of material per month can be stored at each stage. Assuming that the inventory costs are \$12 and \$46 per ton-month, determine the production plan that minimizes total costs over the planning horizon and what the total production costs will be, and what the total inventory costs will be.

Note: The maximum inventory storage constraint can be modeled as an upper bound on y.

In [5]:
u = [40 40 40 40 0 40;             # capacity of each stage (ton)
     80 80  0 80 0 80]            
d = [25, 15, 10, 50, 25, 15]       # demand in each period (ton)
cᵖ = [200, 800]                    # variable production cost for each stage ($/ton)
cⁱ = [12, 46]                      # inventory cost for each stage ($/ton)
y⁰ = [5, 25]                       # initial inventory of each stage (ton)
yᵗ⁺¹ = [10, 15]                    # final inventory of each stage (ton)
m = length(cᵖ)                     # number of production stages
t = length(d)                      # number of periods of production
LBʸ = [y⁰ zeros(m, t-1) yᵗ⁺¹]      # same LB and UB to fix inital and final inventory levels
UBʸ = [y⁰ fill(30, m, t-1) yᵗ⁺¹]   # max inventory storage level
Out[5]:
2×7 Matrix{Int64}:
  5  30  30  30  30  30  10
 25  30  30  30  30  30  15
In [6]:
using JuMP, Cbc

model = Model(Cbc.Optimizer)
@variable(model, 0 <= x[i=1:m, j=1:t] <= u[i,j] )
@variable(model, LBʸ[i,j] <= y[i=1:m, j=1:t+1] <= UBʸ[i,j] )
@objective(model, Min, 
    sum(cᵖ[i]*x[i,j] for i=1:m, j=1:t ) + 
    sum(cⁱ[i]*y[i,j] for i=1:m, j=2:t+1 ) )
# (a) flow-balance all but last stage
@constraint(model, [i=1:m-1, j=1:t], -x[i,j] + x[i+1,j] - y[i,j] + y[i,j+1] == 0 )
# (b) flow-balance for last stage
@constraint(model, [j=1:t], -x[m,j] - y[m,j] + y[m,j+1] == -d[j] )
set_attribute(model, "logLevel", 0)
optimize!(model)
println(solution_summary(model).termination_status)
TCᵒ, xᵒ, yᵒ = objective_value(model), Array(value.(x)), Array(value.(y))

OPTIMAL
Presolve 5 (-7) rows, 12 (-14) columns and 18 (-24) elements
0  Obj 96109.6 Primal inf 80.099997 (4) Dual inf 992 (3)
8  Obj 134300
Optimal - objective value 134300
After Postsolve, objective 134300, infeasibilities - dual 0 (0), primal 0 (0)
Optimal objective 134300 - 8 iterations time 0.002, Presolve 0.00
Out[6]:
(134300.0, [0.0 25.0 … 0.0 40.0; -0.0 30.0 … 0.0 30.0], [5.0 5.0 … 0.0 10.0; 25.0 0.0 … 0.0 15.0])

Optimal production flows (and required demands):

In [7]:
[xᵒ; d']
Out[7]:
3×6 Matrix{Float64}:
  0.0  25.0  30.0  40.0   0.0  40.0
 -0.0  30.0   0.0  70.0   0.0  30.0
 25.0  15.0  10.0  50.0  25.0  15.0

Optimal inventory flows:

In [8]:
yᵒ
Out[8]:
2×7 Matrix{Float64}:
  5.0  5.0   0.0  30.0  -0.0  0.0  10.0
 25.0  0.0  15.0   5.0  25.0  0.0  15.0

Calculate total production cost and total inventory cost:

In [9]:
TCᵖ, TCⁱ = sum(cᵖ.*xᵒ), sum(cⁱ.*yᵒ)
Out[9]:
(131000.0, 4510.0)

Sum of component costs exceeds total cost from solver because it includes first-period (time 0) inventory costs:

In [10]:
TCᵖ + TCⁱ, TCᵒ
Out[10]:
(135510.0, 134300.0)

Exclude first period from inventory costs:

In [11]:
TCⁱ = sum(cⁱ.*yᵒ[:,2:end])
Out[11]:
3300.0

Total component costs now match total cost from solver:

In [12]:
TCᵖ + TCⁱ, TCᵒ
Out[12]:
(134300.0, 134300.0)

2. Single-Product with Fixed Cost: MILP¶

In the production-inventory planning model presented above, unit production costs $c^p$ are constant and the production capacity available at each stage for each period was given; as a result, the model could be solved as an LP. Now, due to scale economies associated with batch production, production costs are separated into a fixed cost, $c^k$, and a variable cost, $c^p$, and the model will include making a decision at each stage for each period regarding how much, if any, capacity is available. The capacity is assumed to be available in discrete units that are either at 100% of the unit's capacity or not available at all. Binary decision variables are used to make these decisions. This internal capacity decision requires the model to be solved as a MILP.

MILP Model¶

$ \begin{eqnarray*} \quad \mbox{Minimize} \quad \sum_{i=1}^{m} \sum_{j=1}^{t} c_i^p x_{ij} + \sum_{i=1}^{m} \sum_{j=2}^{t+1} c_i^i y_{ij} &+& \sum_{i=1}^{m} \sum_{j=1}^{t} c_i^k z_{ij} \\ \quad \mbox{subject to} \quad -x_{ij} + x_{(i+1)j} - y_{ij} + y_{i(j+1)} &=& 0 \mbox{,} &\quad i = 1, \dots, m-1; j = 1, \dots, t &\quad (a) \\ -x_{mj} - y_{mj} + y_{m(j+1)} &=& d_j \mbox{,} &\quad j = 1, \dots, t &\quad (b) \\ x_{ij} &\le& u_{i} z_{ij} \mbox{,} &\quad i = 1, \dots, m; j = 1, \dots, t &\quad (c) \\ y_{i,1} &=& y_i^0 \mbox{,} &\quad i = 1, \dots, m \\ y_{i(t+1)} &=& y_i^{t+1} \mbox{,} &\quad i = 1, \dots, m \\ x, y \ge 0; z &\in& \bigl\{ 0, 1 \bigr\} \mbox{,} \\ \end{eqnarray*} $

where,

$ \begin{eqnarray*} \quad m &=& \mbox{number of production stages} \\ t &=& \mbox{number of periods of production} \\ c_i^p &=& \mbox{variable production cost (dollar/ton) in stage } i \\ x_{ij} &=& \mbox{production (ton) at stage } i \mbox{ in period } j \\ c_i^i &=& \mbox{inventory cost (dollar/ton) in stage } i \\ y_{ij} &=& \mbox{stage-} i \mbox{ inventory (ton) from period } j-1 \mbox{ to } j \\ c_i^k &=& \mbox{fixed production cost (dollar/period) in stage } i \\ z_{ij} &=& \mbox{production indicator at stage } i \mbox{ in period } j \\ d_j &=& \mbox{demand (ton) in period } j \\ u_{i} &=& \mbox{production capacity (ton) of stage } i \mbox{ in period } j \\ y_i^0 &=& \mbox{initial inventory (ton) of stage } i \\ y_i^{t+1} &=& \mbox{final inventory (ton) of stage } i \mbox{.} \\ \end{eqnarray*} $

In the objective function, inventory costs are summed from 2 to t + 1 because the period 1 initial inventory cost was accounted for in period 0 prior to period 1 (rolling horizon). Constraints (a) and (b) in the model are the flow-balance equations for each node. Constraint (c) represents available production capacity and uses the binary variables z to indicate whether production occurs.

Ex: Single-Product with Fixed Cost¶

A single product is produced in a two-stage production process. Stage one has the capacity to produce 40 tons of the product per month, and there is a fixed cost of \$6000 per month and a variable cost of \$200 for each ton of product produced. For stage two, capacity is 80 tons per month, and there is a fixed cost of \$18,000 and a variable cost of \$800. The fixed costs are only incurred when any amount of production occurs at that stage for that month. Demand for the next six months is estimated to be 25, 15, 10, 50, 25, and 15 tons per month. Currently, five tons of material are in storage at stage one and 35 tons at stage two; at the end of the planning period, ten tons of material should be in storage at stage one and 15 at stage two. Up to 30 tons of material per month can be stored at each stage. Assuming that the annual inventory carrying rate is $0.4 \mbox{ \$/\$-yr }$, determine the production plan that minimizes total costs over the planning horizon and what those total costs will be.

The model in Section 1 above can be extended by including a production indicator binary variable to turn capacity on or off at each stage for each period. One modeling issue is determining inventory cost per ton if the production cost includes a fixed component. An approximation is to allocate the fixed cost based on capacity:

$ \quad c_i^i = h\sum_{j=1}^{i} \left( c_j^p + \dfrac{c_j^k}{u_j} \right) $

The maximum inventory storage constraint can be modeled as an upper bound on y.

In [13]:
u = [40, 80]                       # capacity of each stage (ton)
d = [25, 15, 10, 50, 25, 15]       # demand in each period (ton)
cᵖ = [200, 800]                    # variable production cost for each stage ($/ton)
cᵏ = [6000, 18000]                 # fixed production cost for each stage ($/month)
h = 0.4/12                         # monthly inventory carrying rate (1/month)
cⁱ = cumsum(cᵖ + cᵏ./u)*h          # inventory cost for each stage ($/ton)
y⁰ = [5, 25]                       # initial inventory of each stage (ton)
yᵗ⁺¹ = [10, 15]                    # final inventory of each stage (ton)
m = length(u)                      # number of production stages
t = length(d)                      # number of periods of production
LBʸ = [y⁰ zeros(m, t-1) yᵗ⁺¹]      # same LB and UB to fix initial and final inventory levels
UBʸ = [y⁰ fill(30, m, t-1) yᵗ⁺¹]   # max inventory storage level
Out[13]:
2×7 Matrix{Int64}:
  5  30  30  30  30  30  10
 25  30  30  30  30  30  15
In [14]:
using JuMP, Cbc

model = Model(Cbc.Optimizer)
@variable(model, 0 <= x[i=1:m, j=1:t] )
@variable(model, LBʸ[i,j] <= y[i=1:m, j=1:t+1] <= UBʸ[i,j] )
@variable(model, z[i=1:m, j=1:t], Bin )
@objective(model, Min, 
    sum(cᵖ[i]*x[i,j] for i=1:m, j=1:t ) + 
    sum(cⁱ[i]*y[i,j] for i=1:m, j=2:t+1 ) + 
    sum(cᵏ[i]*z[i,j] for i=1:m, j=1:t ) )
# (a) flow-balance all but last stage
@constraint(model, [i=1:m-1, j=1:t], -x[i,j] + x[i+1,j] - y[i,j] + y[i,j+1] == 0 )  
# (b) flow-balance for last stage
@constraint(model, [j=1:t], -x[m,j] - y[m,j] + y[m,j+1] == -d[j] )
# (c) production capacity
@constraint(model, [i=1:m, j=1:t], x[i,j] <= u[i]*z[i,j] )
set_attribute(model, "logLevel", 0)
set_attribute(model, "seconds", 30.0)   # Solution timeout limit
optimize!(model)
println(solution_summary(model).termination_status)
TCᵒ, xᵒ, yᵒ, zᵒ = objective_value(model), 
    Array(value.(x)), Array(value.(y)), Array(round.(value.(z)))

OPTIMAL
Out[14]:
(212275.0, [0.0 24.999999999999986 … 0.0 40.0; 0.0 29.999999999999986 … 0.0 29.999999999999996], [5.0 5.0 … 0.0 10.0; 25.0 0.0 … 0.0 15.0], [0.0 1.0 … 0.0 1.0; 0.0 1.0 … 0.0 1.0])

Optimal production flows (and required demands):

In [15]:
[xᵒ; d']
Out[15]:
3×6 Matrix{Float64}:
  0.0  25.0  30.0  40.0   0.0  40.0
  0.0  30.0   0.0  70.0   0.0  30.0
 25.0  15.0  10.0  50.0  25.0  15.0

Optimal inventory flows:

In [16]:
yᵒ
Out[16]:
2×7 Matrix{Float64}:
  5.0  5.0   0.0  30.0   0.0  0.0  10.0
 25.0  0.0  15.0   5.0  25.0  0.0  15.0

Optimal production indicators:

In [17]:
zᵒ
Out[17]:
2×6 Matrix{Float64}:
 0.0  1.0  1.0  1.0  0.0  1.0
 0.0  1.0  0.0  1.0  0.0  1.0

Calculate total variable production, inventory, and fixed production cost components:

In [18]:
TCᵖ, TCⁱ, TCᵏ = sum(cᵖ.*xᵒ), sum(cⁱ.*yᵒ[:,2:end]), sum(cᵏ.*zᵒ)
Out[18]:
(130999.99999999999, 3274.9999999999986, 78000.0)

Total component costs match total cost from solver:

In [19]:
TCᵖ + TCⁱ + TCᵏ, TCᵒ
Out[19]:
(212274.99999999997, 212275.0)
In [ ]: