Python小白的數學建模課-06 固定費用問題

youcans發表於2021-06-10

Python

Python 例項介紹固定費用問題的建模與求解。

學習 PuLP工具包中處理複雜問題的快捷使用方式。

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

前文講到幾種典型的 0-1 規劃問題，給出了 PuLP 求解的案例。由於 0-1 規劃問題種類很多，又是數模競賽熱點，有必要再結合幾個例項進行介紹。

1. 固定費用問題案例解析

1.1 固定費用問題（Fixed cost problem）

固定費用問題，是指求解生產成本最小問題時，總成本包括固定成本和變動成本，而選擇不同生產方式會有不同的固定成本，因此總成本與選擇的生產方式有關。

固定費用問題，實際上是互斥的目標函式問題，對於不同的生產方式具有多個互斥的目標函式，但只有一個起作用。固定費用問題不能用一般的線性規劃模型求解。

一般地，設有 m 種生產方式可供選擇，採用第 j 種方式時的固定成本為 \(K_j\)、變動成本為 \(c_j\)、產量為 \(x_j\)，則採用各種生產方式的總成本分別為：

\[min\;P_j = \begin{cases} k_j + c_j x_j，&x_j \geq 0\\ 0，&x_j = 0, j=1,...m \end{cases} \]

該類問題的建模方法，為了構造統一的目標函式，可以引入 m 個 0-1 變數 y_j 表示是否採用第 j 種生產方式：

\[y_j = \begin{cases} 0，不採用第\;j\;種生產方式\\ 1，採用第\;j\; 種生產方式 \end{cases} \]

於是可以構造新的目標函式和約束條件：

\[min\;f(x) = \sum_{j=1} ^m (k_j y_j + c_j x_j)\\ s.t.:\;x_j \leq y_j M，j=1,...m \]

M 是一個充分大的常數。

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

1.2 案例問題描述

例題 1：
某服裝廠可以生產 A、B、C 三種服裝，生產不同種類服裝需要租用不同裝置，裝置租金、生產成本、銷售價格等指標如下表所示。

服裝種類	裝置租金	材料成本	銷售價格	人工工時	裝置工時	裝置可用工時
單位	（元）	（元/件）	（元/件）	（小時/件）	（小時/件）	（小時）
A	5000	280	400	5	3	300
B	2000	30	40	1	0.5	300
C	2000	200	300	4	2	300

如果各類服裝的市場需求都足夠大，服裝廠每月可用人工時為 2000h，那麼應該如何安排生產計劃使利潤最大？

1.3 建模過程分析

首先要理解生產某種服裝就會發生裝置租金，租金只與是否生產該產品有關，而與生產數量無關，這就是固定成本。因此本題屬於固定費用問題。

有些同學下意識地認為是從 3 種產品中選擇一種，但題目中並沒有限定必須或只能生產一種產品，因此決策結果可以是都不生產、選擇 1 種或 2 種產品、3 種都生產。

決策結果會是什麼都不生產嗎？有可能的。

每種產品的利潤：（銷售價格 - 材料成本）× 生產數量 - 裝置租金

本題中如果裝置租金很高，決策結果就可能是什麼都不做時利潤最大，這是利潤為 0，至少不虧。

現在可以用固定費用問題的數學模型來描述問題了：

設 \(x_i\) 為是否生產第 \(i\) 種服裝，\(x_i\) 是 0/1變數：

\[x_i = \begin{cases} 0，不生產第\;i\;種服裝\\ 1，生產第\;i\;種服裝，i=1,2,3 \end{cases} \]

設 \(y_i\) 為生產第 \(i\) 種服裝的數量， \(y_i\) 是整數型別。說 \(y_i\) 是實數變數的同學，你經常穿半條褲子嗎？

根據條件確定決策變數的取值範圍。例如，本例中的產量 \(y_i\) 顯然要大於等於 0。進一步地，題目並沒有直接給出 \(y_i\) 的取值上限，但可以從裝置單件工時與裝置可用工時的關係推匯出取值上限為 [100, 600, 150]，也可以從單位人工工時與人工可用工時的關係推匯出上限 [400, 2000, 500]，最後取較小者為 [100, 600, 150]。

數學模型就可以表達為：

\[max\; z = 120 y_1 + 10 y_2 + 100 y_3 - 5000 x_1 - 2000 x_2 - 2000 x_3\\ s.t.:\;\begin{cases} 5 y_1 + y_2 + 4 y_3 \leq 2000\\ 3 y_1 \leq 300 x_1\\ 0.5 y_2 \leq 300 x_2\\ 2 y_3 \leq 300 x_3 \\ 0 \leq y_1 \leq 100\\ 0 \leq y_2 \leq 600\\ 0 \leq y_3 \leq 150\\ \end{cases} \]

1.4 PuLP 求解固定費用問題的程式設計

程式設計求解建立的數學模型，用標準模型的優化演算法對模型求解，得到優化結果。

模型求解的程式設計步驟與之前的線性規劃、整數規劃問題並沒有什麼區別，這就是 PuLP工具包的優勢。

（0）匯入 PuLP庫函式

    import pulp

（1）定義一個規劃問題

    FixedCostP1 = pulp.LpProblem("Fixed_cost_problem", sense=pulp.LpMaximize)  # 定義問題，求最大值

pulp.LpProblem 用來定義問題的建構函式。"FixedCostP1"是使用者定義的問題名。
引數 sense 指定問題求目標函式的最小值/最大值。本例求最大值，選擇 “pulp.LpMaximize” 。

（2）定義決策變數

    x1 = pulp.LpVariable('A', cat='Binary')  # 定義 x1，0-1變數，是否生產 A 產品
    x2 = pulp.LpVariable('B', cat='Binary')  # 定義 x2，0-1變數，是否生產 B 產品
    x3 = pulp.LpVariable('C', cat='Binary')  # 定義 x3，0-1變數，是否生產 C 產品
    y1 = pulp.LpVariable('yieldA', lowBound=0, upBound=100, cat='Integer')  # 定義 y1，整型變數
    y2 = pulp.LpVariable('yieldB', lowBound=0, upBound=600, cat='Integer')  # 定義 y2，整型變數
    y3 = pulp.LpVariable('youCans', lowBound=0, upBound=150, cat='Integer')  # 定義 y3，整型變數

pulp.LpVariable 用來定義決策變數的函式。引數 cat 用來設定變數型別，' Binary ' 表示0/1變數（用於0/1規劃問題），' Integer ' 表示整數變數。'lowBound'、'upBound' 分別表示變數取值範圍的下限和上限。

（3）新增目標函式

    FixedCostP1 += pulp.lpSum(-5000*x1-2000*x2-2000*x3+120*y1+10*y2+100*y3)  # 設定目標函式 f(x)

（4）新增約束條件

    FixedCostP1 += (5*y1 + y2 + 4*y3 <= 2000)  # 不等式約束
    FixedCostP1 += (3*y1 - 300*x1 <= 0)  # 不等式約束
    FixedCostP1 += (0.5*y2 - 300*x2 <= 0)  # 不等式約束
    FixedCostP1 += (2*y3 - 300*x3 <= 0)  # 不等式約束

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

（5）求解

    FixedCostP1.solve()

　　solve() 是求解函式，可以對求解器、求解精度進行設定。

1.5 Python 例程：固定費用問題

# mathmodel07_v1.py
# Demo05 of mathematical modeling algorithm
# Solving assignment problem with PuLP.
# Copyright 2021 Youcans, XUPT
# Crated：2021-06-04
# Python小白的數學建模課 @ Youcans

import pulp      # 匯入 pulp 庫

# 主程式
def main():
    # 固定費用問題(Fixed cost problem)
    print("固定費用問題(Fixed cost problem)")
    # 問題建模：
    """
        決策變數：
            y(i) = 0, 不生產第 i 種產品
            y(i) = 1, 生產第 i 種產品            
            x(i), 生產第 i 種產品的數量, i>=0 整數
            i=1,2,3
        目標函式：
            min profit = 120x1 + 10x2+ 100x3 - 5000y1 - 2000y2 - 2000y3
        約束條件：
            5x1 + x2 + 4x3 <= 2000
            3x1 <= 300y1
            0.5x2 <= 300y2
            2x3 <= 300y3
        變數取值範圍：Youcans XUPT
            0<=x1<=100, 0<=x2<=600, 0<=x3<=150, 整數變數
            y1, y2 ,y3 為 0/1 變數 
    """
    # 1. 固定費用問題(Fixed cost problem), 使用 PuLP 工具包求解
    # (1) 建立優化問題 FixedCostP1: 求最大值(LpMaximize)
    FixedCostP1 = pulp.LpProblem("Fixed_cost_problem_1", sense=pulp.LpMaximize)  # 定義問題，求最大值
    # (2) 建立變數
    x1 = pulp.LpVariable('A', cat='Binary')  # 定義 x1，0-1變數，是否生產 A 產品
    x2 = pulp.LpVariable('B', cat='Binary')  # 定義 x2，0-1變數，是否生產 B 產品
    x3 = pulp.LpVariable('C', cat='Binary')  # 定義 x3，0-1變數，是否生產 C 產品
    y1 = pulp.LpVariable('yieldA', lowBound=0, upBound=100, cat='Integer')  # 定義 y1，整型變數
    y2 = pulp.LpVariable('yieldB', lowBound=0, upBound=600, cat='Integer')  # 定義 y2，整型變數
    y3 = pulp.LpVariable('yieldC', lowBound=0, upBound=150, cat='Integer')  # 定義 y3，整型變數
    # (3) 設定目標函式
    FixedCostP1 += pulp.lpSum(-5000*x1-2000*x2-2000*x3+120*y1+10*y2+100*y3)  # 設定目標函式 f(x)
    # (4) 設定約束條件
    FixedCostP1 += (5*y1 + y2 + 4*y3 <= 2000)  # 不等式約束
    FixedCostP1 += (3*y1 - 300*x1 <= 0)  # 不等式約束
    FixedCostP1 += (0.5*y2 - 300*x2 <= 0)  # 不等式約束
    FixedCostP1 += (2*y3 - 300*x3 <= 0)  # 不等式約束
    # (5) 求解 youcans
    FixedCostP1.solve()
    # (6) 列印結果
    print(FixedCostP1.name)
    if pulp.LpStatus[FixedCostP1.status] == "Optimal":  # 獲得最優解
        for v in FixedCostP1.variables():  # youcans
            print(v.name, "=", v.varValue)  # 輸出每個變數的最優值
        print("Youcans F(x) = ", pulp.value(FixedCostP1.objective))  # 輸出最優解的目標函式值
    return

if __name__ == '__main__':  # Copyright 2021 YouCans, XUPT
    main()  # Python小白的數學建模課 @ Youcans

1.6 Python 例程執行結果

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

Result - Optimal solution found

Fixed_cost_problem_1
A = 1.0
B = 1.0
C = 1.0
yieldA = 100.0
yieldB = 600.0
yieldC = 150.0
Max F(x) =  24000.0

　　從固定費用問題模型的求解結果可知，A、B、C 三種服裝都生產，產量分別為 A/100、B/600、C/150 時獲得最大利潤為：24000。

2. PuLP 求解規劃問題的快捷方法

2.1 PuLP 求解固定費用問題的程式設計

通過從線性規劃、整數規劃、0-1規劃到上例中的混合0-1規劃問題，我們已經充分體會到 PuLP 使用相同的步驟和引數處理不同問題所帶來的便利。

但是，如果問題非常複雜，例如變數數量很多，約束條件複雜，逐個定義變數、逐項編寫目標函式與約束條件的表示式，不僅顯得重複冗長，不方便修改對變數和引數的定義，而且在輸入過程中容易發生錯誤。因此，我們希望用字典、列表、迴圈等快捷方法來進行變數定義、目標函式和約束條件設定。

PuLP 提供了快捷建模的程式設計方案，下面我們仍以上節中的固定費用問題為例進行介紹。本例中的問題、條件和引數都與上節完全相同，以便讀者進行對照比較快捷方法的具體內容。

（0）匯入 PuLP 庫函式

    import pulp

（1）定義一個規劃問題

    FixedCostP2 = pulp.LpProblem("Fixed_cost_problem", sense=pulp.LpMaximize)  # 定義問題，求最大值

（2）定義決策變數

    types = ['A', 'B', 'C']  # 定義產品種類
    status = pulp.LpVariable.dicts("生產決策", types, cat='Binary')  # 定義 0/1 變數，是否生產該產品
    yields = pulp.LpVariable.dicts("生產數量", types, lowBound=0, upBound=600, cat='Integer')  # 定義整型變數

本例中的快捷方法使用列表 types 定義 0/1 變數 status 和整型變數 yields，不論產品的品種有多少，都只有以上幾句，從而使程式大為簡化。

（3）新增目標函式

    fixedCost = {'A':5000, 'B':2000, 'C':2000}  # 各產品的 固定費用
    unitProfit = {'A':120, 'B':10, 'C':100}  # 各產品的 單位利潤
    FixedCostP2 += pulp.lpSum([(yields[i]*unitProfit[i]- status[i]*fixedCost[i]) for i in types])

雖然看起來本例中定義目標函式的程式語句較長，但由於使用字典定義引數、使用 for 迴圈定義目標函式，因此程式更加清晰、簡明、便於修改引數、不容易輸入錯誤。

（4）新增約束條件

    humanHours = {'A':5, 'B':1, 'C':4}  # 各產品的 單位人工工時
    machineHours = {'A':3.0, 'B':0.5, 'C':2.0}  # 各產品的 單位裝置工時
    maxHours = {'A':300, 'B':300, 'C':300}  # 各產品的 最大裝置工時
    FixedCostP2 += pulp.lpSum([humanHours[i] * yields[i] for i in types]) <= 2000  # 不等式約束
    for i in types:
        FixedCostP2 += (yields[i]*machineHours[i] - status[i]*maxHours[i] <= 0)  # 不等式約束

快捷方法對於約束條件的定義與對目標函式的定義相似，使用字典定義引數，使用迴圈定義約束條件，使程式簡單、結構清楚。

注意本例使用了兩種不同的迴圈表達方式：語句內使用 for 迴圈遍歷列表實現所有變數的線性組合，標準的 for 迴圈結構實現多組具有相似結構的約束條件。讀者可以對照數學模型及上例的例程，理解這兩種定義約束條件的快捷方法。

（5）求解和結果的輸出

    # (5) 求解
    FixedCostP2.solve()
    # (6) 列印結果
    print(FixedCostP2.name)
    temple = "品種 %(type)s 的決策是：%(status)s，生產數量為：%(yields)d"
    if pulp.LpStatus[FixedCostP2.status] == "Optimal":  # 獲得最優解
        for i in types:
            output = {'type': i,
                      'status': '同意' if status[i].varValue else '否決',
                      'yields': yields[i].varValue}
            print(temple % output) # youcans@qq.com
        print("最大利潤 = ", pulp.value(FixedCostP2.objective))  # 輸出最優解的目標函式值

由於快捷方法使用列表或字典定義變數，對求解的優化結果也便於實現結構化的輸出。

2.2 Python 例程：PuLP 快捷方法

# mathmodel07_v1.py
# Demo05 of mathematical modeling algorithm
# Solving assignment problem with PuLP.
# Copyright 2021 Youcans, XUPT
# Crated：2021-06-04
# Python小白的數學建模課 @ Youcans

import pulp      # 匯入 pulp 庫


# 主程式
def main():
    # 2. 問題同上，PuLP 快捷方法示例
    # (1) 建立優化問題 FixedCostP2: 求最大值(LpMaximize)
    FixedCostP2 = pulp.LpProblem("Fixed_cost_problem_2", sense=pulp.LpMaximize)  # 定義問題，求最大值
    # (2) 建立變數
    types = ['A', 'B', 'C']  # 定義產品種類
    status = pulp.LpVariable.dicts("生產決策", types, cat='Binary')  # 定義 0/1 變數，是否生產該產品
    yields = pulp.LpVariable.dicts("生產數量", types, lowBound=0, upBound=600, cat='Integer')  # 定義整型變數
    # (3) 設定目標函式
    fixedCost = {'A':5000, 'B':2000, 'C':2000}  # 各產品的 固定費用
    unitProfit = {'A':120, 'B':10, 'C':100}  # 各產品的 單位利潤
    FixedCostP2 += pulp.lpSum([(yields[i]*unitProfit[i]- status[i]*fixedCost[i]) for i in types])
    # (4) 設定約束條件
    humanHours = {'A':5, 'B':1, 'C':4}  # 各產品的 單位人工工時
    machineHours = {'A':3.0, 'B':0.5, 'C':2.0}  # 各產品的 單位裝置工時
    maxHours = {'A':300, 'B':300, 'C':300}  # 各產品的 最大裝置工時
    FixedCostP2 += pulp.lpSum([humanHours[i] * yields[i] for i in types]) <= 2000  # 不等式約束
    for i in types:
        FixedCostP2 += (yields[i]*machineHours[i] - status[i]*maxHours[i] <= 0)  # 不等式約束
    # (5) 求解 youcans
    FixedCostP2.solve()
    # (6) 列印結果
    print(FixedCostP2.name)
    temple = "品種 %(type)s 的決策是：%(status)s，生產數量為：%(yields)d"
    if pulp.LpStatus[FixedCostP2.status] == "Optimal":  # 獲得最優解
        for i in types:
            output = {'type': i,
                      'status': '同意' if status[i].varValue else '否決',
                      'yields': yields[i].varValue}
            print(temple % output)
        print("最大利潤 = ", pulp.value(FixedCostP2.objective))  # 輸出最優解的目標函式值

    return

if __name__ == '__main__':  # Copyright 2021 YouCans, XUPT
    main()  # Python小白的數學建模課 @ Youcans

2.3 Python 例程執行結果

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

Result - Optimal solution found

Fixed_cost_problem_2
品種 A 的決策是：同意，生產數量為：100
品種 B 的決策是：同意，生產數量為：600
品種 C 的決策是：同意，生產數量為：150
最大利潤 =  24000.0

　　本例的問題、條件和引數都與上節完全相同，只是採用 PuLP 提供的快捷建模的程式設計方案，優化結果也與 PuLP 標準方法完全相同，但本例使用了結構化的輸出顯示，使輸出結果更為直觀。

3. 課後練習

修改生產某種服裝的裝置租金，例如將 A 產品租金調整為 10000、20000元，觀察求解結果有何差異？
將各種裝置租金都調整為 20000元，觀察求解結果有何差異？該結果有何現實意義？
如果希望找到影響是否生產某種服裝決策的裝置租金的大小，即租金低於該值就可以生產、高於該值則不能生產，應該如何處理？

【本節完】

版權宣告：

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

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

Crated：2021-06-04

歡迎關注『Python小白的數學建模課 @ Youcans』，每週更新數模筆記
Python小白的數學建模課-01.新手必讀
 Python小白的數學建模課-02.資料匯入
 Python小白的數學建模課-03.線性規劃
 Python小白的數學建模課-04.整數規劃
 Python小白的數學建模課-05.0-1規劃
 Python小白的數學建模課-A1.國賽賽題型別分析
 Python小白的數學建模課-A2.2021年數維杯C題探討
 Python小白的數學建模課-A3.12個新冠疫情數模競賽賽題及短評
 Python數模筆記-PuLP庫
 Python數模筆記-StatsModels統計迴歸
 Python數模筆記-Sklearn
Python數模筆記-NetworkX
Python數模筆記-模擬退火演算法

Python小白的數學建模課-07 選址問題
2021-06-13
Python
Python小白的數學建模課-19.網路流優化問題
2021-08-24
Python優化
Python小白的數學建模課-04.整數規劃
2021-06-03
Python
Python小白的數學建模課-05.0-1規劃
2021-06-06
Python
Python小白的數學建模課-09 微分方程模型
2021-06-21
Python模型
Python小白的數學建模課-03.線性規劃
2021-06-01
Python
Python小白的數學建模課-B5. 新冠疫情 SEIR模型
2021-07-10
Python模型
Python小白的數學建模課-16.最短路徑演算法
2021-08-06
Python演算法
小白學python遇到問題了.......
2019-09-28
Python
小白能學好Python嗎?學費是多少?
2021-01-07
Python
python 來做數學問題
2019-08-04
Python
LeetCode數學問題（Python）
2020-10-08
LeetCodePython
Python數學建模-01.新手必讀
2021-05-29
Python
華中杯數學建模 A題簡單覆盤（附Python原始碼）
2022-08-09
Python原始碼
Python數學建模-02.資料匯入
2021-05-30
Python
python 小白求教一個字串顯示問題
2020-06-18
Python字串
想要學習Python課程，這些問題你必須知道!
2019-12-16
Python
Python小白的爬蟲問題與解決（含程式碼）
2018-05-10
Python爬蟲
用Python解決一個等差數列的求和問題
2020-02-08
Python
小白學python系列-(3)基礎數量型別
2018-04-05
Python型別
Python指令碼應用及學習方法-CSDN公開課-專題視訊課程
2016-01-28
Python指令碼
CSDN學院：24門精品課程免費領取（含Python+人工智慧等專題課程）
2020-04-04
Python人工智慧
2020年數學建模國賽B題解題思路
2020-09-10
06小白學nodejs fs模組其它方法
2021-05-25
NodeJS
Python中的生產者消費者問題
2013-12-04
Python
[全程建模]UML設計類中的實現與方法數量問題
2013-09-14
遊戲數學建模工程手冊
2020-05-18
遊戲
nginx下wordpress 固定連結問題
2014-08-02
Nginx
初學python的常見問題
2020-10-07
Python
小白學 Python（1）：開篇
2019-10-18
Python
python適合小白學嗎？
2021-09-11
Python
小白學python系列-(4)list
2018-04-06
Python
小白學python系列-(5) tuple
2018-04-06
Python
小白學python系列-(8)dict
2018-04-06
Python
資料建模學習筆記 -- 類別不平衡問題
2020-12-09
筆記
[轉載]2016美國數學建模MCM F題（政策）翻譯:難民移入政策建模
2018-02-03
2022全國大學生數學建模A題的思路與解法
2022-10-18
Python入門_給小白的學習指南
2020-02-06
Python