This notebook shows how to solve the diet problem example shown in class, which was taken from R. Larson, Elementary Linear Algebra, Section 9.2:, Example 5.

In [1]:
from pyomo.environ import *
from pyomo.opt import *
opt = solvers.SolverFactory("glpk")

The index sets are lists, and all coefficient vectors/matrices are Python dictionaries with members of the index set as keys:

In [2]:
F = [1,2]
N = ['Calories', 'VitA', 'VitC']
c = {1:12, 2:15}
a = {(1,'Calories'):60, (1,'VitA'):12, (1,'VitC'):10,
     (2,'Calories'):65, (2,'VitA'):6,  (2,'VitC'):30}
r = {'Calories':300, 'VitA':36, 'VitC':90}

model = ConcreteModel()

model.x = Var(F, within=NonNegativeReals)

model.z = Objective(expr = sum(c[i]*model.x[i] for i in F), sense=minimize)

To define a model which works unchanged for any number of nutrient constraints, we write a constraint rule which gets called for each member of the nutrient index set.

In [3]:
def nutrition_rule (model, n):
    return sum(a[i,n]*model.x[i] for i in F) >= r[n]

model.c = Constraint(N, rule=nutrition_rule)

results = opt.solve(model)
In [4]:
model.x.get_values()
Out[4]:
{1: 2.73913043478261, 2: 2.08695652173913}
In [5]:
model.z.expr()
Out[5]:
64.17391304347828