# Python module for simulated annealing - anneal.py - v1.0 - 2 Sep 2009

# Copyright (c) 2009, Richard J. Wagner <wagnerr@umich.edu>

# Permission to use, copy, modify, and/or distribute this software for any

# purpose with or without fee is hereby granted, provided that the above

# copyright notice and this permission notice appear in all copies.

# THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES

# WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF

# MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR

salesman problem to find the shortest route to visit the twenty largest

cities in the United States.

- Changed to slicing lists instead of deepcopy-ing them.

- e.g. state = prevState[:] instead of state = deepcopy(prevState)

- Huge performance enhancement (~5-10x faster)

- Should be identical behavior if the items in the state list are immutable.

- (immutable objects include integers and strings so should be safe)

# How to optimize a system with simulated annealing:

# 1) Define a format for describing the state of the system.

# 2) Define a function to calculate the energy of a state.

# 3) Define a function to make a random change to a state.

# 4) Choose a maximum temperature, minimum temperature, and number of steps.

# 5) Set the annealer to work with your state and functions.

# 6) Study the variation in energy with temperature and duration to find a

# productive annealing schedule.

# 4) Run the automatic annealer which will attempt to choose reasonable values

# for maximum and minimum temperatures and then anneal for the allotted time.

def __init__(self, energy, move):

self.energy = energy # function to calculate energy of a state

self.move = move # function to make a random change to a state

def anneal(self, state, Tmax, Tmin, steps, updates=0):

"""Minimizes the energy of a system by simulated annealing.

state -- an initial arrangement of the system

Tmax -- maximum temperature (in units of energy)

Tmin -- minimum temperature (must be greater than zero)

steps -- the number of steps requested

updates -- the number of updates to print during annealing

Returns the best state and energy found."""

def update(T, E, acceptance, improvement):

"""Prints the current temperature, energy, acceptance rate,

improvement rate, elapsed time, and remaining time.

The acceptance rate indicates the percentage of moves since the last

update that were accepted by the Metropolis algorithm. It includes

moves that decreased the energy, moves that left the energy

unchanged, and moves that increased the energy yet were reached by

The improvement rate indicates the percentage of moves since the

last update that strictly decreased the energy. At high

temperatures it will include both moves that improved the overall

it will tend toward zero as the moves that can decrease the energy

are exhausted and moves that would increase the energy are no longer

elapsed = time.time() - start

print ' Temperature Energy Accept Improve Elapsed Remaining'

print '%12.2f %12.2f %7.2f%% %7.2f%% %s %s' % \

(T, E, 100.0*acceptance, 100.0*improvement,

time_string(elapsed), time_string(remain))

# Precompute factor for exponential cooling from Tmax to Tmin

print 'Exponential cooling requires a minimum temperature greater than zero.'

Tfactor = -math.log( float(Tmax) / Tmin )

updateWavelength = float(steps) / updates

# Attempt moves to new states

if step // updateWavelength > (step-1) // updateWavelength:

update(T, E, float(accepts)/trials, float(improves)/trials)

trials, accepts, improves = 0, 0, 0

# Return best state and energy

return bestState, bestEnergy

def auto(self, state, seconds, steps=2000):

"""Minimizes the energy of a system by simulated annealing with

automatic selection of the temperature schedule.

state -- an initial arrangement of the system

seconds -- time to spend annealing (after exploring temperatures)

steps -- number of steps to spend on each stage of exploration

Returns the best state and energy found."""

def run(state, T, steps):

"""Anneals a system at constant temperature and returns the state,

energy, rate of acceptance, and rate of improvement."""

prevState = copy.deepcopy(state)

return state, E, float(accepts)/steps, float(improves)/steps

print 'Attempting automatic simulated anneal...'

# Find an initial guess for temperature

# This can get in an infinite loop if the energy function always returns

T = abs( self.energy(state) - E )

raise Exception("Couldn't find initial temperature")

print 'Exploring temperature landscape:'

print ' Temperature Energy Accept Improve Elapsed'

def update(T, E, acceptance, improvement):

elapsed = time.time() - start

print '%12.2f %12.2f %7.2f%% %7.2f%% %s' % \

(T, E, 100.0*acceptance, 100.0*improvement, time_string(elapsed))

# Search for Tmax - a temperature that gives 98% acceptance

state, E, acceptance, improvement = run(state, T, steps)

update(T, E, acceptance, improvement)

# Search for Tmin - a temperature that gives 0% improvement

T = round_figures(T/1.5, 2)

update(T, E, acceptance, improvement)

# Calculate anneal duration

elapsed = time.time() - start

duration = round_figures(int(seconds * step / elapsed), 2)

# MP: Don't perform anneal, just return params

#return self.anneal(state, Tmax, Tmin, duration, 20)

return {'tmax': Tmax, 'tmin': Tmin, 'steps': duration}

if __name__ == '__main__':

"""Test annealer with a traveling salesman problem."""

# List latitude and longitude (degrees) for the twenty largest U.S. cities

cities = { 'New York City': (40.72,74.00), 'Los Angeles': (34.05,118.25),

'Chicago': (41.88,87.63), 'Houston': (29.77,95.38),

'Austin': (30.27,97.77), 'Columbus': (39.98,82.98),

'Fort Worth': (32.75,97.33), 'Charlotte': (35.23,80.85),

'Memphis': (35.12,89.97), 'Baltimore': (39.28,76.62) }

"""Calculates distance between two latitude-longitude coordinates."""

R = 3963 # radius of Earth (miles)

lat2, lon2 = math.radians(b[0]), math.radians(b[1])

return math.acos( math.sin(lat1)*math.sin(lat2) +

math.cos(lat1)*math.cos(lat2)*math.cos(lon1-lon2) ) * R

"""Swaps two cities in the route."""

a = random.randint( 0, len(state)-1 )

b = random.randint( 0, len(state)-1 )

state[a], state[b] = state[b], state[a]

"""Calculates the length of the route."""

for i in range(len(state)):

e += distance( cities[state[i-1]], cities[state[i]] )

# Start with the cities listed in random order

# Minimize the distance to be traveled by simulated annealing with a

# manually chosen temperature schedule

annealer = Annealer(route_energy, route_move)

print "%i mile route:" % route_energy(state)

# Minimize the distance to be traveled by simulated annealing with an

# automatically chosen temperature schedule

state, e = annealer.auto(state, 4)

print "%i mile route:" % route_energy(state)