how to implement objective function for unconstrained binary quadratic problem
Steve Reinhardt
entities = ['CA','OR','WA','MT','ID','NV','WY','UT','AZ']
model.entities = Set(initialize=entities, doc='Geographic entities')
colors = ['red','blue','green','yellow']
model.colors = Set(initialize=colors, doc='Colors')
entityColors = [it0+it1 for it0 in entities for it1 in colors]
model.entityColors = Set(initialize=entityColors, ordered=True, doc='entityColors')
model.qubo = Param(model.entityColors,model.entityColors, initialize=quboInit, doc='QUBO')
def objective_rule(model):
vm = [0]*len(model.entityColors)
for (i,key) in zip(range(len(model.entityColors)),model.entityColors.keys()):
vm[i] = sum([it0*it1 for it0,it1 in zip(model.assignedColors.values(),model.qubo[key,:])])
ret = sum([it0*it1 for it0,it1 in zip(vm, model.assignedColors.values())])
print("in objective_rule, returning, ret="+str(ret))
return ret
Siirola, John D
The Abstract/Concrete model thing is a stylistic decision. I can make anything work in either mode, although there are some tasks that are trickier in an Abstract model than they are in a Concrete model. The two modes are really to ease the transition to Pyomo for two different groups of users (folks coming from AMPL tend to resonate with AbstractModels, whereas folks coming from GAMS or AIMMS tend to resonate more with ConcreteModels).
I don’t understand what you mean by the objective being fixed. If there are no variables in your objective and the problem is unconstrained, what is there to optimize?
The best thing I can recommend is to re-read the online documentation and take a look at some of the examples in the PyomoGallery, paying particular attention to how Sets and Params are used in the examples.
The easy way to write your model is:
model.ENTITIES = Set(initialize=['CA','OR','WA','MT','ID','NV','WY','UT','AZ'])
model.COLORS = Set(initialize=['red','blue','green','yellow'])
model.ENTITY_COLORS = model.ENTITIES * model.COLORS
model.qubo = Param(model.ENTITY_COLORS, model.ENTITY_COLORS, initialize=quboInit)
model.assignedColors = Var(model.ENTITY_COLORS, within=Binary)
def objective_rule(m):
return sum( m.qubo[i,i]*m.assignedColors[i] for i in m.ENTITY_COLORS ) \
+ sum( m.qubo[i,j]*m.assignedColors[i]*m.assignedColors[j] \
for i in m.ENTITY_COLORS for j in m.ENTITY_COLORS if i != j)
Note that you will have to modify your quboInit function to take 5 arguments (model, e1, c1, e2, c2).
john
--
You received this message because you are subscribed to the Google Groups "Pyomo Forum" group.
To unsubscribe from this group and stop receiving emails from it, send an email to
pyomo-forum...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
Steve Reinhardt
On Tuesday, April 11, 2017 at 10:45:49 AM UTC-5, John wrote:
The Abstract/Concrete model thing is a stylistic decision. I can make anything work in either mode, although there are some tasks that are trickier in an Abstract model than they are in a Concrete model. The two modes are really to ease the transition to Pyomo for two different groups of users (folks coming from AMPL tend to resonate with AbstractModels, whereas folks coming from GAMS or AIMMS tend to resonate more with ConcreteModels).
I don’t understand what you mean by the objective being fixed. If there are no variables in your objective and the problem is unconstrained, what is there to optimize?
The best thing I can recommend is to re-read the online documentation and take a look at some of the examples in the PyomoGallery, paying particular attention to how Sets and Params are used in the examples.
The easy way to write your model is:
model.ENTITIES = Set(initialize=['CA','OR','WA','MT','ID','NV','WY','UT','AZ'])
model.COLORS = Set(initialize=['red','blue','green','yellow'])
model.ENTITY_COLORS = model.ENTITIES * model.COLORS
model.qubo = Param(model.ENTITY_COLORS, model.ENTITY_COLORS, initialize=quboInit)
model.assignedColors = Var(model.ENTITY_COLORS, within=Binary)
def objective_rule(m):
return sum( m.qubo[i,i]*m.assignedColors[i] for i in m.ENTITY_COLORS ) \
+ sum( m.qubo[i,j]*m.assignedColors[i]*m.assignedColors[j] \
for i in m.ENTITY_COLORS for j in m.ENTITY_COLORS if i != j)
Note that you will have to modify your quboInit function to take 5 arguments (model, e1, c1, e2, c2).
KeyError: "Error accessing indexed component: Index '('WA', 'blue', 'WA', 'blue', 'WA', 'blue', 'WA', 'blue')' is not valid for array component 'qubo'"
ipdb> m.ENTITY_COLORS.keys()[:10]
[('WA', 'blue', 'WA', 'blue'), ('WA', 'blue', 'WA', 'green'), ('WA', 'blue', 'WA', 'yellow'), ('WA', 'blue', 'WA', 'red'), ('WA', 'blue', 'UT', 'blue'), ('WA', 'blue', 'UT', 'green'), ('WA', 'blue', 'UT', 'yellow'), ('WA', 'blue', 'UT', 'red'), ('WA', 'blue', 'CA', 'blue'), ('WA', 'blue', 'CA', 'green')]
I expect ENTITY_COLORS' keys to be (e.g.) ('WA', 'blue'), i.e., 2-dimensional, not 4-dimensional. Am I missing something?
Thanks much! Steve
Siirola, John D
Steve,
> My understanding is that only AbstractModels support loading data from a file to reuse the same model. True?
Sort of, but not really. It is true that AbstractModel supports the “model.create_instance(‘foo.dat’)” syntax. But it is not true that it reuses the model (It actually copies the entire abstract model, then switches the flag on the new model to be concrete, then adds all the data). For ConcreteModel, you would first load the data from the file (in any way you want), then create the ConcreteModel using the data. If you were going to create several concrete models within a single script / invocation, the easiest thing is to declare your model within a function that takes the necessary data as argument(s) and returns the model. The exception to all of this is if your data is in an AMPL-style DAT file. It turns out that DAT files cannot be reliably parsed without a model declaration (you actually need information from the model structure to resolve certain ambiguities with the
> OK, taking the ConcreteModel path, I changed the code […] I expect ENTITY_COLORS' keys to be (e.g.) ('WA', 'blue'), i.e., 2-dimensional, not 4-dimensional
As I write it, it should be 2D. You might be running into something weird with ipython or jupyter where old state is being persisted? Or do you have two models floating around (one named “model” and one named “m”)? We would need to see the whole model to say anything more definitive.
john
From: pyomo...@googlegroups.com [mailto:pyomo...@googlegroups.com] On Behalf Of Steve Reinhardt
Sent: Tuesday, April 11, 2017 10:39 AM
To: Pyomo Forum <pyomo...@googlegroups.com>
--
Steve Reinhardt
On Tuesday, April 11, 2017 at 2:27:10 PM UTC-5, John wrote:
Steve,
> My understanding is that only AbstractModels support loading data from a file to reuse the same model. True?
Sort of, but not really. It is true that AbstractModel supports the “model.create_instance(‘foo.dat’)” syntax. But it is not true that it reuses the model (It actually copies the entire abstract model, then switches the flag on the new model to be concrete, then adds all the data). For ConcreteModel, you would first load the data from the file (in any way you want), then create the ConcreteModel using the data. If you were going to create several concrete models within a single script / invocation, the easiest thing is to declare your model within a function that takes the necessary data as argument(s) and returns the model. The exception to all of this is if your data is in an AMPL-style DAT file. It turns out that DAT files cannot be reliably parsed without a model declaration (you actually need information from the model structure to resolve certain ambiguities with the
> OK, taking the ConcreteModel path, I changed the code […] I expect ENTITY_COLORS' keys to be (e.g.) ('WA', 'blue'), i.e., 2-dimensional, not 4-dimensional
As I write it, it should be 2D. You might be running into something weird with ipython or jupyter where old state is being persisted? Or do you have two models floating around (one named “model” and one named “m”)? We would need to see the whole model to say anything more definitive.
ipdb> m.qubo.keys()[:5]
[('WA', 'red', 'AZ', 'blue'), ('OR', 'red', 'CA', 'red'), ('ID', 'yellow', 'MT', 'blue'), ('CA', 'blue', 'MT', 'green'), ('MT', 'blue', 'CA', 'blue')]
ipdb> m.ENTITY_COLORS.keys()[:5]
[('WA', 'blue', 'WA', 'blue'), ('WA', 'blue', 'WA', 'green'), ('WA', 'blue', 'WA', 'yellow'), ('WA', 'blue', 'WA', 'red'), ('WA', 'blue', 'UT', 'blue')]
#!/usr/bin/env python
# -*- coding: utf-8 -*-
from __future__ import division
from pyomo.environ import *
PENALTY_DEFAULT = -1
PENALTY_MULTIPLE_COLORS = 2
PENALTY_SAME_COLOR = 1
entitiesLut = {
'AZ' : 0,
'CA' : 1,
'ID' : 2,
'MT' : 3,
'NV' : 4,
'OR' : 5,
'UT' : 6,
'WA' : 7,
'WY' : 8,
}
edgesDict = {
('AZ', 'CA') : True,
('AZ', 'NV') : True,
('AZ', 'UT') : True,
('CA', 'AZ') : True,
('CA', 'NV') : True,
('CA', 'OR') : True,
('ID', 'NV') : True,
('ID', 'WY') : True,
('ID', 'UT') : True,
('ID', 'WA') : True,
('ID', 'MT') : True,
('MT', 'ID') : True,
('MT', 'WY') : True,
('NV', 'AZ') : True,
('NV', 'CA') : True,
('NV', 'ID') : True,
('NV', 'UT') : True,
('OR', 'CA') : True,
('OR', 'WA') : True,
('OR', 'ID') : True,
('OR', 'UT') : True,
('UT', 'AZ') : True,
('UT', 'ID') : True,
('UT', 'NV') : True,
('UT', 'OR') : True,
('UT', 'WY') : True,
('WA', 'OR') : True,
('WA', 'ID') : True,
('WY', 'ID') : True,
('WY', 'MT') : True,
('WY', 'UT') : True,
}
model = ConcreteModel(name='mapColoring')
entities = ['CA','OR','WA','MT','ID','NV','WY','UT','AZ']
nEntities = len(entities)
model.ENTITIES = Set(initialize=entities, doc='Geographic entities')
colors = ['red','blue','green','yellow']
nColors = len(colors)
model.COLORS = Set(initialize=colors, doc='Colors')
model.ENTITY_COLORS = model.ENTITIES * model.COLORS
def quboInit(model,iEntity, iColor, jEntity, jColor):
ret = 0
if iEntity==jEntity:
if iColor!=jColor:
ret = PENALTY_MULTIPLE_COLORS # same entity, different color
else:
ret = -1 # same entity, same color
else:
# different entities, same color, adjacent
if iColor==jColor and (iEntity,jEntity) in edgesDict:
ret = PENALTY_SAME_COLOR
return ret
model.qubo = Param(model.ENTITY_COLORS,model.ENTITY_COLORS, initialize=quboInit, doc='QUBO')
def aC_init(model,ign0,ign1,ign2,ign3):
return False
model.assignedColors = Var(model.ENTITY_COLORS,initialize=aC_init,within=Boolean,doc='assigned colors')
def objective_rule(model):
ret = sum( model.qubo[i,i]*model.assignedColors[i] for i in model.ENTITY_COLORS ) \
+ sum( model.qubo[i,j]*model.assignedColors[i]*model.assignedColors[j] \
for i in model.ENTITY_COLORS for j in model.ENTITY_COLORS if i != j)
print("in objective_rule, returning, ret="+str(ret))
return ret
model.objective = Objective(rule=objective_rule, sense=minimize, doc='Define objective function')
def pyomo_postprocess(options=None, instance=None, results=None):
model.assignedColors.display()
# This is an optional code path that allows the script to be run outside of
# pyomo command-line. For example: python transport.py
if __name__ == '__main__':
# This emulates what the pyomo command-line tools does
from pyomo.opt import SolverFactory
import pyomo.environ
opt = SolverFactory("glpk")
results = opt.solve(model) #,solver_io='nl')
#sends results to stdout
results.write()
print("\nDisplaying Solution\n" + '-'*60)
pyomo_postprocess(None, None, results)
Siirola, John D
Ok. You have found a bug with Virtual Sets. A really nasty bug … I vaguely remember something like this popping up a year or two ago, and I thought we squashed it. Until it is fixed, the workaround is to not use a virtual set for ENTITY_COLORS. Instead, declare the ENTITY_COLORS Set as a non-virtual set (i.e.., a set that explicitly contains members):
model.ENTITY_COLORS = Set(initialize=model.ENTITIES * model.COLORS)
…then things will work … up to the ac_init method, which should only take 3 arguments and not 5 … and the fact that glpk doesn’t handle MIQPs (but Gurobi or Cplex do).
The bug has been reported as #134 (https://github.com/Pyomo/pyomo/issues/134)
john
From: pyomo...@googlegroups.com [mailto:pyomo...@googlegroups.com]
On Behalf Of Steve Reinhardt
Sent: Tuesday, April 11, 2017 3:18 PM
To: Pyomo Forum <pyomo...@googlegroups.com>
Subject: Re: [EXTERNAL] how to implement objective function for unconstrained binary quadratic problem
Hi John,
--
Steve Reinhardt
Re: quadratic solver: ipopt is not the simplest to build/install for macOS. Any others you'd suggest?
Thanks much! Steve
Santiago rodriguez
Siirola, John D
Steve,
There are also exact linearizations of QUBO that can be solved by any mixed integer linear solver (e.g., Gurobi, Cplex, SCIP, CBC, or gplk): e.g., replace all bilinear terms x_i*x_j in the objective with an ancillary variable y_ij \in [0,1], and then add the following constraints:
y_ij <= x_i
y_ij <= x_j
y_ij >= x_i + x_j - 1
john