上述问题本质为一个线性规划问题。
而线性规划有以下性质。
- 比例性:每个决策变量对目标函数的贡献,与该决策变量的取值成正比;每个决策变量对每个约束条件右端项的贡献,与该决策变量的取值成正比
- 可加性:各个决策变量对目标函数的贡献,与其他决策变量的取值无关;各个决策变量对每个约束条件右端项的贡献,与其他决策变量的取值无关
- 连续性:每个决策变量的取值都是连续的
最终我们可作出如下假设
模型假设
- 各种蔬菜产品的单价均是与购买量无关的常数,各种蔬菜产品提供的营养价值也是与当前摄入量无关的常数
- 各个蔬菜的购买单价均是与他们相互间购买量无关的常数(一种蔬菜的购买不会影响另一种蔬菜的价格);各种蔬菜餐品提供的营养价值也是与相互的摄入量无关的常数(一种蔬菜的摄入并不会影响到另一种蔬菜的营养价值)
- 各种蔬菜的购买件数(杯数)可以是任意常数
在满足约束满足学生对基本营养的需求情况下,目标函数为投入最小费用。那么我们可以建立数学模型如下
假设苹果、香蕉、胡萝卜、枣汁、鸡蛋的购买量分别为
x
1
,
x
2
,
x
3
,
x
4
,
x
5
x_1,x_2,x_3,x_4,x_5
x1,x2,x3,x4,x5假设每天食材支出为
y
y
y元。那么我们可以建立目标函数为
M
i
n
y
=
10
x
1
+
15
x
2
+
5
x
3
+
60
x
4
+
8
x
5
Min \quad y=10x_1+15x_2+5x_3+60x_4+8x_5
Miny=10x1+15x2+5x3+60x4+8x5
同样我们也可以建立约束条件
S
.
T
.
0.3
∗
x
1
+
1.2
∗
x
2
+
0.7
∗
x
3
+
3.5
∗
x
4
+
5.5
∗
x
5
>
=
50
73
∗
x
1
+
96
∗
x
2
+
20253
∗
x
3
+
890
∗
x
4
+
279
∗
x
5
>
=
4000
9.6
∗
x
1
+
7
∗
x
2
+
19
∗
x
3
+
57
∗
x
4
+
22
∗
x
5
>
=
1000
x
i
>
=
0
,
i
=
1
,
2
,
3
,
4
,
5
S.T.\\ 0.3*x_1+1.2*x_2+0.7*x_3+3.5*x_4+5.5*x_5>=50\\ 73*x_1+ 96*x_2+ 20253*x_3+ 890*x_4+ 279*x_5>=4000\\ 9.6*x_1+ 7*x_2+ 19*x_3+ 57 *x_4+22*x_5>=1000\\ x_i>=0\qquad,i=1,2,3,4,5
S.T.0.3∗x1+1.2∗x2+0.7∗x3+3.5∗x4+5.5∗x5>=5073∗x1+96∗x2+20253∗x3+890∗x4+279∗x5>=40009.6∗x1+7∗x2+19∗x3+57∗x4+22∗x5>=1000xi>=0,i=1,2,3,4,5
我们利用python求解得到结果。
代码
#程序文件Pex5_4.py
import numpy as np
from cvxopt import matrix, solvers
c=matrix([10.,15,5,60,8]); #系数向量
A=matrix([[-0.3,-1.2,-0.7,-3.5,-5.5],
[-73,-96,-20253,-890,-279],
[-9.6,-7,-19,-57,-22],
[-1,0,0,0,0],
[0,-1,0,0,0],
[0,0,-1,0,0],
[0,0,0,-1,0],
[0,0,0,0,-1],
]).T#约束条件左端系数
b=matrix([-50.,-4000,-1000,0,0,0,0,0]);
sol=solvers.lp(c,A,b)
print("最优解为:\n",sol['x'])
print("最优值为:",sol['primal objective'])
'''
运行结果
pcost dcost gap pres dres k/t
0: 1.0965e+03 2.1410e+04 2e+02 2e-03 1e+03 1e+00
1: 4.1684e+02 1.1224e+05 1e+04 1e-02 7e+03 7e+01
2: 5.4657e+02 7.8504e+03 4e+02 8e-04 5e+02 1e+02
3: 4.5291e+02 3.3434e+03 2e+02 3e-04 2e+02 5e+01
4: 1.2653e+03 1.3518e+04 6e+03 1e-03 8e+02 4e+02
5: 2.8744e+03 2.1860e+04 1e+05 1e-03 9e+02 6e+03
6: 2.0199e+03 7.8726e+03 3e+04 6e-04 4e+02 3e+01
7: 1.3657e+03 2.6481e+03 4e+03 1e-04 8e+01 7e+01
8: 8.6825e+02 1.6685e+03 3e+03 8e-05 5e+01 8e+01
9: 4.3787e+02 5.4587e+02 3e+02 1e-05 6e+00 1e+01
10: 3.6601e+02 3.6929e+02 8e+00 3e-07 2e-01 3e-01
11: 3.6390e+02 3.6470e+02 2e+00 7e-08 4e-02 1e-01
12: 3.5529e+02 3.5833e+02 6e+01 1e-07 7e-02 2e+00
13: 3.1848e+02 3.2539e+02 1e+02 7e-08 4e-02 6e+00
14: 3.0232e+02 3.0603e+02 6e+01 4e-08 2e-02 3e+00
15: 2.7112e+02 2.7126e+02 3e+00 2e-09 1e-03 1e-01
16: 2.6938e+02 2.6938e+02 3e-02 2e-11 1e-05 1e-03
17: 2.6936e+02 2.6936e+02 3e-04 2e-13 1e-07 1e-05
18: 2.6936e+02 2.6936e+02 3e-06 1e-13 1e-09 1e-07
Optimal solution found.
最优解为:
[ 7.15e-09]
[ 8.94e-09]
[ 4.94e+01]
[ 1.67e-09]
[ 2.81e+00]
最优值为: 269.3602711693341
'''
或者:
import cvxpy as cp
from numpy import array
c=array([10.,15,5,60,8]) #定义目标向量
a=array([[-0.3,-1.2,-0.7,-3.5,-5.5],
[-73,-96,-20253,-890,-279],
[-9.6,-7,-19,-57,-22],
[-1,0,0,0,0],
[0,-1,0,0,0],
[0,0,-1,0,0],
[0,0,0,-1,0],
[0,0,0,0,-1]]) #定义约束矩阵
b=array([-50.,-4000,-1000,0,0,0,0,0]) #定义约束条件的右边向量
x=cp.Variable(5,integer=False) #定义5个决策变量,True为整数规划,False为实数
obj=cp.Minimize(c*x) #构造目标函数
cons=[a*x<=b, x>=0] #构造约束条件
prob=cp.Problem(obj, cons) #构建问题模型
prob.solve(solver='GLPK_MI',verbose =True) #求解问题
print("最优值为:",prob.value)
print("最优解为:\n",x.value)
'''
运行结果:
===============================================================================
CVXPY
v1.1.12
===============================================================================
(CVXPY) Aug 11 05:13:02 PM: Your problem has 5 variables, 2 constraints, and 0 parameters.
(CVXPY) Aug 11 05:13:02 PM: It is compliant with the following grammars: DCP, DQCP
(CVXPY) Aug 11 05:13:02 PM: (If you need to solve this problem multiple times, but with different data, consider using parameters.)
(CVXPY) Aug 11 05:13:02 PM: CVXPY will first compile your problem; then, it will invoke a numerical solver to obtain a solution.
-------------------------------------------------------------------------------
Compilation
-------------------------------------------------------------------------------
(CVXPY) Aug 11 05:13:02 PM: Compiling problem (target solver=GLPK_MI).
(CVXPY) Aug 11 05:13:02 PM: Reduction chain: Dcp2Cone -> CvxAttr2Constr -> ConeMatrixStuffing -> GLPK_MI
(CVXPY) Aug 11 05:13:02 PM: Applying reduction Dcp2Cone
(CVXPY) Aug 11 05:13:02 PM: Applying reduction CvxAttr2Constr
(CVXPY) Aug 11 05:13:02 PM: Applying reduction ConeMatrixStuffing
(CVXPY) Aug 11 05:13:02 PM: Applying reduction GLPK_MI
(CVXPY) Aug 11 05:13:02 PM: Finished problem compilation (took 2.800e-02 seconds).
-------------------------------------------------------------------------------
Numerical solver
-------------------------------------------------------------------------------
(CVXPY) Aug 11 05:13:02 PM: Invoking solver GLPK_MI to obtain a solution.
-------------------------------------------------------------------------------
Summary
-------------------------------------------------------------------------------
(CVXPY) Aug 11 05:13:02 PM: Problem status: optimal
(CVXPY) Aug 11 05:13:02 PM: Optimal value: 2.694e+02
(CVXPY) Aug 11 05:13:02 PM: Compilation took 2.800e-02 seconds
(CVXPY) Aug 11 05:13:02 PM: Solver (including time spent in interface) took 1.800e-02 seconds
最优值为: 269.36026936026735
最优解为:
[-0. -0. 49.38271605 -0. 2.80583614]
'''
最终得到结果上述蔬菜分别购买:
0. -0. 49.38271605 0. 2.80583614
最终支出为
269.36026936026735
269.36026936026735
269.36026936026735(角)
灵敏性分析
接下来我们介绍几个概念:
- 影子价格:影子价格通俗理解
接下来我们求解下面两问,最终我们根据上述求解结果分别购买49.38的胡萝卜和2.80的鸡蛋。我们可以算出若按照这个方案安排,营养含量如下图
由于python灵敏度性不太完善,影子价格的计算不太方便,因此对于后两问需重新调整模型进行计算(后续可考虑Scipy Linprog中的影子价格)