Access dual variables

79 views
Skip to first unread message

Pancho

unread,
May 14, 2012, 1:56:05 PM5/14/12
to coopr forum
Hi again,

I'm trying to implement Benders Decomposition, based on the example
included with the installation. Performing a sanity check on the code,
I found a really odd behavior of dual variables when solving
subproblems iteratively.
For example, here I'm not even doing Benders, I'm just giving a
subproblem different values of the complicating variables (variables
from master) and checking the duals.

In the subproblem I'm imposing

X=Xparam (dual)

and getting the dual of this constraint using
(inst.Imposeconstraints[r].dual).

The main code is:

max_iterations = 10

# the main benders loop.
for i in xsequence(max_iterations):

print "\nIteration=",i

for instance in sub_insts:
for r in mstr_inst.R:
if i is 5:
instance.Xparam[r] = 0
else:
instance.Xparam[r] = 10000
instance.preprocess() # re-generate constraints based on new
parameter data

# solve the subproblems.
solve_all_instances(solver_manager, solver, sub_insts,
suffixes=['rc', 'dual'])

For X=0, I get the right duals, however, for X=10000 I get the right
duals (all None, or zero) only until I impose X=0 (for i=5). For the
next iterations (i=6..10), I impose X=10000 again, but I keep getting
the duals from imposing X=0 (nonzero), its like they get "stuck" in
those values.

Any suggestions why this could be happening?

Xparam is a parameter of the subproblems declared as:

model.Xparam = Param(model.R)

The constraint where I impose the X and from which I'm getting the
duals is defined as:

#Impose solution from master
def imposeconstraints_rule(model, r):
return model.Xsubproblem[r]-model.Xparam[r] == 0.0
model.Imposeconstraints = Constraint(model.R,
rule=imposeconstraints_rule)

Thanks for your help
Francisco

David Woodruff

unread,
May 14, 2012, 2:01:01 PM5/14/12
to coopr...@googlegroups.com
A good guess is that solve failed. Can you check the solver status?
Dave

Pancho

unread,
May 14, 2012, 2:29:47 PM5/14/12
to coopr forum
Sorry, but how can you check the solver status when using the command
solve_all_instances(solver_manager, solver, sub_insts,
suffixes=['rc', 'dual'])?

I know how to do it when solving a regular pyomo model, by just
using .write()....but here?

Thanks for your quick response

Francisco

David Woodruff

unread,
May 14, 2012, 4:55:55 PM5/14/12
to coopr...@googlegroups.com
So is "solve_all_instances" your code? If so, this text might help you:

After a solve, the results object has a member Solution.Status that contains the
solver status. The following snippet shows an example of access via a
print statement:
----
instance = model.create()
results = opt.solve(instance)
print "The solver returned a status of:"+str(results.Solution.Status)
----

The use of the Python str function to cast the value to a be string makes it
easy to test it. In particular, the value 'optimal' indicates that the
solver succeeded.

Pancho

unread,
May 14, 2012, 5:32:19 PM5/14/12
to coopr forum
I can see how to do it by using opt.solve, but it seems that
solve_all_instances doesn't return anything.
This is the definition of this function:

def solve_all_instances(solver_manager, solver, instances, **kwds):
action_handles = []
instance_map = {}
kwds['opt'] = solver
for instance in instances:
action_handle = solver_manager.queue(instance, **kwds)
instance_map[action_handle] = instance
action_handles.append(action_handle)
solver_manager.wait_all(action_handles)
for action_handle in action_handles:
results = solver_manager.get_results(action_handle)
instance_map[action_handle].load(results)

I can see how you could get the results by defining :

results = solver_manager.get_results(action_handle)

But outside of function solve_all_instances, can I define
action_handle?
solver_manager.get_results(inst) doesn't work

Watson, Jean-paul

unread,
May 14, 2012, 5:57:13 PM5/14/12
to coopr...@googlegroups.com
One thing to try here is to add the keyword argument "tee=True" (at least
I think that's it) to your invocation for solve_all_instances. This will
spit a lot of output to the screen, which will at least be able to tell us
if the sub-problems are actually feasible.

We could modify solve_all_instances to return a list of results - Dave and
Pancho: Should I do this?

jpw

Watson, Jean-paul

unread,
May 14, 2012, 5:58:51 PM5/14/12
to coopr...@googlegroups.com
The "solve_all_instances" was a quick+dirty routine I introduced a while
back, to basically implement a barrier synchronization after a bunch of
sub-problem solves. The function signature and return value can and should
be improved (see another e-mail from me in this thread, sent a few minutes
ago).

jpw

Pancho

unread,
May 14, 2012, 6:14:36 PM5/14/12
to coopr forum
Thanks David and JP. It seems that all the problems are feasible, but
the duals get "stuck" after the 5th iteration.

Iteration 4 (Correct)

Variables : 7
Objective nonzeros : 3
Linear constraints : 10 [Less: 6, Equal: 4]
Nonzeros : 19
RHS nonzeros : 7

Variables : Min LB: 0.000000 Max UB: all
infinite
Objective nonzeros : Min : 4.000000 Max :
24.00000
Linear constraints :
Nonzeros : Min : 0.5000000 Max :
8.000000
RHS nonzeros : Min : 1.000000 Max :
10000.00
CPLEX> Tried aggregator 1 time.
LP Presolve eliminated 10 rows and 7 columns.
All rows and columns eliminated.
Presolve time = 0.00 sec.

Dual simplex - Optimal: Objective = -6.2600000000e+03
Solution time = 0.00 sec. Iterations = 0 (0)
Deterministic time = 0.01 ticks (7.66 ticks/sec)


duals= None
X[Finishing]= 10000.0

duals= None
X[Lumber]= 10000.0

duals= None
X[Carpentry]= 10000.0
#################################################################################

Iteration 5 (Correct)

duals= None
X[Finishing]= 0.0

duals= None
X[Lumber]= 0.0

duals= -12.0
X[Carpentry]= 0.0

Variables : 7
Objective nonzeros : 3
Linear constraints : 10 [Less: 6, Equal: 4]
Nonzeros : 19
RHS nonzeros : 4

Variables : Min LB: 0.000000 Max UB: all
infinite
Objective nonzeros : Min : 4.000000 Max :
24.00000
Linear constraints :
Nonzeros : Min : 0.5000000 Max :
8.000000
RHS nonzeros : Min : 1.000000 Max :
225.0000
CPLEX> Tried aggregator 1 time.
LP Presolve eliminated 10 rows and 7 columns.
All rows and columns eliminated.
Presolve time = 0.00 sec.

Dual simplex - Optimal: Objective = 0.0000000000e+00
Solution time = 0.00 sec. Iterations = 0 (0)
Deterministic time = 0.01 ticks (7.54 ticks/sec)

#################################################################################

Iterations >= 6 (Incorrect!)

duals= None
X[Finishing]= 10000.0

duals= None
X[Lumber]= 10000.0

duals= -12.0 (this dual should be zero, just like in iteration 4!)
X[Carpentry]= 10000.0

Variables : 7
Objective nonzeros : 3
Linear constraints : 10 [Less: 6, Equal: 4]
Nonzeros : 19
RHS nonzeros : 7

Variables : Min LB: 0.000000 Max UB: all
infinite
Objective nonzeros : Min : 4.000000 Max :
24.00000
Linear constraints :
Nonzeros : Min : 0.5000000 Max :
8.000000
RHS nonzeros : Min : 1.000000 Max :
10000.00
CPLEX> Tried aggregator 1 time.
LP Presolve eliminated 10 rows and 7 columns.
All rows and columns eliminated.
Presolve time = 0.00 sec.

Dual simplex - Optimal: Objective = -6.2600000000e+03
Solution time = 0.00 sec. Iterations = 0 (0)
Deterministic time = 0.01 ticks (7.66 ticks/sec)


Should I reset something between iterations? I could not warm start
the solver then though...

Pancho



On May 14, 5:58 pm, "Watson, Jean-paul" <jwat...@sandia.gov> wrote:
> The "solve_all_instances" was a quick+dirty routine I introduced a while
> back, to basically implement a barrier synchronization after a bunch of
> sub-problem solves. The function signature and return value can and should
> be improved (see another e-mail from me in this thread, sent a few minutes
> ago).
>
> jpw
>

Pancho

unread,
May 14, 2012, 10:11:19 PM5/14/12
to coopr forum
It seems that the problem is when the value of
inst.Imposeconstraints[r].dual=None, right after an iteration where
inst.Imposeconstraints[r].dual=(number). Then, in the new iteration
inst.Imposeconstraints[r].dual returns (number) instead of None.

Francisco

Pancho

unread,
May 15, 2012, 12:46:34 AM5/15/12
to coopr forum
Just to add more info, I extracted this from the file CPLEX.py.

Are the default duals=None?
Why constraint.dual = float(field_value) only if float(field_value) !=
0.0?


elif (tokens[0] == "constraint") and ((extract_duals is
True) or (extract_slacks is True)):
constraint_name = None
constraint_dual = None
constaint = None # cache the solution constraint
reference, as the getattr is expensive.
for i in xrange(1,len(tokens)):
field_name = tokens[i].split('=')[0]
field_value = tokens[i].split('=')
[1].lstrip("\"").rstrip("\"")
if field_name == "name":
constraint_name = field_value
constraint = soln.constraint[constraint_name]
elif (extract_duals is True) and (field_name ==
"dual"): # for LPs
# assumes the name field is first.
if float(field_value) != 0.0:
constraint.dual = float(field_value)
elif (extract_slacks is True) and (field_name ==
"slack"): # for MIPs
# assumes the name field is first.
if float(field_value) != 0.0:
constraint.slack = float(field_value)

Watson, Jean-paul

unread,
May 15, 2012, 8:03:10 AM5/15/12
to coopr...@googlegroups.com
That's a great question. I'll take a look at this tonight (I'm traveling
today).

The short answer is that we haven't had to much use in terms of the suffix
information, so you may be encountering something stupid-ish in the
implementation.

jpw

Watson, Jean-paul

unread,
May 15, 2012, 12:40:41 PM5/15/12
to coopr...@googlegroups.com
Hi Pancho,

I just went looking through old notes and the SVN revision log - and there
is no real good reason why the "!= 0" test is in there. The short version
is that it is both historical, and now a bug. A very long time ago, we
were trying to save space in the results object, so we were not storing
any values that were equal to 0.0. But we also 0-initialized all values
before loading a result. We culled the latter, but missed the former when
dealing with duals.

I just committed a fix to repository. You can update this file using "svn
update", if you installed from the trunk. Or do the obvious patch until
you re-install from the next official Coopr release.

Thanks for pointing this out. I have added a TODO on my end to make our
suffix tests more rigorous. The only place we currently have regression
tests for suffixes is in the benders example, which is clearly
insufficient.

jpw

On 5/14/12 10:46 PM, "Pancho" <elpanc...@gmail.com> wrote:

David Woodruff

unread,
May 15, 2012, 12:45:23 PM5/15/12
to coopr...@googlegroups.com
Pancho,
If you did not install trunk in the beginning, then I suggest that
you do so now rather than patching. See the download for VOTD or trunk
at the bottom of
of the downloads page.
For unix use:
python coopr_install --trunk
Dave

Watson, Jean-paul

unread,
May 15, 2012, 4:38:07 PM5/15/12
to coopr...@googlegroups.com
Thanks for pointing out this (now obsolete) "feature"! :)

Jpw



----- Original Message -----
From: Pancho [mailto:elpanc...@gmail.com]
Sent: Tuesday, May 15, 2012 02:34 PM
To: coopr forum <coopr...@googlegroups.com>
Subject: [EXTERNAL] Re: Access dual variables

It works now! Thanks a lot for your help.

Pancho

On May 15, 12:45 pm, David Woodruff <david.l.woodr...@gmail.com>
wrote:
> ...
>
> read more »


Pancho

unread,
May 15, 2012, 4:34:56 PM5/15/12
to coopr forum
It works now! Thanks a lot for your help.

Pancho

On May 15, 12:45 pm, David Woodruff <david.l.woodr...@gmail.com>
wrote:
> ...
>
> read more »
Reply all
Reply to author
Forward
0 new messages