Changeset 637


Ignore:
Timestamp:
01/07/13 15:18:17 (3 years ago)
Author:
mmckerns
Message:

add _map_solver attribute, merge trialPop to trialSolution (2D for map solver);
refactor diffev* solvers to allow stepping over a single iteration;
account for 2D trialSolution in termination and strategy;

Location:
mystic
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • mystic/mystic/abstract_map_solver.py

    r631 r637  
    107107       #self.signal_handler   = None 
    108108       #self._handle_sigint   = False 
     109        trialPop = [[0.0 for i in range(dim)] for j in range(self.nPop)] 
     110        self.trialSolution    = trialPop 
     111        self._map_solver      = True 
    109112 
    110113        # import 'map' defaults 

  • mystic/mystic/abstract_solver.py

    r636 r637  
    113113        self.population       = [[0.0 for i in range(dim)] for j in range(NP)] 
    114114        self.trialSolution    = [0.0] * dim 
     115        self._map_solver      = False 
    115116        self._bestEnergy      = None 
    116117        self._bestSolution    = None 

  • mystic/mystic/differential_evolution.py

    r636 r637  
    191191        return 
    192192 
     193    def _SetEvaluationLimits(self, iterscale=None, evalscale=None): 
     194        """set the evaluation limits""" 
     195        if iterscale is None: iterscale = 10 
     196        if evalscale is None: evalscale = 1000 
     197        # if SetEvaluationLimits not applied, use the solver default 
     198        if self._maxiter is None: 
     199            self._maxiter = self.nDim * self.nPop * iterscale 
     200        if self._maxfun is None: 
     201            self._maxfun = self.nDim * self.nPop * evalscale 
     202        return 
     203 
     204    def _Decorate(self, costfunction, ExtraArgs=None): 
     205        """decorate cost function with bounds, penalties, monitors, etc""" 
     206        self._fcalls, costfunction = wrap_function(costfunction, ExtraArgs, self._evalmon) 
     207        if self._useStrictRange: 
     208            for i in range(self.nPop): 
     209                self.population[i] = self._clipGuessWithinRangeBoundary(self.population[i]) 
     210            costfunction = wrap_bounds(costfunction, self._strictMin, self._strictMax) 
     211        costfunction = wrap_penalty(costfunction, self._penalty) 
     212        return costfunction 
     213 
     214    def Step(self, costfunction, strategy=None): 
     215        """perform a single optimization iteration""" 
     216        if not self.generations: # this is the first iteration 
     217            self.population[0] = asfarray(self.population[0]) 
     218            # decouple bestSolution from population and bestEnergy from popEnergy 
     219            self.bestSolution = self.population[0] 
     220            self.bestEnergy = self.popEnergy[0] 
     221 
     222        for candidate in range(self.nPop): 
     223            if not self.generations: 
     224                # generate trialSolution (within valid range) 
     225                self.trialSolution[:] = self.population[candidate] 
     226            if strategy: 
     227                # generate trialSolution (within valid range) 
     228                strategy(self, candidate) 
     229            # apply constraints 
     230            self.trialSolution[:] = self._constraints(self.trialSolution) 
     231            # apply penalty 
     232           #trialEnergy = self._penalty(self.trialSolution) 
     233            # calculate cost 
     234            trialEnergy = costfunction(self.trialSolution) 
     235 
     236            if trialEnergy < self.popEnergy[candidate]: 
     237                # New low for this candidate 
     238                self.popEnergy[candidate] = trialEnergy 
     239                self.population[candidate][:] = self.trialSolution 
     240                self.UpdateGenealogyRecords(candidate, self.trialSolution[:]) 
     241 
     242                # Check if all-time low 
     243                if trialEnergy < self.bestEnergy: 
     244                    self.bestEnergy = trialEnergy 
     245                    self.bestSolution[:] = self.trialSolution 
     246 
     247        # log bestSolution and bestEnergy (includes penalty) 
     248        self._stepmon(self.bestSolution[:], self.bestEnergy, self.id) 
     249        return 
     250 
     251    def _process_inputs(self, kwds): 
     252        """process and activate input settings""" 
     253        #allow for inputs that don't conform to AbstractSolver interface 
     254        from mystic.strategy import Best1Bin 
     255        strategy=Best1Bin    #mutation strategy (see mystic.strategy) 
     256        callback=None        #user-supplied function, called after each step 
     257        disp=0               #non-zero to print convergence messages 
     258        probability=0.9      #potential for parameter cross-mutation 
     259        scale=0.8            #multiplier for mutation impact 
     260        if not kwds.has_key('strategy'): kwds['strategy'] = strategy 
     261        if not kwds.has_key('callback'): kwds['callback'] = callback 
     262        if not kwds.has_key('disp'): kwds['disp'] = disp 
     263        self.probability = kwds.pop('CrossProbability', probability) 
     264        self.scale = kwds.pop('ScalingFactor', scale) 
     265        # backward compatibility 
     266        if kwds.has_key('EvaluationMonitor'): \ 
     267           self.SetEvaluationMonitor(kwds.pop('EvaluationMonitor')) 
     268        if kwds.has_key('StepMonitor'): \ 
     269           self.SetGenerationMonitor(kwds.pop('StepMonitor')) 
     270        if kwds.has_key('penalty'): \ 
     271           self.SetPenalty(kwds.pop('penalty')) 
     272        if kwds.has_key('constraints'): \ 
     273           self.SetConstraints(kwds.pop('constraints')) 
     274        return kwds 
     275 
    193276    def Solve(self, costfunction, termination, sigint_callback=None, 
    194277                                               ExtraArgs=(), **kwds): 
    … …  
    223306    disp -- non-zero to print convergence messages. 
    224307        """ 
    225         #allow for inputs that don't conform to AbstractSolver interface 
    226         self.population[0] = asfarray(self.population[0]) 
    227         from mystic.strategy import Best1Bin 
    228         strategy=Best1Bin    #mutation strategy (see mystic.strategy) 
    229         CrossProbability=0.9 #potential for parameter cross-mutation 
    230         ScalingFactor=0.8    #multiplier for mutation impact 
    231         callback=None        #user-supplied function, called after each step 
    232         disp=0               #non-zero to print convergence messages 
    233         if kwds.has_key('strategy'): strategy = kwds['strategy'] 
    234         if kwds.has_key('CrossProbability'): \ 
    235            CrossProbability = kwds['CrossProbability'] 
    236         if kwds.has_key('ScalingFactor'): ScalingFactor = kwds['ScalingFactor'] 
    237         if kwds.has_key('callback'): callback = kwds['callback'] 
    238         if kwds.has_key('disp'): disp = kwds['disp'] 
    239         # backward compatibility 
    240         if kwds.has_key('EvaluationMonitor'): \ 
    241            self.SetEvaluationMonitor(kwds['EvaluationMonitor']) 
    242         if kwds.has_key('StepMonitor'): \ 
    243            self.SetGenerationMonitor(kwds['StepMonitor']) 
    244         if kwds.has_key('penalty'): \ 
    245            self.SetPenalty(kwds['penalty']) 
    246         if kwds.has_key('constraints'): \ 
    247            self.SetConstraints(kwds['constraints']) 
    248         #------------------------------------------------------------- 
    249  
     308        # process and activate input settings 
     309        settings = {} 
     310        settings.update(self._process_inputs(kwds)) 
     311        disp = settings['disp'] 
     312        callback = settings['callback'] 
     313        strategy = settings['strategy'] #XXX: specific to diffev* ? 
     314 
     315        # set up signal handler 
    250316        import signal 
    251317        self._EARLYEXIT = False 
    252  
    253         self._fcalls, costfunction = wrap_function(costfunction, ExtraArgs, self._evalmon) 
    254         if self._useStrictRange: 
    255             for i in range(self.nPop): 
    256                 self.population[i] = self._clipGuessWithinRangeBoundary(self.population[i]) 
    257             costfunction = wrap_bounds(costfunction, self._strictMin, self._strictMax) 
    258         costfunction = wrap_penalty(costfunction, self._penalty) 
    259  
    260         #generate signal_handler 
    261318        self._generateHandler(sigint_callback)  
    262319        if self._handle_sigint: signal.signal(signal.SIGINT, self.signal_handler) 
    263320 
    264         id = self.id 
    265         self.probability = CrossProbability 
    266         self.scale = ScalingFactor 
    267  
    268         # decouple bestSolution from population and bestEnergy from popEnergy 
    269         self.bestSolution = self.population[0] 
    270         self.bestEnergy = self.popEnergy[0] 
    271  
    272         #set initial solution and energy by running a single iteration 
    273         for candidate in range(self.nPop): 
    274             # generate trialSolution (within valid range) 
    275             self.trialSolution[:] = self.population[candidate] 
    276             # apply constraints 
    277             self.trialSolution[:] = self._constraints(self.trialSolution) 
    278             # apply penalty 
    279            #trialEnergy = self._penalty(self.trialSolution) 
    280             # calculate cost 
    281             trialEnergy = costfunction(self.trialSolution) 
    282  
    283             if trialEnergy < self.popEnergy[candidate]: 
    284                 # New low for this candidate 
    285                 self.popEnergy[candidate] = trialEnergy 
    286                 self.population[candidate][:] = self.trialSolution 
    287                 self.UpdateGenealogyRecords(candidate, self.trialSolution[:]) 
    288  
    289                 # Check if all-time low 
    290                 if trialEnergy < self.bestEnergy: 
    291                     self.bestEnergy = trialEnergy 
    292                     self.bestSolution[:] = self.trialSolution 
    293  
    294         # log bestSolution and bestEnergy (includes penalty) 
    295         self._stepmon(self.bestSolution[:], self.bestEnergy, id) 
    296         termination(self) #XXX: initialize termination conditions, if needed 
     321        # decorate cost function with bounds, penalties, monitors, etc 
     322        costfunction = self._Decorate(costfunction, ExtraArgs) 
     323 
     324        # the initital optimization iteration 
     325        self.Step(costfunction) 
     326        # initialize termination conditions, if needed 
     327        termination(self) 
    297328        if callback is not None: 
    298329            callback(self.bestSolution) 
    299330          
    300         # if SetEvaluationLimits not applied, use the solver default 
    301         if self._maxiter is None: 
    302             self._maxiter = self.nDim * self.nPop * 10  #XXX: better defaults? 
    303         if self._maxfun is None: 
    304             self._maxfun = self.nDim * self.nPop * 1000 #XXX: better defaults? 
    305  
     331        # impose the evaluation limits 
     332        self._SetEvaluationLimits() 
     333 
     334        # the main optimization loop 
    306335        while not self._terminated(termination) and not self._EARLYEXIT: 
    307             for candidate in range(self.nPop): 
    308                 # generate trialSolution (within valid range) 
    309                 strategy(self, candidate) 
    310                 # apply constraints 
    311                 self.trialSolution[:] = self._constraints(self.trialSolution) 
    312                 # apply penalty 
    313                #trialEnergy = self._penalty(self.trialSolution) 
    314                 # calculate cost 
    315                 trialEnergy = costfunction(self.trialSolution) 
    316  
    317                 if trialEnergy < self.popEnergy[candidate]: 
    318                     # New low for this candidate 
    319                     self.popEnergy[candidate] = trialEnergy 
    320                     self.population[candidate][:] = self.trialSolution 
    321                     self.UpdateGenealogyRecords(candidate, self.trialSolution[:]) 
    322  
    323                     # Check if all-time low 
    324                     if trialEnergy < self.bestEnergy: 
    325                         self.bestEnergy = trialEnergy 
    326                         self.bestSolution[:] = self.trialSolution 
    327  
    328             # log bestSolution and bestEnergy (includes penalty) 
    329             self._stepmon(self.bestSolution[:], self.bestEnergy, id) 
     336            self.Step(costfunction, strategy=strategy) 
    330337            if callback is not None: 
    331338                callback(self.bestSolution) 
    332339 
     340        # handle signal interrupts 
    333341        signal.signal(signal.SIGINT,signal.default_int_handler) 
    334342 
    … …  
    371379        return 
    372380 
    373     def Solve(self, costfunction, termination, sigint_callback=None, 
    374                                                ExtraArgs=(), **kwds): 
    375         """Minimize a function using differential evolution. 
    376  
    377 Description: 
    378  
    379     Uses a differential evolution algorith to find the minimum of 
    380     a function of one or more variables. This implementation holds 
    381     the current generation invariant until the end of each iteration. 
    382  
    383 Inputs: 
    384  
    385     costfunction -- the Python function or method to be minimized. 
    386     termination -- callable object providing termination conditions. 
    387  
    388 Additional Inputs: 
    389  
    390     sigint_callback -- callback function for signal handler. 
    391     ExtraArgs -- extra arguments for func. 
    392  
    393 Further Inputs: 
    394  
    395     strategy -- the mutation strategy for generating new trial 
    396         solutions [default = Best1Bin] 
    397     CrossProbability -- the probability of cross-parameter mutations 
    398         [default = 0.9] 
    399     ScalingFactor -- multiplier for the impact of mutations on the 
    400         trial solution [default = 0.8] 
    401     callback -- an optional user-supplied function to call after each 
    402         iteration.  It is called as callback(xk), where xk is 
    403         the current parameter vector.  [default = None] 
    404     disp -- non-zero to print convergence messages. 
    405         """ 
    406         #allow for inputs that don't conform to AbstractSolver interface 
    407         self.population[0] = asfarray(self.population[0]) 
    408         from mystic.strategy import Best1Bin 
    409         strategy=Best1Bin    #mutation strategy (see mystic.strategy) 
    410         CrossProbability=0.9 #potential for parameter cross-mutation 
    411         ScalingFactor=0.8    #multiplier for mutation impact 
    412         callback=None        #user-supplied function, called after each step 
    413         disp=0               #non-zero to print convergence messages 
    414         if kwds.has_key('strategy'): strategy = kwds['strategy'] 
    415         if kwds.has_key('CrossProbability'): \ 
    416            CrossProbability = kwds['CrossProbability'] 
    417         if kwds.has_key('ScalingFactor'): ScalingFactor = kwds['ScalingFactor'] 
    418         if kwds.has_key('callback'): callback = kwds['callback'] 
    419         if kwds.has_key('disp'): disp = kwds['disp'] 
    420         # backward compatibility 
    421         if kwds.has_key('EvaluationMonitor'): \ 
    422            self.SetEvaluationMonitor(kwds['EvaluationMonitor']) 
    423         if kwds.has_key('StepMonitor'): \ 
    424            self.SetGenerationMonitor(kwds['StepMonitor']) 
    425         if kwds.has_key('penalty'): \ 
    426            self.SetPenalty(kwds['penalty']) 
    427         if kwds.has_key('constraints'): \ 
    428            self.SetConstraints(kwds['constraints']) 
    429         #------------------------------------------------------------- 
    430  
    431         import signal 
    432         self._EARLYEXIT = False 
    433  
     381    def _SetEvaluationLimits(self, iterscale=None, evalscale=None): 
     382        """set the evaluation limits""" 
     383        if iterscale is None: iterscale = 10 
     384        if evalscale is None: evalscale = 1000 
     385        # if SetEvaluationLimits not applied, use the solver default 
     386        if self._maxiter is None: 
     387            self._maxiter = self.nDim * self.nPop * iterscale 
     388        if self._maxfun is None: 
     389            self._maxfun = self.nDim * self.nPop * evalscale 
     390        return 
     391 
     392    def _Decorate(self, costfunction, ExtraArgs=None): 
     393        """decorate cost function with bounds, penalties, monitors, etc""" 
    434394       #FIXME: EvaluationMonitor fails for MPI, throws error for 'pp' 
    435395        from python_map import python_map 
    … …  
    443403            costfunction = wrap_bounds(costfunction, self._strictMin, self._strictMax) 
    444404        costfunction = wrap_penalty(costfunction, self._penalty) 
    445  
    446         #generate signal_handler 
    447         self._generateHandler(sigint_callback)  
    448         if self._handle_sigint: signal.signal(signal.SIGINT, self.signal_handler) 
    449  
    450         id = self.id 
    451         self.probability = CrossProbability 
    452         self.scale = ScalingFactor 
    453  
    454         # break link between population and popEnergy 
    455         self.bestSolution = self.population[0] 
    456         self.bestEnergy = self.popEnergy[0] 
    457  
    458         trialPop = [[0.0 for i in range(self.nDim)] for j in range(self.nPop)] 
    459         #set initial solution and energy by running a single iteration 
     405        return costfunction 
     406 
     407    def Step(self, costfunction, strategy=None): 
     408        """perform a single optimization iteration""" 
     409        if not self.generations: # this is the first iteration 
     410            self.population[0] = asfarray(self.population[0]) 
     411            # decouple bestSolution from population and bestEnergy from popEnergy 
     412            self.bestSolution = self.population[0] 
     413            self.bestEnergy = self.popEnergy[0] 
     414 
    460415        for candidate in range(self.nPop): 
    461             # generate trialSolution (within valid range) 
    462             self.trialSolution[:] = self.population[candidate] 
     416            if not self.generations: 
     417                # generate trialSolution (within valid range) 
     418                self.trialSolution[candidate][:] = self.population[candidate] 
     419            if strategy: 
     420                # generate trialSolution (within valid range) 
     421                strategy(self, candidate) 
    463422            # apply constraints 
    464             self.trialSolution[:] = self._constraints(self.trialSolution) 
    465             trialPop[candidate][:] = self.trialSolution 
     423            self.trialSolution[candidate][:] = self._constraints(self.trialSolution[candidate]) 
    466424 
    467425        mapconfig = dict(nnodes=self._nnodes, launcher=self._launcher, \ 
    … …  
    469427                         timelimit=self._timelimit, scheduler=self._scheduler, \ 
    470428                         ncpus=self._ncpus, servers=self._servers) 
     429 
    471430        # apply penalty 
    472        #trialEnergy = map(self._penalty, trialPop)#, **mapconfig) 
     431       #trialEnergy = map(self._penalty, self.trialSolution)#, **mapconfig) 
    473432        # calculate cost 
    474         trialEnergy = self._map(costfunction, trialPop, **mapconfig) 
     433        trialEnergy = self._map(costfunction, self.trialSolution, **mapconfig) 
    475434 
    476435        for candidate in range(self.nPop): 
    … …  
    478437                # New low for this candidate 
    479438                self.popEnergy[candidate] = trialEnergy[candidate] 
    480                 self.population[candidate][:] = trialPop[candidate] 
    481                 self.UpdateGenealogyRecords(candidate, trialPop[candidate][:]) 
     439                self.population[candidate][:] = self.trialSolution[candidate] 
     440                self.UpdateGenealogyRecords(candidate, self.trialSolution[candidate][:]) 
    482441 
    483442                # Check if all-time low 
    484443                if trialEnergy[candidate] < self.bestEnergy: 
    485444                    self.bestEnergy = trialEnergy[candidate] 
    486                     self.bestSolution[:] = trialPop[candidate] 
     445                    self.bestSolution[:] = self.trialSolution[candidate] 
    487446 
    488447        # log bestSolution and bestEnergy (includes penalty) 
    489448       #FIXME: StepMonitor works for 'pp'? 
    490         self._stepmon(self.bestSolution[:], self.bestEnergy, id) 
    491         termination(self) #XXX: initialize termination conditions, if needed 
     449        self._stepmon(self.bestSolution[:], self.bestEnergy, self.id) 
     450        return 
     451 
     452    def _process_inputs(self, kwds): 
     453        """process and activate input settings""" 
     454        #allow for inputs that don't conform to AbstractSolver interface 
     455        from mystic.strategy import Best1Bin 
     456        strategy=Best1Bin    #mutation strategy (see mystic.strategy) 
     457        callback=None        #user-supplied function, called after each step 
     458        disp=0               #non-zero to print convergence messages 
     459        probability=0.9      #potential for parameter cross-mutation 
     460        scale=0.8            #multiplier for mutation impact 
     461        if not kwds.has_key('strategy'): kwds['strategy'] = strategy 
     462        if not kwds.has_key('callback'): kwds['callback'] = callback 
     463        if not kwds.has_key('disp'): kwds['disp'] = disp 
     464        self.probability = kwds.pop('CrossProbability', probability) 
     465        self.scale = kwds.pop('ScalingFactor', scale) 
     466        # backward compatibility 
     467        if kwds.has_key('EvaluationMonitor'): \ 
     468           self.SetEvaluationMonitor(kwds.pop('EvaluationMonitor')) 
     469        if kwds.has_key('StepMonitor'): \ 
     470           self.SetGenerationMonitor(kwds.pop('StepMonitor')) 
     471        if kwds.has_key('penalty'): \ 
     472           self.SetPenalty(kwds.pop('penalty')) 
     473        if kwds.has_key('constraints'): \ 
     474           self.SetConstraints(kwds.pop('constraints')) 
     475        return kwds 
     476 
     477    def Solve(self, costfunction, termination, sigint_callback=None, 
     478                                               ExtraArgs=(), **kwds): 
     479        """Minimize a function using differential evolution. 
     480 
     481Description: 
     482 
     483    Uses a differential evolution algorith to find the minimum of 
     484    a function of one or more variables. This implementation holds 
     485    the current generation invariant until the end of each iteration. 
     486 
     487Inputs: 
     488 
     489    costfunction -- the Python function or method to be minimized. 
     490    termination -- callable object providing termination conditions. 
     491 
     492Additional Inputs: 
     493 
     494    sigint_callback -- callback function for signal handler. 
     495    ExtraArgs -- extra arguments for func. 
     496 
     497Further Inputs: 
     498 
     499    strategy -- the mutation strategy for generating new trial 
     500        solutions [default = Best1Bin] 
     501    CrossProbability -- the probability of cross-parameter mutations 
     502        [default = 0.9] 
     503    ScalingFactor -- multiplier for the impact of mutations on the 
     504        trial solution [default = 0.8] 
     505    callback -- an optional user-supplied function to call after each 
     506        iteration.  It is called as callback(xk), where xk is 
     507        the current parameter vector.  [default = None] 
     508    disp -- non-zero to print convergence messages. 
     509        """ 
     510        # process and activate input settings 
     511        settings = {} 
     512        settings.update(self._process_inputs(kwds)) 
     513        disp = settings['disp'] 
     514        callback = settings['callback'] 
     515        strategy = settings['strategy'] #XXX: specific to diffev* ? 
     516 
     517        # set up signal handler 
     518        import signal 
     519        self._EARLYEXIT = False 
     520        self._generateHandler(sigint_callback)  
     521        if self._handle_sigint: signal.signal(signal.SIGINT, self.signal_handler) 
     522 
     523        # decorate cost function with bounds, penalties, monitors, etc 
     524        costfunction = self._Decorate(costfunction, ExtraArgs) 
     525 
     526        # the initital optimization iteration 
     527        self.Step(costfunction) 
     528        # initialize termination conditions, if needed 
     529        termination(self) 
    492530        if callback is not None: 
    493531            callback(self.bestSolution) 
    494532          
    495         # if SetEvaluationLimits not applied, use the solver default 
    496         if self._maxiter is None: 
    497             self._maxiter = self.nDim * self.nPop * 10  #XXX: better defaults? 
    498         if self._maxfun is None: 
    499             self._maxfun = self.nDim * self.nPop * 1000 #XXX: better defaults? 
    500  
     533        # impose the evaluation limits 
     534        self._SetEvaluationLimits() 
     535 
     536        # the main optimization loop 
    501537        while not self._terminated(termination) and not self._EARLYEXIT: 
    502             for candidate in range(self.nPop): 
    503                 # generate trialSolution (within valid range) 
    504                 strategy(self, candidate) 
    505                 # apply constraints 
    506                 self.trialSolution[:] = self._constraints(self.trialSolution) 
    507                 trialPop[candidate][:] = self.trialSolution 
    508  
    509             # apply penalty 
    510            #trialEnergy = map(self._penalty, trialPop)#, **mapconfig) 
    511             # calculate cost 
    512             trialEnergy = self._map(costfunction, trialPop, **mapconfig) 
    513            ##XXX:[HACK] return trialPop b/c may be altered by constraints 
    514            #trials = map(costfunction, trialPop) 
    515            #trialEnergy = [i for i,j in trials] 
    516            #trialPop = [j for i,j in trials] #FIXME: breaks on inf (out of bounds) 
    517            ##XXX:[HACK] -end- 
    518  
    519             for candidate in range(self.nPop): 
    520                 if trialEnergy[candidate] < self.popEnergy[candidate]: 
    521                     # New low for this candidate 
    522                     self.popEnergy[candidate] = trialEnergy[candidate] 
    523                     self.population[candidate][:] = trialPop[candidate] 
    524                     self.UpdateGenealogyRecords(candidate, trialPop[candidate][:]) 
    525  
    526                     # Check if all-time low 
    527                     if trialEnergy[candidate] < self.bestEnergy: 
    528                         self.bestEnergy = trialEnergy[candidate] 
    529                         self.bestSolution[:] = trialPop[candidate] 
    530  
    531             # log bestSolution and bestEnergy (includes penalty) 
    532            #FIXME: StepMonitor works for 'pp'? 
    533             self._stepmon(self.bestSolution[:], self.bestEnergy, id) 
     538            self.Step(costfunction, strategy=strategy) 
    534539            if callback is not None: 
    535540                callback(self.bestSolution) 
    536541 
     542        # handle signal interrupts 
    537543        signal.signal(signal.SIGINT,signal.default_int_handler) 
    538544 

  • mystic/mystic/strategy.py

    r124 r637  
    3333    r1,r2 = get_random_candidates(inst.nPop, candidate, 2)  
    3434    n = random.randrange(inst.nDim) 
    35     
    36     inst.trialSolution[:] = inst.population[candidate][:] 
    37  
    38     i = 0 
    39     while 1: 
    40         if random.random() >= inst.probability or i == inst.nDim: 
    41             break 
    42         inst.trialSolution[n] = inst.bestSolution[n] + \ 
    43                                 inst.scale * (inst.population[r1][n] - \ 
    44                                               inst.population[r2][n]) 
     35 
     36    if inst._map_solver: 
     37        trialSolution = inst.trialSolution[candidate] 
     38    else: 
     39        trialSolution = inst.trialSolution 
     40    trialSolution[:] = inst.population[candidate] 
     41 
     42    i = 0 
     43    while 1: 
     44        if random.random() >= inst.probability or i == inst.nDim: 
     45            break 
     46        trialSolution[n] = inst.bestSolution[n] + \ 
     47                           inst.scale * (inst.population[r1][n] - \ 
     48                                         inst.population[r2][n]) 
    4549        n = (n + 1) % inst.nDim 
    4650        i += 1 
    … …  
    5862    r1,r2 = get_random_candidates(inst.nPop, candidate, 2)  
    5963 
    60     inst.trialSolution[:] = inst.population[candidate][:] 
     64    if inst._map_solver: 
     65        trialSolution = inst.trialSolution[candidate] 
     66    else: 
     67        trialSolution = inst.trialSolution 
     68    trialSolution[:] = inst.population[candidate] 
    6169 
    6270    # Randomly chosen index between [0, ND-1] (See Eq.4 of [1] ) 
    … …  
    6775        if i==n or cross < inst.probability: 
    6876            # this component of trial vector will come from vector v 
    69             inst.trialSolution[i] = inst.bestSolution[i] + \ 
    70                                     inst.scale * (inst.population[r1][i] - \ 
    71                                                   inst.population[r2][i]) 
     77            trialSolution[i] = inst.bestSolution[i] + \ 
     78                               inst.scale * (inst.population[r1][i] - \ 
     79                                             inst.population[r2][i]) 
    7280 
    7381#   inst._keepSolutionWithinRangeBoundary(inst.bestSolution) 
    … …  
    8290    n = random.randrange(inst.nDim) 
    8391 
    84     inst.trialSolution[:] = inst.population[candidate][:] 
    85  
    86     i = 0 
    87     while 1: 
    88         if random.random() >= inst.probability or i == inst.nDim: 
    89             break 
    90         inst.trialSolution[n] = inst.population[r1][n] + \ 
    91                                 inst.scale * (inst.population[r2][n] - \ 
    92                                               inst.population[r3][n]) 
     92    if inst._map_solver: 
     93        trialSolution = inst.trialSolution[candidate] 
     94    else: 
     95        trialSolution = inst.trialSolution 
     96    trialSolution[:] = inst.population[candidate] 
     97 
     98    i = 0 
     99    while 1: 
     100        if random.random() >= inst.probability or i == inst.nDim: 
     101            break 
     102        trialSolution[n] = inst.population[r1][n] + \ 
     103                           inst.scale * (inst.population[r2][n] - \ 
     104                                         inst.population[r3][n]) 
    93105        n = (n + 1) % inst.nDim 
    94106        i += 1 
    … …  
    108120    n = random.randrange(inst.nDim) 
    109121 
    110     inst.trialSolution[:] = inst.population[candidate][:] 
    111     i = 0 
    112     while 1: 
    113         if random.random() >= inst.probability or i == inst.nDim: 
    114             break 
    115         inst.trialSolution[n] += inst.scale * (inst.bestSolution[n] - \ 
    116                                                inst.trialSolution[n]) + \ 
    117                                  inst.scale * (inst.population[r1][n] - \ 
    118                                                inst.population[r2][n]) 
     122    if inst._map_solver: 
     123        trialSolution = inst.trialSolution[candidate] 
     124    else: 
     125        trialSolution = inst.trialSolution 
     126    trialSolution[:] = inst.population[candidate][:] 
     127    i = 0 
     128    while 1: 
     129        if random.random() >= inst.probability or i == inst.nDim: 
     130            break 
     131        trialSolution[n] += inst.scale * (inst.bestSolution[n] - \ 
     132                                          trialSolution[n]) + \ 
     133                            inst.scale * (inst.population[r1][n] - \ 
     134                                          inst.population[r2][n]) 
    119135        n = (n + 1) % inst.nDim 
    120136        i += 1 
    … …  
    131147    n = random.randrange(inst.nDim) 
    132148 
    133     inst.trialSolution[:] = inst.population[candidate][:] 
    134     i = 0 
    135     while 1: 
    136         if random.random() >= inst.probability or i == inst.nDim: 
    137             break 
    138         inst.trialSolution[n] = inst.bestSolution[n] + \ 
    139                                 inst.scale * (inst.population[r1][n] + \ 
    140                                               inst.population[r2][n] - \ 
    141                                               inst.population[r3][n] - \ 
    142                                               inst.population[r4][n]) 
     149    if inst._map_solver: 
     150        trialSolution = inst.trialSolution[candidate] 
     151    else: 
     152        trialSolution = inst.trialSolution 
     153    trialSolution[:] = inst.population[candidate] 
     154    i = 0 
     155    while 1: 
     156        if random.random() >= inst.probability or i == inst.nDim: 
     157            break 
     158        trialSolution[n] = inst.bestSolution[n] + \ 
     159                           inst.scale * (inst.population[r1][n] + \ 
     160                                         inst.population[r2][n] - \ 
     161                                         inst.population[r3][n] - \ 
     162                                         inst.population[r4][n]) 
    143163        n = (n + 1) % inst.nDim 
    144164        i += 1 
    … …  
    155175    n = random.randrange(inst.nDim) 
    156176 
    157     inst.trialSolution[:] = inst.population[candidate][:] 
    158     i = 0 
    159     while 1: 
    160         if random.random() >= inst.probability or i == inst.nDim: 
    161             break 
    162         inst.trialSolution[n] = inst.population[r1][n] + \ 
    163                                 inst.scale * (inst.population[r2][n] + \ 
    164                                               inst.population[r3][n] - \ 
    165                                               inst.population[r4][n] - \ 
    166                                               inst.population[r5][n]) 
     177    if inst._map_solver: 
     178        trialSolution = inst.trialSolution[candidate] 
     179    else: 
     180        trialSolution = inst.trialSolution 
     181    trialSolution[:] = inst.population[candidate] 
     182    i = 0 
     183    while 1: 
     184        if random.random() >= inst.probability or i == inst.nDim: 
     185            break 
     186        trialSolution[n] = inst.population[r1][n] + \ 
     187                           inst.scale * (inst.population[r2][n] + \ 
     188                                         inst.population[r3][n] - \ 
     189                                         inst.population[r4][n] - \ 
     190                                         inst.population[r5][n]) 
    167191        n = (n + 1) % inst.nDim 
    168192        i += 1 
    … …  
    179203    n = random.randrange(inst.nDim) 
    180204 
    181     inst.trialSolution[:] = inst.population[candidate][:] 
    182     i = 0 
    183     while 1: 
    184         if random.random() >= inst.probability or i == inst.nDim: 
    185             break 
    186         inst.trialSolution[n] = inst.population[r1][n] + \ 
    187                                 inst.scale * (inst.population[r2][n] -\ 
    188                                               inst.population[r3][n]) 
     205    if inst._map_solver: 
     206        trialSolution = inst.trialSolution[candidate] 
     207    else: 
     208        trialSolution = inst.trialSolution 
     209    trialSolution[:] = inst.population[candidate] 
     210    i = 0 
     211    while 1: 
     212        if random.random() >= inst.probability or i == inst.nDim: 
     213            break 
     214        trialSolution[n] = inst.population[r1][n] + \ 
     215                           inst.scale * (inst.population[r2][n] -\ 
     216                                         inst.population[r3][n]) 
    189217        n = (n + 1) % inst.nDim 
    190218        i += 1 
    … …  
    202230    n = random.randrange(inst.nDim) 
    203231 
    204     inst.trialSolution[:] = inst.population[candidate][:] 
    205     i = 0 
    206     while 1: 
    207         if random.random() >= inst.probability or i == inst.nDim: 
    208             break 
    209         inst.trialSolution[n] += inst.scale * (inst.bestSolution[n] - \ 
    210                                                inst.trialSolution[n]) + \ 
    211                                  inst.scale * (inst.population[r1][n] - \ 
    212                                                inst.population[r2][n]) 
     232    if inst._map_solver: 
     233        trialSolution = inst.trialSolution[candidate] 
     234    else: 
     235        trialSolution = inst.trialSolution 
     236    trialSolution[:] = inst.population[candidate] 
     237    i = 0 
     238    while 1: 
     239        if random.random() >= inst.probability or i == inst.nDim: 
     240            break 
     241        trialSolution[n] += inst.scale * (inst.bestSolution[n] - \ 
     242                                          trialSolution[n]) + \ 
     243                            inst.scale * (inst.population[r1][n] - \ 
     244                                          inst.population[r2][n]) 
    213245        n = (n + 1) % inst.nDim 
    214246        i += 1 
    … …  
    225257    n = random.randrange(inst.nDim) 
    226258 
    227     inst.trialSolution[:] = inst.population[candidate][:] 
    228     i = 0 
    229     while 1: 
    230         if random.random() >= inst.probability or i == inst.nDim: 
    231             break 
    232         inst.trialSolution[n] = inst.bestSolution[n] + \ 
    233                                 inst.scale * (inst.population[r1][n] + \ 
    234                                               inst.population[r2][n] - \ 
    235                                               inst.population[r3][n] - \ 
    236                                               inst.population[r4][n]) 
     259    if inst._map_solver: 
     260        trialSolution = inst.trialSolution[candidate] 
     261    else: 
     262        trialSolution = inst.trialSolution 
     263    trialSolution[:] = inst.population[candidate] 
     264    i = 0 
     265    while 1: 
     266        if random.random() >= inst.probability or i == inst.nDim: 
     267            break 
     268        trialSolution[n] = inst.bestSolution[n] + \ 
     269                           inst.scale * (inst.population[r1][n] + \ 
     270                                         inst.population[r2][n] - \ 
     271                                         inst.population[r3][n] - \ 
     272                                         inst.population[r4][n]) 
    237273        n = (n + 1) % inst.nDim 
    238274        i += 1 
    … …  
    249285    n = random.randrange(inst.nDim) 
    250286 
    251     inst.trialSolution[:] = inst.population[candidate][:] 
    252     i = 0 
    253     while 1: 
    254         if random.random() >= inst.probability or i == inst.nDim: 
    255             break 
    256         inst.trialSolution[n] = inst.population[r1][n] + \ 
    257                                 inst.scale * (inst.population[r2][n] + \ 
    258                                               inst.population[r3][n] - \ 
    259                                               inst.population[r4][n] - \ 
    260                                               inst.population[r5][n]) 
     287    if inst._map_solver: 
     288        trialSolution = inst.trialSolution[candidate] 
     289    else: 
     290        trialSolution = inst.trialSolution 
     291    trialSolution[:] = inst.population[candidate] 
     292    i = 0 
     293    while 1: 
     294        if random.random() >= inst.probability or i == inst.nDim: 
     295            break 
     296        trialSolution[n] = inst.population[r1][n] + \ 
     297                           inst.scale * (inst.population[r2][n] + \ 
     298                                         inst.population[r3][n] - \ 
     299                                         inst.population[r4][n] - \ 
     300                                         inst.population[r5][n]) 
    261301        n = (n + 1) % inst.nDim 
    262302        i += 1 

  • mystic/mystic/termination.py

    r632 r637  
    211211        best = numpy.array(inst.bestSolution) 
    212212        trial = numpy.array(inst.trialSolution) 
    213         update = best - trial #XXX: if inf - inf ? 
    214         answer = numpy.add.reduce(abs(update)) <= tolerance 
     213        update = abs(best - trial) #XXX: if inf - inf ? 
     214        answer = numpy.add.reduce(update.T) 
     215        if isinstance(answer, numpy.ndarray): # if trialPop, take 'best' answer 
     216            answer = max(answer)              #XXX: is this 'best' or 'worst'? 
     217        answer = answer <= tolerance 
    215218        if answer: return info(doc) 
    216219        return info(null) 

  • mystic/tests/solver_test_compare.py

    r635 r637  
    6969  print "comparing against known results" 
    7070  sol = solvers.diffev(rosen, x0, npop=40, disp=0, full_output=True) 
    71   assert almostEqual(sol[1], 0.0020640145337293249, tol=1e-3) 
     71  assert almostEqual(sol[1], 0.0020640145337293249, tol=3e-3) 
    7272  sol = solvers.diffev2(rosen, x0, npop=40, disp=0, full_output=True) 
    73   assert almostEqual(sol[1], 0.0017516784703663288) 
     73  assert almostEqual(sol[1], 0.0017516784703663288, tol=3e-3) 
    7474  sol = solvers.fmin_powell(rosen, x0, disp=0, full_output=True) 
    7575  assert almostEqual(sol[1], 8.3173488898295291e-23) 
Note: See TracChangeset for help on using the changeset viewer.