source: coopr.pysp/trunk/coopr/pysp/scenariotree.py @ 3217

Last change on this file since 3217 was 3217, checked in by jwatson, 11 years ago

Various updates to support heteogeneous index sets in PH for different nodes in the scenario tree - more work / testing remains.

File size: 41.5 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
11import sys
12import types
13import copy
14import os.path
15import traceback
16
17from math import fabs
18
19from coopr.pyomo import *
20from phutils import *
21
22class ScenarioTreeNode(object):
23
24   #
25   # update the _solutions attribute of a tree node, given a specific
26   # variable/match-template/variable-index triple.
27   #
28   def _update_solution_map(self, variable, match_template):
29
30      variable_indices = self._variable_indices[variable.name]
31
32      # don't bother copying bounds for variables, as the values stored
33      # here are computed elsewhere - and that entity is responsible for
34      # ensuring feasibility. this also leaves room for specifying infeasible
35      # or partial solutions.
36      new_variable_index = variable._index
37      new_variable_name = variable.name
38      new_variable = None
39      if (len(new_variable_index) is 1) and (None in new_variable_index):
40         new_variable = Var(name=new_variable_name)
41      else:
42         new_variable = Var(new_variable_index, name=new_variable_name)
43      new_variable.construct()
44
45      # by default, deactive all variable values - we're copying the
46      # full index set of a variable, which is necessarily wasteful and
47      # incorrect. then, do a second pass and activate those indicies
48      # that are actually associated with this tree node.
49      for index in new_variable_index:
50         new_variable[index].deactivate()
51
52      for index in variable_indices:
53         new_variable[index].activate()
54
55      self._solutions[new_variable_name] = new_variable
56
57   #
58   # initialize the _solutions attribute of a tree node, from scratch.
59   #
60
61   def _initialize_solution_map(self):
62
63      # clear whatever was there before.
64      self._solutions = {}
65
66      # NOTE: Given the expense, this should really be optional - don't
67      #       construct unless we're actually going to use!
68      # for each variable referenced in the stage, clone the variable
69      # for purposes of storing solutions. we are being wasteful in
70      # terms copying indices that may not be referenced in the stage.
71      # this is something that we might revisit if space/performance
72      # is an issue (space is the most likely issue)
73      for variable, match_template in self._stage._variables:
74         self._update_solution_map(variable, match_template)
75
76      self._solution_map_initialized = True
77
78
79   """ Constructor
80   """
81   def __init__(self, name, conditional_probability, stage, initialize_solution=False):
82
83      self._name = name
84      self._stage = stage
85      self._parent = None
86      self._children = [] # a collection of ScenarioTreeNodes
87      self._conditional_probability = conditional_probability # conditional on parent
88      self._scenarios = [] # a collection of all Scenarios passing through this node in the tree
89      self._variable_indices = {} # a map from a variable name to the indices blended at this node.
90
91      # general use statistics for the variables at each node.
92      # each attribute is a map between the variable name and a
93      # parameter (over the same index set) encoding the corresponding
94      # statistic computed over all scenarios for that node. the
95      # parameters are named as the source variable name suffixed
96      # by one of: "NODEMIN", "NODEAVG", and "NODEMAX".
97      # NOTE: the averages are probability_weighted - the min/max
98      #       values are not.
99      # NOTE: the parameter names are basically irrelevant, and the
100      #       convention is assumed to be enforced by whoever populates
101      #       these parameters.
102      self._averages = {}
103      self._minimums = {}
104      self._maximums = {}
105
106      # a flag indicating whether the _solutions attribute has been properly initialized.
107      self._solution_map_initialized = False
108
109      # solution (variable) values for this node. assumed to be distinct
110      # from self._averages, as the latter are not necessarily feasible.
111      # objects in the map are actual pyomo Var instances; keys are
112      # variable names.
113      self._solutions = {}
114
115      if initialize_solution is True:
116         self._initialize_solution_map()
117
118   #
119   # given a set of scenario instances, compute the set of indices being blended
120   # for each variable at this node.
121   #
122
123   def defineVariableIndexSets(self, scenario_instances):
124
125      # find a representative scenario instance belonging to this node. the
126      # first scenario is as good as any.
127      scenario_instance = scenario_instances[self._scenarios[0]._name]
128
129      for reference_variable, match_template in self._stage._variables:
130
131         # the stage variable simply references the variable object in the
132         # reference scenario instance - we need to grab the variable in the
133         # scenario instance, as the index set might be different.
134         variable_name = reference_variable.name
135
136         instance_variable = getattr(scenario_instance, variable_name)
137
138         match_indices = extractVariableIndices(instance_variable, match_template)
139
140         self._variable_indices[variable_name] = match_indices
141
142   #
143   # copies the parameter values values from the _averages attribute
144   # into the _solutions attribute - only for active variable values.
145   #
146
147   def snapshotSolutionFromAverages(self):
148
149      if self._solution_map_initialized is False:
150         self._initialize_solution_map()
151
152      for variable_name, variable in self._solutions.items():
153
154         # try and grab the corresponding averages parameter - if it
155         # doesn't exist, throw an exception.
156         average_parameter = None
157         try:
158            average_parameter = self._averages[variable_name]
159         except:
160            raise RuntimeError, "No averages parameter present on tree node="+self._name+" for variable="+variable_name
161
162         for index in variable._index:
163            if variable[index].active is True:
164               variable[index] = value(average_parameter[index])
165
166   #
167   # computes the solution values from the composite scenario instances at this tree node.
168   # the input scenario_instance_map is a map from scenario name to instance objects.
169   #
170
171   def snapshotSolutionFromInstances(self, scenario_instance_map):
172
173      if self._solution_map_initialized is False:
174         self._initialize_solution_map()
175
176      for variable_name, variable in self._solutions.items():
177         for index in variable._index:
178            if variable[index].active is True:
179               node_probability = 0.0
180               avg = 0.0
181               num_scenarios_with_index = 0
182               for scenario in self._scenarios:
183                  scenario_instance = scenario_instance_map[scenario._name]
184                  node_probability += scenario._probability
185                  scenario_variable = getattr(scenario_instance, variable.name)
186                  if index in scenario_variable:
187                     num_scenarios_with_index = num_scenarios_with_index + 1
188                     var_value = value(getattr(scenario_instance, variable.name)[index])
189                     avg += (scenario._probability * var_value)
190               if num_scenarios_with_index > 0:
191                  variable[index] = avg / node_probability
192
193   #
194   # a utility to compute the cost of the current node plus the expected costs of child nodes.
195   #
196
197   def computeExpectedNodeCost(self, scenario_instance_map):
198
199      # IMPT: This implicitly assumes convergence across the scenarios - if not, garbage results.
200      instance = scenario_instance_map[self._scenarios[0]._name]
201      stage_cost_variable = instance.active_components(Var)[self._stage._cost_variable[0].name][self._stage._cost_variable[1]]
202      if stage_cost_variable.value is not None:
203         my_cost = stage_cost_variable.value
204      else:
205         # depending on the situation (and the model), the stage cost variables might not have values.
206         my_cost = 0.0
207      child_cost = 0.0
208      for child in self._children:
209         child_cost += (child._conditional_probability * child.computeExpectedNodeCost(scenario_instance_map))
210      return my_cost + child_cost
211
212
213class Stage(object):
214
215   """ Constructor
216   """
217   def __init__(self, *args, **kwds):
218
219      self._name = ""
220
221      # a collection of ScenarioTreeNode objects.
222      self._tree_nodes = []
223
224      # a collection of pairs consisting of (1) a reference to a Pyomo model Var object (in the reference scenario instance) and
225      # (2) the original match template string (for output purposes). the specific indices that match belong to the tree node.
226      self._variables = []
227
228      # a tuple consisting of (1) a reference to a pyomo model Var that computes the stage-specific cost and (2) the corresponding
229      # index. the index *is* the sole index in the cost variable, as the cost variable refers to a single variable index.
230      self._cost_variable = (None, None)
231
232   #
233   # add a new variable to the stage, which will include updating the solution maps for each associated ScenarioTreeNode.
234   #
235   def add_variable(self, variable, match_template):
236
237      self._variables.append((variable, match_template))
238
239      for tree_node in self._tree_nodes:
240         tree_node._update_solution_map(variable, match_template)
241
242class Scenario(object):
243
244   """ Constructor
245   """
246   def __init__(self, *args, **kwds):
247
248      self._name = None
249      self._leaf_node = None  # allows for construction of node list
250      self._node_list = []    # sequence from parent to leaf of ScenarioTreeNodes
251      self._probability = 0.0 # the unconditional probability for this scenario, computed from the node list
252
253class ScenarioTree(object):
254
255   # a utility to construct the stage objects for this scenario tree.
256   # operates strictly by side effects, initializing the self
257   # _stages and _stage_map attributes.
258   def _construct_stages(self, stage_ids, stage_variable_ids, stage_cost_variable_ids):
259
260      # construct the stage objects, which will leave them
261      # largely uninitialized - no variable information, in particular.
262      for stage_name in stage_ids:
263         new_stage = Stage()
264         new_stage._name = stage_name
265         self._stages.append(new_stage)
266         self._stage_map[stage_name] = new_stage
267
268      # initialize the variable collections (blended and cost) for each stage.
269      # these are references, not copies.
270      for stage_id in stage_variable_ids.keys():
271
272         if stage_id not in self._stage_map.keys():
273            raise ValueError, "Unknown stage=" + stage_id + " specified in scenario tree constructor (stage->variable map)"
274
275         stage = self._stage_map[stage_id]
276         variable_ids = stage_variable_ids[stage_id]
277
278         for variable_string in variable_ids:
279
280            if isVariableNameIndexed(variable_string) is True:
281
282               variable_name, index_template = extractVariableNameAndIndex(variable_string)
283
284               # validate that the variable exists and extract the reference.
285               if variable_name not in self._reference_instance.active_components(Var):
286                  raise ValueError, "Variable=" + variable_name + " associated with stage=" + stage_id + " is not present in model=" + self._reference_instance.name
287               variable = self._reference_instance.active_components(Var)[variable_name]
288
289               stage._variables.append((variable, index_template))
290
291            else:
292
293               # verify that the variable exists.
294               if variable_string not in self._reference_instance.active_components(Var).keys():
295                  raise RuntimeError, "Unknown variable=" + variable_string + " associated with stage=" + stage_id + " is not present in model=" + self._reference_instance.name
296
297               variable = self._reference_instance.active_components(Var)[variable_string]
298
299               # 9/14/2009 - now forcing the user to explicit specify the full
300               # match template (e.g., "foo[*,*]") instead of just the variable
301               # name (e.g., "foo") to represent the set of all indices.
302
303               # if the variable is a singleton - that is, non-indexed - no brackets is fine.
304               # we'll just tag the var[None] variable value with the (suffix,value) pair.
305               if None not in variable._index:
306                  raise RuntimeError, "Variable="+variable_string+" is an indexed variable, and templates must specify an index match; encountered in scenario tree specification for model="+self._reference_instance.name
307
308               match_indices = []
309               match_indices.append(None)
310
311               stage._variables.append((variable, ""))
312
313      for stage_id in stage_cost_variable_ids.keys():
314
315         if stage_id not in self._stage_map.keys():
316            raise ValueError, "Unknown stage=" + stage_id + " specified in scenario tree constructor (stage->cost variable map)"
317         stage = self._stage_map[stage_id]
318
319         cost_variable_string = value(stage_cost_variable_ids[stage_id]) # de-reference is required to access the parameter value
320
321         # to be extracted from the string.
322         cost_variable_name = None
323         cost_variable = None
324         cost_variable_index = None
325
326         # do the extraction.
327         if isVariableNameIndexed(cost_variable_string) is True:
328
329            cost_variable_name, index_template = extractVariableNameAndIndex(cost_variable_string)
330
331            # validate that the variable exists and extract the reference.
332            if cost_variable_name not in self._reference_instance.active_components(Var):
333               raise ValueError, "Variable=" + cost_variable_name + " associated with stage=" + stage_id + " is not present in model=" + self._reference_instance.name
334            cost_variable = self._reference_instance.active_components(Var)[cost_variable_name]
335
336            # extract all "real", i.e., fully specified, indices matching the index template.
337            match_indices = extractVariableIndices(cost_variable, index_template)
338
339            # only one index can be supplied for a stage cost variable.
340            if len(match_indices) > 1:
341                msg = 'Only one index can be specified for a stage cost '     \
342                      'variable - %s match template "%s" for variable "%s" ;' \
343                      ' encountered in scenario tree specification for model' \
344                      ' "%s"'
345                raise RuntimeError, msg % (
346                  len(match_indices),
347                  index_template,
348                  cost_variable_name,
349                  self._reference_instance.name
350                )
351            elif len(match_indices) == 0:
352                msg = 'Stage cost index not found: %s[%s]\n'                  \
353                      'Do you have an off-by-one miscalculation, or did you ' \
354                      'forget to specify it in ReferenceModel.dat?'
355                raise RuntimeError, msg % ( cost_variable_name, index_template )
356
357            cost_variable_index = match_indices[0]
358
359         else:
360
361            cost_variable_name = cost_variable_string
362
363            # validate that the variable exists and extract the reference
364            if cost_variable_name not in self._reference_instance.active_components(Var):
365               raise ValueError, "Cost variable=" + cost_variable_name + " associated with stage=" + stage_id + " is not present in model=" + self._reference_instance.name
366            cost_variable = self._reference_instance.active_components(Var)[cost_variable_name]
367
368         # store the validated info.
369         stage._cost_variable = (cost_variable, cost_variable_index)
370
371   """ Constructor
372       Arguments:
373           scenarioinstance     - the reference (deterministic) scenario instance.
374           scenariotreeinstance - the pyomo model specifying all scenario tree (text) data.
375           scenariobundlelist   - a list of scenario names to retain, i.e., cull the rest to create a reduced tree!
376   """
377   def __init__(self, *args, **kwds):
378
379      self._name = None # some arbitrary identifier
380      self._reference_instance = kwds.pop( 'scenarioinstance', None ) # the reference (deterministic) base model
381
382      # the core objects defining the scenario tree.
383      self._tree_nodes = [] # collection of ScenarioTreeNodes
384      self._stages = [] # collection of Stages - assumed to be in time-order. the set (provided by the user) itself *must* be ordered.
385      self._scenarios = [] # collection of Scenarios
386
387      # dictionaries for the above.
388      self._tree_node_map = {}
389      self._stage_map = {}
390      self._scenario_map = {}
391
392      # a boolean indicating how data for scenario instances is specified.
393      # possibly belongs elsewhere, e.g., in the PH algorithm.
394      self._scenario_based_data = None
395
396      scenario_tree_instance = kwds.pop( 'scenariotreeinstance', None )
397      scenario_bundle_list = kwds.pop( 'scenariobundlelist', None )
398
399      # process the keyword options
400      for key in kwds.keys():
401         print >>sys.stderr, "Unknown option '%s' specified in call to ScenarioTree constructor" % key
402
403      if self._reference_instance is None:
404         raise ValueError, "A reference scenario instance must be supplied in the ScenarioTree constructor"
405
406      if scenario_tree_instance is None:
407         raise ValueError, "A scenario tree instance must be supplied in the ScenarioTree constructor"
408
409      node_ids = scenario_tree_instance.Nodes
410      node_child_ids = scenario_tree_instance.Children
411      node_stage_ids = scenario_tree_instance.NodeStage
412      node_probability_map = scenario_tree_instance.ConditionalProbability
413      stage_ids = scenario_tree_instance.Stages
414      stage_variable_ids = scenario_tree_instance.StageVariables
415      stage_cost_variable_ids = scenario_tree_instance.StageCostVariable
416      scenario_ids = scenario_tree_instance.Scenarios
417      scenario_leaf_ids = scenario_tree_instance.ScenarioLeafNode
418      scenario_based_data = scenario_tree_instance.ScenarioBasedData
419
420      # save the method for instance data storage.
421      self._scenario_based_data = scenario_based_data()
422
423      # the input stages must be ordered, for both output purposes and knowledge of the final stage.
424      if stage_ids.ordered is False:
425         raise ValueError, "An ordered set of stage IDs must be supplied in the ScenarioTree constructor"
426
427      #
428      # construct the actual tree objects
429      #
430
431      # construct the stage objects w/o any linkages first; link them up
432      # with tree nodes after these have been fully constructed.
433      self._construct_stages(stage_ids, stage_variable_ids, stage_cost_variable_ids)
434
435      # construct the tree node objects themselves in a first pass,
436      # and then link them up in a second pass to form the tree.
437      # can't do a single pass because the objects may not exist.
438      for tree_node_name in node_ids:
439
440         if tree_node_name not in node_stage_ids:
441            raise ValueError, "No stage is assigned to tree node=" + tree_node._name
442
443         stage_name = value(node_stage_ids[tree_node_name])
444         if stage_name not in self._stage_map.keys():
445            raise ValueError, "Unknown stage=" + stage_name + " assigned to tree node=" + tree_node._name
446
447         new_tree_node = ScenarioTreeNode(tree_node_name,
448                                          value(node_probability_map[tree_node_name]),
449                                          self._stage_map[stage_name],
450                                          self._reference_instance)
451
452         self._tree_nodes.append(new_tree_node)
453         self._tree_node_map[tree_node_name] = new_tree_node
454         self._stage_map[stage_name]._tree_nodes.append(new_tree_node)
455
456      # link up the tree nodes objects based on the child id sets.
457      for this_node in self._tree_nodes:
458         this_node._children = []
459         if this_node._name in node_child_ids: # otherwise, you're at a leaf and all is well.
460            child_ids = node_child_ids[this_node._name]
461            for child_id in child_ids:
462               if child_id in self._tree_node_map.keys():
463                  child_node = self._tree_node_map[child_id]
464                  this_node._children.append(self._tree_node_map[child_id])
465                  if child_node._parent is None:
466                     child_node._parent = this_node
467                  else:
468                     raise ValueError, "Multiple parents specified for tree node="+child_id+"; existing parent node="+child_node._parent._name+"; conflicting parent node="+this_node._name
469               else:
470                  raise ValueError, "Unknown child tree node=" + child_id + " specified for tree node=" + this_node._name
471
472      # at this point, the scenario tree nodes and the stages are set - no
473      # two-pass logic necessary when constructing scenarios.
474      for scenario_name in scenario_ids:
475
476         # IMPT: the name of the scenario is assumed to have no '_' (underscore) characters in the identifier.
477         #       this is required when writing the extensive form, e.g., the scenario is used in the extensive
478         #       form as a prefix on variable and constraint names. this convention simplifies parsing on the
479         #       back end; if the underscore isn't used as a reserved character, then some other separator
480         #       symbol would be required, or we would have to engage in some complex prefix matching with
481         #       all possible scenario names.
482         if string.find(scenario_name, "_") != -1:
483            raise ValueError, "By convention, scenario names in PySP cannot contain underscore (_) characters; the scenario in violation="+scenario_name
484
485         new_scenario = Scenario()
486         new_scenario._name=scenario_name
487
488         if scenario_name not in scenario_leaf_ids.keys():
489            raise ValueError, "No leaf tree node specified for scenario=" + scenario_name
490         else:
491            scenario_leaf_node_name = value(scenario_leaf_ids[scenario_name])
492            if scenario_leaf_node_name not in self._tree_node_map.keys():
493               raise ValueError, "Uknown tree node=" + scenario_leaf_node_name + " specified as leaf of scenario=" + scenario_name
494            else:
495               new_scenario._leaf_node = self._tree_node_map[scenario_leaf_node_name]
496
497         current_node = new_scenario._leaf_node
498         probability = 1.0
499         while current_node is not None:
500            new_scenario._node_list.append(current_node)
501            current_node._scenarios.append(new_scenario) # links the scenarios to the nodes to enforce necessary non-anticipativity
502            probability *= current_node._conditional_probability
503            current_node = current_node._parent
504         new_scenario._node_list.reverse()
505         new_scenario._probability = probability
506
507         self._scenarios.append(new_scenario)
508         self._scenario_map[scenario_name] = new_scenario
509
510      # for output purposes, it is useful to known the maximal length of identifiers
511      # in the scenario tree for any particular category. I'm building these up
512      # incrementally, as they are needed. 0 indicates unassigned.
513      self._max_scenario_id_length = 0
514
515      # does the actual traversal to populate the members.
516      self.computeIdentifierMaxLengths()
517
518      # if a sub-bundle of scenarios has been specified, mark the
519      # active scenario tree components and compress the tree.
520      if scenario_bundle_list is not None:
521         self.compress(scenario_bundle_list)
522
523   #
524   # given a set of scenario instances, compute the set of variable indices being blended at each node.
525   # this can't be done until the scenario instances are available, as different scenarios can have
526   # different index sets.
527   #
528
529   def defineVariableIndexSets(self, scenario_instances):
530
531      for tree_node in self._tree_nodes:
532
533         tree_node.defineVariableIndexSets(scenario_instances)
534
535   #
536   # is the indicated scenario in the tree?
537   #
538
539   def contains_scenario(self, name):
540      return name in self._scenario_map.keys()
541
542   #
543   # get the scenario object from the tree.
544   #
545
546   def get_scenario(self, name):
547      return self._scenario_map[name]
548
549   #
550   # compute the scenario cost for the input instance, i.e.,
551   # the sum of all stage cost variables.
552   #
553
554   def compute_scenario_cost(self, instance):
555      aggregate_cost = 0.0
556      for stage in self._stages:
557            instance_cost_variable = instance.active_components(Var)[stage._cost_variable[0].name][stage._cost_variable[1]]
558            if instance_cost_variable.value is not None:
559                # depending on the situation (and the model), the stage cost variables might not have values.
560               aggregate_cost += instance_cost_variable
561      return aggregate_cost
562
563   #
564   # utility for compressing or culling a scenario tree based on
565   # a provided list of scenarios (specified by name) to retain -
566   # all non-referenced components are eliminated. this particular
567   # method compresses *in-place*, i.e., via direct modification
568   # of the scenario tree structure.
569   #
570
571   def compress(self, scenario_bundle_list):
572
573      # scan for and mark all referenced scenarios and
574      # tree nodes in the bundle list - all stages will
575      # obviously remain.
576      for scenario_name in scenario_bundle_list:
577         if scenario_name not in self._scenario_map:
578            raise ValueError, "Scenario="+scenario_name+" selected for bundling not present in scenario tree"
579         scenario = self._scenario_map[scenario_name]
580         scenario.retain = True
581
582         # chase all nodes comprising this scenario,
583         # marking them for retention.
584         for node in scenario._node_list:
585            node.retain = True
586
587      # scan for any non-retained scenarios and tree nodes.
588      scenarios_to_delete = []
589      tree_nodes_to_delete = []
590      for scenario in self._scenarios:
591         if hasattr(scenario, "retain") is True:
592            delattr(scenario, "retain")
593            pass
594         else:
595            scenarios_to_delete.append(scenario)
596            del self._scenario_map[scenario._name]
597
598      for tree_node in self._tree_nodes:
599         if hasattr(tree_node, "retain") is True:
600            delattr(tree_node, "retain")
601            pass
602         else:
603            tree_nodes_to_delete.append(tree_node)
604            del self._tree_node_map[tree_node._name]
605
606      # JPW does not claim the following routines are
607      # the most efficient. rather, they get the job
608      # done while avoiding serious issues with
609      # attempting to remove elements from a list that
610      # you are iterating over.
611
612      # delete all references to unmarked scenarios
613      # and child tree nodes in the scenario tree node
614      # structures.
615      for tree_node in self._tree_nodes:
616         for scenario in scenarios_to_delete:
617            if scenario in tree_node._scenarios:
618               tree_node._scenarios.remove(scenario)
619         for node_to_delete in tree_nodes_to_delete:
620            if node_to_delete in tree_node._children:
621               tree_node._children.remove(node_to_delete)
622
623      # delete all references to unmarked tree nodes
624      # in the scenario tree stage structures.
625      for stage in self._stages:
626         for tree_node in tree_nodes_to_delete:
627            if tree_node in stage._tree_nodes:
628               stage._tree_nodes.remove(tree_node)
629
630      # delete all unreferenced entries from the core scenario
631      # tree data structures.
632      for scenario in scenarios_to_delete:
633         self._scenarios.remove(scenario)
634      for tree_node in tree_nodes_to_delete:
635         self._tree_nodes.remove(tree_node)
636
637
638      # re-normalize the conditional probabilities of the
639      # children at each tree node.
640      for tree_node in self._tree_nodes:
641         sum_child_probabilities = 0.0
642         for child_node in tree_node._children:
643            sum_child_probabilities += child_node._conditional_probability
644         for child_node in tree_node._children:
645            child_node._conditional_probability = child_node._conditional_probability / sum_child_probabilities
646
647      # re-compute the absolute scenario probabilities based
648      # on the re-normalized conditional node probabilities.
649      for scenario in self._scenarios:
650         probability = 1.0
651         for tree_node in scenario._node_list:
652            probability = probability * tree_node._conditional_probability
653         scenario._probability = probability
654
655   #
656   # returns the root node of the scenario tree
657   #
658
659   def findRootNode(self):
660
661      for tree_node in self._tree_nodes:
662         if tree_node._parent is None:
663            return tree_node
664      return None
665
666   #
667   # a utility function to compute, based on the current scenario tree content,
668   # the maximal length of identifiers in various categories.
669   #
670
671   def computeIdentifierMaxLengths(self):
672
673      self._max_scenario_id_length = 0
674      for scenario in self._scenarios:
675         if len(scenario._name) > self._max_scenario_id_length:
676            self._max_scenario_id_length = len(scenario._name)
677
678   #
679   # a utility function to (partially, at the moment) validate a scenario tree
680   #
681
682   def validate(self):
683
684      # for any node, the sum of conditional probabilities of the children should sum to 1.
685      for tree_node in self._tree_nodes:
686         sum_probabilities = 0.0
687         if len(tree_node._children) > 0:
688            for child in tree_node._children:
689               sum_probabilities += child._conditional_probability
690            if abs(1.0 - sum_probabilities) > 0.000001:
691               print "The child conditional probabilities for tree node=" + tree_node._name + " sum to " + `sum_probabilities`
692               return False
693
694      # ensure that there is only one root node in the tree
695      num_roots = 0
696      root_ids = []
697      for tree_node in self._tree_nodes:
698         if tree_node._parent is None:
699            num_roots += 1
700            root_ids.append(tree_node._name)
701
702      if num_roots != 1:
703         print "Illegal set of root nodes detected: " + `root_ids`
704         return False
705
706      # there must be at least one scenario passing through each tree node.
707      for tree_node in self._tree_nodes:
708         if len(tree_node._scenarios) == 0:
709            print "There are no scenarios associated with tree node=" + tree_node._name
710            return False
711
712      return True
713
714   #
715   # copies the parameter values stored in any tree node _averages attribute
716   # into any tree node _solutions attribute - only for active variable values.
717   #
718
719   def snapshotSolutionFromAverages(self):
720
721      for tree_node in self._tree_nodes:
722
723         tree_node.snapshotSolutionFromAverages()
724
725   #
726   # computes the variable values at each tree node from the input scenario instances.
727   #
728
729   def snapshotSolutionFromInstances(self, scenario_instance_map):
730
731      for tree_node in self._tree_nodes:
732
733         tree_node.snapshotSolutionFromInstances(scenario_instance_map)
734
735   #
736   # a utility to determine the stage to which the input variable belongs.
737   # this is horribly inefficient, lacking an inverse map. fortunately,
738   # it isn't really called that often (yet). stage membership is determined
739   # by comparing the input variable name with the reference instance
740   # variable name (which is what the scenario tree refers to) and the
741   # associated indices.
742   #
743
744   def variableStage(self, variable, index):
745
746      for stage in self._stages:
747         for (stage_var, match_template) in stage._variables:
748            if (variable.name == stage_var.name) and (index in match_indices):
749               return stage
750
751         # IMPT: this is a temporary hack - the real fix is to force users to
752         # have every variable assigned to some stage in the model, either
753         # automatically or explicitly.
754         if (variable.name == stage._cost_variable[0].name):
755            return stage
756
757      raise RuntimeError, "The variable="+str(variable.name)+", index="+str(index)+" does not belong to any stage in the scenario tree"
758
759   #
760   # a utility to determine the stage to which the input constraint "belongs".
761   # a constraint belongs to the latest stage in which referenced variables
762   # in the constraint appears in that stage.
763   # input is a constraint is of type "Constraint", and an index of that
764   # constraint - which might be None in the case of singleton constraints.
765   # currently doesn't deal with SOS constraints, for no real good reason.
766   # returns an instance of a Stage object.
767   # IMPT: this method works on the canonical representation ("repn" attribute)
768   #       of a constraint. this implies that pre-processing of the instance
769   #       has been performed.
770   # NOTE: there is still the issue of whether the contained variables really
771   #       belong to the same model, but that is a different issue we won't
772   #       address right now (e.g., what does it mean for a constraint in an
773   #       extensive form binding instance to belong to a stage?).
774   #
775
776   def constraintStage(self, constraint, index):
777
778      largest_stage_index = -1
779      largest_stage = None
780
781      canonical_repn = constraint[index].repn
782      for degree, terms in canonical_repn.items():
783         if (degree != -1) and (degree != 0): # ignore constant terms and the variable definitions themselves.
784            for var_key, coefficient in terms.items():
785               var_value = canonical_repn[-1][var_key.keys()[0]]
786               var_stage = self.variableStage(var_value.var, var_value.index)
787               var_stage_index = self._stages.index(var_stage)
788
789               if var_stage_index > largest_stage_index:
790                  largest_stage_index = var_stage_index
791                  largest_stage = var_stage
792
793      return largest_stage
794
795   #
796   # a utility function to pretty-print the static/non-cost information associated with a scenario tree
797   #
798
799   def pprint(self):
800
801      print "Scenario Tree Detail"
802
803      print "----------------------------------------------------"
804      if self._reference_instance is not None:
805         print "Model=" + self._reference_instance.name
806      else:
807         print "Model=" + "Unassigned"
808      print "----------------------------------------------------"
809      print "Tree Nodes:"
810      print ""
811      for tree_node_name in sorted(self._tree_node_map.keys()):
812         tree_node = self._tree_node_map[tree_node_name]
813         print "\tName=" + tree_node_name
814         if tree_node._stage is not None:
815            print "\tStage=" + str(tree_node._stage._name)
816         else:
817            print "\t Stage=None"
818         if tree_node._parent is not None:
819            print "\tParent=" + tree_node._parent._name
820         else:
821            print "\tParent=" + "None"
822         if tree_node._conditional_probability is not None:
823            print "\tConditional probability=%4.4f" % tree_node._conditional_probability
824         else:
825            print "\tConditional probability=" + "***Undefined***"
826         print "\tChildren:"
827         if len(tree_node._children) > 0:
828            for child_node in sorted(tree_node._children, cmp=lambda x,y: cmp(x._name, y._name)):
829               print "\t\t" + child_node._name
830         else:
831            print "\t\tNone"
832         print "\tScenarios:"
833         if len(tree_node._scenarios) == 0:
834            print "\t\tNone"
835         else:
836            for scenario in sorted(tree_node._scenarios, cmp=lambda x,y: cmp(x._name, y._name)):
837               print "\t\t" + scenario._name
838         print ""
839      print "----------------------------------------------------"
840      print "Stages:"
841      for stage_name in sorted(self._stage_map.keys()):
842         stage = self._stage_map[stage_name]
843         print "\tName=" + str(stage_name)
844         print "\tTree Nodes: "
845         for tree_node in sorted(stage._tree_nodes, cmp=lambda x,y: cmp(x._name, y._name)):
846            print "\t\t" + tree_node._name
847         print "\tVariables: "
848         for (variable, index_template) in stage._variables:
849            print "\t\t",variable.name,":",index_template
850         print "\tCost Variable: "
851         if stage._cost_variable[1] is None:
852            print "\t\t" + stage._cost_variable[0].name
853         else:
854            print "\t\t" + stage._cost_variable[0].name + indexToString(stage._cost_variable[1])
855         print ""
856      print "----------------------------------------------------"
857      print "Scenarios:"
858      for scenario_name in sorted(self._scenario_map.keys()):
859         scenario = self._scenario_map[scenario_name]
860         print "\tName=" + scenario_name
861         print "\tProbability=%4.4f" % scenario._probability
862         if scenario._leaf_node is None:
863            print "\tLeaf node=None"
864         else:
865            print "\tLeaf node=" + scenario._leaf_node._name
866         print "\tTree node sequence:"
867         for tree_node in scenario._node_list:
868            print "\t\t" + tree_node._name
869         print ""
870      print "----------------------------------------------------"
871
872   #
873   # a utility function to pretty-print the solution associated with a scenario tree
874   #
875
876   def pprintSolution(self, epsilon=1.0e-5):
877
878      print "----------------------------------------------------"
879      print "Tree Nodes:"
880      print ""
881      for tree_node_name in sorted(self._tree_node_map.keys()):
882         tree_node = self._tree_node_map[tree_node_name]
883         print "\tName=" + tree_node_name
884         if tree_node._stage is not None:
885            print "\tStage=" + tree_node._stage._name
886         else:
887            print "\t Stage=None"
888         if tree_node._parent is not None:
889            print "\tParent=" + tree_node._parent._name
890         else:
891            print "\tParent=" + "None"
892         print "\tVariables: "
893         for (variable, index_template) in tree_node._stage._variables:
894            indices = tree_node._variable_indices[variable.name]
895            solution_variable = tree_node._solutions[variable.name]
896            if (len(indices) == 1) and (indices[0] == None):
897               # if this is a singleton variable, then it should necessarily be active -
898               # otherwise, it wouldn't be referenced in the stage!!!
899               value = solution_variable[None].value
900               if fabs(value) > epsilon:
901                  print "\t\t"+variable.name+"="+str(value)
902            else:
903               for index in indices:
904                  if (solution_variable[index].active is True) and (index in solution_variable):
905                     value = solution_variable[index].value
906                     if (value is not None) and (fabs(value) > epsilon):
907                        print "\t\t"+variable.name+indexToString(index)+"="+str(value)
908         print ""
909
910   #
911   # a utility function to pretty-print the cost information associated with a scenario tree
912   #
913
914   def pprintCosts(self, scenario_instance_map):
915
916      print "Scenario Tree Costs"
917      print "***WARNING***: Assumes full (or nearly so) convergence of scenario solutions at each node in the scenario tree - computed costs are invalid otherwise"
918
919      print "----------------------------------------------------"
920      if self._reference_instance is not None:
921         print "Model=" + self._reference_instance.name
922      else:
923         print "Model=" + "Unassigned"
924
925      print "----------------------------------------------------"
926      print "Tree Nodes:"
927      print ""
928      for tree_node_name in sorted(self._tree_node_map.keys()):
929         tree_node = self._tree_node_map[tree_node_name]
930         print "\tName=" + tree_node_name
931         if tree_node._stage is not None:
932            print "\tStage=" + tree_node._stage._name
933         else:
934            print "\t Stage=None"
935         if tree_node._parent is not None:
936            print "\tParent=" + tree_node._parent._name
937         else:
938            print "\tParent=" + "None"
939         if tree_node._conditional_probability is not None:
940            print "\tConditional probability=%4.4f" % tree_node._conditional_probability
941         else:
942            print "\tConditional probability=" + "***Undefined***"
943         print "\tChildren:"
944         if len(tree_node._children) > 0:
945            for child_node in sorted(tree_node._children, cmp=lambda x,y: cmp(x._name, y._name)):
946               print "\t\t" + child_node._name
947         else:
948            print "\t\tNone"
949         print "\tScenarios:"
950         if len(tree_node._scenarios) == 0:
951            print "\t\tNone"
952         else:
953            for scenario in sorted(tree_node._scenarios, cmp=lambda x,y: cmp(x._name, y._name)):
954               print "\t\t" + scenario._name
955         print "\tExpected node cost=%10.4f" % tree_node.computeExpectedNodeCost(scenario_instance_map)
956         print ""
957
958      print "----------------------------------------------------"
959      print "Scenarios:"
960      print ""
961      for scenario_name, scenario in sorted(self._scenario_map.iteritems()):
962         instance = scenario_instance_map[scenario_name]
963
964         print "\tName=" + scenario_name
965         print "\tProbability=%4.4f" % scenario._probability
966
967         if scenario._leaf_node is None:
968            print "\tLeaf Node=None"
969         else:
970            print "\tLeaf Node=" + scenario._leaf_node._name
971
972         print "\tTree node sequence:"
973         for tree_node in scenario._node_list:
974            print "\t\t" + tree_node._name
975
976         aggregate_cost = 0.0
977         for stage in self._stages:
978            instance_cost_variable = instance.active_components(Var)[stage._cost_variable[0].name][stage._cost_variable[1]]
979            if instance_cost_variable.value is not None:
980               print "\tStage=%20s     Cost=%10.4f" % (stage._name, instance_cost_variable.value)
981               cost = instance_cost_variable.value
982            else:
983               print "\tStage=%20s     Cost=%10s" % (stage._name, "Not Rprted.")
984               cost = 0.0
985            aggregate_cost += cost
986         print "\tTotal scenario cost=%10.4f" % aggregate_cost
987         print ""
988      print "----------------------------------------------------"
989
Note: See TracBrowser for help on using the repository browser.