Python之建模规划篇--整数规划

  • 基本介绍
  • 整数规划的分类
  • 整数规划的特点
  • 求解方法分类
    • 0 - 1 型整数规划
    • 蒙特卡洛法 (随机取样法)
    • 整数线性规划的计算机求解
    • 分枝定界法
    • Python 实现 (分支定界代码)

基本介绍

规划中的变量(部分或全部)限制为整数时,称为整数规划。若在线性规划模型中,变量限制为整数,则称为整数线性规划。目前所流行的求解整数规划的方法,往往只适用于整数线性规划。目前还没有一种方法能有效地求解一切整数规划。

整数规划的分类

如不加特殊说明,一般指整数线性规划。对于整数线性规划模型大致可分为两类:

  • 变量全限制为整数时,称纯(完全)整数规划。
  • 变量部分限制为整数的,称混合整数规划。

整数规划的特点

  • 原线性规划有最优解,当自变量限制为整数后,其整数规划解出现下述情况:
    ①原线性规划最优解全是整数,则整数规划最优解与线性规划最优解一致。
    ②整数规划无可行解
  • 整数规划最优解不能按照实数最优解简单取整而获得。

求解方法分类

  • 分枝定界法—可求纯或混合整数线性规划。
  • 割平面法—可求纯或混合整数线性规划。
  • 隐枚举法—求解“0-1”整数规划:
    ①过滤隐枚举法;
    ②分枝隐枚举法。
  • 匈牙利法—解决指派问题(“0-1”规划特殊情形)。
  • 蒙特卡洛法—求解各种类型规划。

0 - 1 型整数规划

0 −1型整数规划是整数规划中的特殊情形,它的变量 xj 仅取值0或1。这时xj 称为0−1变量,或称二进制变量。xj 仅取值0 或1 这个条件可由下述约束条件:

所代替,是和一般整数规划的约束条件形式一致的。在实际问题中,如果引入 0 −1变量,就可以把有各种情况需要分别讨论的线性规划问题统一在一个问题中讨论了。我们先介绍引入0 −1变量的实际问题,再研究解法。
比如有一些相互排斥的约束条件,就是一种0-1问题,如运输方式只能选择一种,用车或者用船等类似的
除此之外,还有关于固定费用的问题,在讨论线性规划时,有些问题是要求使成本为最小。那时总设固定成本为常数,并在线性规划的模型中不必明显列出。但有些固定费用(固定成本)的问题不能用一般线性规划来描述,但可改变为混合整数规划来解决

蒙特卡洛法 (随机取样法)

蒙特卡洛方法也称为计算机随机模拟方法,它源于世界著名的赌城一摩纳哥的Monte Carlo(蒙特卡洛)。它是基于对大量事件的统计结果来实现–些确定性问题的计算。使用蒙特卡洛方法必须使用计算机生成相关分布的随机数,Matlab和python等各种编程语言都给出了生成各种随机数的命令。
如,给个例子

前面介绍的常用的整数规划求解方法,主要是针对线性整数规划而言,而对于非线性整数规划目前尚未有一种成熟而准确的求解方法,因为非线性规划本身的通用有效解法尚未找到,更何况是非线性整数规划。
然而,尽管整数规划由于限制变量为整数而增加了难度;然而又由于整数解是有限个,于是为枚举法提供了方便。当然,当自变量维数很大和取值范围很宽情况下,企图用显枚举法(即穷举法)计算出最优值是不现实的,但是应用概率理论可以证明,在一定的计算量的情况下,完全可以得出一个满意解。


如果用显枚举法试探,共需计算(100)5 = 1010个点,其计算量非常之大。然而应用蒙特卡洛去随机计算106个点,便可找到满意解,那么这种方法的可信度究竟怎样
呢?
下面就分析随机取样采集106个点计算时,应用概率理论来估计一下可信度。
不失一般性,假定一个整数规划的最优点不是孤立的奇点。
假设目标函数落在高值区的概率分别为 0.01,0.00001,则当计算106个点后,有
任一个点能落在高值区的概率分别为

首先编写M 文件mente.m 定义目标函数f 和约束向量函数g,程序如下:

function [f,g]=mengte(x);
f=x(1)^2+x(2)^2+3*x(3)^2+4*x(4)^2+2*x(5)-8*x(1)-2*x(2)-3*x(3)-...
x(4)-2*x(5);
g=[sum(x)-400
x(1)+2*x(2)+2*x(3)+x(4)+6*x(5)-800
2*x(1)+x(2)+6*x(3)-200
x(3)+x(4)+5*x(5)-200]

编写M文件mainint.m如下求问题的解

rand('state',sum(clock));  %初始化随机数发生器
p0=0;
tic    %计时开始
for i=1:10^6x=randi([0,99],1,5); %产生一行五列的区间[0,99]上的随机整数[f,g]=mengte(x);if all(g<=0)if p0<fx0=x; p0=f; %记录下当前较好的解endend
end
x0,p0
toc    %计时结束

由于是随机模拟,所以最后的答案都是不同的

整数线性规划的计算机求解

整数规划问题的求解使用Lingo等专用软件比较方便。对于整数线性规划问题,也可以使用Matlab的intlinprog函数求解,但使用Matlab软件求解数学规划问题有–个缺陷,即必须把所有的决策变量化成一-维决策向量,实际上对于多维变量的数学规划问题,用Matlab软件求解,需要做–个变量替换,把多维决策变量化成–维决策向量,变量替换后,约束条件很难写出;而使用Lingo软件求解数学规划问题是不需要做变换的,使用起来相对比较容易。


这里用一个指派问题作为例子进行解决

若用Matlab解决,代码如下

clc, clear
c=[3 8 2 10 3;8 7 2 9 7;6 4 2 7 58 4 2 3 5;9 10 6 9 10];
c=c(:); a=zeros(10,25); intcon=1:25;
for i=1:5a(i,(i-1)*5+1:5*i)=1;a(5+i,i:5:25)=1;
end
b=ones(10,1); lb=zeros(25,1); ub=ones(25,1);
x=intlinprog(c,intcon,[],[],a,b,lb,ub);
x=reshape(x,[5,5])

这样可以得到最优的指派结果

分枝定界法

对有约束条件的最优化问题(其可行解为有限数)的所有可行解空间恰当地进行系统搜索,这就是分枝与定界内容。通常,把全部可行解空间反复地分割为越来越小的子集,称为分枝;并且对每个子集内的解集计算一个目标下界(对于最小值问题),这称为定界。在每次分枝后,凡是界限超出已知可行解集目标值的那些子集不再进一步分枝,这样,许多子集可不予考虑,这称剪枝。这就是分枝定界法的主要思路。

分枝定界法可用于解纯整数或混合的整数规划问题。在本世纪六十年代初由 Land Doig 和Dakin 等人提出的。由于这方法灵活且便于用计算机求解,所以现在它已是解整数规划的重要方法。目前已成功地应用于求解生产进度问题、旅行推销员问题、工厂选址问题、背包问题及分配问题等。
设有最大化的整数规划问题 A ,与它相应的线性规划为问题B ,从解问题B 开始,若其最优解不符合 A的整数条件,那么B的最优目标函数必是 A的最优目标函数z的上界,记作z1 ;而 A的任意可行解的目标函数值将是z的一个下界z2 。分枝定界法就是将B的可行域分成子区域的方法。逐步减小z 1和增大z2 ,最终求到z*。现用下例来说明:



从以上解题过程可得用分枝定界法求解整数规划(最大化)问题的步骤为:
开始,将要求解的整数规划问题称为问题 A ,将与它相应的线性规划问题称为问题B 。
(i)解问题B 可能得到以下情况之一:
(a) B 没有可行解,这时A 也没有可行解,则停止.
(b) B 有最优解,并符合问题A 的整数条件, B 的最优解即为A 的最优解,则停止。
(c) B 有最优解,但不符合问题A 的整数条件,记它的目标函数值为z1 。
(ii)用观察法找问题 A的一个整数可行解,一般可取x j , j = 0,1.2,…,n,试探,求得其目标函数值,并记作z2。以z*表示问题 A的最优目标函数值;这时有

z2 ≤ z* ≤ z1 进行迭代。 第一步:分枝,在 B 的最优解中任选一个不符合整数条件的变量x j ,其值为b j,以b[ j ] 表示小于b j的最大整数。构造两个约束条件 x j ≤ [b j] 和 x j ≥ [b j] + 1 将这两个约束条件,分别加入问题B ,求两个后继规划问题B 1 和B 2。不考虑整数条件求解这两个后继问题。
定界,以每个后继问题为一分枝标明求解的结果,与其它问题的解的结果中,找出 最优目标函数值最大者作为新的上界z 1。从已符合整数条件的各分支中,找出目标函数 值为最大者作为新的下界z2,若无作用z 不变。
第二步:比较与剪枝,各分枝的最优目标函数中若有小于z2 者,则剪掉这枝,即 以后不再考虑了。若大于z2 ,且不符合整数条件,则重复第一步骤。一直到最后得到 z* = z2 为止。得最优整数解x j* , j = 1,2,.....,n

Python 实现 (分支定界代码)

  • 整数规划的模型与线性规划基本相同,只是额外增加了部分变量为整数的约束
  • 整数规划求解的基本框架是分支定界法,首先去除整数约束得到“松弛模型”,使用线性规划的方法求解。
  • 若有某个变量不是整数,在松弛模型.上分别添加约束: x≤floor(A)和x≥ceil(A),然后再分别求解,这个过程叫做分支。当节点求解结果中所有变量都是整数时,停止分支。这样不断迭代,形成了一颗树。
  • 所谓定界,指的是叶子节点产生后,相当于给问题定了一个下界。之后在求解过程中一旦某个节点的目标函数值小于这个下界,那就直接pass,不再进行分支了;每次新产生叶子节点,则更新下界。
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 test1():""" 此测试的真实最优解为 [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 1's result:", solver.opt_val, solver.opt_x)print("Test 1's true optimal x: [4, 2]\n")def test2():""" 此测试的真实最优解为 [2, 4] """c = np.array([3, 13])A = np.array([[2, 9], [11, -8]])b = np.array([40, 82])Aeq = Nonebeq = Nonebounds = [(0, None), (0, None)]solver = ILP(c, A, b, Aeq, beq, bounds)solver.solve()print("Test 2's result:", solver.opt_val, solver.opt_x)print("Test 2's true optimal x: [2, 4]\n")if __name__ == '__main__':test1()test2()

可以得到结果,就可以求出整数规划了(个人感觉整数规划没有matlab那么好,matlab有直接的函数)

每日一句
Never underestimate your power to change yourself!(永远不要低估你改变自我的能力! )

Python之建模规划篇--整数规划相关推荐

  1. Python之建模规划篇--线性规划

    Python之建模规划篇--线性规划 基本介绍 线性规划的实例与定义 线性规划问题的解的概念 求解线性规划的Matlab 解法 Python解法 Python Scipy库实现 Python plup ...

  2. python非线性规划求解_Python之建模规划篇--非线性规划

    Python之建模规划篇--非线性规划 基本介绍 如果目标函数或约束条件中包含非线性函数,就称这种规划问题为非线性规划问 题.一般说来,解非线性规划要比解线性规划问题困难得多.而且,也不象线性规划有 ...

  3. Python之建模数值逼近篇–最小二乘拟合

    Python之建模数值逼近篇–最小二乘拟合 介绍 系数ak的确定 函数rk(x)r_k(x)rk​(x)的选取 理解和区别 样例 介绍 曲线拟合问题的提法是,已知一组(二维)数据,即平面上的n个点(x ...

  4. Python之建模数值逼近篇--一维插值

    Python之建模数值逼近篇--一维插值 基本介绍 拉格朗日插值 分段插值 样条插值 概念 二次样条函数 三次样条函数 线性插值与样条插值 样例 1 高阶样条插值 样例2: 基本介绍 首先介绍一下插值 ...

  5. python数学建模

    1.python之建模规划篇 1.1 scipy库 例如:from scipy import 想导入的模块名 scipy具体模块如下: scipy.cluster向量量化 scipy.constant ...

  6. Python小白的数学建模课-04.整数规划

    整数规划与线性规划的差别只是变量的整数约束. 问题区别一点点,难度相差千万里. 选择简单通用的编程方案,让求解器去处理吧. 『Python小白的数学建模课 @ Youcans』带你从数模小白成为国赛达 ...

  7. 数学建模第四天:数学建模算法篇之整数规划、指派问题及其求解方法

    目录 一.前言 二.整数规划模型 1.整数规划特征 2.分枝定界法 ①分枝定界法的步骤 ②实际解题 三.0-1整数规划 1.隐枚举法 ①隐枚举法的步骤: ②案例 2.匈牙利法 ①指派问题 ②匈牙利法步 ...

  8. 数学建模——规划问题

    文章目录 线性规划 整数规划 一般的整数线性规划问题 0-1整数规划 广义指派问题 非线性规划 二次规划 线性规划  运筹学对于线性规划问题直接使用图解法,单纯形法利用求解.在python中可以直接使 ...

  9. Python数学建模系列(五):微分方程

    文章目录 前言 往期文章 1.微分方程分类 2.微分方程解析解 3.微分方程数值解 3.1 场线图与数值解 3.2 洛伦兹曲线与数值解 4.传染病模型 模型一:SI-Model 模型二:SIS mod ...

最新文章

  1. ctimespan 获取毫秒_VC++中通过CTime类获取日期差
  2. ML之分类预测之ElasticNet之PLoR:在二分类数据集上调用Glmnet库训练PLoR模型(T2)
  3. turnitin时间
  4. rpg人物制作软件_新机制和随机性的完美结合!新RPG《元素梦境》参上
  5. Java|达梦工作笔记-达梦数据库同步工具(JDBC)
  6. MS SQL数据库备份和恢复存储过程
  7. 【比赛分享】互联网新闻情感分析复赛top8(8/2745)解决方案及总结
  8. (六)6.6 Neurons Networks PCA
  9. 面试题:原型Bean在一个线程多次获取是否一样?
  10. java json 变量所有的属性
  11. 英特尔nuc能代替主机吗_拆了拆了!Intel NUC装机!小机箱退烧器啊!主机显示器合体...
  12. VS2010 快捷键
  13. 解决ThinkPad S3-440无法睡眠问题
  14. Python中对列表中的字典元素进行排序
  15. APP Launch 优化
  16. 如何一键关闭win安全中心(Windows Defender )
  17. 什么是CDN,简单了解CDN
  18. 本地调试获取微信code网页授权,免部署(前端+开发者工具)
  19. 医疗行业容灾备份解决方案
  20. 原生js源码之Array数组的every方法

热门文章

  1. VC 开机自动启动程序代码
  2. TOF/结构光camera区别、TOF同时成像深度图、IR图原理?
  3. 怎样才能设计一个“易用性”好的网站?
  4. WPF布局控件与子控件的HorizontalAlignment/VerticalAlignment属性之间的关系
  5. 团队成员的角色类型及其给我们的启示
  6. 试卷: 【2022】小米秋招笔试-软件开发-卷2
  7. Wpf依赖属性和附加属性在样式中的应用
  8. php研究所 百科_松下幸之助_PHP研究所
  9. 北京市重点区域5G网络实测分析
  10. Docker 制作带有中文字体的镜像