本文章为上篇建模学习打卡第二天的续
文章目录
一、本次问题
二、本题理解
三、问题求解
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 += objective for cons in constraints : prob += cons print(prob) status = prob.solve() if status != 1 : return None else : 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 pulp def solve_ilp(objective, constraints): print(objective) print(constraints) prob = pulp.LpProblem('LP1', pulp.LpMaximize) prob += objective for cons in constraints: prob += cons print(prob) status = prob.solve() if status != 1: # print 'status' # print status return None else: # 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.maxsize self.UPPER_BOUND = sys.maxsize self.opt_val = None self.opt_x = None self.Q = Queue() # 这些参数在每轮计算中都不会改变 self.c = -c self.A_eq = A_eq self.b_eq = b_eq self.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.fun if self.opt_val is None or self.opt_val < -res.fun: self.opt_val = -res.fun self.opt_x = res.x continue # 进行分枝 else: # 寻找 x 中第一个不是整数的,取其下标 idx idx = 0 for i, x in enumerate(res.x): if not x.is_integer(): break idx += 1 # 构建新的约束条件(分割 new_con1 = np.zeros(A_ub.shape[1]) new_con1[idx] = -1 new_con2 = np.zeros(A_ub.shape[1]) new_con2[idx] = 1 new_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 = None beq = None bounds = [(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 Queue class ILP(): def __init__(self, c, A_ub, b_ub, A_eq, b_eq, bounds): # 全局参数 self.LOWER_BOUND = -sys.maxsize self.UPPER_BOUND = sys.maxsize self.opt_val = None self.opt_x = None self.Q = Queue() # 这些参数在每轮计算中都不会改变 self.c = -c self.A_eq = A_eq self.b_eq = b_eq self.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.fun if self.opt_val is None or self.opt_val < -res.fun: self.opt_val = -res.fun self.opt_x = res.x continue # 进行分枝 else: # 寻找 x 中第一个不是整数的,取其下标 idx idx = 0 for i, x in enumerate(res.x): if not x.is_integer(): break idx += 1 # 构建新的约束条件(分割 new_con1 = np.zeros(A_ub.shape[1]) new_con1[idx] = -1 new_con2 = np.zeros(A_ub.shape[1]) new_con2[idx] = 1 new_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 = None beq = None bounds = [(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;符合题意
欢迎分享,转载请注明来源:内存溢出
评论列表(0条)