3D expression in JuMP

Specific Domains Optimization (Mathematical)
,
1

I want to build the following expression in the optimisation problem defined below.
A_{i0s} = x_i
A_{its} = A_{i,t-1,s} + y_{its} - z_{its} \quad \forall s, \forall i, \forall t > 0
I have given the code that I tried to build this expression, but I am getting it wrong. Any ideas? The variables y and z index t (time) from 1:T, but the expression need to index time from 0:T, where t =0 refers to the initial position.

Update: Just to add, I am able to achieve the desired outcome by defining A as variable and then defining a set of constraints. But I think these constraints can be avoided by defining this as expression.

model = Model(Cbc.Optimizer)
n = 10
T = 5
S= collect(1:3)
@variable(model, 0 <= x[i=1:n] <= 1)
@variable(model, 0 <= y[i=1:n, t=1:T, s in S] <= 1)
@variable(model, 0 <= z[i=1:n, t=1:T, s in S] <= 1)

Code for building expression (but incorrect)

# For t = 0
for i=1:n
    for s in S
        add_to_expression!(A[i, 0, s], x[i])
    end
end

# For t = 1:T
for t=1:T
    for i=1:n
        for s in S
            add_to_expression!(A[i, t, s], A[i, t-1, s] + y[i, t, s] - z[i, t, s])
        end
    end
end
# Define it as variable instead 
# Define an auxiliary variable for A
@variable(model, A[1:n, 0:T, S])

# Constraints for the initialization at t=0
for i = 1:n
    for s in S
        @constraint(model, A[i, 0, s] == x[i])
    end
end

# Constraints for the iteration relationship from t=1 to T
for t = 1:T
    for i = 1:n
        for s in S
            @constraint(model, A[i, t, s] == A[i, t-1, s] + y[i, t, s] - z[i, t, s])
        end
    end
end
1 Like
2

Here’s how I would write it:

using JuMP
n = 2
T = 2
S = ["a", "b"]
model = Model()
@variable(model, x[1:n])
@variable(model, y[1:n, 1:T, S])
@variable(model, z[1:n, 1:T, S])
@expression(model, A[i in 1:n, t in 0:T, s in S], AffExpr(0.0))
for i in 1:n, s in S, t in 0:T
    A[i, t, s] = if t == 0
        x[i]
    else
        @expression(model, A[i, t-1, s] + y[i, t, s] - z[i, t, s])
    end
end

As a note of caution: if you (or future readers) are generalizing this, the A[i, t-1, s] term must have already been created in the for-loop. The A[i, t, s] will not update if you later modify A[i, t-1, s].

1 Like
3

Two additional comments:
First, I would recommend you use HiGHS rather than Cbc. The former is a lot faster and numerically robust than Cbc. It’s open-source, and you can access it via the HiGHS.jl package.

Second, since x, y and z are decision variables, I would expect that defining A as a variable and adding constraints as in your initial message would be more efficient.
This is because, if you expand the A[i,t,s] expression, it is a the sum of x[i] and t y and z variables, which is 2*t+1 coefficients. That means that overall, your problem has O(n \times T^2 \times S). In contrast, when you declare A as a variable, you get O(n\times T \times S) coefficients. If T is large, the gain can be substantial.
LP and MILP solvers like sparse problems :slight_smile:

2 Likes
4

Thank you @odow and @mtanneau very useful

1 Like