source: coopr.pysp/trunk/coopr/pysp/phinit.py @ 2446

Last change on this file since 2446 was 2446, checked in by jwatson, 9 years ago

Added --output-scenario-tree-solution option to PySP runph script, which will:
1) Create a solution from the node averages in a scenario tree -and-
2) Output the full solution in scenario tree format (which includes the leaf nodes).

  • Property svn:executable set to *
File size: 27.6 KB
Line 
1#  _________________________________________________________________________
2#
3#  Coopr: A COmmon Optimization Python Repository
4#  Copyright (c) 2009 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
12import sys
13import os
14from optparse import OptionParser
15
16import pyutilib.services
17import pyutilib.misc
18import textwrap
19import traceback
20
21# garbage collection control.
22import gc
23
24# for profiling
25import cProfile
26import pstats
27
28# for serializing
29import pickle
30
31from coopr.pysp.convergence import *
32from coopr.pysp.scenariotree import *
33from coopr.pysp.ph import *
34from coopr.pysp.ef import *
35from coopr.opt.base import SolverFactory
36from coopr.opt.parallel import SolverManagerFactory
37
38#
39# utility method to construct an option parser for ph arguments, to be
40# supplied as an argument to the runph method.
41#
42
43def construct_ph_options_parser(usage_string):
44
45   parser = OptionParser()
46   parser.add_option("--verbose",
47                     help="Generate verbose output for both initialization and execution. Default is False.",
48                     action="store_true",
49                     dest="verbose",
50                     default=False)
51   parser.add_option("--report-solutions",
52                     help="Always report PH solutions after each iteration. Enabled if --verbose is enabled. Default is False.",
53                     action="store_true",
54                     dest="report_solutions",
55                     default=False)
56   parser.add_option("--output-scenario-tree-solution",
57                     help="Report the full solution (even leaves) in scenario tree format upon termination. Values represent averages, so convergence is not an issue. Default is False.",
58                     action="store_true",
59                     dest="output_scenario_tree_solution",
60                     default=False)
61   parser.add_option("--report-weights",
62                     help="Always report PH weights prior to each iteration. Enabled if --verbose is enabled. Default is False.",
63                     action="store_true",
64                     dest="report_weights",
65                     default=False)
66   parser.add_option("--model-directory",
67                     help="The directory in which all model (reference and scenario) definitions are stored. Default is \".\".",
68                     action="store",
69                     dest="model_directory",
70                     type="string",
71                     default=".")
72   parser.add_option("--instance-directory",
73                     help="The directory in which all instance (reference and scenario) definitions are stored. Default is \".\".",
74                     action="store",
75                     dest="instance_directory",
76                     type="string",
77                     default=".")
78   parser.add_option("--solver",
79                     help="The type of solver used to solve scenario sub-problems. Default is cplex.",
80                     action="store",
81                     dest="solver_type",
82                     type="string",
83                     default="cplex")
84   parser.add_option("--solver-manager",
85                     help="The type of solver manager used to coordinate scenario sub-problem solves. Default is serial.",
86                     action="store",
87                     dest="solver_manager_type",
88                     type="string",
89                     default="serial")
90   parser.add_option("--scenario-solver-options",
91                     help="Solver options for all PH scenario sub-problems",
92                     action="append",
93                     dest="scenario_solver_options",
94                     type="string",
95                     default=[])
96   parser.add_option("--scenario-mipgap",
97                     help="Specifies the mipgap for all PH scenario sub-problems",
98                     action="store",
99                     dest="scenario_mipgap",
100                     type="float",
101                     default=None)
102   parser.add_option("--ef-solver-options",
103                     help="Solver options for the extension form problem",
104                     action="append",
105                     dest="ef_solver_options",
106                     type="string",
107                     default=[])
108   parser.add_option("--ef-mipgap",
109                     help="Specifies the mipgap for the EF solve",
110                     action="store",
111                     dest="ef_mipgap",
112                     type="float",
113                     default=None)
114   parser.add_option("--max-iterations",
115                     help="The maximal number of PH iterations. Default is 100.",
116                     action="store",
117                     dest="max_iterations",
118                     type="int",
119                     default=100)
120   parser.add_option("--default-rho",
121                     help="The default (global) rho for all blended variables. Default is 1.",
122                     action="store",
123                     dest="default_rho",
124                     type="float",
125                     default=1.0)
126   parser.add_option("--rho-cfgfile",
127                     help="The name of a configuration script to compute PH rho values. Default is None.",
128                     action="store",
129                     dest="rho_cfgfile",
130                     type="string",
131                     default=None)
132   parser.add_option("--bounds-cfgfile",
133                     help="The name of a configuration script to set variable bound values. Default is None.",
134                     action="store",
135                     dest="bounds_cfgfile",
136                     default=None)
137   parser.add_option("--enable-termdiff-convergence",
138                     help="Terminate PH based on the termdiff convergence metric. Default is True.",
139                     action="store_true",
140                     dest="enable_termdiff_convergence",
141                     default=True)
142   parser.add_option("--enable-normalized-termdiff-convergence",
143                     help="Terminate PH based on the normalized termdiff convergence metric. Default is True.",
144                     action="store_true",
145                     dest="enable_normalized_termdiff_convergence",
146                     default=False)
147   parser.add_option("--termdiff-threshold",
148                     help="The convergence threshold used in the term-diff and normalized term-diff convergence criteria. Default is 0.01.",
149                     action="store",
150                     dest="termdiff_threshold",
151                     type="float",
152                     default=0.01)
153   parser.add_option("--enable-free-discrete-count-convergence",
154                     help="Terminate PH based on the free discrete variable count convergence metric. Default is False.",
155                     action="store_true",
156                     dest="enable_free_discrete_count_convergence",
157                     default=False)
158   parser.add_option("--free-discrete-count-threshold",
159                     help="The convergence threshold used in the criterion based on when the free discrete variable count convergence criterion. Default is 20.",
160                     action="store",
161                     dest="free_discrete_count_threshold",
162                     type="float",
163                     default=20)
164   parser.add_option("--enable-ww-extensions",
165                     help="Enable the Watson-Woodruff PH extensions plugin. Default is False.",
166                     action="store_true",
167                     dest="enable_ww_extensions",
168                     default=False)
169   parser.add_option("--ww-extension-cfgfile",
170                     help="The name of a configuration file for the Watson-Woodruff PH extensions plugin. Default is wwph.cfg.",
171                     action="store",
172                     dest="ww_extension_cfgfile",
173                     type="string",
174                     default="")
175   parser.add_option("--ww-extension-suffixfile",
176                     help="The name of a variable suffix file for the Watson-Woodruff PH extensions plugin. Default is wwph.suffixes.",
177                     action="store",
178                     dest="ww_extension_suffixfile",
179                     type="string",
180                     default="")
181   parser.add_option("--user-defined-extension",
182                     help="The name of a python module specifying a user-defined PH extension plugin.",
183                     action="store",
184                     dest="user_defined_extension",
185                     type="string",
186                     default=None)
187   parser.add_option("--write-ef",
188                     help="Upon termination, write the extensive form of the model - accounting for all fixed variables.",
189                     action="store_true",
190                     dest="write_ef",
191                     default=False)
192   parser.add_option("--solve-ef",
193                     help="Following write of the extensive form model, solve it.",
194                     action="store_true",
195                     dest="solve_ef",
196                     default=False)
197   parser.add_option("--ef-output-file",
198                     help="The name of the extensive form output file (currently only LP format is supported), if writing of the extensive form is enabled. Default is efout.lp.",
199                     action="store",
200                     dest="ef_output_file",
201                     type="string",
202                     default="efout.lp")
203   parser.add_option("--suppress-continuous-variable-output",
204                     help="Eliminate PH-related output involving continuous variables.",
205                     action="store_true",
206                     dest="suppress_continuous_variable_output",
207                     default=False)
208   parser.add_option("--keep-solver-files",
209                     help="Retain temporary input and output files for scenario sub-problem solves",
210                     action="store_true",
211                     dest="keep_solver_files",
212                     default=False)
213   parser.add_option("--output-solver-logs",
214                     help="Output solver logs during scenario sub-problem solves",
215                     action="store_true",
216                     dest="output_solver_logs",
217                     default=False)
218   parser.add_option("--output-ef-solver-log",
219                     help="Output solver log during the extensive form solve",
220                     action="store_true",
221                     dest="output_ef_solver_log",
222                     default=False)
223   parser.add_option("--output-solver-results",
224                     help="Output solutions obtained after each scenario sub-problem solve",
225                     action="store_true",
226                     dest="output_solver_results",
227                     default=False)
228   parser.add_option("--output-times",
229                     help="Output timing statistics for various PH components",
230                     action="store_true",
231                     dest="output_times",
232                     default=False)
233   parser.add_option("--disable-warmstarts",
234                     help="Disable warm-start of scenario sub-problem solves in PH iterations >= 1. Default is False.",
235                     action="store_true",
236                     dest="disable_warmstarts",
237                     default=False)
238   parser.add_option("--drop-proximal-terms",
239                     help="Eliminate proximal terms (i.e., the quadratic penalty terms) from the weighted PH objective. Default is False.",
240                     action="store_true",
241                     dest="drop_proximal_terms",
242                     default=False)
243   parser.add_option("--retain-quadratic-binary-terms",
244                     help="Do not linearize PH objective terms involving binary decision variables",
245                     action="store_true",
246                     dest="retain_quadratic_binary_terms",
247                     default=False)
248   parser.add_option("--linearize-nonbinary-penalty-terms",
249                     help="Approximate the PH quadratic term for non-binary variables with a piece-wise linear function, using the supplied number of equal-length pieces from each bound to the average",
250                     action="store",
251                     dest="linearize_nonbinary_penalty_terms",
252                     type="int",
253                     default=0)
254   parser.add_option("--breakpoint-strategy",
255                     help="Specify the strategy to distribute breakpoints on the [lb, ub] interval of each variable when linearizing. 0 indicates uniform distribution. 1 indicates breakpoints at the node min and max, uniformly in-between. 2 indicates more aggressive concentration of breakpoints near the observed node min/max.",
256                     action="store",
257                     dest="breakpoint_strategy",
258                     type="int",
259                     default=0)
260   parser.add_option("--checkpoint-interval",
261                     help="The number of iterations between writing of a checkpoint file. Default is 0, indicating never.",
262                     action="store",
263                     dest="checkpoint_interval",
264                     type="int",
265                     default=0)
266   parser.add_option("--restore-from-checkpoint",
267                     help="The name of the checkpoint file from which PH should be initialized. Default is \"\", indicating no checkpoint restoration",
268                     action="store",
269                     dest="restore_from_checkpoint",
270                     type="string",
271                     default="")
272   parser.add_option("--profile",
273                     help="Enable profiling of Python code.  The value of this option is the number of functions that are summarized.",
274                     action="store",
275                     dest="profile",
276                     type="int",
277                     default=0)
278   parser.add_option("--disable-gc",
279                     help="Disable the python garbage collecter. Default is False.",
280                     action="store_true",
281                     dest="disable_gc",
282                     default=False)
283
284   parser.usage=usage_string
285
286   return parser
287
288#
289# Create the reference model / instance and scenario tree instance for PH.
290#
291
292def load_reference_and_scenario_models(options):
293
294   #
295   # create and populate the reference model/instance pair.
296   #
297
298   reference_model = None
299   reference_instance = None
300
301   try:
302      reference_model_filename = os.path.expanduser(options.model_directory)+os.sep+"ReferenceModel.py"
303      if options.verbose is True:
304         print "Scenario reference model filename="+reference_model_filename
305      model_import = pyutilib.misc.import_file(reference_model_filename)
306      if "model" not in dir(model_import):
307         print ""
308         print "***ERROR: Exiting test driver: No 'model' object created in module "+reference_model_filename
309         return None, None, None, None
310
311      if model_import.model is None:
312         print ""
313         print "***ERROR: Exiting test driver: 'model' object equals 'None' in module "+reference_model_filename
314         return None, None, None, None
315 
316      reference_model = model_import.model
317   except IOError:
318      print "***ERROR: Failed to load scenario reference model from file="+reference_model_filename
319      return None, None, None, None
320
321   try:
322      reference_instance_filename = os.path.expanduser(options.instance_directory)+os.sep+"ReferenceModel.dat"
323      if options.verbose is True:
324         print "Scenario reference instance filename="+reference_instance_filename
325      reference_instance = reference_model.create(reference_instance_filename)
326   except IOError:
327      print "***ERROR: Failed to load scenario reference instance data from file="+reference_instance_filename
328      return None, None, None, None
329
330   #
331   # create and populate the scenario tree model
332   #
333
334   from coopr.pysp.util.scenariomodels import scenario_tree_model
335   scenario_tree_instance = None
336
337   try:
338      scenario_tree_instance_filename = os.path.expanduser(options.instance_directory)+os.sep+"ScenarioStructure.dat"
339      if options.verbose is True:
340         print "Scenario tree instance filename="+scenario_tree_instance_filename
341      scenario_tree_instance = scenario_tree_model.create(scenario_tree_instance_filename)
342   except IOError:
343      print "***ERROR: Failed to load scenario tree reference instance data from file="+scenario_tree_instance_filename
344      return None, None, None, None
345   
346   #
347   # construct the scenario tree
348   #
349   scenario_tree = ScenarioTree(scenarioinstance=reference_instance,
350                                scenariotreeinstance=scenario_tree_instance)
351
352   return reference_model, reference_instance, scenario_tree, scenario_tree_instance
353
354#
355# Create a PH object from a checkpoint.
356#
357def create_ph_from_checkpoint(options):
358
359   # we need to load the reference model, as pickle doesn't save contents of .py files!
360   try:
361      reference_model_filename = os.path.expanduser(options.model_directory)+os.sep+"ReferenceModel.py"
362      if options.verbose is True:
363         print "Scenario reference model filename="+reference_model_filename
364      model_import = pyutilib.misc.import_file(reference_model_filename)
365      if "model" not in dir(model_import):
366         print "***ERROR: Exiting test driver: No 'model' object created in module "+reference_model_filename
367         return
368
369      if model_import.model is None:
370         print "***ERROR: Exiting test driver: 'model' object equals 'None' in module "+reference_model_filename
371         return None
372 
373      reference_model = model_import.model
374   except IOError:
375      print "***ERROR: Failed to load scenario reference model from file="+reference_model_filename
376      return None
377
378   # import the saved state
379     
380   try:
381      checkpoint_file = open(options.restore_from_checkpoint,"r")
382      ph = pickle.load(checkpoint_file)
383      checkpoint_file.close()
384     
385   except IOError, msg:
386      raise RuntimeError, msg
387
388   # tell PH to build the right solver manager and solver TBD - AND PLUGINS, BUT LATER
389
390   raise RuntimeError, "Checkpoint restoration is not fully supported/tested yet!"
391
392   return ph
393
394#
395# Create a PH object from a checkpoint.
396#
397def create_ph_from_scratch(options, reference_model, reference_instance, scenario_tree):
398
399   #
400   # print the input tree for validation/information purposes.
401   #
402   if options.verbose is True:
403      scenario_tree.pprint()
404
405   #
406   # validate the tree prior to doing anything serious
407   #
408   print ""
409   if scenario_tree.validate() is False:
410      print "***ERROR: Scenario tree is invalid****"
411      return None
412   else:
413      if options.verbose is True:
414         print "Scenario tree is valid!"
415   print ""
416
417   #
418   # if any of the ww extension configuration options are specified without the
419   # ww extension itself being enabled, halt and warn the user - this has led
420   # to confusion in the past, and will save user support time.
421   #
422   if len(options.ww_extension_cfgfile) > 0 and options.enable_ww_extensions is False:
423      print "***ERROR: A configuration file was specified for the WW extension module, but the WW extensions are not enabled!"
424      return None
425
426   if len(options.ww_extension_suffixfile) > 0 and options.enable_ww_extensions is False:
427      print "***ERROR: A suffix file was specified for the WW extension module, but the WW extensions are not enabled!"         
428      return None
429
430   #
431   # if a breakpoint strategy is specified without linearization eanbled, halt and warn the user.
432   #
433   if (options.breakpoint_strategy > 0) and (options.linearize_nonbinary_penalty_terms == 0):
434      print "***ERROR: A breakpoint distribution strategy was specified, but linearization is not enabled!"
435      return None
436
437   #
438   # deal with any plugins. ww extension comes first currently, followed by an option user-defined plugin.
439   # order only matters if both are specified.
440   #
441   if options.enable_ww_extensions is True:
442
443      from coopr.pysp import wwphextension
444
445      plugin = ExtensionPoint(IPHExtension)
446      if len(options.ww_extension_cfgfile) > 0:
447         plugin.service()._configuration_filename = options.ww_extension_cfgfile
448      if len(options.ww_extension_suffixfile) > 0:
449         plugin.service()._suffix_filename = options.ww_extension_suffixfile
450
451   if options.user_defined_extension is not None:
452      print "Trying to import user-defined PH extension module="+options.user_defined_extension
453      # JPW removed the exception handling logic, as the module importer
454      # can raise a broad array of exceptions.
455      __import__(options.user_defined_extension)
456      print "Module successfully loaded"
457
458   #
459   # construct the convergence "computer" class.
460   #
461   converger = None
462   # go with the non-defaults first, and then with the default.
463   if options.enable_free_discrete_count_convergence is True:
464      converger = NumFixedDiscreteVarConvergence(convergence_threshold=options.free_discrete_count_threshold)
465   elif options.enable_normalized_termdiff_convergence is True:
466      converger = NormalizedTermDiffConvergence(convergence_threshold=options.termdiff_threshold)     
467   else:
468      converger = TermDiffConvergence(convergence_threshold=options.termdiff_threshold)     
469
470   
471   #
472   # construct and initialize PH
473   #
474   ph = ProgressiveHedging(max_iterations=options.max_iterations, \
475                           rho=options.default_rho, \
476                           rho_setter=options.rho_cfgfile, \
477                           bounds_setter=options.bounds_cfgfile, \
478                           solver=options.solver_type, \
479                           solver_manager=options.solver_manager_type, \
480                           output_scenario_tree_solution=options.output_scenario_tree_solution, \
481                           scenario_solver_options=options.scenario_solver_options, \
482                           scenario_mipgap=options.scenario_mipgap, \
483                           keep_solver_files=options.keep_solver_files, \
484                           output_solver_log=options.output_solver_logs, \
485                           output_solver_results=options.output_solver_results, \
486                           verbose=options.verbose, \
487                           report_solutions=options.report_solutions, \
488                           report_weights=options.report_weights, \
489                           output_times=options.output_times, \
490                           disable_warmstarts=options.disable_warmstarts,
491                           drop_proximal_terms=options.drop_proximal_terms,
492                           retain_quadratic_binary_terms=options.retain_quadratic_binary_terms, \
493                           linearize_nonbinary_penalty_terms=options.linearize_nonbinary_penalty_terms, \
494                           breakpoint_strategy=options.breakpoint_strategy, \
495                           checkpoint_interval=options.checkpoint_interval)
496
497   ph.initialize(scenario_data_directory_name=os.path.expanduser(options.instance_directory), \
498                 model=reference_model, \
499                 model_instance=reference_instance, \
500                 scenario_tree=scenario_tree, \
501                 converger=converger)
502
503   if options.suppress_continuous_variable_output is True:
504      ph._output_continuous_variable_stats = False # clutters up the screen, when we really only care about the binaries.
505
506   return ph
507
508#
509# Given a PH object, execute it and optionally solve the EF at the end.
510#
511
512def run_ph(options, ph):
513
514   #
515   # at this point, we have an initialized PH object by some means.
516   #
517   start_time = time.time()
518
519   #
520   # kick off the solve
521   #
522   ph.solve()
523
524   end_time = time.time()
525
526   print ""
527   print "Total PH execution time=%8.2f seconds" %(end_time - start_time)
528   print ""
529   if options.output_times is True:
530      ph.print_time_stats()
531
532   #
533   # write the extensive form, accounting for any fixed variables.
534   #
535   if (options.write_ef is True) or (options.solve_ef is True):
536      print ""
537      print "Writing EF for remainder problem"
538      print ""
539      create_and_write_ef(ph._scenario_tree, ph._instances, os.path.expanduser(options.ef_output_file))
540
541   #
542   # solve the extensive form.
543   #
544   if options.solve_ef is True:
545      print ""
546      print "Solving extensive form written to file="+os.path.expanduser(options.ef_output_file)
547      print ""
548
549      ef_solver = SolverFactory(options.solver_type)
550      if ef_solver is None:
551         raise ValueError, "Failed to create solver of type="+options.solver_type+" for use in extensive form solve"
552      if len(options.ef_solver_options) > 0:
553         print "Initializing ef solver with options="+str(options.ef_solver_options)         
554         ef_solver.set_options("".join(options.ef_solver_options))
555      if options.ef_mipgap is not None:
556         if (options.ef_mipgap < 0.0) or (options.ef_mipgap > 1.0):
557            raise ValueError, "Value of the mipgap parameter for the EF solve must be on the unit interval; value specified=" + `options.ef_mipgap`
558         else:
559            ef_solver.mipgap = options.ef_mipgap
560
561      ef_solver_manager = SolverManagerFactory(options.solver_manager_type)
562      if ef_solver is None:
563         raise ValueError, "Failed to create solver manager of type="+options.solver_type+" for use in extensive form solve"
564
565      print "Queuing extensive form solve"
566      ef_action_handle = ef_solver_manager.queue(os.path.expanduser(options.ef_output_file), opt=ef_solver, warmstart=False, tee=options.output_ef_solver_log)
567      print "Waiting for extensive form solve"
568      ef_results = ef_solver_manager.wait_for(ef_action_handle)
569      print "Extensive form solve results:"
570      ef_results.write(num=1)
571   
572#
573# The main PH initialization / runner routine.
574#
575
576def exec_ph(options):
577
578   ph = None
579
580   # if we are restoring from a checkpoint file, do so - otherwise, construct PH from scratch.
581   if len(options.restore_from_checkpoint) > 0:
582
583      ph = create_ph_from_checkpoint(options)
584     
585   else:
586
587      reference_model, reference_instance, scenario_tree, scenario_tree_instance = load_reference_and_scenario_models(options)
588      if reference_model is None or reference_instance is None or scenario_tree is None:
589         return
590      ph = create_ph_from_scratch(options, reference_model, reference_instance, scenario_tree)
591
592   if ph is None:
593      print "***FAILED TO CREATE PH OBJECT"
594      return
595
596   run_ph(options, ph)
597
598def run(args=None):
599
600    #
601    # Top-level command that executes the extensive form writer.
602    # This is segregated from run_ef_writer to enable profiling.
603    #
604
605    #
606    # Parse command-line options.
607    #
608    try:
609       ph_options_parser = construct_ph_options_parser("runph [options]")
610       (options, args) = ph_options_parser.parse_args(args=args)
611    except SystemExit:
612       # the parser throws a system exit if "-h" is specified - catch
613       # it to exit gracefully.
614       return
615
616    if options.disable_gc is True:
617       gc.disable()
618    else:
619       gc.enable()
620
621    if options.profile > 0:
622        #
623        # Call the main PH routine with profiling.
624        #
625        tfile = pyutilib.services.TempfileManager.create_tempfile(suffix=".profile")
626        tmp = cProfile.runctx('exec_ph(options)',globals(),locals(),tfile)
627        p = pstats.Stats(tfile).strip_dirs()
628        p.sort_stats('time', 'cum')
629        p = p.print_stats(options.profile)
630        p.print_callers(options.profile)
631        p.print_callees(options.profile)
632        p = p.sort_stats('cum','calls')
633        p.print_stats(options.profile)
634        p.print_callers(options.profile)
635        p.print_callees(options.profile)
636        p = p.sort_stats('calls')
637        p.print_stats(options.profile)
638        p.print_callers(options.profile)
639        p.print_callees(options.profile)
640        pyutilib.services.TempfileManager.clear_tempfiles()
641        ans = [tmp, None]
642    else:
643        #
644        # Call the main PH routine without profiling.
645        #
646        ans = exec_ph(options)
647
648    gc.enable()
649   
650    return ans
Note: See TracBrowser for help on using the repository browser.