import pulp
def dietlp(price,nutrition,requirements, nmin = 0, nmax=100):
diet = pulp.LpProblem('Diet Problem', pulp.LpMinimize)
x = pulp.LpVariable.dicts("x", list(price.index), nmin, nmax, pulp.LpInteger)
diet += pulp.lpSum([price.ix[idx,'price']*x[idx] for idx in list(price.index)])
for n in nutrition.columns:
diet += pulp.lpSum([nutrition.ix[idx,n]*x[idx] for idx in list(price.index)]) >= requirements.ix[n,'n_min'], n + '_min'
diet += pulp.lpSum([nutrition.ix[idx,n]*x[idx] for idx in list(price.index)]) <= requirements.ix[n,'n_max'], n + '_max'
return diet
diet = dietlp(price,nutrition,requirements,0,10)
def solvelp(lp):
lp.solve();
obj = pd.DataFrame(
[[lp.objective.name,
pulp.LpSenses[lp.sense],
pulp.LpStatus[lp.status],
pulp.value(lp.objective)]],
columns=["Objective", "Sense", "Status", "Value"])
var = pd.DataFrame(
[[x.name, x.cat, x.lowBound, x.upBound, x.varValue] for x in lp.variables()],
columns=["Name","Category","Lower Bound","Upper Bound","Value"])
con = pd.DataFrame(
[[c.name,
pulp.LpConstraintSenses[c.sense],
c.getLb(),
c.getUb(),
c.slack,
pulp.value(c) - c.constant] for (k,c) in lp.constraints.items()],
columns=["Name", 'Sense', "Lower Bound", "Upper Bound", "Slack", "Value"])
return {'Objective':obj, 'Variables':var, 'Constraints':con}
soln = solvelp(diet)
from IPython.display import display,HTML
display(HTML('<h3>Model</h3>'))
display(diet)
for k in soln.keys():
display(HTML('<h3>' + k + '</h3>'))
display(soln[k])
Diet Problem: MINIMIZE 3.19*x_BEEF + 2.59*x_CHK + 2.29*x_FSH + 2.89*x_HAM + 1.89*x_MCH + 1.99*x_MTL + 1.99*x_SPG + 2.49*x_TUR + 0.0 SUBJECT TO A_min: 60 x_BEEF + 8 x_CHK + 8 x_FSH + 40 x_HAM + 15 x_MCH + 70 x_MTL + 25 x_SPG + 60 x_TUR >= 700 A_max: 60 x_BEEF + 8 x_CHK + 8 x_FSH + 40 x_HAM + 15 x_MCH + 70 x_MTL + 25 x_SPG + 60 x_TUR <= 10000 B1_min: 10 x_BEEF + 20 x_CHK + 15 x_FSH + 35 x_HAM + 15 x_MCH + 15 x_MTL + 25 x_SPG + 15 x_TUR >= 700 B1_max: 10 x_BEEF + 20 x_CHK + 15 x_FSH + 35 x_HAM + 15 x_MCH + 15 x_MTL + 25 x_SPG + 15 x_TUR <= 10000 B2_min: 15 x_BEEF + 20 x_CHK + 10 x_FSH + 10 x_HAM + 15 x_MCH + 15 x_MTL + 15 x_SPG + 10 x_TUR >= 700 B2_max: 15 x_BEEF + 20 x_CHK + 10 x_FSH + 10 x_HAM + 15 x_MCH + 15 x_MTL + 15 x_SPG + 10 x_TUR <= 10000 C_min: 20 x_BEEF + 10 x_FSH + 40 x_HAM + 35 x_MCH + 30 x_MTL + 50 x_SPG + 20 x_TUR >= 700 C_max: 20 x_BEEF + 10 x_FSH + 40 x_HAM + 35 x_MCH + 30 x_MTL + 50 x_SPG + 20 x_TUR <= 10000 VARIABLES 0 <= x_BEEF <= 10 Integer 0 <= x_CHK <= 10 Integer 0 <= x_FSH <= 10 Integer 0 <= x_HAM <= 10 Integer 0 <= x_MCH <= 10 Integer 0 <= x_MTL <= 10 Integer 0 <= x_SPG <= 10 Integer 0 <= x_TUR <= 10 Integer
Objective | Sense | Status | Value | |
---|---|---|---|---|
0 | OBJ | Minimize | Optimal | 95.56 |
Name | Category | Lower Bound | Upper Bound | Value | |
---|---|---|---|---|---|
0 | x_BEEF | Integer | 0 | 10 | 2 |
1 | x_CHK | Integer | 0 | 10 | 10 |
2 | x_FSH | Integer | 0 | 10 | 2 |
3 | x_HAM | Integer | 0 | 10 | 0 |
4 | x_MCH | Integer | 0 | 10 | 10 |
5 | x_MTL | Integer | 0 | 10 | 10 |
6 | x_SPG | Integer | 0 | 10 | 10 |
7 | x_TUR | Integer | 0 | 10 | 0 |
Name | Sense | Lower Bound | Upper Bound | Slack | Value | |
---|---|---|---|---|---|---|
0 | A_min | >= | 700 | NaN | -616 | 1316 |
1 | A_max | <= | NaN | 10000 | 8684 | 1316 |
2 | B1_min | >= | 700 | NaN | -100 | 800 |
3 | B1_max | <= | NaN | 10000 | 9200 | 800 |
4 | B2_min | >= | 700 | NaN | -0 | 700 |
5 | B2_max | <= | NaN | 10000 | 9300 | 700 |
6 | C_min | >= | 700 | NaN | -510 | 1210 |
7 | C_max | <= | NaN | 10000 | 8790 | 1210 |
ax = soln['Variables'].ix[:,'Value'].plot(kind='bar')
ax.set_xticklabels(list(soln['Variables'].ix[:,'Name']));
nmax = arange(0,60,1)
obj = zeros(size(nmax))
for n in nmax:
soln = solvelp(dietlp(price,nutrition,requirements,0,n));
obj[n-1] = soln['Objective'].ix[0,'Value'] if (soln['Objective'].ix[0,'Status'] == 'Optimal') else NaN
plot(nmax,obj,'r.')
xlabel('n_max')
ylabel('Cost')
<matplotlib.text.Text at 0x10a6dead0>
%%script glpsol -m /dev/stdin -o /dev/stdout --out output
set COMPONENTS;
set PRODUCTS;
set CRUDES;
param feed {COMPONENTS,CRUDES};
param avail {CRUDES};
data;
set COMPONENTS :=
Butane
Naptha
Gasoil
Residue ;
set CRUDES :=
Crude_1
Crude_2 ;
param feed : Crude_1 Crude_2 :=
Butane 0.03 0.05
Naphtha 0.15 0.20
Gasoil 0.40 0.35
Residue 0.15 0.10 ;
param avail :=
Crude_1 250000
Crude_2 400000 ;
end;
print output
GLPSOL: GLPK LP/MIP Solver, v4.52 Parameter(s) specified in the command line: -m /dev/stdin -o /dev/stdout Reading model section from /dev/stdin... Reading data section from /dev/stdin... 30 lines were read Model has been successfully generated GLPK Simplex Optimizer, v4.52 0 rows, 0 columns, 0 non-zeros ~ 0: obj = 0.000000000e+00 infeas = 0.000e+00 OPTIMAL SOLUTION FOUND Time used: 0.0 secs Memory used: 0.1 Mb (73651 bytes) Writing basic solution to `/dev/stdout'... Problem: stdin Rows: 0 Columns: 0 Non-zeros: 0 Status: OPTIMAL Objective: 0 (MINimum) No. Row name St Activity Lower bound Upper bound Marginal ------ ------------ -- ------------- ------------- ------------- ------------- No. Column name St Activity Lower bound Upper bound Marginal ------ ------------ -- ------------- ------------- ------------- ------------- Karush-Kuhn-Tucker optimality conditions: KKT.PE: max.abs.err = 0.00e+00 on row 0 max.rel.err = 0.00e+00 on row 0 High quality KKT.PB: max.abs.err = 0.00e+00 on row 0 max.rel.err = 0.00e+00 on row 0 High quality KKT.DE: max.abs.err = 0.00e+00 on column 0 max.rel.err = 0.00e+00 on column 0 High quality KKT.DB: max.abs.err = 0.00e+00 on row 0 max.rel.err = 0.00e+00 on row 0 High quality End of output