Express logical constraints : one constraint implies another in linear programming (open source solver)

huangapple
go评论80阅读模式
英文:

Express logical constraints : one constraint implies another in linear programming (open source solver)

问题

How to express the following in linear programming with open-source solvers:

  1. 2*x + 3*y >= 15 ==> m + n <= 5

In open-source solvers, you can use the following syntax:

  1. 2*x + 3*y >= 15
  2. m + n <= 5

This represents the same constraint without the need for if-then notation, as open-source solvers typically allow you to directly specify linear constraints like this.

英文:

How can I express the following in linear programming:

  1. 2*x + 3*y &gt;= 15 ==&gt; m + n &lt;= 5

x, y, m and n are decision variables (continuous or integer) with relevant upper and lower bounds

I know with CPLEX (OPL and with python API), we can write if-then constraints like :

  1. 2*x + 3*y &gt;= 15 =&gt; m + n &lt;= 5 (OPL)
  2. model.add_if_then(2*x + 3*y &gt;= 15, m + n &lt;= 5) [python API, docplex]

But how to express above in an open source solver ?

答案1

得分: 1

使用指示约束，您可以将逻辑建模如下：

  1. 2*x + 3*y >= 15 ==> w == 1
  2. w == 1 ==> m + n <= 5

其中，w是一个二进制变量。指示约束通常（至少在Gurobi中）的形式是：如果二进制变量（这里是 w ）为真，则将强制满足线性约束。因此，在Gurobi中，您需要表达以下关系：

  1. 2*x + 3*y >= 15 ==> w == 1

以及其逻辑的反命题：

  1. w == 0 ==> 2*x + 3*y < 15

数学规划框架中不允许严格不等式，但我们可以使用某个小的 epsilon > 0 来近似这个约束。请注意，如果所有变量都是整数，我们可以使用 epsilon = 1。因此最终得到：

  1. w == 0 ==> 2*x + 3*y < 15 - epsilon
  2. w == 1 ==> m + n <= 5

在Gurobi中的实现如下：

  1. model.addGenConstrIndicator(w, 0, 2*x + 3*y <= 15 - epsilon)
  2. model.addGenConstrIndicator(w, 1, m + n <= 5)

对于开源求解器，您可以使用 big M 来表达以下内容：

  1. 2*x + 3*y >= 15 ==> w == 1
  3. 实现：
  4. a*x - (M + e) * w <= b - e
  5. a、x 表示系数和变量
  6. e（读作：epsilon）是一个小的容忍值，超过这个值，我们将视约束为被打破。
  7. 如果变量是整数，我们可以取 e = 1
  8. M（big_M）是 2*x + 3*y - 15 的上界
  9. b 是约束的右手边

以及：

  1. w == 1 ==> m + n <= 5
  3. 实现：
  4. a*x - b <= M * (1 - w)
  5. M（big_M）是 m + n - 5 的上界

在Google OR-Tools中的线性求解器（Python API）中的实现如下：

  1. # 请注意，代码中可能存在格式化问题，需要进一步验证。

这个实现假设了一些变量的上下界，并设置了目标函数以最大化。在解决问题后，您可以通过查询变量的解来获取结果。

英文:

Using indicator constraints, you can model the logic as follows:

  1. 2*x + 3*y &gt;= 15 ==&gt; w == 1
  2. w == 1 ==&gt; m + n &lt;= 5

w is a binary variable

indicator are generally (at least in Gurobi) of the form : if binary variable (here w) is true, then that would force the linear constraint to be satisfied. So in Gurobi, you would have to express the following :
2*x + 3*y &gt;= 15 ==&gt; w == 1 with its logical contra-positive i.e.
w == 0 ==&gt; 2*x + 3*y &lt; 15

Strict inequalities aren't allowed in a math programming framework, but we can closely approximate this constraint with 2*x + 3*y &lt;= 15 - epsilon for some small epsilon &gt; 0. Note that we can use epsilon = 1 if all variables are integers.

So finally:

  1. w == 0 ==&gt; 2*x + 3*y &lt;= 15 - epsilon
  2. w == 1 ==&gt; m + n &lt;= 5

In Gurobi :

  1. model.addGenConstrIndicator(w, 0, 2*x + 3*y &lt;= 15 - epsilon)
  2. model.addGenConstrIndicator(w, 1, m + n &lt;= 5)

I am not sure whether indicator constraints are supported in open source solvers.

For open source solvers, however, you can use big M to express following:

  1. 2*x + 3*y &gt;= 15 ==&gt; w == 1
  3. implementation:
  4. a*x - (M + e) * w &lt;= b - e
  5. a, x represent coefficients and variables
  6. e (read: epsilon) is some small tolerance value beyond which we will 
  7.    regard the constraint as having been broken. 
  8.    incase variables are integers, we can take e = 1
  9. M (big_M) is the upper bound on 2*x + 3*y - 15
  10. b is the RHS of the constraint
  1. w == 1 ==&gt; m + n &lt;= 5
  3. implementation:
  4. a*x - b &lt;= M * (1 - w)
  5. M (big_M) is the upper bound on : m + n - 5

Below is the implementation in Google OR-Tools, linear solver (python API)

  1. from ortools.linear_solver import pywraplp
  3. s = pywraplp.Solver(&quot;&quot;, pywraplp.Solver.SCIP_MIXED_INTEGER_PROGRAMMING)
  5. e = 10e-5 #epsilon
  7. x = s.NumVar(lb = 0, ub = 5, name = &quot;x&quot;)
  8. y = s.NumVar(lb = 0, ub = 5, name = &quot;y&quot;)
  10. m = s.NumVar(lb = 0, ub = 10, name = &quot;x&quot;)
  11. n = s.NumVar(lb = 0, ub = 10, name = &quot;x&quot;)
  13. # 2*x + 3*y &gt;= 15 ==&gt; m + n &lt;= 5
  15. # introduce w as a binary variable, then, enforce following
  16. # 2*x + 3*y &gt;= 15 ==&gt; w == 1
  17. # w == 1 ==&gt; m + n &lt;= 5
  19. # to test the above, lets fix x = 3 and y = 4
  20. # inequality : 2*x + 3*y &gt;= 15 now will be true,
  21. # and m + n &lt;= 5 should hold true post the solve
  22. x.SetLb(3)
  23. x.SetUb(3)
  25. y.SetLb(4)
  26. y.SetUb(4)
  28. # add w as a binary variable:
  29. w = s.BoolVar(&quot;&quot;)
  31. # 2*x + 3*y &gt;= 15 ==&gt; w == 1
  33. # a*x - (M + e) * w &lt;= b - e
  34. # a, x represent coefficients and variables
  35. # e (read: epsilon) is some small tolerance value beyond which we will 
  36. #   regard the constraint as having been broken 
  37. #   is variables are integers, we can take e = 1
  38. # M (big_M) is the upper bound on 2*x + 3*y - 15
  39. # b is the RHS of the constraint
  41. big_M_1 = 2*x.ub() + 3*y.ub() - 15
  42. s.Add(2*x + 3*y - (big_M_1 + e) * w &lt;= 15 - e)
  44. # now we express : w == 1 ==&gt; m + n &lt;= 5
  46. # a*x - b &lt;= M * (1 - w)
  47. # M (big_M) is the upper bound on : m + n - 5
  49. big_M_2 = m.ub() + n.ub() - 5
  50. s.Add(m + n - 5 &lt;= big_M_2 * (1 - w))
  52. s.Maximize(x + y + m + n)
  54. s.Solve()
  56. s.Objective().Value() # objective function value
  58. print(&quot;x == &quot; + str(x.SolutionValue())) # x = 3
  59. print(&quot;y == &quot; + str(y.SolutionValue())) # y = 4
  60. print(&quot;m == &quot; + str(m.SolutionValue())) # m = 5
  61. print(&quot;n == &quot; + str(n.SolutionValue())) # n = 0
  63. print(&quot;w == &quot; + str(w.SolutionValue())) # w = 1

as you see since we fixed x, y to be 3 and 4 respectively, 2x + 3y >= 15 is satisfied, which forces m + n <= 5 (m = 5 and n = 0 in optimal solution). This is despite of the fact that we are maximising x + y + m + n (without the constraint, m, n would take their upper bound values.

答案2

得分: 1

以下是您要翻译的代码部分：

  1. import pulp
  3. m = pulp.LpVariable('m', lowBound=0, upBound=40)
  4. n = pulp.LpVariable('n', lowBound=0, upBound=40)
  5. x = pulp.LpVariable('x', lowBound=2.4)
  6. y = pulp.LpVariable('y', lowBound=3.2)  # Increase lower bound to see constraint changes
  7. b = pulp.LpVariable('b', cat=pulp.LpBinary)
  8. prob = pulp.LpProblem(name='conditional_demo', sense=pulp.LpMinimize)
  9. prob.objective = x + y - m - n
  11. # If 2x + 3y < 15, b == 0.
  12. # If 2x + 3y >= 15, b == 1.
  13. # Epsilon must be used for the first expression, otherwise b is unconstrained at 2x + 3y == 15
  14. # p must be a coefficient sufficiently large that it accounts for the highest possible values of x, y
  15. epsilon = 1e-3
  16. p = 3
  17. prob.addConstraint(b >= ((2*x + 3*y)/15 - 1)/p + epsilon)
  18. prob.addConstraint(b <= (2*x + 3*y)/15)
  20. # If b == 0, then m + n is unconstrained
  21. # If b == 1, then m + n <= 5
  22. # M should be sufficiently large that it will always exceed realistic values of m + n - 5
  23. M = 100
  24. prob.addConstraint(m + n <= 5 + (1-b)*M)
  26. print(prob)
  27. prob.solve()
  28. assert prob.status == pulp.LpStatusOptimal
  29. print('m =', m.value())
  30. print('n =', n.value())
  31. print('x =', x.value())
  32. print('y =', y.value())
  33. print('b =', b.value(), 'because 2x+3y =', 2*x.value() + 3*y.value())
  1. m = 5.0
  2. n = 0.0
  3. x = 2.4
  4. y = 5.2
  5. b = 1.0 because 2x+3y = 20.400000000000002

请注意：由于对m和n的约束是单向的，b也可以是单向的；也就是说，您可以将其从二进制变量放宽为具有单边整数上限的变量。

英文:

The proper solution is strongly dependent on the possible ranges of x, y, a, b and their integrality. The following is a demonstration only, and you'll need to tune the parameters for whatever it is you're actually doing, in particular the values of M and p for the big-M constraints.

  1. import pulp
  3. m = pulp.LpVariable(&#39;m&#39;, lowBound=0, upBound=40)
  4. n = pulp.LpVariable(&#39;n&#39;, lowBound=0, upBound=40)
  5. x = pulp.LpVariable(&#39;x&#39;, lowBound=2.4)
  6. y = pulp.LpVariable(&#39;y&#39;, lowBound=3.2)  # Increase lower bound to see constraint changes
  7. b = pulp.LpVariable(&#39;b&#39;, cat=pulp.LpBinary)
  8. prob = pulp.LpProblem(name=&#39;conditional_demo&#39;, sense=pulp.LpMinimize)
  9. prob.objective = x + y - m - n
  11. # If 2x + 3y &lt; 15, b == 0.
  12. # If 2x + 3y &gt;= 15, b == 1.
  13. # Epsilon must be used for the first expression, otherwise b is unconstrained at 2x + 3y == 15
  14. # p must be a coefficient sufficiently large that it accounts for the highest possible values of x,y
  15. epsilon = 1e-3
  16. p = 3
  17. prob.addConstraint(b &gt;= ((2*x + 3*y)/15 - 1)/p + epsilon)
  18. prob.addConstraint(b &lt;= (2*x + 3*y)/15)
  20. # If b == 0, then m + n is unconstrained
  21. # If b == 1, then m + n &lt;= 5
  22. # M should be sufficiently large that it will always exceed realistic values of m + n - 5
  23. M = 100
  24. prob.addConstraint(m + n &lt;= 5 + (1-b)*M)
  26. print(prob)
  27. prob.solve()
  28. assert prob.status == pulp.LpStatusOptimal
  29. print(&#39;m =&#39;, m.value())
  30. print(&#39;n =&#39;, n.value())
  31. print(&#39;x =&#39;, x.value())
  32. print(&#39;y =&#39;, y.value())
  33. print(&#39;b =&#39;, b.value(), &#39;because 2x+3y =&#39;, 2*x.value() + 3*y.value())
  1. m = 5.0
  2. n = 0.0
  3. x = 2.4
  4. y = 5.2
  5. b = 1.0 because 2x+3y = 20.400000000000002

Note: since constraints on m,n are one-sided, b may also be one-sided; that is, you can loosen it from a binary variable to a singly-bounded integral variable:

  1. import pulp
  3. m = pulp.LpVariable(&#39;m&#39;, lowBound=0, upBound=40)
  4. n = pulp.LpVariable(&#39;n&#39;, lowBound=0, upBound=40)
  5. x = pulp.LpVariable(&#39;x&#39;, lowBound=6)
  6. y = pulp.LpVariable(&#39;y&#39;, lowBound=6)  # Change lower bound to see constraint changes
  7. b = pulp.LpVariable(&#39;b&#39;, cat=pulp.LpInteger, upBound=1)
  8. prob = pulp.LpProblem(name=&#39;conditional_demo&#39;, sense=pulp.LpMinimize)
  9. prob.objective = x + y - m - n
  11. # If 2x + 3y &lt; 15, b == 0.
  12. # If 2x + 3y &gt;= 15, b == 1.
  13. # Epsilon must be used for the first expression, otherwise b is unconstrained at 2x + 3y == 15
  14. # p must be a coefficient sufficiently large that it accounts for the highest possible values of x,y
  15. epsilon = 1e-3
  16. p = 3
  17. prob.addConstraint(b &gt;= ((2*x + 3*y)/15 - 1)/p + epsilon)
  18. prob.addConstraint(b &lt;= (2*x + 3*y)/15)
  20. # If b == 0, then m + n is unconstrained
  21. # If b == 1, then m + n &lt;= 5
  22. # M should be sufficiently large that it will always exceed realistic values of m + n - 5
  23. M = 100
  24. prob.addConstraint(m + n &lt;= 5 + (1-b)*M)
  26. print(prob)
  27. prob.solve()
  28. assert prob.status == pulp.LpStatusOptimal
  29. print(&#39;m =&#39;, m.value())
  30. print(&#39;n =&#39;, n.value())
  31. print(&#39;x =&#39;, x.value())
  32. print(&#39;y =&#39;, y.value())
  33. print(&#39;b =&#39;, b.value(), &#39;because 2x+3y =&#39;, 2*x.value() + 3*y.value())

huangapple
  • 本文由 发表于 2023年7月6日 14:37:55
  • 转载请务必保留本文链接：https://go.coder-hub.com/76626108.html
  • linear-programming
  • mixed-integer-programming
匿名

发表评论

匿名网友

:?: :razz: :sad: :evil: :!: :smile: :oops: :grin: :eek: :shock: :???: :cool: :lol: :mad: :twisted: :roll: :wink: :idea: :arrow: :neutral: :cry: :mrgreen:

确定