Solving linear programming problems using Python
Why solve extreme problems
In practice, very often problems arise for the solution of which optimization methods are used. In everyday life, with multiple choices, for example, gifts for the new year, we intuitively solve the problem of minimum costs with a given quality of purchases.
Unfortunately, one cannot always rely on intuition. Suppose you are an employee of a commercial company and are responsible for advertising. The cost of advertising per month should not exceed 10,000 monetary units (e). A minute of radio advertising costs CU 5, and a television advertisement of 90 CU The company intends to use radio advertising three times more often than television advertising. Practice shows that 1 minute of television advertising provides sales 30 times greater than 1 minute of radio advertising.
Your task is to determine such a distribution of funds between the two mentioned types of advertising at which the company's sales will be maximum. You first select the variables, namely the monthly volume in minutes on television - x1, and on radio - x2. Now it’s not difficult to create the following system:
30x1 + x2 –increase sales from advertising;
90x1 + 5x2 <= 10,000 - funds limit;
x2 = 3x1 - the ratio of the times of the radio and the body of the advertisement.
Now forget about advertising for a while and try to summarize the above. There may be many such tasks as a fixed one, but they have the following common features. Mandatory is the presence of a target function linearly dependent on variables, in our case it is 30x1 + x2which, with the found values of the incoming variables, should have a single maximum value. Moreover, the condition of non-negativity of the incoming variables is satisfied automatically. Then again linear equalities and inequalities follow in an amount depending on the presence of conditions. So we formulated one group of linear programming problems.
We consider another large group of linear programming problems using the so-called transport problem as an example. Suppose you are an employee of a commercial company that provides transportation services. There are suppliers of goods with warehouses in different three cities, and the volumes of homogeneous products in these warehouses are respectively equal to a1, a2, a3. There are consumers in the other three cities who need to bring goods from suppliers in volumes of b1, b2, b3, respectively. Also known are the shipping costs from 1 to 9 goods from suppliers to consumers, according to the table.
If x1 ... xn is the quantity of cargo transported, then the function of the goal will be the total cost of transportation:
F (x) = c1 * x1 + c2 * x2 + c3 * x3 + c4 * x4 + c5 * x5 + c6 * x6 + c7 * x7 + c8 * x8 + c9 * x9.
Conditions to be recorded. in the form of inequalities:
x1 + x2 + x3 <= 20 - you can’t get more than the supplier has
x4 + x5 + x6 <= 45
x7 + x8 + x9 <= 30
Conditions to be recorded. in the form of equalities:
x1 + x4 + x7 = b1– how much is needed and bring
x2 + x5 + x8 = b2
x3 + x6 + x9 = b3
Here we additionally need the conditions for the non-negativity of the variables x since they are meaningfully not negative and the minimum F is sought ( x). These inequalities are not given.
Now you know how to build goal functions and conditions for the basic tasks of linear programming. But when you read in the special literature about the geometric, simplex, artificial basis methods for solving these problems, you abandoned both advertising and logistics. But you can find a simple and understandable solution in Python.
Choosing Python libraries to solve typical linear programming problems
For linear programming in Python, I know of three specialized libraries. Consider the solution to both of these problems in each of the libraries. In addition to the interface and optimization results, we will also evaluate performance. Since we only need a qualitative difference in speed, we will use the simplest listing for this, averaging the results of each program launch.
Optimization with pulp library [1].
Listing of the program for solving the "About Advertising" task
from pulp import *
import time
start = time.time()
x1 = pulp.LpVariable("x1", lowBound=0)
x2 = pulp.LpVariable("x2", lowBound=0)
problem = pulp.LpProblem('0',pulp.LpMaximize)
problem += 30*x1 +x2, "Функция цели"
problem += 90*x1+ 5*x2 <= 10000, "1"
problem +=x2 ==3*x1, "2"
problem.solve()
print ("Результат:")
for variable in problem.variables():
print (variable.name, "=", variable.varValue)
print ("Прибыль:")
print (value(problem.objective))
stop = time.time()
print ("Время :")
print(stop - start)
In the fiscal tenge of the program, the ratios already familiar to us for the maximum advertising profit are 30 * x1 + x2, the conditions for limiting costs, marked “1” for comparison. We have not forgotten about the relationship between the times of using the radio and the advertising body, marked in fox tenge as “2”. The purpose of the other operators is obvious. Details can be found in [1].
The results of solving the optimization problem using pulp .
Result:
x1 = 95.238095
x2 = 285.71429
Profit:
3142.85714
Time:
0.10001182556152344
As a result, we got the values of advertising times at which the expected profit from its use will be maximum.
Optimization with the cvxopt library [2].
Listing of the program to solve the problem "About advertising".
from cvxopt.modeling import variable, op
import time
start = time.time()
x = variable(2, 'x')
z=-(30*x[0] +1*x[1])#Функция цели
mass1 = (90*x[0] + 5*x[1] <= 10000) #"1"
mass2 = (3*x[0] -x[1] == 0) # "2"
x_non_negative = (x >= 0) #"3"
problem =op(z,[mass1,mass2,x_non_negative])
problem.solve(solver='glpk')
problem.status
print ("Прибыль:")
print(abs(problem.objective.value()[0]))
print ("Результат:")
print(x.value)
stop = time.time()
print ("Время :")
print(stop - start)
The structure of the program is similar to the previous one, but there are two significant differences. Firstly, the cvxopt library is configured to find the minimum of the target function, and not the maximum. Therefore, the objective function is taken with a negative minus sign - (30 * x [0] + 1 * x [1]). The resulting negative value is derived in absolute value. Secondly, the restriction on non-negativity of non-negative variables has been introduced. Whether this affected the result, we now see.
Results of solving the optimization problem using cvxopt.
Profit:
3142.857142857143
Result :
[9.52e + 01]
[2.86e + 02]
Time:
0.041656494140625
No significant changes have occurred in comparison with the pulp library, with the exception of the running time of the program.
Optimization with scipy library . optimiz e [3].
Listing of the program to solve the problem "About advertising".
from scipy.optimize import linprog
import time
start = time.time()
c = [-30,-1] #Функция цели
A_ub = [[90,5]] #'1'
b_ub = [10000]#'1'
A_eq = [[3,-1]] #'2'
b_eq = [0] #'2'
print (linprog(c, A_ub, b_ub, A_eq, b_eq))
stop = time.time()
print ("Время :")
print(stop - start)
A quick glance at the listing is enough to understand that we are dealing with a fundamentally different approach to data entry. Although the numbers given in the listings help clarify the principle of data organization, by comparison, I will give explanations. The list c = [-30, -1] contains the coefficients of the target function with the opposite sign, since linprog () is looking for a minimum. The matrix A_ub contains the coefficients of the variables for the conditions in the form of inequalities. For our task, this is 90x1 + 5x2 <= 10000. The values on the right side of inequality-1000 are placed on the b_ub list. The matrix A_eq contains the coefficients of the variables for the conditions in the form of equalities. For our task, 3x1-x2 = 0 , and the zero on the right side is placed in the b_eq list.
The results of solving the optimization problem using scipy. optimize.
Fun: -3142.8571428571431
message: 'Optimization terminated successfully.'
nit: 2
slack: array ([0.])
status: 0
success: True
x: array ([95.23809524, 285.71428571])
Time:
0.03020191192626953
The calculation results are shown in more detail here, but the results themselves are the same as in the previous libraries.
Based on the results of the solution of the “On Advertising” task, an intermediate conclusion can be drawn regarding the use of the scipy library. optimize provides greater speed and streamlined form of source data. However, without the results of solving the transport problem, it is too early to make a final conclusion.
I give a solution to the transport problem, but without detailed explanations, since the main stages of the solution are already described in detail.
Optimization with pulp library .
Listing programs to solve the transport problem.
from pulp import *
import time
start = time.time()
x1 = pulp.LpVariable("x1", lowBound=0)
x2 = pulp.LpVariable("x2", lowBound=0)
x3 = pulp.LpVariable("x3", lowBound=0)
x4 = pulp.LpVariable("x4", lowBound=0)
x5 = pulp.LpVariable("x5", lowBound=0)
x6 = pulp.LpVariable("x6", lowBound=0)
x7 = pulp.LpVariable("x7", lowBound=0)
x8 = pulp.LpVariable("x8", lowBound=0)
x9 = pulp.LpVariable("x9", lowBound=0)
problem = pulp.LpProblem('0',pulp.LpMaximize)
problem += -7*x1 - 3*x2 - 6* x3 - 4*x4 - 8*x5 -2* x6-1*x7- 5*x8-9* x9, "Функция цели"
problem +=x1 + x2 +x3<= 74,"1"
problem +=x4 + x5 +x6 <= 40, "2"
problem +=x7 + x8+ x9 <= 36, "3"
problem +=x1+ x4+ x7 == 20, "4"
problem +=x2+x5+ x8 == 45, "5"
problem +=x3 + x6+x9 == 30, "6"
problem.solve()
print ("Результат:")
for variable in problem.variables():
print (variable.name, "=", variable.varValue)
print ("Стоимость доставки:")
print (abs(value(problem.objective)))
stop = time.time()
print ("Время :")
print(stop - start)
The solution of the transport problem requires minimizing the cost of delivery, so the goal function is introduced with a minus sign, and displayed on the absolute value.
The results of solving the transport problem using pulp.
Result:
x1 = 0.0
x2 = 45.0
x3 = 0.0
x4 = 0.0
x5 = 0.0
x6 = 30.0
x7 = 20.0
x8 = 0.0
x9 = 0.0
Shipping Cost:
215.0
Time:
0.19006609916687012
Optimization with cvxopt library .
Listing programs to solve the transport problem.
from cvxopt.modeling import variable, op
import time
start = time.time()
x = variable(9, 'x')
z=(7*x[0] + 3*x[1] +6* x[2] +4*x[3] + 8*x[4] +2* x[5]+x[6] + 5*x[7] +9* x[8])
mass1 = (x[0] + x[1] +x[2] <= 74)
mass2 = (x[3] + x[4] +x[5] <= 40)
mass3 = (x[6] + x[7] + x[8] <= 36)
mass4 = (x[0] + x[3] + x[6] == 20)
mass5 = (x[1] +x[4] + x[7] == 45)
mass6 = (x[2] + x[5] + x[8] == 30)
x_non_negative = (x >= 0)
problem =op(z,[mass1,mass2,mass3,mass4 ,mass5,mass6, x_non_negative])
problem.solve(solver='glpk')
problem.status
print("Результат:")
print(x.value)
print("Стоимость доставки:")
print(problem.objective.value()[0])
stop = time.time()
print ("Время :")
print(stop - start)
The results of solving the transport problem using cvxopt.
Result:
[0.00e + 00]
[4.50e + 01]
[0.00e + 00]
[0.00e + 00]
[0.00e + 00]
[3.00e + 01]
[2.00e + 01]
[0.00e + 00]
[0.00e + 00]
Delivery cost:
215.0
Time:
0.03001546859741211
Optimization with scipy library. optimize.
Listing programs to solve the transport problem.
from scipy.optimize import linprog
import time
start = time.time()
c = [7, 3, 6,4,8,2,1,5,9]
A_ub = [[1,1,1,0,0,0,0,0,0],
[0,0,0,1,1,1,0,0,0],
[0,0,0,0,0,0,1,1,1]]
b_ub = [74,40,36]
A_eq = [[1,0,0,1,0,0,1,0,0],
[0,1,0,0,1,0,0,1,0],
[0,0,1,0,0,1,0,0,1]]
b_eq = [20,45,30]
print(linprog(c, A_ub, b_ub, A_eq, b_eq))
stop = time.time()
print ("Время :")
print(stop - start)
Results of solving a transport problem using scipy optimize.
fun: 215.0
message: 'Optimization terminated successfully.'
nit: 9
slack: array ([29., 10., 16.])
status: 0
success: True
x: array ([0., 45., 0., 0., 0., 30., 20., 0., 0.])
Time:
0.009982585906982422
An analysis of the solution of two typical linear programming problems using three libraries of similar purpose does not raise doubts about the choice of scipy library. optimize as a leader in data entry compactness and speed.
What's new for using scipy library. optimize when solving linear programming problems
Obtaining from the source matrix, a list for the goal function, as well as filling in the matrices of inequalities A_ub and equalities A_eq programmatically, will simplify the work of data entry and increase the dimension of the original matrix. This will allow solving real optimization problems. Let's see how this can be done with a simple demo that does not claim to be perfect code.
Spoiler heading
import numpy as np
from scipy.optimize import linprog
b_ub = [74,40,36]
b_eq = [20,45,30]
A=np.array([[7, 3,6],[4,8,2],[1,5,9]])
m, n = A.shape
c=list(np.reshape(A,n*m))# Преобразование матрицы A в список c.
A_ub= np.zeros([m,m*n])
for i in np.arange(0,m,1):# Заполнение матрицы условий –неравенств.
for j in np.arange(0,n*m,1):
if i*n<=j<=n+i*n-1:
A_ub [i,j]=1
A_eq= np.zeros([m,m*n])
for i in np.arange(0,m,1):# Заполнение матрицы условий –равенств.
k=0
for j in np.arange(0,n*m,1):
if j==k*n+i:
A_eq [i,j]=1
k=k+1
print(linprog(c, A_ub, b_ub, A_eq, b_eq))
Now only the matrix A itself and the lists of the right-hand sides of b_ub inequalities and b_ub - equalities are introduced.
The result of the program is predictable.
fun: 215.0
message: 'Optimization terminated successfully.'
nit: 9
slack: array ([29., 10., 16.])
status: 0
success: True
x: array ([0., 45., 0., 0., 0., 30., 20., 0., 0.])
The conclusion is private
When solving linear programming problems, it is advisable to use the scipy.optimize library according to the technique described in the article, and fill in the condition matrices programmatically. In this case, you do not need special knowledge about methods for solving optimization problems.
General conclusion
Recently, various Python libraries have appeared that solve similar problems. The decision to use a particular library is often subjective. Therefore, it is advisable to conduct a comparative analysis of them for the area of your tasks.