Linear programming in python by the scipy library

In my first publication, I want to talk about how you can quickly and simply solve the problem of linear programming with the help of the wonderful scipy library. For similar tasks in python there is also a pulp , but for beginners in scipy there is a clearer syntax.

Why might linear programming be needed in practice? As a rule, it solves the problem of minimizing the function f (x) (or the inverse maximization problem for - f (x)).

Here I will not give theoretical calculations (you can see here ), and consider a specific example.
')
So, the task.

We have 8 factories that produce a certain amount of products every week. We need to distribute products to 13 stores so as to maximize the total profit, while it is allowed to close unprofitable stores.

We know:

1. Factory productivity (field of supply in kg of products per week), their coordinates (X, Y)
   
   

   fact x y supply 0 F_01 59.250276 59.871183 389 1 F_02 84.739320 14.179336 409 2 F_03 42.397937 42.474530 124 3 F_04 19.539202 13.714643 70 4 F_05 41.280669 37.860993 386 5 F_06 37.159066 41.353602 196 6 F_07 96.890453 64.420010 394 7 F_08 86.267499 81.662811 365 

   
   
2. Store productivity (field of demand in kg products per week), their coordinates (X, Y)
   
   

     shop x y demand 0 S_01 13.490869 73.269974 200 1 S_02 85.435435 66.637250 20 2 S_03 28.578297 8.997380 320 3 S_04 31.324145 91.839907 360 4 S_05 40.338575 15.487028 360 5 S_06 41.642451 42.121572 120 6 S_07 53.983692 20.950457 360 7 S_08 75.761895 87.067552 60 8 S_09 81.836739 36.799647 80 9 S_10 54.260517 25.920108 100 10 S_11 67.918105 68.108601 340 11 S_12 92.200710 10.898110 360 12 S_13 19.966539 39.046271 60 

   
   
3. We also know that the cost of producing candy: 200 rubles / kg; proceeds from the sale of chocolates: 800 rubles / kg; delivery cost of sweets: 20 rubles / kg / km; weekly expenses for the store: 10 000 rub.
   

Now we impose the first logical constraint:

1. general supply <= general demand, that is, factories cannot ship more products to the store than they can sell.

We introduce the notation:

1. costkg = production costs;
2. revkg = production revenue (sales price per kg);
3. devkg = logistics costs;
4. costweek = store maintenance costs.

To work with the data, we will make a cross join of factories with stores in order to get all the combinations of supplies (factory - store).

We get a table of this type:

  fact supply shop demand distance cost_kg rev_kg del_kg cost_week 99 F_08 365 S_09 80 44.863164 200 800 20 10000 100 F_08 365 S_10 100 55.742703 200 800 20 10000 101 F_08 365 S_11 340 13.554211 200 800 20 10000 102 F_08 365 S_12 360 70.764701 200 800 20 10000 103 F_08 365 S_13 60 42.616540 200 800 20 10000 

Now let's take a look at the function that we will optimize (so far without considering the content of the stores).

$$$profit = volume * (revkg - costkg - devkg * distance)$

Here profit is the profit in rubles from shipments to a particular store, and volume is the volume of shipments of the factory to the store.

In total we will have 104 terms of profit (because we have 104 combinations of factories and stores).

The final form of the function will look like:

$$$Profit = sum (profit_1 ... profit_n) - 13 * 10000$

13 * 10000 is the cost of maintaining stores.

Now go to the most linear programming. Description of the module scipy.optimize you can find here . In general:

 res = linprog(c, A_ub=A, b_ub=b, bounds=(x0_bounds, x1_bounds), options={"disp": True}) 

Here:

1. c are the coefficients at our volume (i.e., costs). We will have them negative, because we solve the maximization problem (minimization-f (x));
2. A_ub is the value at volume for the condition imposed by us, b_ub is the value of the limit;
3. x0_bounds and x1_bounds lower and upper bounds on volume

The limitations are as follows:

1. The volume for each row in our table does not exceed the supply;
2. sum for factory volume in rows does not exceed supply;
3. the amount of volume for the store in rows does not exceed demand;
4. The volume takes values ​​from zero to max (supply) for a specific factory.

Code below:

 ## ff - table with factories and shops coefs = [] for f in ff.iterrows(): coefs.append((f[1]['rev_kg'] - f[1]['cost_kg'] - f[1]['del_kg']*f[1]['distance'])*(-1)) A = [] b = [] for f in ff['fact'].values: A.append((ff['fact'] == f)*1) b.append(ff[ff['fact'] ==f]['supply'].max()) for f in ff['shop'].values: A.append((ff['shop'] == f)*1) b.append(ff[ff['shop'] ==f]['demand'].max()) x0_bounds = [] x1_bounds = [] for f in ff.iterrows(): x0_bounds.append(0) x1_bounds.append(f[1]['demand']) x0_bounds = tuple(x0_bounds) x1_bounds = tuple(x1_bounds) A.append(coefs) b.append(-10000*13) res = linprog(coefs, A_ub=A, b_ub=b, bounds=list(zip(x0_bounds, x1_bounds)), options={"disp": True, 'maxiter' : 50000} ) 

Output:

Optimization terminated successfully.
Current function value: -948302.914122
Iterations: 20

We calculate the net profit and see how much each factory will supply and which stores should be closed:

 ff['supply_best'] = res.x ff['stay_opened'] = (ff['supply_best'] > 0)*1 ff['profit'] = (ff['supply_best']*(ff['rev_kg']- ff['cost_kg'] - ff['distance'] * ff['del_kg']))*ff['stay_opened'] net_profit = ff['profit'].sum() - ff[ff['stay_opened']==1]['shop'].nunique()*10000 grouped = ff.groupby(['fact', 'shop'])['supply_best'].sum().reset_index() f = {'supply_best': 'sum', 'supply': 'max'} ff.groupby('fact')['supply_best', 'supply'].agg(f) 

Reading profit - 828302.9141220199 rubles

Supplies of each factory:

  fact shop supply_best 0 F_01 S_01 166.0 17 F_02 S_05 360.0 24 F_02 S_12 49.0 31 F_03 S_06 120.0 38 F_03 S_13 4.0 50 F_04 S_12 70.0 58 F_05 S_07 346.0 61 F_05 S_10 40.0 73 F_06 S_09 80.0 74 F_06 S_10 60.0 77 F_06 S_13 56.0 78 F_07 S_01 34.0 79 F_07 S_02 20.0 88 F_07 S_11 340.0 94 F_08 S_04 305.0 98 F_08 S_08 60.0 

Deliveries to the store:

 supply_best demand shop S_01 200.0 200 S_02 20.0 20 S_03 0.0 320 S_04 305.0 360 S_05 360.0 360 S_06 120.0 120 S_07 346.0 360 S_08 60.0 60 S_09 80.0 80 S_10 100.0 100 S_11 340.0 340 S_12 119.0 360 S_13 60.0 60 

We see that S_03 should be closed, and in S_12 about 30% of demand will be supplied.

Thus, with the help of simple manipulations, we decided, albeit a basic and highly simplified, but often encountered in practice, task of optimizing the supply chain.

If you have questions and comments, I will be glad to answer.