source: branches/collapse/collapse_code.py @ 855

Revision 855, 7.4 KB checked in by mmckerns, 5 months ago (diff)

updated copyright to 2016
Line 
1#!/usr/bin/env python
2#
3# Author: Lan Huong Nguyen (lanhuong @stanford)
4# Copyright (c) 2012-2016 California Institute of Technology.
5# License: 3-clause BSD.  The full license text is available at:
6#  - http://mmckerns.github.io/project/mystic/browser/mystic/LICENSE
7
8################################################################################
9
10                           ### COLLAPSE CODE ###                   
11
12################################################################################
13MINMAX = -1  ## NOTE: sup = maximize = -1; inf = minimize = 1 
14debug = True #False
15from mystic.termination import Or
16from mystic.termination import ChangeOverGeneration as COG
17from collapse_termination import CLPS_Weight as WCLPS
18from collapse_termination import CLPS_Position as PCLPS
19###################### FUNCTION: TERMINATION CONDITIONS ########################
20
21def termCond(npts, tol, tolW, tolP, ngen, ngcol, indexW_ids, indexP_ids):
22   
23  x = []
24  _COG = COG(tol,ngen) 
25  _WCLPS = []
26  _PCLPS = []
27 
28  for m in range(len(indexW_ids)):
29      for idx in indexW_ids[m]:
30          _WCLPS.append(WCLPS(npts, tol_weight = tolW, generations=ngcol, \
31                              axisW = m, indexW = idx))
32         
33  for m in range(len(indexP_ids)):
34      for pair in indexP_ids[m]:
35          _PCLPS.append( PCLPS(npts, tol_position=tolP, generations=ngcol, \
36                                axisP = m, index1 = pair[0], index2  = pair[1]))
37  x.extend(_WCLPS)
38  x.extend(_PCLPS)
39  x.append(_COG)
40  trm = Or(*x)
41  return trm
42 
43#################### FUNCTIONs: NEW CONSTRAINTS AFTER COLLAPSE #################
44
45def _set_weight_to_zero(rv, npts, axisW, indexW):
46    from mystic.math.discrete import product_measure
47    c = product_measure()
48    c.load(rv, npts)
49    c[axisW][indexW].weight = 0
50    return c.flatten()
51 
52def _set_same_positions(rv, npts, axisP, index1, index2):
53    from mystic.math.discrete import product_measure
54    c = product_measure()
55    c.load(rv, npts)
56    c[axisP][index1].weight += c[axisP][index2].weight
57    c[axisP][index2].weight = 0
58    c[axisP][index2].position = c[axisP][index1].position
59    return c.flatten()
60 
61def new_constraints(rv, npts, _constraints, list_collapsed_pairs, list_vanished_weights):
62    for ele in list_vanished_weights:
63      rv = _set_weight_to_zero(rv, npts, axisW = ele[0], indexW = ele[1])
64    for mem in list_collapsed_pairs:
65      rv = _set_same_positions(rv, npts, axisP = mem[0], \
66              index1 = mem[1][0], index2 = mem[1][1])           
67    return _constraints(rv)   
68 
69################################################################################
70
71                               ### SOLVER  ###
72
73############################## RUN THE SOLVER ##################################
74def collapseSolver(solver, cost, npts, lb, ub, tol, tolW, tolP_range_fr, ngen, ngcol, \
75                   Best1Exp, crossover,  percent_change, stepmon, _constraints): 
76   
77  from mystic.math.measures import _nested_split
78 
79  ############################ TOLERANCE BOUNDS ################################
80  w_lb, x_lb = _nested_split(lb, npts)
81  w_ub, x_ub = _nested_split(ub, npts)
82  range_x = []
83 
84  for i in range(len(x_lb)):
85      range_x.append(abs(x_lb[i][0] - x_ub[i][0])) 
86  tolP = [ele*tolP_range_fr for ele in range_x]
87  if debug: print 'tolP = %s' %tolP 
88 
89  ######################### TERMINATION CHECK LIST #############################
90 
91  indexW_ids = [range(nx) for nx in npts]  # list of indices for CLPS_Weight termination check
92  indexP_ids = []  # list of indices for CLPS_Position termination check 
93 
94  for i in range(len(npts)):
95      c = []
96      for idx in range(npts[i]):
97          for jdx in range(idx,npts[i]):
98              if idx !=jdx:
99                  c.append([idx, jdx])
100      indexP_ids.append(c)   
101     
102  ####################### NEW CONDSTRANIS COLLAPSE LIST ########################
103  list_vanished_weights = [] # list of [measure, index] lists
104  list_collapsed_pairs = [] # list of [measure, pair] lists, where pair = [index1, index2]
105 
106  ##############################################################################
107  collapse_data = {} # dictoranty {iteration number: termination messages}
108  tot_evaluations = 0 # total number of function eveluations
109  _COG = COG(tol,ngen)
110  trm = termCond(npts, tol, tolW, tolP, ngen, ngcol, indexW_ids, indexP_ids)
111 
112  while _COG(solver) != True:         
113    solver.Solve(cost,termination=trm,strategy=Best1Exp, \
114               CrossProbability=crossover,ScalingFactor=percent_change)
115    tot_evaluations += solver.evaluations   
116    collapse_data[stepmon._step] = trm(solver, True)
117   
118    if debug:
119      print 'last_itr_evaluations = %s' %solver.evaluations
120      print "Generation  %d " %stepmon._step
121   
122    ##################### PARSING TERMINATION MESSAGE ##########################
123    int_trm = trm(solver, info = 'self') 
124   
125    for ele in range(len(int_trm)):
126      cond = int_trm[ele]
127      string = cond(solver, info = True)
128      string1 = string.replace(',', ' ')
129      string2 = string1.replace('&', ' ')
130      string3 = string2.replace(':', ' ')
131      lst = string3.split()
132     
133      if lst[0] == 'CLPS_Weight':
134          for i in range(len(lst)):
135              if lst[i] == 'measure':
136                  axisW = int(lst[i+1])
137                  if lst[i+2] == 'index': indexW = int(lst[i+3])
138                 
139          if debug:
140            print 'axisW = %s' %axisW
141            print 'indexW = %s' %indexW
142         
143          ###################### UPDATING COLLAPSE LISTS #######################   
144         
145          list_vanished_weights.append([axisW, indexW])
146          indexW_ids[axisW].remove(indexW)
147                 
148      if lst[0] == 'CLPS_Position':
149          for i in range(len(lst)): 
150              if lst[i] == 'measure':       
151                  axisP = int(lst[i+1])
152                  if lst[i+2] == 'indices':
153                      index1 = int(lst[i+3].lstrip('('))
154                      index2 = int(lst[i+4].rstrip(')'))
155          if debug:
156            print "axisP = %s" %axisP
157            print "index1 = %s" %index1
158            print "index2 = %s" %index2
159   
160          ###################### UPDATING COLLAPSE LISTS #######################   
161          list_collapsed_pairs.append([axisP, [index1, index2]])
162          indexP_ids[axisP].remove([index1,index2])         
163     
164      ##### AFTER COLLAPSE OCCURRED CHANGE THE TERMINATION AND CONSTRAINTS #####
165      if debug:
166              print 'list_vanished_weights = '
167              print  list_vanished_weights
168              print 'list_collapsed_pairs = '
169              print  list_collapsed_pairs                     
170
171      def _new_constraints(rv):
172        return new_constraints(rv, npts, _constraints, list_collapsed_pairs, list_vanished_weights)     
173      solver.SetConstraints(_new_constraints)
174     
175      trm = termCond(npts, tol, tolW, tolP, ngen, ngcol, indexW_ids, indexP_ids)
176                             
177 
178  solved = solver.bestSolution
179 #print "solved: %s" % solver.Solution()
180  func_max = MINMAX * solver.bestEnergy       #NOTE: -solution assumes -Max
181 #func_max = 1.0 + MINMAX*solver.bestEnergy   #NOTE: 1-sol => 1-success = fail
182  func_evals = solver.evaluations
183  from mystic.munge import write_support_file
184  write_support_file(stepmon)
185  return solved, func_max, tot_evaluations, collapse_data
Note: See TracBrowser for help on using the repository browser.