Python數學建模系列(二):規划問題之整數規划


@

前言

Hello!小伙伴!
非常感謝您閱讀海轟的文章,倘若文中有錯誤的地方,歡迎您指出~
 
自我介紹 ଘ(੭ˊᵕˋ)੭
昵稱:海轟
標簽:程序猿|C++選手|學生
簡介:因C語言結識編程,隨后轉入計算機專業,有幸拿過一些國獎、省獎...已保研。目前正在學習C++/Linux/Python
學習經驗:扎實基礎 + 多做筆記 + 多敲代碼 + 多思考 + 學好英語!
 
初學Python 小白階段
文章僅作為自己的學習筆記 用於知識體系建立以及復習
題不在多 學一題 懂一題
知其然 知其所以然!
 
本文僅從Pyhton如何解決建模問題出發
未對建模思路等進行深一步探索

推薦(排版簡潔 適合閱讀)
https://mp.weixin.qq.com/s/O7YIWsJZY6rWfbqsjcv88g

整數規划

整數規划的模型與線性規划基本相同,只是額外增加了部分變量為整數的約束

整數規划求解的基本框架是分支定界法,首先去除整數約束得到"松弛模型"。使用線性規划的方法求解。

若有某個變量不是整數,在松弛模型.上分別添加約束:x≤floor(A)和x≥ceil(A),然后再分別求解,這個過程叫做分支。當節點求解結果中所有變量都是整數時。停止分支。這樣不斷迭代,形成了一顆樹。

所謂定界,指的是葉子節點產生后,相當於給問題定了一個下界。之后在求解過程中一旦某個節點的目標函數值小於這個下界,那就直接pass,不再進行分支了;每次新產生葉子節點,則更新下界。

例題

\(min \quad z = 3x_1 + 4x_2 + x_3\) 的最小值

\[min \quad z = 3x_1 + 4x_2 + x_3\\ s.t. = \begin{cases} x_1 + 6x_2 + 2x_3 >= 5 \\ 2x_1 >= 3\\ x_1,x_2,x_3 >=0(且都為整數) \end{cases} \]

方法一:分支定界法(使用scipy庫)

Demo代碼

# 運行環境:Vs Code
import math
from scipy.optimize import linprog
import sys

def integerPro(c, A, b, Aeq, beq,t=1.0E-8):
    res = linprog(c, A_ub=A, b_ub=b, A_eq=Aeq, b_eq=beq)
    bestVal = sys.maxsize
    bestX = res.x
    if not(type(res.x) is float or res.status != 0): 
        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]
        newCon1 = [0]*len(A[0])
        newCon2 = [0]*len(A[0])
        newCon1[ind] = -1
        newCon2[ind] = 1
        newA1 = A.copy()
        newA2 = A.copy()
        newA1.append(newCon1)
        newA2.append(newCon2)
        newB1 = b.copy()
        newB2 = b.copy()
        newB1.append(-math.ceil(bestX[ind]))
        newB2.append(math.floor(bestX[ind]))
        r1 = integerPro(c, newA1, newB1, Aeq, beq)
        r2 = integerPro(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]]
beq = [0]
print(integerPro(c, A, b, Aeq, beq))

運行結果

(8.000000000001586, array([2., 0., 2.]))
# 或者
(8.000000000001586, array([2.00000000e+00, 1.83247535e-13, 2.00000000e+00]))

在這里插入圖片描述

方法二:使用pulp庫進行求解

只需要在設置變量的時候

設置參數cat='Integer' 即可

  • Continuous:連續
  • Binary:0 或 1
  • Integer:整數

Demo代碼

import pulp as pp

# 參數設置
c = [3,4,1]        #目標函數未知數前的系數

A_gq = [[1,6,2],[2,0,0]]   # 大於等於式子 未知數前的系數集合 二維數組 
b_gq = [5,3]         # 大於等於式子右邊的數值 一維數組


# 確定最大最小化問題,當前確定的是最小化問題
m = pp.LpProblem(sense=pp.LpMinimize)

# 定義三個變量放到列表中 生成x1 x2 x3
x = [pp.LpVariable(f'x{i}',lowBound=0,cat='Integer') for i in [1,2,3]]

# 定義目標函數,並將目標函數加入求解的問題中 
m += pp.lpDot(c,x) # lpDot 用於計算點積 

# 設置比較條件
for i in range(len(A_gq)):# 大於等於
    m += (pp.lpDot(A_gq[i],x) >= b_gq[i])

# 求解
m.solve()

# 輸出結果
print(f'優化結果:{pp.value(m.objective)}')
print(f'參數取值:{[pp.value(var) for var in x]}')

運行結果

優化結果:8.0
參數取值:[2.0, 0.0, 2.0]

結語

學習來源:B站及其課堂PPT,對其中代碼進行了復現

鏈接:https://www.bilibili.com/video/BV12h411d7Dm?from=search&seid=5685064698782810720

文章僅作為學習筆記,記錄從0到1的一個過程

希望對您有所幫助,如有錯誤歡迎小伙伴指正~

我是 海轟ଘ(੭ˊᵕˋ)੭

如果您覺得寫得可以的話,請點個贊吧

謝謝支持 ❤️

在這里插入圖片描述


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM