Solve a system of linear equations and linear inequalities

Question

I have to get the min and max y for a linear expression, restricted by some linear inequalities in python.

You can see the equation and inequalities here that I have entered into Desmos:

3x+12y = 1000
x > 30
x < 160
y < 60
y > 10
x + y > 180

I can solve them by hand by drawing and crossing out the inequalities. But I cannot do that in Python. What I have tried so far in Python is to get y=83.33 when x=0; x=333.33 when y=0; After getting the min and max x,y I then apply the inequalities 1 by 1. But with every inequality I have to add the previous ones, and also do check if x or y has gone over certain range and so far and it is almost certain I will miss a check.

I looked at numpy and sympy, but could not figure out how to solve this using them. Can you suggest what / how can I use in order to get the range which the white arrow shows on the picture?

Dunderhead

This is a problem in linear programming, where your equality and inequalities are the limitations and you want to minimize and maximize the expression y. Do a search on linear programming in Python, and if you have a specific problem implementing it let us know. You could start with scipy.optimize.linprog. Did you look into that when you looked at scipy? – Allanadale

Answer 1

Yours is a problem in linear programming, where your equality and inequalities are the limitations and you want to minimize (then later maximize) the expression y. The equality, inequalities, and expression are all linear, so that makes it linear programming. The scipy package, using the scipy.optimize.linprog function, can do this kind of linear programming.

Here is commented code to do what you want. Note that all the inequalities were slightly changed to include equality, which is necessary to have a maximum or minimum value of y. To find the maximum value of y the code instead finds the minimum value of -y then prints the additive inverse of that, since linprog minimizes the objective function. Finally, the inequality restrictions must be "less than or equal to" in linprog, so I multiplied both sides of your inequality x + y > 180 by -1 to get one, namely -x + -y <= -180. Ask if you have any questions.

from scipy.optimize import linprog

# Set up values relating to both minimum and maximum values of y
coefficients_inequalities = [[-1, -1]]  # require -1*x + -1*y <= -180
constants_inequalities = [-180]
coefficients_equalities = [[3, 12]]  # require 3*x + 12*y = 1000
constants_equalities = [1000]
bounds_x = (30, 160)  # require 30 <= x <= 160
bounds_y = (10, 60)  # require 10 <= y <= 60

# Find and print the minimal value of y
coefficients_min_y = [0, 1]  # minimize 0*x + 1*y
res = linprog(coefficients_min_y,
              A_ub=coefficients_inequalities,
              b_ub=constants_inequalities,
              A_eq=coefficients_equalities,
              b_eq=constants_equalities,
              bounds=(bounds_x, bounds_y))
print('Minimum value of y =', res.fun)

# Find and print the maximal value of y = minimal value of -y
coefficients_max_y = [0, -1]  # minimize 0*x + -1*y
res = linprog(coefficients_max_y,
              A_ub=coefficients_inequalities,
              b_ub=constants_inequalities,
              A_eq=coefficients_equalities,
              b_eq=constants_equalities,
              bounds=(bounds_x, bounds_y))
print('Maximum value of y =', -res.fun)  # opposite of value of -y

The printout from that code is

Minimum value of y = 43.3333333333
Maximum value of y = 51.1111111111

which is correct within the precision of floating point. If you need the corresponding values of x, see the value of res.x which is an array that gives the values of both x and y at the desired point--x is res.x[0] and y is res.x[1].

Allanadale

Thanks alot I got it to work with the actual problem i had (it had more inequalities). But unfortunately I am looping this code about 1500 times, which is done in about 8s. Is there a way to just check if there is root for x in y in the equation given the inequalities (i think it would be faster without minimizing and maximizing it). I am marking your answer as the solution btw. – Dunderhead

I love COPY & RUN 100% working examples. And this is one of them. Great work. – Defilade

Answer 2

Following up on the hint provided in the comment, you can solve your problem using scipy.optimize.linprog as follows...

In [1]: import numpy as np

In [2]: from scipy import optimize

In [3]: c = np.zeros(2)

In [4]: A_ub = np.array([[-1, -1]])

In [5]: b_ub = np.array([-180])

In [6]: A_eq = np.array([[3, 12]])

In [7]: b_eq = np.array([1000])

In [8]: bounds = [(30, 160), (10, 60)]

In [9]: optimize.linprog(c, A_ub, b_ub, A_eq, b_eq, bounds, method="simplex")
Out[11]: 
     fun: -0.0
 message: 'Optimization terminated successfully.'
     nit: 5
   slack: array([ 31.11111111,  98.88888889,   0.        ,  41.11111111,   8.88888889])
  status: 0
 success: True
       x: array([ 128.88888889,   51.11111111])

The trick is to note that solving a system of equations can be expressed as an optimization problem with a constant zero objective function. This is why I set c=np.zeros(2) above.

Explication

It seems the OP wanted to "get the min and max y", not just one feasible point. – Aaronson

Answer 3

You could try cvxpy. It's a solver for convex problems, so it may be overkill for your needs. I like how easy it is to encode your problem with it. Most other libraries require you to encode your problem into several matrices.

Chili

(1) cvxpy has no own solver at all, except unconstrained LS. All the others are external solvers. (2) cvxpy is a modelling-tool but can't express all convex-problems. It only can build problem compliant with Disciplined Convex Programming (DCP < Convex). (3) I love cvxpy, but i don't see the point to use it here. Approaching LPs with general conic-programming will lose performance & robustness (i ignored the alpha-version 1.0 with it's reductions; but even with those: there is not a concept of variable-bounds, which are helping in the real-world). – Filigree

I agree with you on your first point, I should have been clearer there. I was able to solve the problem faster with cvxpy both writing it and running the solver than the scipy solution. – Chili