Python小白的數學建模課-03.線性規劃

youcans發表於2021-06-01
Python

  • 線性規劃是很多數模培訓講的第一個演算法，演算法很簡單，思想很深刻。

  • 要通過線性規劃問題，理解如何學習數學建模、如何選擇程式設計演算法。

  • 『Python小白的數學建模課 @ Youcans』帶你從數模小白成為國賽達人。


1. 求解方法、演算法和程式設計方案

線性規劃 （Linear Programming，LP） 是很多數模培訓講的第一個演算法，演算法很簡單，思想很深刻。

線性規劃問題是中學數學的內容，雞兔同籠就是一個線性規劃問題。數學規劃的題目在高考中也經常出現，有直接給出線性約束條件求線性目標函式極值，有間接給出約束條件求線性目標函式極值，還有已知約束條件求非線性目標函式極值問題。

因此，線性規劃在數學建模各類問題和演算法中確實是比較簡單的問題。下面我們通過這個比較簡單、也比較熟悉的問題，分析一下數學模型問題的方法、演算法和學習方案。探討這些容易混淆的概念，也便於大家理解本系列教程的初衷和特色。


歡迎關注 『Python小白的數學建模課 @ Youcans』，每週更新數模筆記
Python小白的數學建模課-01.新手必讀
Python小白的數學建模課-02.資料匯入
Python小白的數學建模課-03.線性規劃
Python數模筆記-PuLP庫
Python數模筆記-StatsModels統計迴歸
Python數模筆記-Sklearn
Python數模筆記-NetworkX
Python數模筆記-模擬退火演算法


1.1 線性規劃問題的求解方法

解決線性規劃問題有很多數學方法，例如：

  • 圖解法， 用幾何作圖的方法並求出其最優解，中學就講過這種方法，在經濟學研究中十分常用；
  • 矩陣法， 引進鬆弛變數將線性規劃問題轉換成增廣矩陣形式後逐次求解， 是單純性法之前的典型方法；
  • 單純性法， 利用多面體在可行域內逐步構造新的頂點來不斷逼近最優解，是線性規劃研究的里程碑，至今仍然是最重要的方法之一；
  • 內點法，通過選取可行域內部點沿下降方向不斷迭代來達到最優解，是目前理論上最好的線性規劃問題求解方法；
  • 啟發式方法，依靠經驗準則不斷迭代改進來搜尋最優解 ，如貪心法、模擬退火、遺傳演算法、神經網路。

雖然不同的求解方法都是面對線性規劃問題，也就必然會殊途同歸，但它們在思想上就存在著本質區別，在求解方法和步驟上也就完全不同。

不誇張地說，對於很多小白，學沒學過單純性法，對於學習啟發式方法可能完全沒有區別。

這意味著什麼呢？這就是說，對於非數學專業的同學，對於學習數學建模的同學，針對每一類問題，完全沒必要學習各種解決方法。即便是數學專業的同學，也不可能在數模學習期間把各種方法都學會。

對於小白，本文推薦選擇較為通用、相對簡單（思路簡單、程式簡單）的方法來進行學習，沒必要貪多求新。


1.2 線性規劃的最快演算法

演算法，跟方法有什麼不同呢？

演算法的定義是“解題方案的準確而完整的描述”，是一系列解決問題的清晰指令，演算法代表著用系統的方法描述解決問題的策略機制。

我對“方法”的理解是思想方法，是求解問題總體框架，而“演算法”是具體和明確的實現步驟，在計算機程式設計中相當於詳細的流程圖。

在每一種方法的基本思想和方案提出後，往往都會有很多變形、改進和發展的演算法。極少的改進演算法具有實質貢獻而成為主流的經典演算法，即便如此往往也只是效能、效率上的提升，對於求解數模競賽中的問題基本沒有影響。

而絕大多數改進演算法只是針對某些特殊情況、特殊問題（自稱）有效，常用於大量的灌水論文。對於數學建模來說，學習基本演算法或者目前的經典演算法就足夠了，不需要聽信改進演算法中自稱的優點，那都是莆田系的廣告。

有一種例外情況，就是一些演算法是有適用範圍和限制條件的。舉個例子，內點法的基本演算法不能處理等式約束，最短路徑問題中 Dijkstar演算法不能處理負權邊。這種情況下如果選錯演算法，問題是無法求解的。所以對我們來說，搞清楚演算法的適用範圍，比理解演算法本身更重要。

回到本節的標題，對於線性規劃問題，什麼演算法是最快的呢？答案是：猜。不是讓你猜，而是說求解線性規劃問題，猜起來比較快。不是開玩笑，我是認真的。

佐治亞理工學院彭泱教授在 2021年計算機理論頂會 SODA2021 獲得最佳論文（Best paper award at ACM-SIAM symposium on discrete algorithms 2021），正是研究線性規劃問題的求解——“Solving sparse linear systems faster than matrix multiplication”，所用的全新思路是：猜，反覆猜，迭代猜。



當然，猜起來比較快只是在某些特殊條件下才有效的，至於在什麼條件下猜，怎麼猜，這不是我們所要關心，所能理解的問題了。只是以此說明，簡單的問題也有複雜的情況，每個問題都有很多求解的思路、方法和演算法。


1.3 選擇適合自己的程式設計方案

程式設計方案是我杜撰的術語。我所要表達意思是，在選擇了求解方法和演算法以後，是自己按照演算法步驟一步步程式設計實現，或者找到例程除錯使用，還是呼叫第三方工具包/庫函式來完成呢？

首先，對於學習數學建模、參加數模競賽，不建議自己按照演算法步驟去程式設計。我們在《01.新手必讀》中討論過這個問題，對於數學小白兼計算機小白，這樣做既不可行也沒必要；即使你願意挑戰自我去試試，那其實已經是走在學習另一門計算機或演算法課程的路上了。

其次，要不要找到例程自己除錯、使用？很多數模培訓就是這麼說，這麼做的，而且把這些收集的例程當作核心機密吸引同學。我不反對這樣做，這種學習方法對於理解演算法、提高程式設計能力很有幫助；但是並不推薦這樣做，原因是：（1）我認為學習數學建模、參加數模競賽，重點應該放在識別問題、分析問題、解決問題，能使用演算法和程式設計就足夠了；（2）第三方庫與例程沒有本質區別，第三方庫就是經典的、規範的、標準化的例程，既然選擇例程為什麼不選擇優秀的例程——第三方庫呢？（3）大部分例程都存在很多問題，即使除錯通過仍然有很多坑，而且新手難以識別。

所以我是明確推薦優選直接使用第三方庫來解決問題，這也是 Python 語言“不要重複造輪子”的思想。

進一步地，很多工具包/庫函式都能實現常用的演算法，應該如何選擇呢？

如果你對某個工具包已經很熟悉，又能實現所要的演算法，這當然是理想的選擇。如果你是小白，就跟著我走吧。

本系列選擇第三方工具包的原則是：（1）優選常用的工具包；（2）優選通用功能的工具包和函式（例如，最好既能實現線性規劃，又能實現整數規劃、非線性規劃）；（3）優選安裝簡單、使用簡單、配置靈活的工具包；（4）優選兼模型檢驗、圖形繪製的工具包。


2. PuLP庫求解線性規劃問題

2.1 線性規劃問題的描述

線性規劃是研究線性等式或不等式約束條件下求解線性目標函式的極值問題，常用於解決資源分配、生產排程和混合問題。

一般線性規劃問題的標準形式為：

\[max\;f(x) = \sum_{j=1} ^n c_j x_j\\ s.t.:\begin{cases} \sum_{j=1} ^n a_{ij} x_j = b_i, \\ x_j \geq 0 \end{cases} \]

滿足所有約束條件的解，稱為線性規劃問題的可行解；所有可行解構成的集合，稱為可行域。

使目標函式達到最小值的解，稱為最優解。

線性規劃問題的建模和求解，通常按照以下步驟進行：

  1. 問題定義，確定決策變數、目標函式和約束條件；
  2. 模型構建，由問題描述建立數學方程，並轉化為標準形式的數學模型；
  3. 模型求解，用標準模型的優化演算法對模型求解，得到優化結果。

很多 Python 的第三方包，都提供求解線性規劃問題的演算法，有的工具包還提供整數規劃、非線性規劃的演算法。例如：

  • Scipy 庫提供瞭解簡單線性或非線性規劃問題，但是不能求解如揹包問題的0-1規劃問題，或整數規劃問題，混合整數規劃問題。
  • PuLP 可以求解線性規劃、整數規劃、0-1規劃、混合整數規劃問題，提供多種針對不同型別問題的求解器。
  • Cvxpy 是一種凸優化工具包，可以求解線性規劃、整數規劃、0-1規劃、混合整數規劃、二次規劃和幾何規劃問題。

此外，SKlearn、DOcplex、Pymprog 等很多第三方工具包也都能求解線性規劃問題。


2.2 PuLP 求解線性規劃問題的步驟

例題 1：

max fx = 2x1 + 3x2 - 5x3
s.t. x1 + 3x2 + x3 <= 12
2x1 - 5x2 + x3 >= 10
x1 + x2 + x3 = 7
x1, x2, x3 >=0

下面以該題為例講解 PuLP 求解線性規劃問題的步驟：

（0）匯入 PuLP庫函式

    import pulp

（1）定義一個規劃問題

    MyProbLP = pulp.LpProblem("LPProbDemo1", sense=pulp.LpMaximize)

pulp.LpProblem 是定義問題的建構函式。
"LPProbDemo1"是使用者定義的問題名（用於輸出資訊）。
引數 sense 用來指定求最小值/最大值問題，可選引數值：LpMinimize、LpMaximize 。本例 “sense=pulp.LpMaximize” 表示求目標函式的最大值。

（2）定義決策變數

    x1 = pulp.LpVariable('x1', lowBound=0, upBound=7, cat='Continuous') 
    x2 = pulp.LpVariable('x2', lowBound=0, upBound=7, cat='Continuous')
    x3 = pulp.LpVariable('x3', lowBound=0, upBound=7, cat='Continuous')

pulp.LpVariable 是定義決策變數的函式。
'x1' 是使用者定義的變數名。
引數 lowBound、upBound 用來設定決策變數的下界、上界；可以不定義下界/上界，預設的下界/上界是負無窮/正無窮。本例中 x1,x2,x3 的取值區間為 [0,7]。
引數 cat 用來設定變數型別，可選引數值：'Continuous' 表示連續變數（預設值）、' Integer ' 表示離散變數（用於整數規劃問題）、' Binary ' 表示0/1變數（用於0/1規劃問題）。

（3）新增目標函式

    MyProbLP += 2*x1 + 3*x2 - 5*x3  	# 設定目標函式

新增目標函式使用 "問題名 += 目標函式式" 格式。

（4）新增約束條件

    MyProbLP += (2*x1 - 5*x2 + x3 >= 10)  # 不等式約束
    MyProbLP += (x1 + 3*x2 + x3 <= 12)  # 不等式約束
    MyProbLP += (x1 + x2 + x3 == 7)  # 等式約束

　　新增約束條件使用 "問題名 += 約束條件表示式" 格式。
　　約束條件可以是等式約束或不等式約束，不等式約束可以是 小於等於 或 大於等於，分別使用關鍵字">="、"<="和"=="。

（5）求解

    MyProbLP.solve()
    print("Status:", pulp.LpStatus[MyProbLP.status]) # 輸出求解狀態
    for v in MyProbLP.variables():
        print(v.name, "=", v.varValue)  # 輸出每個變數的最優值
    print("F(x) = ", pulp.value(MyProbLP.objective))  #輸出最優解的目標函式值

solve() 是求解函式。PuLP預設採用 CBC 求解器來求解優化問題，也可以呼叫其它的優化器來求解，如：GLPK，COIN CLP/CBC，CPLEX，和GUROBI，但需要另外安裝。　


2.3 Python例程：線性規劃問題

例程 1：求解線性規劃問題

# mathmodel04_v1.py
# Demo01 of mathematical modeling algorithm
# Solving linear programming with PuLP.
# Copyright 2021 Youcans, XUPT
# Crated：2021-05-28

import pulp
MyProbLP = pulp.LpProblem("LPProbDemo1", sense=pulp.LpMaximize)  # 求最大值
x1 = pulp.LpVariable('x1', lowBound=0, upBound=7, cat='Continuous') 
x2 = pulp.LpVariable('x2', lowBound=0, upBound=7, cat='Continuous') 
x3 = pulp.LpVariable('x3', lowBound=0, upBound=7, cat='Continuous') 
MyProbLP += 2*x1 + 3*x2 - 5*x3  	# 設定目標函式
MyProbLP += (2*x1 - 5*x2 + x3 >= 10)  # 不等式約束
MyProbLP += (x1 + 3*x2 + x3 <= 12)  # 不等式約束
MyProbLP += (x1 + x2 + x3 == 7)  # 等式約束
MyProbLP.solve()  # youcans@xupt
print("Status:", pulp.LpStatus[MyProbLP.status]) # 輸出求解狀態
for v in MyProbLP.variables():  # youcans
    print(v.name, "=", v.varValue)  # 輸出每個變數的最優值
print("Max F(x) = ", pulp.value(MyProbLP.objective))  #輸出最優解的目標函式值

例程 1 執行結果：

Welcome to the CBC MILP Solver 
Version: 2.9.0 
Build Date: Feb 12 2015 

Status: Optimal
x1 = 6.4285714
x2 = 0.57142857
x3 = 0.0
Max F(x) =  14.57142851

例程01 程式說明：

  1. 用 PuLP 庫求解線性規劃問題，可以選擇求最大值或最小值，可以按照問題的數學描述，直接輸入目標函式、等式約束和不等式約束，不等式約束可以選擇 <= 或 >=，不需要進行轉換。這中方式簡單直觀，非常適合初學者掌握。
  2. 對於較大規模線性規劃問題， PuLP 庫支援用字典型別（dict）建立多個變數，設定目標函式和約束條件。詳見拙文《Python數模筆記-PuLP庫（2）線性規劃進階》。

3. 小結

求解線性規劃問題的方法非常簡單，本文實際上並未講解具體的演算法。

希望通過對求解方法、演算法和程式設計方案的講解，闡明作者對於數學建模學什麼、怎麼學的理解，也使讀者能瞭解本系列教程的特點：本教程不打算詳細講解各種演算法的具體方法，重點介紹如何使用第三方包實現演算法、解決問題。

【本節完】


版權說明：

歡迎關注『Python小白的數學建模課 @ Youcans』 原創作品

原創作品，轉載必須標註原文連結。

Copyright 2021 Youcans, XUPT

Crated：2021-05-29


歡迎關注 『Python小白的數學建模課 @ Youcans』，每週更新數模筆記
Python小白的數學建模課-01.新手必讀
Python小白的數學建模課-02.資料匯入
Python小白的數學建模課-03.線性規劃
Python數模筆記-PuLP庫（1）線性規劃入門
Python數模筆記-PuLP庫（2）線性規劃進階
Python數模筆記-PuLP庫（3）線性規劃例項
Python數模筆記-StatsModels 統計迴歸（1）簡介
Python數模筆記-StatsModels 統計迴歸（2）線性迴歸
Python數模筆記-StatsModels 統計迴歸（3）模型資料的準備
Python數模筆記-StatsModels 統計迴歸（4）視覺化
Python數模筆記-Sklearn （1）介紹
Python數模筆記-Sklearn （2）聚類分析
Python數模筆記-Sklearn （3）主成分分析
Python數模筆記-Sklearn （4）線性迴歸
Python數模筆記-Sklearn （5）支援向量機
Python數模筆記-NetworkX（1）圖的操作
Python數模筆記-NetworkX（2）最短路徑
Python數模筆記-NetworkX（3）條件最短路徑
Python數模筆記-NetworkX（4）最小生成樹
Python數模筆記-NetworkX（5）關鍵路徑法
Python數模筆記-模擬退火演算法（1）多變數函式優化
Python數模筆記-模擬退火演算法（2）約束條件的處理
Python數模筆記-模擬退火演算法（3）整數規劃問題
Python數模筆記-模擬退火演算法（4）旅行商問題


相關文章