Optimization problem with Julia - Launch sequence

It looks like this:

using JuMP
using GLPK
using MathOptInterface
using GLPKMathProgInterface
using Gurobi

n = 12

env = Gurobi.Env()


m = Model(Gurobi.Optimizer)
set_optimizer_attribute(m, "NonConvex", 2)

@variable(m, x[i = 1:n, j = 1:n] >= 0)
@variable(m, y[i = 1:n, j = 1:n], Bin)
@variable(m, 1 <= nr[i = 1:n, j = 1:n] <= n, Int)


start = [18, 15, 17, 17, 20, 13, 15, 17, 17, 21, 18, 15]
pop = [100, 80, 130, 90, 1500, 1200, 200, 250, 800, 1000, 1100, 750]
growth = [0.2, 0.15, 0.22, 0.24, 0.12, 0.42, 0.1, 0.11, 0.14, 0.2, 0.05, 0.05]
mgrowth = zeros(n,n)

for i = 1:n
    for j = 1:n
        mgrowth[i,j] = (1 + growth[i])^(n-j)
    end
end

@objective(m, Max, sum(x[i,j]*y[i,j] * pop[j] * mgrowth[i,j] for i =1:n, j = 1:n) )

for j = 1:n
    @constraint(m, sum(y[i,j] for i = 1:n) == 1 )
end

for i = 1:n
    @constraint(m, sum(y[i,j] for j = 1:n) == 1 )
end

for i = 1:n
    @constraint(m, x[i,1]*y[i,1] <= start[i] )
end

for j = 2:n
    @constraint(m, nr[1,j] >= sum(y[2,i] + y[3,i] for i = 1:j-1) )
    @constraint(m, nr[1,j] * x[1,j] == sum(y[2,i]*x[2,i] + y[3,i]*x[3,i] for i = 1:j-1)  )

    @constraint(m, nr[2,j] >= sum(y[3,i] + y[1,i] + y[12,i] for i = 1:j-1) )
    @constraint(m, nr[2,j] * x[2,j] == sum(y[3,i]*x[3,i] + y[1,i]*x[1,i] + y[12,i]*x[12,i]  for i = 1:j-1) )

    @constraint(m, nr[3,j] >= sum(y[2,i] + y[1,i] + y[4,i] + y[5,i] for i = 1:j-1) )
    @constraint(m, nr[3,j] * x[3,j] == sum(y[2,i]*x[2,i] + y[1,i]*x[1,i] + y[4,i]*x[4,i] + y[5,i]*x[5,i] for i = 1:j-1)  )

    @constraint(m, nr[4,j] >= sum(y[1,i] + y[2,i] + y[3,i] for i = 1:j-1) )
    @constraint(m, nr[4,j] * x[4,j] == sum(y[1,i]*x[1,i] + y[2,i]*x[2,i] + y[3,i]*x[3,i]  for i = 1:j-1) )

    @constraint(m, nr[5,j] >= sum(y[1,i] + y[2,i] + y[3,i] + y[4,i] + y[10,i] + y[11,i] + y[12,i] for i = 1:j-1) )
    @constraint(m, nr[5,j] * x[5,j] == sum(y[1,i]*x[1,i] + y[2,i]*x[2,i] + y[3,i]*x[3,i] + y[4,i]*x[4,i]
     + y[10,i]*x[10,i] + y[11,i]*x[11,i] + y[12,i]*x[12,i] for i = 1:j-1)  )

    @constraint(m, nr[6,j] >= sum(y[2,i] + y[12,i] for i = 1:j-1) )
    @constraint(m, nr[6,j] * x[6,j] == sum(y[2,i]*x[2,i] + y[12,i]*x[12,i] for i = 1:j-1)  )

    @constraint(m, nr[7,j] >= sum(y[8,i] for i = 1:j-1) )
    @constraint(m, nr[7,j] * x[7,j] == sum(y[8,i]*x[8,i] for i = 1:j-1)  )

    @constraint(m, nr[8,j] >= sum(y[7,i] for i = 1:j-1) )
    @constraint(m, nr[8,j] * x[8,j] == sum(y[7,i]*x[7,i] for i = 1:j-1)  )

    @constraint(m, nr[9,j] >= sum(y[10,i] + y[11,i] + y[12,i] for i = 1:j-1) )
    @constraint(m, nr[9,j] * x[9,j] == sum(y[10,i]*x[10,i] + y[11,i]*x[11,i] + y[12,i]*x[12,i]  for i = 1:j-1) )

    @constraint(m, nr[10,j] >= sum(y[9,i] + y[11,i] for i = 1:j-1) )
    @constraint(m, nr[10,j] * x[10,j] == sum(y[9,i]*x[9,i] + y[11,i]*x[11,i] for i = 1:j-1)  )

    @constraint(m, nr[11,j] >= sum(y[7,i] + y[8,i] + y[9,i] + y[10,i] + y[12,i] for i = 1:j-1) )
    @constraint(m, nr[11,j] * x[11,j] == sum(y[7,i]*x[7,i] + y[8,i]*x[8,i] + y[9,i]*x[9,i] + y[10,i]*x[10,i] + y[12,i]*x[12,i]  for i = 1:j-1) )

    @constraint(m, nr[12,j] >= sum(y[1,i] + y[2,i] + y[3,i] + y[4,i] + y[5,i] + y[6,i] + y[7,i] + y[8,i] + y[9,i] for i = 1:j-1) )
    @constraint(m, nr[12,j] * x[12,j] == sum(y[1,i]*x[1,i] + y[2,i]*x[2,i] + y[3,i]*x[3,i] + y[4,i]*x[4,i]
     + y[5,i]*x[5,i] + y[6,i]*x[6,i] + y[7,i]*x[7,i] + y[8,i]*x[8,i] + y[9,i]*x[9,i] for i = 1:j-1)  )
end

optimize!(m)
println(JuMP.value.(x))
println(JuMP.value.(y))
println(JuMP.value.(nr))
println(JuMP.objective_value(m))
println(JuMP.termination_status(m))

As you can see, I added a mgrowth to the objective which is a shop-specific monthly growth in the customer size.