This tutorial includes everything you need to set up Decision Optimization engines, build a mathematical programming model, then use the progress listeners to monitor progress, capture intermediate solutions and stop the solve on your own criteria.
When you finish this tutorial, you'll have a foundational knowledge of Prescriptive Analytics.
This notebook is part of Prescriptive Analytics for Python
It requires either an installation of CPLEX Optimizers or it can be run on IBM Watson Studio Cloud (Sign up for a free IBM Cloud account and you can start using Watson Studio Cloud right away).
Table of contents:
Prescriptive analytics (Decision Optimization) technology recommends actions that are based on desired outcomes. It takes into account specific scenarios, resources, and knowledge of past and current events. With this insight, your organization can make better decisions and have greater control of business outcomes.
Prescriptive analytics is the next step on the path to insight-based actions. It creates value through synergy with predictive analytics, which analyzes data to predict future outcomes.
Prescriptive analytics takes that insight to the next level by suggesting the optimal way to handle that future situation. Organizations that can act fast in dynamic conditions and make superior decisions in uncertain environments gain a strong competitive advantage.
With prescriptive analytics, you can:
In [1]:
from docplex.mp.environment import Environment
Environment().print_information()
You start with a very simple infeasible model: you have three variables, each of which is greater than the previous one, in a cyclic fashion. Of course this leads to an infeasible model
In [2]:
from docplex.mp.model import Model
def build_infeasible_cyclic_model3():
m = Model(name='cyclic3')
x,y,z = m.continuous_var_list(keys=['x', 'y', 'z'], name=str)
m.add( y >= x+1, name="y_gt_x")
m.add( z >= y+1, name="z_gt_y")
m.add( x >= z+1, name="x_gt_z")
# add another constraint, should noever appear in conflicts
m.add(x + y + z <= 33)
return m
cycle3 = build_infeasible_cyclic_model3()
cycle3.print_information()
As expected, the model is infeasible.
In [3]:
s = cycle3.solve(log_output=True)
assert s is None
print("the model is infeasible")
First, you can use the Conflict refiner on this model. The conflict refiner computes a minimal cluster of constraints, which causes the infeasibility.
Using the conflict refiner requires the following steps:
ConflictRefiner instanceThe output is an object of type ConflictRefinerResults which holds all information about the minimal conflict.
Displaying this result object lists all modeling objects which belong to the minimal conflict.
In [4]:
from docplex.mp.conflict_refiner import ConflictRefiner
cr = ConflictRefiner()
crr = cr.refine_conflict(cycle3, display=True)
Another way to handle infeasibilities is to use the relaxer (class docplex.mp.relaxer.Relaxer). The relaxer tries to find a minimal feasible relaxation of the model, by relaxing certain constraints.
For example, a constraint x == 1 can be relaxed with a slack of 1 to accept a value of 2 for x.
the relaxer tries to minimize the total value of slack in finding a feasible relaxation.
In [5]:
from docplex.mp.relaxer import Relaxer
rx = Relaxer()
rs = rx.relax(cycle3)
rx.print_information()
rs.display()
The relaxer has relaxed one constraint ( x>= z+1) by 3, and found a solution wiht x=0, x=1, z=2, breaking the cyclic chain of constraints.
Unlike the conflict refiner, the relaxer provides a relaxed solution to the initial model, with minimal slack. But there's more: the relaxer can also search for the best business objective, once it has found the minimal slack.
To illustrate, add an objective to your model, try to minimize z, and see what happens:
In [6]:
# retrieve the z variable using Model.get_var_by_name()
z = cycle3.get_var_by_name('z')
assert z
cycle3.minimize(z)
rs = rx.relax(cycle3)
rx.print_information()
rs.display()
Th relaxed solution has changed, finding a minimum objective of 0 for z, and the relaxations have also changed: now two constraints are relaxed, but the total absolute slack remains unchanged, equal to 3.
To summarize, the relaxer finds a relaxed solution in two steps:
In this case, the minimal conflict contains the three cyclic constraints, but not the fourth (x+y+z <= 30).
The model aims at minimizing the production cost for a number of products while satisfying customer demand. Production is constrained by the company's resources.
The model first declares the products and the resources. The data consists of the description of the products (the demand, the inside and outside costs, and the resource consumption) and the capacity of the various resources.
The variables for this problem are the quantities produced for each products.
Of course, this model is naive, un-realistic and not robust, as it fails as soon as demand is not satisfied. But this is excatly why it is well adapted to show you how to repair infeasible models.
First define the data:
In [7]:
# costs are stored in a dict from product names to cost
COSTS = {"kluski": 0.6, "capellini": 0.8, "fettucine": 0.3}
# demands are stored in a dict from product names to demands
DEMANDS = {"kluski": 10, "capellini": 20, "fettucine": 30}
# resources are stored as a dict of resource name to resource capacity
RESOURCES = {"flour": 40, "eggs": 60}
CONSUMPTIONS = {("kluski", "flour"): 0.5,
("kluski", "eggs"): 0.2,
("capellini", "flour"): 0.4,
("capellini", "eggs"): 0.4,
("fettucine", "flour"): 0.3,
("fettucine", "eggs"): 0.6}
In [8]:
from docplex.mp.model import Model
import six
def build_production_problem(costs, resources, consumptions, demands, **kwargs):
products = [p for p, _ in six.iteritems(costs)]
mdl = Model(name='pasta_production', **kwargs)
# --- decision variables ---
mdl.q_vars = mdl.continuous_var_dict(products, name="q")
# --- constraints ---
# demand satisfaction
mdl.add_constraints((mdl.q_vars[p] >= demands[p], 'ct_demand_%s' % p) for p in products)
# --- resource capacity ---
mdl.add_constraints((mdl.sum(mdl.q_vars[p] * consumptions[p, res] for p in products) <= cap,
'ct_res_%s' % res) for res, cap in six.iteritems(resources))
# --- objective ---
mdl.minimize(mdl.dotf(mdl.q_vars, lambda p: costs[p]))
return mdl
In [9]:
pasta1 = build_production_problem(COSTS, RESOURCES, CONSUMPTIONS, DEMANDS)
pasta1.print_information()
This default model is feasible, solve the model
In [10]:
s1 = pasta1.solve()
s1.display()
In [11]:
demands2 = {p: 2*d for p, d in six.iteritems(DEMANDS)}
pasta2 = build_production_problem(COSTS, RESOURCES, CONSUMPTIONS, demands2)
s2 = pasta2.solve()
if s2 is None:
print("!! Pasta production with double demand is impossible")
Now double the demand.
In [23]:
from docplex.mp.conflict_refiner import ConflictRefiner
crr = ConflictRefiner().refine_conflict(pasta2, display=True)
crr.display()
Not surprisingly, you can see that the conflict involves all three demands but also the flour resource constraint, but we have no idea which quantity demands cannot be satisfied.
In [13]:
from docplex.mp.relaxer import Relaxer
# create an instance of relaxer
rx = Relaxer()
rs = rx.relax(pasta2)
rx.print_information()
rs.display()
The relaxer managed to satisfy all demands by relaxing the flour constraint by an amount of 4. What does this mean?
To explain, first remember what this flour constraint was all about. You can use the Model.get_constraint_by_name method to retrieve a constraint from its name.
In [14]:
# get back the constraint from its name
ctf = pasta2.get_constraint_by_name("ct_res_flour")
assert ctf is not None
print(str(ctf))
Now you can see what the left hand side evaluates in the relaxed solution
In [15]:
ctf.lhs.solution_value
Out[15]:
This explains the relaxation of 4 for the flour resource constraint
It might well happen that the relaxation found by the relaxer does not make sense in real world. For example, in our production example, resource constraints can be impossible to relax, but demands could be.
This is where priorities enter the game. By setting priorities, users can control how the relaxer chooses constraints to relax. In the following code, you can set a HIGH priority to resource constraints (you could even make them mandatory) and a LOW priority to demand constraints.
In [16]:
from docplex.mp.basic import Priority
for ctr in pasta2.find_matching_linear_constraints("ct_res"):
ctr.priority = Priority.HIGH
for ctd in pasta2.find_matching_linear_constraints("ct_dem"):
ctd.priority = Priority.LOW
rx2 = Relaxer()
rs2 = rx2.relax(pasta2)
rx2.print_information()
rs2.display()
In this new relaxed solution, all resource constraints are satisfied, but one demand is not: kluski demand has an unfilled quantity of 8.
Setting constraint priorities explicitly is the most basic way to control relaxation, but there are others. A function can be used: the relaxer will call the function for each constraint to determine its priority. Possible values are:
Constraints with higher priority are less likely to be relaxed than constraints with lower priorities. Still, relaxation of a high-priority constraint cannot be ruled out, if it is the only way to provide a relaxed solution.
In this section, you can see how to use a function to compute the priority of a constraint. The function must take a constraint and return a priority (an enumerated type, see docplex.mp.basic.Priority
First, you reset all priorities to None (the default).
In [17]:
# reset all priorities
for c in pasta2.iter_constraints():
c.priority = None
In [18]:
# define the constraint -> priority function
def map_priority(ct):
ctname = ct.name
if not ctname:
return Priority.MANDATORY
elif "ct_res" in ctname:
return Priority.HIGH
elif "ct_dem" in ctname:
return Priority.LOW
else:
# will not be relaxed
return Priority.MANDATORY
# create a new instance of Relaxer with this function.
rx3 = Relaxer(prioritizer=map_priority)
# use it to relax pasta2 model
rx3.relax(pasta2)
# display relaxation.
rx3.print_information()
As expected, you get the same result as with the explicit priorities: an unsatisfied demand of 8 for kluski.
Note that relaxer can also accept a dictionary of constraints to priorities.
Now that you know about setting priorities, you can understand the default behavior of the relaxer class: for each constraint, you use either its explicit priority (if set) or the defaultMEDIUM priority.
If no priority has been set, all constraints are considered equally relaxable.
In [19]:
m4 = Model(name='m4')
ijs = m4.integer_var_list(keys=["i", "j", "k", "l"], name =str, lb = 1)
m4.add(m4.sum(ijs) <= 2)
s4 = m4.solve()
assert s4 is None
ConflictRefiner().refine_conflict(m4, display=True);
The resulting conflict contains the sum constraint and the three lower bounds.
In [20]:
r4 = Relaxer()
r4.relax(m4)
r4.print_information()
Changing variable lower bounds to constraints allows the relaxer to take them into account. You can set a LOW priority to lower bounds, so we expect the default relaxer to relax them before the sum constraint, which will be considered with the default MEDIUM priority.
In [21]:
for v in m4.iter_variables():
v.lb = 0
clb = m4.add(v >= 1, "{0}_ge_1".format(v.name))
clb.priority = Priority.LOW
In [22]:
r4 = Relaxer()
r4.relax(m4)
r4.print_information()
As expected, the relaxer has relaxed two lower bound constraints, but not the sum constraint.
You have learned how to use both the conflict refiner and the relaxer, and the differences between them
In constrast, the relaxer provides a relaxed solution, and indicates which constraints are relaxed, and with what quantities. It does not consider variables bounds. It requires a mapping from constraints to Priority objects, which can take many forms: a function, a dictionary,...
Copyright © 2017-2019 IBM. Sample Materials.
In [ ]: