数学建模代码实现

1、线性规划

下面是代码实现

 导入包并把约束转化成标准格式

from scipy import optimize
import numpy as np
 
 
c = [2, 3, -5]
A = [[-2, 5, -1], [1, 3, 1]]
b = [-10, 12]
Aeq = [[1, 1, 1]]
beq = [7]
x1 = (0, None)
x2 = (0, None)
x3 = (0, None)

LP求解函数

 
#-*- coding:utf-8 -*-
#导入包
from scipy import optimize
import numpy as np
def LP(m='',clist=[],Alist=[],blist=[],Aeqlist=[],beqlist=[],all_x=()):
    #c,A,b,Aeq,beq,LB,UB,X0,OPTIONS
    c = np.array(clist)
    A = np.array(Alist)
    b = np.array(blist)
    Aeq = np.array(Aeqlist)
    beq = np.array(beqlist)
    #求解
    if m == 'min':
        res = optimize.linprog(c, A, b, Aeq, beq, bounds=all_x)
        fun = res.fun
        x = res.x
    else:
        res = optimize.linprog(-c, A, b, Aeq, beq, bounds=all_x)
        fun = -(res.fun)
        x = res.x
    return fun,x

main函数，方便其它调用：（这个是方便写成python文件直接用命令行执行）

 
#-*- coding:utf-8 -*-
import LP
import sys
if __name__ == '__main__':
    m = sys.argv[1]
    clist = list(eval(sys.argv[2]))
    Alist = list(eval(sys.argv[3]))
    blist = list(eval(sys.argv[4]))
    Aeqlist = list(eval(sys.argv[5]))
    beqlist =list(eval(sys.argv[6]))
    all_x = tuple(eval(sys.argv[7]))
    r=LP.LP(m=m,clist=clist,Alist=Alist,blist=blist,Aeqlist=Aeqlist,beqlist=beqlist,all_x=all_x)
    print(r)

2、分支定界算法

2.1 概念

分支定界算法（Branch and bound）始终围绕着一棵搜索树进行。分支的过程就是不断给树增加子节点的过程，而定界就是在分支的过程中检查子问题的上下界，如果子问题不能产生一比当前最优解还要优的解，那么砍掉这一枝。知道所有子问题都不能产生一个更优的解时，算法结束。

 具体到整数线性规划问题，这个用来对比“叶子”是否更优，能不能被砍掉的标准就是目前最优的一套整数解，例如最小化问题，记录这套整数解的目标函数的值，这个值是最小的。任何一个枝叶解出来的目标函数如果比这个还大，下面继续分支增加约束了之后更加不能得到更优的解，所以这些所有不小于最优整数解目标函数的分支应该全部砍掉。

但是一开始很可能没有一套整数解作为标准来比较，那么对于最小化问题，这个用来参照的目标函数就是无穷大，任何一个解出来的分支不管是不是整数解都不会被砍掉。这个参照更新发生在每次出现更优的整数解时，所以现在出现任何一个整数解都要优于这个无穷大，这个解就会变成参照。

如上图，4节点得到最优解Z=10（最大问题），其他3,7,8,6节点都不能优于这个整数解，往下面继续分支添加约束只会让结果越来越差，所以找到全局最优解4节点。

 下面上代码：

import numpy
from scipy.optimize import linprog
import sys
import math
 
#判断是否为整数（当一个值与他的下取整或上取整相差小于1e-6时，为整数）
def judge(L, t=1e-6):
    for i in L:
        if (i-math.floor(i)) < t or (math.ceil(i)-i) < t:
            return False
    return True
 
 
def intergerpro(c, A, B, Aeq, Beq, t=1.0e-6):
    res = linprog(c, A, B, Aeq, Beq)
    if not res.success:#整数规划问题不可解时终止分支（定界）
        return [sys.maxsize]
    if judge(res.x):
        bestX = [10000]*len(c)#抄的视频上的代码，意义不明
    else:
        bestX = res.x
    bestVal = sum(x*y for x, y in zip(c, bestX))
    if all(((x-math.floor(x)) < t or (math.ceil(x)-x) < t) for x in bestX):
        return (bestVal, bestX)#整数解
    else:#加约束条件
        ind = [i for i, x in enumerate(bestX) if (x-math.floor(x) > t) and (math.ceil(x)-x) > t][0]
        newcol1 = [0] * len(A[0])
        newcol2 = [0] * len(A[0])
        newcol1[ind] = -1#新不等式约束的系数
        newcol2[ind] = 1#新不等式约束的系数
        newA1 = A.copy()
        newA2 = A.copy()
        newA1.append(newcol1)#追加新不等式约束系数
        newA2.append(newcol2)#追加新不等式约束系数
        newB1 = B.copy()
        newB2 = B.copy()
        newB1.append(-math.ceil(bestX[ind]))#新不等式约束的值
        newB2.append(math.floor(bestX[ind]))#新不等式约束的值
        r1 = intergerpro(c, newA1, newB1, Aeq, Beq)#分支
        r2 = intergerpro(c, newA2, newB2, Aeq, Beq)#分支
        if r1[0] < r2[0]:
            return r1
        else:
            return r2
 
 
c = [3, 4, 1]#问题方程
A = [[-1, -6, -2],[-2, 0, 0]]#不等式约束的系数
B = [-5, -3]#不等式约束右边的值
Aeq = [[0, 0, 0]]#等式约束（aeq与beq都为0，代表0*x1+0*x2+0*x3=0，因此这里不写这个参数也行）
Beq = [0]#等式约束值
print(intergerpro(c, A, B, Aeq, Beq))

上面的代码其实没有剪枝，bound只有整数解，把所有分支往下死里搜索知道碰到整数解。用迭代算法深度有限，抓住一个分支往下死里搜索，搜不到就死循环。最下面的分支如果是整数解就由上面return这个整数解，然后上面一层一层迭代返回那些非整数解，当然，这是它这个节点的兄弟节点也往下算清楚的情况下，返回这两个兄弟中比较牛的那个（但其实意义不大了）。结果没有整理，直接一股脑迭代输出，没时间修改了。

有机会再把代码完善了 

2、非线性规划

非线性规划可简单分为两种，目标函数为凸函数or非凸函数

凸函数的非线性规划，比如，有很多常用的库完成，比如cvxpy

非凸函数的非线性规划（求极值），可尝试一下方法：

• 纯数学方法，求导求极值
• 神经网络、深度学习（反向传播算法中链式求导过程）
• scipy.optimize.minimize(接近matlab）

scipy.optimize.minimize(fun,x0,args=(),method=None,jac=None,hess=None,hessp=None,bounds=None,constraints=(),tol=None,callback=None,options=None)
 
fun:求最小值的目标函数
args:常数值
method:求极值方法，一般默认
constraints:约束条件
x0:变量的初始猜测值，注意minimize是局部最优

x0初始值不同，结果可能不同，最后要经过调参的过程。可以写一个循环调参，靠算力提升结果。

昨天的代码其实是matlab，用matlab深度学习确实挺奇怪的（