source: coopr.pyomo/stable/2.2/coopr/pyomo/io/cpxlp.py @ 2230

Last change on this file since 2230 was 2230, checked in by wehart, 11 years ago

Merged revisions 2208-2229 via svnmerge from
https://software.sandia.gov/svn/public/coopr/coopr.pyomo/trunk

........

r2209 | jwatson | 2010-01-30 17:01:12 -0700 (Sat, 30 Jan 2010) | 3 lines


Fix to new LP writer to output variable label, instead of name. The name has embedded quotes, which hoses all kinds of things. The label doesn't.

........

r2210 | jwatson | 2010-01-30 20:11:21 -0700 (Sat, 30 Jan 2010) | 3 lines


Additional minor fixings to make new LP writer work with PH - which it now does.

........

r2211 | jwatson | 2010-01-30 20:17:19 -0700 (Sat, 30 Jan 2010) | 3 lines


Suppressing warning for "constant constraint" when writing LP files - this happens all the time in PH and other contexts, where variable fixing is active.

........

r2213 | wehart | 2010-01-30 23:03:26 -0700 (Sat, 30 Jan 2010) | 2 lines


Adding the ability to process multiple *.dat files on the Pyomo command line

........

r2214 | wehart | 2010-01-31 14:44:12 -0700 (Sun, 31 Jan 2010) | 3 lines


Reworking the canonical expression definition to include a
dictionary from variable ID -> _VarValue object.

........

r2217 | wehart | 2010-01-31 21:48:34 -0700 (Sun, 31 Jan 2010) | 2 lines


Update to NL file writer, after recent changes in the canonical representation.

........

r2218 | wehart | 2010-01-31 22:11:59 -0700 (Sun, 31 Jan 2010) | 7 lines


Reworking the symbol_map generated for the NL file.


Fixing the use of this within PyomoModel? when loading data.


Adding logic to specify the 'solver' option when an optimization interface
is specified (e.g. asl:cplexamp).

........

r2226 | wehart | 2010-01-31 22:27:19 -0700 (Sun, 31 Jan 2010) | 2 lines


Removing some debugging IO.

........

r2227 | wehart | 2010-02-01 14:40:08 -0700 (Mon, 01 Feb 2010) | 7 lines


This is a major rework of the ModelData? logic. I've extracted the
actual operations for data loading into different plugins. These
plugins share a common code base for processing the data into a form that
Pyomo can use, but that logic is now separate from the ModelData? class itself.
This will provide a path-forward for cleaning up this code, as well as
modularizing the parsing of AMPL *.dat files.

........

r2228 | jwatson | 2010-02-01 14:49:40 -0700 (Mon, 01 Feb 2010) | 3 lines


Final fix to LP writer to (optionally) account for the models in which a variable originates - the PySP extensive form writer now passes all tests.

........

File size: 19.9 KB
Line 
1#  _________________________________________________________________________
2#
3#  Coopr: A COmmon Optimization Python Repository
4#  Copyright (c) 2008 Sandia Corporation.
5#  This software is distributed under the BSD License.
6#  Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
7#  the U.S. Government retains certain rights in this software.
8#  For more information, see the Coopr README.txt file.
9#  _________________________________________________________________________
10
11#
12# Problem Writer for CPLEX LP Format Files
13#
14
15from coopr.opt import ProblemFormat
16from coopr.pyomo.base import expr, Var, Constraint, Objective
17from coopr.pyomo.base.var import _VarValue, _VarBase
18from coopr.pyomo.base.param import _ParamValue
19from coopr.pyomo.base.numtypes import minimize, maximize
20from coopr.pyomo.base.expr import *
21from coopr.pyomo.base import *
22from coopr.opt.base import problem, AbstractProblemWriter
23from coopr.pyomo.expr import *
24import math
25
26def convert_name(namestr):
27    namestr = namestr.replace('[','(')
28    namestr = namestr.replace(']',')')
29    return namestr
30
31class ProblemWriter_cpxlp(AbstractProblemWriter):
32
33    def __init__(self):
34        AbstractProblemWriter.__init__(self,ProblemFormat.cpxlp)
35
36        # temporarily placing attributes used in extensive-form writing
37        # here, for eventual migration to the base class.
38        self._output_objectives = True
39        self._output_constraints = True       
40        self._output_variables = True # means we're outputting *some* variables - types determined by following flags.
41
42        # building on the above, partition out the variables that should
43        # be written if the _output_variables attribute is set to True.
44        # this is useful when writing components of multiple models to a
45        # single LP file. unfortunately, the CPLEX LP format requires all
46        # integer and binary (and continuous) variables to appear in a
47        # single continuous block.
48        self._output_continuous_variables = True
49        self._output_integer_variables = True
50        self._output_binary_variables = True
51
52        # do I pre-pend the model name to all identifiers I output?
53        # useful in cases where you are writing components of multiple
54        # models to a single output LP file.
55        self._output_prefixes = False 
56
57    def __call__(self, model, filename):
58        if filename is None:
59           filename = model.name + ".lp"
60        OUTPUT=open(filename,"w")
61        self._print_model_LP(model,OUTPUT)
62        OUTPUT.close()
63        return filename, None
64
65    def _get_bound(self, exp):
66        if isinstance(exp,expr._IdentityExpression):
67           return self._get_bound(exp._args[0])
68        elif isinstance(exp,NumericConstant):
69           return exp.value
70        else:
71           raise ValueError, "ERROR: nonconstant bound: " + str(exp)
72           return None
73
74    def _print_bound(self, exp, OUTPUT, offset=0.0):
75        if isinstance(exp,expr._IdentityExpression):
76           self._print_bound(exp._args[0], OUTPUT, offset)
77        elif isinstance(exp,NumericConstant) or isinstance(exp,_ParamValue):
78           print >>OUTPUT, exp.value+offset
79        else:
80           # Idea: Bundle the following into a single string and pass with thrown exception
81           print "ERROR: nonconstant constraint bound"
82           print "Expression type=" + exp.__class__.__name__
83           if not exp is None:
84                print "Offending expression:"
85                exp.pprint()
86           raise ValueError
87
88    def _collect_key(self, a):
89        return a[0]
90
91    def _print_expr(self, x, OUTPUT, print_offset=False):
92        max_terms_on_line=5 # this is the line-limit hack
93        terms_output = 0
94        #
95        # Linear
96        #
97        if 1 in x:
98            keys = x[1].keys()
99            keys.sort()
100            for id in keys:
101                coef = x[1][id]
102                prefix = ""
103                if self._output_prefixes is True:
104                   parent_var = x[-1][id.keys()[0]].var
105                   prefix = convert_name(parent_var.model.name) + "_"
106                name = prefix + convert_name(x[-1][id.keys()[0]].label)
107                print >>OUTPUT, ('- %f' % math.fabs(coef) if coef < 0.0 else '+ %f' % coef), name,
108                terms_output += 1
109                if terms_output >= max_terms_on_line:
110                    print >>OUTPUT, ""
111                    print >>OUTPUT, "     ",
112                    terms_output = 0
113        #
114        # Quadratic
115        #
116        if 2 in x:
117            keys = x[2].keys()
118            keys.sort()
119            for id in keys:
120                coef = x[2][id]
121                print >>OUTPUT, ('- %f' % math.fabs(coef) if coef < 0.0 else '+ %f' % coef),
122                for var in id:
123                    prefix = ""
124                    if self._output_prefixes is True:
125                       parent_var = x[-1][var].var
126                       prefix = convert_name(parent_var.model.name) + "_"
127                    name = prefix + convert_name(x[-1][var].label)
128                    if id[var] > 1:
129                        print "%s^%d" % (name,id[var]),
130                    else:
131                        print name,
132                terms_output += 1
133                if terms_output >= max_terms_on_line:
134                    print >>OUTPUT, ""
135                    print >>OUTPUT, "     ",
136                    terms_output = 0
137        #
138        # Constant offset
139        #
140        offset=0.0
141        if 0 in x:
142            offset = x[0][None]
143        if print_offset and offset != 0.0:
144            print >>OUTPUT, ('- %f' % math.fabs(offset) if offset < 0.0 else '+ %f' % offset), "ONE_VAR_CONSTANT"
145            terms_output += 1
146        #
147        # Return constant offset
148        #
149        return offset
150
151    def _print_quadterm(self, x, is_minimizing, OUTPUT):
152        # the LP format doesn't allow for expression of constant terms in the objective.
153        # a work-around involves tracking the sum of constant terms in the quadratic
154        # quadratic terms, and then writing that out with a dummy variable forced equal
155        # to one.
156        print >>OUTPUT, ""
157        for arg in x._args:
158            if isinstance(arg,expr._ProductExpression):
159                # NOTE: we need to handle quadratics defined with the 'pow' term here.
160                # WARNING: The following code is very specific to progressive hedging, with
161                #          a very specific format assumed. we do need to handle more general
162                #          expressions, but we'll worry about that at a latter time.
163                blend = arg._numerator[0]()
164
165                if blend is 1:
166
167                   rho = arg._numerator[1]()
168
169                   pow_expression = arg._numerator[2]
170               
171                   base = pow_expression._args[0]
172                   exponent = pow_expression._args[1]
173
174                   if not isinstance(base,expr._SumExpression):
175                       raise ValueError, "Quadratic term base must be a _SumExpression"
176                   if not isinstance(exponent,numvalue.NumericConstant):
177                       raise ValueError, "Quadratic term exponent must be a NumericConstant"
178                   variable = base._args[0]
179                   offset = base._args[1]
180                   if variable.status is not VarStatus.unused:
181
182                      if is_minimizing is True:
183                         print >>OUTPUT, " + [" + str(rho) + " " + convert_name(variable.label) + "^2]/2",
184                      else:
185                         print >>OUTPUT, " - [" + str(rho) + " " + convert_name(variable.label) + "^2]/2",
186
187                      if (is_minimizing is True):
188                         if offset.value < 0.0:
189                            print >>OUTPUT, " + " + str(abs(rho*offset.value)) + " " + convert_name(variable.label),                   
190                         else:
191                            print >>OUTPUT, " - " + str(rho*offset.value) + " " + convert_name(variable.label),
192                      else:
193                         if offset.value < 0.0:
194                            print >>OUTPUT, " - " + str(abs(rho*offset.value)) + " " + convert_name(variable.label),                   
195                         else:
196                            print >>OUTPUT, " + " + str(rho*offset.value) + " " + convert_name(variable.label),
197
198                      objective_offset = (rho * offset.value * offset.value / 2.0)
199                      if is_minimizing is True:
200                         print >>OUTPUT, " + " + str(objective_offset) + " ONE_VAR_CONSTANT"
201                      else:
202                         print >>OUTPUT, " - " + str(objective_offset) + " ONE_VAR_CONSTANT"
203                     
204            elif isinstance(arg,numvalue.NumericConstant):
205                # this is the "0.0" element that forms the initial expression - the
206                # quadratic sub-expressions aren't known to the presolve routines.
207                # ideally unnecessary - hacked in for now.
208                pass
209
210            else:
211                print `arg`
212                raise ValueError, "Unknown expression sub-type found in quadratic objective expression"   
213   
214
215    def _print_model_LP(self, model, OUTPUT):
216
217        _obj = model.active_components(Objective)
218       
219        #
220        # Objective
221        #
222        if self._output_objectives is True:
223           printed_quadterm = False
224           if len(_obj) == 0:
225              raise ValueError, "ERROR: No objectives defined for input model=" + str(model.name) + "; cannot write legal LP file"
226           if _obj[ _obj.keys()[0] ].sense == maximize:
227              print >>OUTPUT, "max "
228           else:
229              print >>OUTPUT, "min "
230           print >>OUTPUT, "obj: ",
231           obj = _obj[ _obj.keys()[0] ]
232           for key in obj:
233                if is_constant(obj[key].repn):
234                    print "WARNING: ignoring objective %s[%s] which is constant" % (str(obj),str(key))
235                    continue
236                if is_nonlinear(obj[key].repn) and not is_quadratic(obj[key].repn):
237                    raise ValueError, "Cannot write legal LP file.  Objective %s[%s] has nonlinear terms that are not quadratic." % (str(obj),str(key))
238                self._print_expr(obj[key].repn, OUTPUT, print_offset=True)
239           if obj._quad_subexpr is not None:
240               self._print_quadterm(obj._quad_subexpr, (_obj[ _obj.keys()[0] ].sense == minimize), OUTPUT)
241               printed_quadterm = True               
242           print >>OUTPUT, ""
243       
244        #
245        # Constraints
246        #
247        # if there are no non-trivial constraints, you'll end up with an empty constraint block. CPLEX is OK with
248        # this, but GLPK isn't. And eliminating the constraint block (i.e., the "s.t." line) causes GLPK to whine
249        # elsewhere. output a warning if the constraint block is empty, so users can quickly determine the cause
250        # of the solve failure.
251
252        if self._output_constraints is True:
253
254           # for now, if this routine isn't writing everything, then assume a meta-level handler is
255           # dealing with the writing of the transitional elements - these should probably be done in
256           # the form of "end_objective()" and "end_constraint()" helper methods.                   
257           if (self._output_objectives is True) and (self._output_variables is True):
258              print >>OUTPUT, "s.t."
259           
260           CON = model.active_components(Constraint)
261           have_nontrivial=False
262           for key in CON:
263             if CON[key].trivial:
264                continue
265             have_nontrivial=True
266             i=0
267             C = CON[key]
268             for cndx in C:
269               if not C[cndx].active:
270                    continue
271               try:
272                    # there are conditions, e.g., when fixing variables, under which a constraint block might be empty.
273                    # ignore these, for both practical reasons and the fact that the CPLEX LP format requires a variable
274                    # in the constraint body. it is also possible that the body of the constraint consists of only a
275                    # constant, in which case the "variable" of
276                    if is_constant(C[cndx].repn):
277                        # this happens *all* the time in many applications, including PH - so suppress the warning.
278                        #print "WARNING: ignoring constraint %s[%s] which is constant" % (str(C),str(cndx))
279                        continue
280                    if is_nonlinear(C[cndx].repn):
281                        raise ValueError, "Cannot write legal LP file.  Constraint %s[%s] has a body with nonlinear terms." % (str(C),str(cndx))
282                    prefix = ""
283                    if self._output_prefixes is True:
284                        if C.model is None:
285                           raise RuntimeError, "Constraint="+C._data[cndx].label+" has no model attribute - no label prefix can be assigned"
286                        prefix = C.model.name+"_"
287
288                    if C._data[cndx]._equality:
289                         #
290                         # Equality constraint
291                         #
292                         print >>OUTPUT, prefix + "c_e_" + convert_name(C._data[cndx].label) + "_: ",
293                         offset = self._print_expr(C[cndx].repn, OUTPUT)
294                         print >>OUTPUT, "=",
295                         self._print_bound(C._data[cndx].lower, OUTPUT, -offset)
296                         print >>OUTPUT, ""
297                         i=i+1
298                    else:
299                         #
300                         # Inequality constraint
301                         #
302                         if C._data[cndx].lower is not None:
303                            print >>OUTPUT, prefix + "c_l_" + convert_name(C._data[cndx].label) + "_: ",
304                            offset = self._print_expr(C[cndx].repn, OUTPUT)
305                            print >>OUTPUT, ">=",
306                            self._print_bound(C._data[cndx].lower, OUTPUT, -offset)
307                            print >>OUTPUT, ""
308                            i=i+1
309                         if C._data[cndx].upper is not None:
310                            print >>OUTPUT, prefix + "c_u_" + convert_name(C._data[cndx].label) + "_: ",
311                            offset = self._print_expr(C[cndx].repn, OUTPUT)
312                            print >>OUTPUT, "<=",
313                            self._print_bound(C._data[cndx].upper, OUTPUT, -offset)
314                            print >>OUTPUT, ""
315                            i=i+1
316               except ValueError, msg:
317                  print msg
318                  raise ValueError, "ERROR: Failed to output constraint="+C.name+", index="+`cndx`+" in generation of LP format file"
319           if not have_nontrivial:
320             print "WARNING: Empty constraint block written in LP format - solver may error"
321         
322           # the CPLEX LP format doesn't allow constants in the objective (or constraint body), which is a bit silly.
323           # to avoid painful book-keeping, we introduce the following "variable", constrained to the value 1. this is
324           # used when quadratic terms are present. worst-case, if not used, is that CPLEX easily pre-processes it out.
325           print >>OUTPUT, "c_e_ONE_VAR_CONSTANT: ONE_VAR_CONSTANT = 1.0"
326           print >>OUTPUT, ""
327           
328        #
329        # Bounds
330        #
331
332        if self._output_variables is True:
333
334           # for now, if this routine isn't writing everything, then assume a meta-level handler is
335           # dealing with the writing of the transitional elements - these should probably be done in
336           # the form of "end_objective()" and "end_constraint()" helper methods.                   
337           if (self._output_objectives is True) and (self._output_constraints is True):
338              print >>OUTPUT, "bounds "
339
340           # scan all variables even if we're only writing a subset of them. required
341           # because we don't store maps by variable type currently.
342           
343           # track the number of integer and binary variables, so you can output their status later.
344           niv = nbv = 0
345           VAR = model.active_components(Var)
346           for var in VAR.values():
347              if isinstance(var.domain, IntegerSet):
348                 niv += 1
349              elif isinstance(var.domain, BooleanSet):
350                 nbv += 1
351
352              if self._output_continuous_variables is True:
353                 for ndx in var._varval.keys():
354                    if not var._varval[ndx].active:
355                        continue
356                    prefix = ""
357                    if self._output_prefixes is True:
358                       prefix = convert_name(var.model.name)+"_"
359
360                    if var[ndx].id != -1: # if the variable isn't referenced in the model, don't output bounds...
361                       # in the CPLEX LP file format, the default variable bounds are 0 and +inf.
362                       # these bounds are in conflict with Pyomo, which assumes -inf and inf (which
363                       # we would argue is more rational).
364                       print >>OUTPUT,"   ",
365                       if var[ndx].lb is not None:
366                          print >>OUTPUT, str(value(var[ndx].lb())) + " <= ",
367                       else:
368                          print >>OUTPUT, " -inf <= ",
369                       name_to_output = prefix+convert_name(var[ndx].label)
370                       if name_to_output == "e":
371                          raise ValueError, "Attempting to write variable with name=e in a CPLEX LP formatted file - will cause a parse failure due to confusion with numeric values expressed in scientific notation"
372                       print >>OUTPUT, name_to_output,
373                       if var[ndx].ub is not None:
374                          print >>OUTPUT, " <= " + str(value(var[ndx].ub())),
375                          print >>OUTPUT, ""
376                       else:
377                          print >>OUTPUT, " <= +inf"
378                       
379           if (niv > 0) and (self._output_integer_variables is True): 
380
381              # if we're outputting the whole model, then assume we can output the "integer"
382              # header. if not, then assume a meta-level process is taking care of it.
383              if (self._output_objectives is True) and (self._output_constraints is True):             
384                 print >>OUTPUT, "integer"
385
386              prefix = ""
387              if self._output_prefixes is True:
388                 prefix = convert_name(var.model.name)+"_"
389
390              for var in VAR.values():
391                 if isinstance(var.domain, IntegerSet):
392                    for ndx in var.keys():
393                       if not var[ndx].active:
394                            continue
395                       if var[ndx].id != -1: # if the variable isn't referenced, skip.
396                          print >>OUTPUT, "   ", prefix+convert_name(var[ndx].label)
397
398           if (nbv > 0) and (self._output_binary_variables is True):
399
400              # if we're outputting the whole model, then assume we can output the "binary"
401              # header. if not, then assume a meta-level process is taking care of it.
402              if (self._output_objectives is True) and (self._output_constraints is True):             
403                 print >>OUTPUT, "binary"
404
405              prefix = ""
406              if self._output_prefixes is True:
407                 prefix = convert_name(var.model.name)+"_"
408
409              for var in VAR.values():
410                 if isinstance(var.domain, BooleanSet):
411                    for ndx in var.keys():
412                       if not var[ndx].active:
413                            continue
414                       if var[ndx].id != -1: # if the variable isn't referenced, skip.
415                          print >>OUTPUT, "   ", prefix+convert_name(var[ndx].label)
416
417        #
418        # wrap-up
419        #
420
421        if (self._output_objectives is True) and (self._output_constraints is True) and (self._output_variables is True):
422           #
423           # End
424           #
425           print >>OUTPUT, "end "
426       
427problem.WriterRegistration(str(ProblemFormat.cpxlp), ProblemWriter_cpxlp)
Note: See TracBrowser for help on using the repository browser.