Introduction to Linear Programming (LP) with Python PuLP
ตัวอย่างการแก้ปัญหา LP ด้วย PuLP. พิจารณาปัญหา Product Mix อย่างง่าย: เราต้องการตัดสินใจว่าจะผลิตสินค้าแต่ละประเภทเป็นจำนวนเท่าใด เพื่อให้ได้กำไรรวมสูงสุด ภายใต้ข้อจำกัดเรื่องทรัพยากร. ทรัพยากรสามารถเป็นวัตถุดิบ แรงงานคนหรือเครื่องจักรเป็นต้น.
โจทย์: บริษัทผลิตและขายเครื่องดื่มหลักสองชนิด คือ
- Mojito
- Singapore Sling
ในการผลิตเครื่องดื่มสองชนิดนี้ต้องใช้ทรัพยากร (วัตถุดิบ) สามอย่าง คือ
- Rum
- Gin
- Lemon Juice
โดยมีอัตราการใช้เป็นดังนี้
- หนึ่งขวดของ Mojito ใช้ rum 1.0 หน่วย และใช้ lemon juice 0.75 หน่วย
- หนึ่งขวดของ Singapore Sling ใช้ gin 1.5 หน่วย และใช้ lemon juice 2.0 หน่วย
ในสัปดาห์นี้ เรามีวัตถุดิบ (gin, rum, lemon juice) อยู่อย่างจำกัด คือ
gin มีอยู่ 60 หน่วย
rum มีอยู่ 80 หน่วย
lemon juice มีอยู่ 120 หน่วย
ต้นทุนในการผลิต Mojito และ Singapore Sling เท่ากับ 50 และ 80 บาทต่อขวดตามลำดับ. ราคาขายต่อขวดของ Mojito และ Singapore Sling เท่ากับ 200 และ 280 บาทต่อขวดตามลำดับ. จะได้
- กำไรต่อขวด Mojito คือ 200–50 = 150 บาทต่อขวด
- กำไรต่อขวดของ Singapore Sling คือ 280–80 = 200 บาทต่อขวด
สังเกตว่า กำไรต่อขวดของ Singapore Sling ถึงแม้จะสูงกว่าของ Mojito, แต่ในการผลิตใช้ lemon juice เป็นจำนวนมากกว่า. จึงไม่ชัดเจนว่า แผนการผลิตควรเป็นอย่างไร.
เราต้องการทราบว่า ในสัปดาห์นี้ ควรผลิต Mojito และ Singapore Sling อย่างเท่าใด เพื่อให้ได้กำไรรวมสูงสุด ภายใต้ทรัพยากรที่จำกัด.
ปัญหาข้างต้น สามารถเขียนเป็น mathematical programming model ได้ดังนี้
Decision variables:
Xm = จำนวน Mojito ที่จะผลิตในสัปดาห์นี้ (หน่วย: ขวด)
Xs = จำนวน Singapore Sling ที่จะผลิตในสัปดาห์นี้ (หน่วย: ขวด)Objective function:
Maximize z = 150Xm + 200Xs
Subject to:
1.5Xs ≤ 60 (gin constraint)
1.0Xm ≤80 (rum constraint)
0.75Xm + 2.0Xs ≤ 120 (lemon juice constraint)
Xm, Xs ≥ 0 (non-negativity constraint)
สังเกตว่า ฟังก์ชันทั้งหมดใน objective และ constraint เป็น linear, ตัวแบบข้างต้นจึงเป็น linear programming (LP).
Simple model
เราจะแสดงการใช้ PuLP ใน Python เพื่อแก้ปัญหา LP ข้างต้น
- import functions ทั้งหมดของ PuLP เข้ามา
- สร้าง model เก็บไว้ใน prob โดยใช้ฟังก์ชัน
LpProblem
ซึ่งมี 2 parameters, ตัวแรกเป็นชื่ออธิบายปัญหา (ในที่นี้คือ string “ProductMixSimple”), ตัวหลังเป็นLpMinimize
หรือLpMaximize
ขึ้นอยู่กับปัญหา - ระบุ decision variables ใช้ฟังก์ชัน
LpVariable
ซึ่งมีได้ 4 parameters ดังนี้ 1) ตัวแรกเป็นชื่ออธิบายว่าตัวแปรนี้คืออะไร เช่น ตัวแปร Xm แทนจำนวน “Mojito”; 2) ตัวสองเป็น lower bound ซึ่งในที่นี้ทั้ง Xm และ Xs ≥ 0; 3) ตัวสามเป็น upper bound ซึ่งในที่นี้ไม่ได้ระบุไว้; 4) ตัวสุดท้ายเป็นLpContinuous
(ซึ่งเป็น default) หรือเป็นLpInteger
หากต้องการให้คำตอบเป็นจำนวนเต็ม. - เพิ่ม objective function ด้วย operator += และ , ตามหลังด้วยชื่อวัตถุประสงค์ที่ต้องการจะ optimize ซึ่งในที่นี้คือกำไรรวม “TotalProfit”. (สัญลักษณ์
+=
ใช้ใน code ภาษาอื่นด้วย เช่น ใน C ในการวน loop, สามารถเขียนเป็นt = t+ 1
หรือt += 1
ก็ได้.) - เพิ่ม constraints สำหรับวัตถุดิบทั้งสามอย่าง. สังเกตว่าเราไม่จำเป็นต้องใส่ non-negativity constraint เข้าไปเนื่องจากได้กำหนดไว้แล้วตั้งแต่ตอนประกาศตัวแปร Xm,Xs. เป็นอันเสร็จสิ้นการสร้าง model.
- หาคำตอบที่ดีสุด ใช้ attribute
solve()
บน model prob ที่สร้างไว้ในบรรทัด 2. - แสดงสถานะที่ได้จากการ solve ซึ่งจะเป็น Not Solved, Infeasible, Unbounded, Undefined, Optimal อย่างใดอย่างหนึ่ง
- แสดง optimal solution ที่ได้โดยใช้ for loop วนไปบน
variables
ซึ่งเป็น attribute ของ prob. ในโจทย์นี้ แผนการผลิตที่ดีสุด คือ ผลิต Mojito 80 ขวด และผลิต Sling 30 ขวด. หากเราต้องการแสดงเฉพาะจำนวนการผลิตของ Mojito สามารถทำได้โดยใช้ attributevarValue
ของตัวแปรนั้น นั่นคือprint(Xm.varValue)
จะแสดงค่า 80. - แสดง optimal objective value โดยใช้ attribute
objective
ของ prob. จากแผนการผลิตนี้ ทำให้ได้กำไรสูงสุดคือ 18,000 บาทในสัปดาห์นี้.
สังเกต constraints ของ gin และ rum สามารถเขียนได้เป็น Xs ≤ 40 และ Xm ≤ 80. ดังนั้น เราไม่จำเป็นต้องใส่ constraints ทั้งสองนี้ถ้าหากเราระบุ upper bound ของ decision variable ตั้งแต่ต้น
Xm = LpVariable("Mojito", 0, 80)
Xs = LpVariable("Sling", 0, 40)
Better model
ตัวแบบข้างต้น สามารถเขียนให้อ่านง่ายขึ้น โดยแยกเป็นส่วน input data และส่วนของ model.
ในส่วนของ data มี list PROD เก็บชื่อประเภทเครื่องดื่ม. มี dictionary profit เก็บกำไรต่อขวดของเครื่องดื่มแต่ละประเภท เป็นต้น.
ในส่วนของ model ระบุ decision variable ที่ชื่อ quant (เก็บเป็น dictionary) โดยใช้ฟังก์ชัน LpVariable.dicts
. ฟังก์ชัน lpSum
ทำให้เขียน objective และ constraint functions ได้สั้นลง.
เราสามารถ optimize ได้เช่นเดิม
สังเกตชื่อ decision variable ขึ้นต้นด้วย “Drink” ตามด้วยที่ชื่อระบุไว้ใน PROD เช่น Drink_Mojito.
หากต้องการนำผลลัพธ์บางตัวไปใช้ต่อในส่วนอื่นๆ ของโปรแกรม สามารถเรียกใช้ตัวแปร quant (ซึ่งได้ประกาศไว้ด้วย LpVariable.dicts
) เช่น
print(quant["Sling"].varValue)
แสดงค่า จำนวนผลิตที่เหมาะสมสุดของ Sling (ซึ่งคือ 30.0)
print(value(prob3.objective))
แสดงค่า optimal objective value กำไรสูงสุดที่ได้ (ซึ่งคือ 18000.0)
หาก decision variables มีทั้ง lower bounds และ upper bounds เช่น
20 ≤ Xm ≤ 80
10 ≤ Xs ≤ 40
เราสามารถสร้าง dictionary quantUB, quantLB
quantUB = {"Mojito": 80,
"Sling": 40}
quantLB = {"Mojito": 20,
"Sling": 10}
และหลังจากได้ประกาศตัวแปร quant เป็น decision variables แล้วเราสามารถใช้ attribute bounds
for i in PROD:
quant[i].bounds(quantLB[i], quantUB[i])
เพื่อใส่ bounds ให้กับ decision variables.
Model with many constraints
หากมี constraints เป็นจำนวนมาก สามารถใช้ for loop เพิ่ม constraints เข้าไปดังนี้
ในส่วน data, สังเกตการประกาศตัวแปร rateDict จะเป็น dictionary ที่มี dictionary ซ้อนกันอยู่ เช่น อัตราการใช้ lemon ใน Mojito คือ 0.75 หน่วยต่อขวด และใน Sling คือ 2.0 หน่วยต่อขวด.
ในส่วน model, สังเกตการระบุ resource constraint ใช้ for loop โดยมี iteration variable i วนบน RES. สำหรับแต่ละวัตถุดิบ i,ใน lpSum
เป็นการบวกวัตถุดิบที่ใช้ในการผลิตเครื่องดื่มทุกชนิด (ทุกเครื่องดื่มชนิด j).
อนึ่ง หากเรามีหลายสายการผลิตซึ่งทำสินค้าประเภทแตกต่างกัน, การแยก data และ model นั้นทำให้เราสามารถเปลี่ยนเฉพาะตรงส่วน data เท่านั้น และยังคงสามารถใช้ส่วนของ model ตามเดิมได้.
การหา optimal solution ทำได้เช่นเดิม.
เราสามารถตรวจสอบตัวแบบที่สร้างขึ้นได้ (ที่อยู่ใน prob3) โดย writeLP ออกมาเป็น file
prob3.writeLP('ProductMix3.lp')
เมื่อเปิด file ProfuctMix3.lp โดยใช้ text editor จะได้
สังเกตชื่อ constraint ขึ้นต้นด้วย “ConRes” ตามด้วยที่ชื่อระบุไว้ใน RES เช่น ConResLemon. ค่า parameters ใน input data ได้ถูกใส่เข้าไปแล้ว เช่น อัตราการใช้ lemon ในการผลิต Mojito
print(rateDict["Lemon"]["Mojito"])
จะได้ 0.75 ซึ่งอยู่ใน constraint ConResLemon.