Note
You can download this example as a Jupyter notebook or start it in interactive mode.
Solve a basic optimization model
In this example, we explain the basic functions of the linopy Model class. First, we are setting up a very simple linear optimization model, given by
Minimize:
subject to:
[1]:
from linopy import Model
[2]:
m = Model()
The Model class serves at a container for all the relevant data.
Adding variables
Let’s add the two variables. Note that a variable can always be assigned with a lower and an upper bound. In this case, both x and y have a lower bound of zero (coming from the first two constraints). Note, the default for lower and upper bounds are minus and plus infinity.
[3]:
x = m.add_variables(lower=0, name='x')
y = m.add_variables(lower=0, name='y');
x and y are linopy variables. Each of them consist of an array with the variable references in the model, that is the names of variables that will finally be used when solving the model.
[4]:
x
[4]:
Variable
--------
x ∈ [0, inf]
Since both x and y are scalar variables, so their arrays also contain just one variable reference. The variable x points to the optimisation variable 1 and the variable y points to the optimisation variable 2. Later we will see the benefit of this behaviour.
Adding Constraints
Constraints consist of the left hand side (lhs) and the righ hand side (rhs). The lhs must contain all the variables with the corresponding coefficients, while the rhs is just constant. The first lhs of our two remaining constraints is \(3x + 7y\), which we can write just exactly in this way
[5]:
3*x + 7*y
[5]:
LinearExpression:
-----------------
3.0 x + 7.0 y
When assigning to the model, we call the function m.add_constraints.
[6]:
m.add_constraints(3*x + 7*y >= 10)
m.add_constraints(5*x + 2*y >= 3);
Adding the Objective
We do the same for defining the objective while the objective function only consists of a linear expression.
[7]:
m.add_objective(x + 2*y)
[7]:
LinearExpression:
-----------------
1.0 x + 2.0 y
[8]:
m.solve()
Restricted license - for non-production use only - expires 2024-10-28
Read LP format model from file /tmp/linopy-problem-_gb0qr9s.lp
Reading time = 0.00 seconds
obj: 2 rows, 2 columns, 4 nonzeros
Gurobi Optimizer version 10.0.1 build v10.0.1rc0 (linux64)
CPU model: Intel(R) Xeon(R) Platinum 8259CL CPU @ 2.50GHz, instruction set [SSE2|AVX|AVX2|AVX512]
Thread count: 1 physical cores, 2 logical processors, using up to 2 threads
Optimize a model with 2 rows, 2 columns and 4 nonzeros
Model fingerprint: 0x4c5ee539
Coefficient statistics:
Matrix range [2e+00, 7e+00]
Objective range [1e+00, 2e+00]
Bounds range [0e+00, 0e+00]
RHS range [3e+00, 1e+01]
Presolve time: 0.00s
Presolved: 2 rows, 2 columns, 4 nonzeros
Iteration Objective Primal Inf. Dual Inf. Time
0 0.0000000e+00 1.625000e+00 0.000000e+00 0s
2 2.8620690e+00 0.000000e+00 0.000000e+00 0s
Solved in 2 iterations and 0.01 seconds (0.00 work units)
Optimal objective 2.862068966e+00
[8]:
('ok', 'optimal')
The solution of the linear problem is stored in m.solution in fom of a xarray.Dataset.
[9]:
m.solution
[9]:
<xarray.Dataset>
Dimensions: ()
Data variables:
x float64 0.03448
y float64 1.414Expanding the dimensionality
Now comes the interesting part. Suppose the two variables x and y are a function of time t and a modified problem setup like
Minimize:
subject to:
whereas t spans all the range from 0 to 10.
First, we define a new model (we just overwrite the old m).
[10]:
m = Model()
Again, we define x and y using the add_variables function, but now we are adding a coords argument. This automatically creates optimization variables for all coordinates, in this case time-steps.
[11]:
import pandas as pd
time = pd.Index(range(10), name='time')
x = m.add_variables(lower=0, coords=[time], name='x', )
y = m.add_variables(lower=0, coords=[time], name='y')
We again write the constraints out using the syntax from above, while multiplying the rhs with t. Note that the coordinates from the lhs and the rhs have to match.
Note: In case lhs and rhs have different sets of coordinates, the constraint creation is broadcasted over all combinations of coordinates.
[12]:
factor = pd.Series(time, index=time)
m.add_constraints(3*x + 7*y >= 10*factor, name='Constraint1')
m.add_constraints(5*x + 2*y >= 3*factor, name='Constraint2')
m
[12]:
Linopy model
============
Variables:
----------
Dimensions: (time: 10)
Coordinates:
* time (time) int64 0 1 2 3 4 5 6 7 8 9
Data:
x (time) int64 0 1 2 3 4 5 6 7 8 9
y (time) int64 10 11 12 13 14 15 16 17 18 19
Constraints:
------------
Dimensions: (time: 10)
Coordinates:
* time (time) int64 0 1 2 3 4 5 6 7 8 9
Data:
Constraint1 (time) int64 0 1 2 3 4 5 6 7 8 9
Constraint2 (time) int64 10 11 12 13 14 15 16 17 18 19
Status:
-------
initialized
When we add the objective, we use the sum function of linopy.LinearExpression. This stacks the term_ dimension, such that all terms of the time dimension are rewritten into one big expression.
[13]:
obj = (x + 2*y).sum()
m.add_objective(obj)
[13]:
LinearExpression:
-----------------
1.0 x[0] + 2.0 y[0] + 1.0 x[1] ... 2.0 y[8] + 1.0 x[9] + 2.0 y[9]
[14]:
m.solve()
Read LP format model from file /tmp/linopy-problem-t6bdcgkw.lp
Reading time = 0.00 seconds
obj: 20 rows, 20 columns, 40 nonzeros
Gurobi Optimizer version 10.0.1 build v10.0.1rc0 (linux64)
CPU model: Intel(R) Xeon(R) Platinum 8259CL CPU @ 2.50GHz, instruction set [SSE2|AVX|AVX2|AVX512]
Thread count: 1 physical cores, 2 logical processors, using up to 2 threads
Optimize a model with 20 rows, 20 columns and 40 nonzeros
Model fingerprint: 0x89bb3e20
Coefficient statistics:
Matrix range [2e+00, 7e+00]
Objective range [1e+00, 2e+00]
Bounds range [0e+00, 0e+00]
RHS range [3e+00, 9e+01]
Presolve removed 2 rows and 2 columns
Presolve time: 0.01s
Presolved: 18 rows, 18 columns, 36 nonzeros
Iteration Objective Primal Inf. Dual Inf. Time
0 0.0000000e+00 7.312500e+01 0.000000e+00 0s
18 1.2879310e+02 0.000000e+00 0.000000e+00 0s
Solved in 18 iterations and 0.01 seconds (0.00 work units)
Optimal objective 1.287931034e+02
[14]:
('ok', 'optimal')
[15]:
m.solution.to_dataframe().plot(grid=True, ylabel='Optimal Value')
[15]:
<Axes: xlabel='time', ylabel='Optimal Value'>
You can easily add dimensions in this way. Since the underlying code is highly relying on xarray broadcasting on new dimensions is automatically supported, so even the rhs can introduce new dimensions. Just be careful that all variable dimensions have the correct dimension name.
We limit the slope of the y variable by adding a constraint in the form of
[16]:
lhs = (y - y.shift(time=1)).sel(time=time[1:])
m.add_constraints(lhs, '<=', 0.5, name='Limited growth y')
[16]:
Constraint `Limited growth y` (time: 10) - 1 masked entries
-----------------------------------------------------------
[1]: 1.0 y[1] - 1.0 y[0] <= 0.5
[2]: 1.0 y[2] - 1.0 y[1] <= 0.5
[3]: 1.0 y[3] - 1.0 y[2] <= 0.5
[4]: 1.0 y[4] - 1.0 y[3] <= 0.5
[5]: 1.0 y[5] - 1.0 y[4] <= 0.5
[6]: 1.0 y[6] - 1.0 y[5] <= 0.5
[7]: 1.0 y[7] - 1.0 y[6] <= 0.5
[8]: 1.0 y[8] - 1.0 y[7] <= 0.5
[9]: 1.0 y[9] - 1.0 y[8] <= 0.5
[17]:
m.solve()
m.solution.to_dataframe().plot(grid=True, ylabel='Optimal Value');
Read LP format model from file /tmp/linopy-problem-wwt1ohtb.lp
Reading time = 0.00 seconds
obj: 29 rows, 20 columns, 58 nonzeros
Gurobi Optimizer version 10.0.1 build v10.0.1rc0 (linux64)
CPU model: Intel(R) Xeon(R) Platinum 8259CL CPU @ 2.50GHz, instruction set [SSE2|AVX|AVX2|AVX512]
Thread count: 1 physical cores, 2 logical processors, using up to 2 threads
Optimize a model with 29 rows, 20 columns and 58 nonzeros
Model fingerprint: 0xbd65942b
Coefficient statistics:
Matrix range [1e+00, 7e+00]
Objective range [1e+00, 2e+00]
Bounds range [0e+00, 0e+00]
RHS range [5e-01, 9e+01]
Presolve removed 2 rows and 1 columns
Presolve time: 0.00s
Presolved: 27 rows, 19 columns, 54 nonzeros
Iteration Objective Primal Inf. Dual Inf. Time
0 0.0000000e+00 7.312500e+01 0.000000e+00 0s
19 1.4158621e+02 0.000000e+00 0.000000e+00 0s
Solved in 19 iterations and 0.01 seconds (0.00 work units)
Optimal objective 1.415862069e+02
