本文章為上篇建模學習打卡第二天的續

文章目錄

一、本次問題

?二、本題理解

三、問題求解

1.lingo實現

（1）先拋除整數約束條件對問題求解

?（2）加入整數約束條件求解

2.python實現求解

（1）先拋除整數約束條件對問題求解

?（2）加入整數約束條件求解實現? ?通過 pulp 庫求解

?（3）加入整數約束條件求解實現? 分枝界定法求解

---

一、本次問題

?二、本題理解

目標函數：

max = 40x1+90x2

一級約束條件：

9x1+7x2<=56 7x1+20x2<=70 x1,x2 >= 0

二級約束條件：

x1，x2全為整數

三、問題求解

1.lingo實現

lingo編寫代碼時，每行代碼結束后必須以 ‘ ; ’ 結束，否則無法運行。

（1）先拋除整數約束條件對問題求解

基礎線性規劃實現（matlab，lingo）_菜菜笨小孩的博客-CSDN博客

? ? ? ? lingo代碼實現：（l無其他條件下，ingo中默認變量大于等于0）

max = 40*x1+90*x2; 9*x1+7*x2<=56; 7*x1+20*x2<=70;

????????結果：最優解 x1=4.80916 , x2 = 1.816794 ; 最優值為355.8779；顯然不符合題意

?（2）加入整數約束條件求解

首先，需要引出lingo的變量界定函數 @gin(x) --- 將x限制為整數條件

? ? ? ?ingo代碼實現:通過變量界定函數將x1，x2限制為整數約束

max = 40*x1+90*x2; 9*x1+7*x2<=56; 7*x1+20*x2<=70; @gin(x1); @gin(x2);

????????結果：最優解 x1=4 , x2 = 2 ; 最優值為340；符合題意

?lingo實現求解到此結束。

2.python實現求解

（1）先拋除整數約束條件對問題求解

基礎線性規劃實現---python_菜菜笨小孩的博客-CSDN博客

????????python代碼實現如下：詳解請看上方python基礎線性規劃的文章

#導入包 from scipy import optimize as opt import numpy as np#確定c,A,b,Aeq,beq c = np.array([40,90]) #目標函數變量系數 A = np.array([[9,7],[7,20]]) #不等式變量系數 b = np.array([56,70]) #不等式變量值 Aeq = np.array([[0,0]]) #等式變量系數 beq = np.array([0]) #等式變量值 #限制 lim1=(0,None) lim2=(0,None) #求解 res = opt.linprog(-c,A,b,Aeq,beq,bounds=(lim1,lim2)) #輸出結果 print(res)

????????結果：最優解 x1=4.80916 , x2 = 1.816794 ; 最優值為355.8779；顯然不符合題意

?（2）加入整數約束條件求解實現 ? 通過 pulp 庫求解

安裝庫：我用python3.8安裝成功了

pip install pulp

python代碼實現：

1.導入庫

import pulp as pulp

2.定義解決問題的函數

def solve_ilp(objective , constraints) :print(objective)print(constraints)prob = pulp.LpProblem('LP1' , pulp.LpMaximize)prob += objectivefor cons in constraints :prob += consprint(prob)status = prob.solve()if status != 1 :return Noneelse :return [v.varValue.real for v in prob.variables()]

3.設置目標函數和約束條件等

V_NUM = 2 #本題變量個數 # 變量，直接設置下限 variables = [pulp.LpVariable('X%d' % i, lowBound=0, cat=pulp.LpInteger) for i in range(0, V_NUM)] # 目標函數 c = [40, 90] objective = sum([c[i] * variables[i] for i in range(0, V_NUM)]) # 約束條件 constraints = [] a1 = [9, 7] constraints.append(sum([a1[i] * variables[i] for i in range(0, V_NUM)]) <= 56) a2 = [7, 20] constraints.append(sum([a2[i] * variables[i] for i in range(0, V_NUM)]) <= 70)

4.求解：

res = solve_ilp(objective, constraints) print(res) #輸出結果

完整代碼如下：

# -*- coding: utf-8 -*- import pulp as pulpdef solve_ilp(objective, constraints):print(objective)print(constraints)prob = pulp.LpProblem('LP1', pulp.LpMaximize)prob += objectivefor cons in constraints:prob += consprint(prob)status = prob.solve()if status != 1:# print 'status'# print statusreturn Noneelse:# return [v.varValue.real for v in prob.variables()]return [v.varValue.real for v in prob.variables()]V_NUM = 2 # 變量，直接設置下限 variables = [pulp.LpVariable('X%d' % i, lowBound=0, cat=pulp.LpInteger) for i in range(0, V_NUM)] # 目標函數 c = [40, 90] objective = sum([c[i] * variables[i] for i in range(0, V_NUM)]) # 約束條件 constraints = [] a1 = [9, 7] constraints.append(sum([a1[i] * variables[i] for i in range(0, V_NUM)]) <= 56) a2 = [7, 20] constraints.append(sum([a2[i] * variables[i] for i in range(0, V_NUM)]) <= 70) # print (constraints)res = solve_ilp(objective, constraints) print(res) #輸出解

結果：最優解 x1=4 , x2 = 2 ; 最優值為340；符合題意

?（3）加入整數約束條件求解實現? 分枝界定法求解

何為分枝界定法，請看詳解https://blog.csdn.net/qq_25990967/article/details/121211474

python代碼實現：

1.導入庫

from scipy.optimize import linprog import numpy as np import math import sys from queue import Queue

2.定義整數線性規劃類

class ILP()

3.定義分枝界定法函數

def __init__(self, c, A_ub, b_ub, A_eq, b_eq, bounds):# 全局參數self.LOWER_BOUND = -sys.maxsizeself.UPPER_BOUND = sys.maxsizeself.opt_val = Noneself.opt_x = Noneself.Q = Queue()# 這些參數在每輪計算中都不會改變self.c = -cself.A_eq = A_eqself.b_eq = b_eqself.bounds = bounds# 首先計算一下初始問題r = linprog(-c, A_ub, b_ub, A_eq, b_eq, bounds)# 若最初問題線性不可解if not r.success:raise ValueError('Not a feasible problem!')# 將解和約束參數放入隊列self.Q.put((r, A_ub, b_ub))def solve(self):while not self.Q.empty():# 取出當前問題res, A_ub, b_ub = self.Q.get(block=False)# 當前最優值小于總下界，則排除此區域if -res.fun < self.LOWER_BOUND:continue# 若結果 x 中全為整數，則嘗試更新全局下界、全局最優值和最優解if all(list(map(lambda f: f.is_integer(), res.x))):if self.LOWER_BOUND < -res.fun:self.LOWER_BOUND = -res.funif self.opt_val is None or self.opt_val < -res.fun:self.opt_val = -res.funself.opt_x = res.xcontinue# 進行分枝else:# 尋找 x 中第一個不是整數的，取其下標 idxidx = 0for i, x in enumerate(res.x):if not x.is_integer():breakidx += 1# 構建新的約束條件（分割new_con1 = np.zeros(A_ub.shape[1])new_con1[idx] = -1new_con2 = np.zeros(A_ub.shape[1])new_con2[idx] = 1new_A_ub1 = np.insert(A_ub, A_ub.shape[0], new_con1, axis=0)new_A_ub2 = np.insert(A_ub, A_ub.shape[0], new_con2, axis=0)new_b_ub1 = np.insert(b_ub, b_ub.shape[0], -math.ceil(res.x[idx]), axis=0)new_b_ub2 = np.insert(b_ub, b_ub.shape[0], math.floor(res.x[idx]), axis=0)# 將新約束條件加入隊列，先加最優值大的那一支r1 = linprog(self.c, new_A_ub1, new_b_ub1, self.A_eq,self.b_eq, self.bounds)r2 = linprog(self.c, new_A_ub2, new_b_ub2, self.A_eq,self.b_eq, self.bounds)if not r1.success and r2.success:self.Q.put((r2, new_A_ub2, new_b_ub2))elif not r2.success and r1.success:self.Q.put((r1, new_A_ub1, new_b_ub1))elif r1.success and r2.success:if -r1.fun > -r2.fun:self.Q.put((r1, new_A_ub1, new_b_ub1))self.Q.put((r2, new_A_ub2, new_b_ub2))else:self.Q.put((r2, new_A_ub2, new_b_ub2))self.Q.put((r1, new_A_ub1, new_b_ub1))

4.定義求解問題中的變量級約束條件

def test():""" 此測試的真實最優解為 [4, 2] """c = np.array([40, 90])A = np.array([[9, 7], [7, 20]])b = np.array([56, 70])Aeq = Nonebeq = Nonebounds = [(0, None), (0, None)]solver = ILP(c, A, b, Aeq, beq, bounds)solver.solve()print("Test 's result:", solver.opt_val, solver.opt_x)print("Test 's true optimal x: [4, 2]\n")

5.求解并輸出

if __name__ == '__main__':test()

完整代碼如下：

from scipy.optimize import linprog import numpy as np import math import sys from queue import Queueclass ILP():def __init__(self, c, A_ub, b_ub, A_eq, b_eq, bounds):# 全局參數self.LOWER_BOUND = -sys.maxsizeself.UPPER_BOUND = sys.maxsizeself.opt_val = Noneself.opt_x = Noneself.Q = Queue()# 這些參數在每輪計算中都不會改變self.c = -cself.A_eq = A_eqself.b_eq = b_eqself.bounds = bounds# 首先計算一下初始問題r = linprog(-c, A_ub, b_ub, A_eq, b_eq, bounds)# 若最初問題線性不可解if not r.success:raise ValueError('Not a feasible problem!')# 將解和約束參數放入隊列self.Q.put((r, A_ub, b_ub))def solve(self):while not self.Q.empty():# 取出當前問題res, A_ub, b_ub = self.Q.get(block=False)# 當前最優值小于總下界，則排除此區域if -res.fun < self.LOWER_BOUND:continue# 若結果 x 中全為整數，則嘗試更新全局下界、全局最優值和最優解if all(list(map(lambda f: f.is_integer(), res.x))):if self.LOWER_BOUND < -res.fun:self.LOWER_BOUND = -res.funif self.opt_val is None or self.opt_val < -res.fun:self.opt_val = -res.funself.opt_x = res.xcontinue# 進行分枝else:# 尋找 x 中第一個不是整數的，取其下標 idxidx = 0for i, x in enumerate(res.x):if not x.is_integer():breakidx += 1# 構建新的約束條件（分割new_con1 = np.zeros(A_ub.shape[1])new_con1[idx] = -1new_con2 = np.zeros(A_ub.shape[1])new_con2[idx] = 1new_A_ub1 = np.insert(A_ub, A_ub.shape[0], new_con1, axis=0)new_A_ub2 = np.insert(A_ub, A_ub.shape[0], new_con2, axis=0)new_b_ub1 = np.insert(b_ub, b_ub.shape[0], -math.ceil(res.x[idx]), axis=0)new_b_ub2 = np.insert(b_ub, b_ub.shape[0], math.floor(res.x[idx]), axis=0)# 將新約束條件加入隊列，先加最優值大的那一支r1 = linprog(self.c, new_A_ub1, new_b_ub1, self.A_eq,self.b_eq, self.bounds)r2 = linprog(self.c, new_A_ub2, new_b_ub2, self.A_eq,self.b_eq, self.bounds)if not r1.success and r2.success:self.Q.put((r2, new_A_ub2, new_b_ub2))elif not r2.success and r1.success:self.Q.put((r1, new_A_ub1, new_b_ub1))elif r1.success and r2.success:if -r1.fun > -r2.fun:self.Q.put((r1, new_A_ub1, new_b_ub1))self.Q.put((r2, new_A_ub2, new_b_ub2))else:self.Q.put((r2, new_A_ub2, new_b_ub2))self.Q.put((r1, new_A_ub1, new_b_ub1))def test():""" 此測試的真實最優解為 [4, 2] """c = np.array([40, 90])A = np.array([[9, 7], [7, 20]])b = np.array([56, 70])Aeq = Nonebeq = Nonebounds = [(0, None), (0, None)]solver = ILP(c, A, b, Aeq, beq, bounds)solver.solve()print("Test 's result:", solver.opt_val, solver.opt_x)print("Test 's true optimal x: [4, 2]\n")if __name__ == '__main__':test()

結果：最優解 x1=4 , x2 = 2 ; 最優值為340；符合題意

總結

以上是生活随笔為你收集整理的整数线性规划实现（lingo，python分枝界定法）的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。