This page was automatically generated by NetLogo 5.0.
The applet requires Java 5 or higher. Java must be enabled in your browser settings. Mac users must have Mac OS X 10.4 or higher. Windows and Linux users may obtain the latest Java from Sun's Java site.
powered by NetLogo
view/download model file: apo.nlogo
This simulation examines the application of artificial physics to the task of noisy function optimization. We call this “artificial physics optimization” (APO).
Multiple particles use F = ma and a “split Newtonian” force law to self-organize into a triangular lattice. However, the algorithm has been extended to perform function optimization as follows.
There is a function f(x, y) defined for all (x, y) points in the environment. Furthermore, this function is noisy, i.e., if you compute f(x, y) multiple times, you will get multiple function values for the same values of x and y. This occurs in the real world because the function values can only be measured via noisy sensors, or because the function is itself changing over time.
Suppose two neighboring particles in the lattice are labeled ‘a’ and ‘b.’ These two particles are at coordinate positions (x_a, y_a) and (x_b, y_b) respectively. Let f_a = f(x_a, y_a) and f_b = f(x_b, y_b) be the (noisy) fitness values associated with particles ‘a’ and ‘b.’
In addition to the normal split Newtonian force that creates triangular lattices, one more force is applied to neighboring pairs of particles that depends on the function values. Suppose that we want to find the minimum of the function. Then if f_a < f_b, particle ‘b’ is attracted to ‘a,’ while particle ‘a’ is repelled from ‘b.’ Otherwise, particle ‘a’ is attracted to ‘b,’ while particle ‘b’ is repelled from ‘a.’ If we want to find the maximum of the function, this logic is reversed. Note that this is a deliberate breaking of the Newtonian assumption that forces are “equal and opposite.”
In this simulation we assume that our particles are robots that are trying to find the minimum amount of chemical in some region. Chapter 19 extends this work to multiple dimensions and considers the situation when the optimum moves, thereby making the task much more difficult.
Press SETUP to initialize the particles and the patches. The color of each patch is determined by the value of the function at that patch. You will notice a small delay, because it takes a while to compute the value of the function at every patch.
Once SETUP is complete, click MOVE ROBOTS. They will move towards the optimum of the function. In this simulation we are minimizing, and the minimal point is shown as a blue disk in the simulation.
The CLEAR PATCHES button will erase the patches, turning them black.
The NUMBER_OF_ROBOTS slider allows you to control the number of robots created at initialization. Changing this slider while the simulation is running will have no effect.
All other sliders will affect the simulation when it is running.
The most important slider is the NOISE slider. This controls the amount of noise in the function. If you move the NOISE slider all the way to the right, the amount of noise is 5,000,000. This is a very high level of noise.
The DISTANCE TO OPTIMUM monitor indicates the distance of the center of mass of the formation from the minimum at the current moment. The DISTANCE TO OPTIMUM graph displays the distance as a function of time.
This simulation has two other nice features. First, you can “zoom” into the function, allowing you to see fine-level structure embedded in the function. Each time you click ZOOM IN you “zoom” closer by a factor of 10. The RESET button brings you back to the default view. The ZOOM IN button should be clicked after SETUP has been clicked and before you click MOVE ROBOTS.
Also, when the simulation is running, you can mouse click in the graphics pane. This moves the optimum to the location of your mouse. This allows you to move the optimum in real time, so you can see whether the formation of robots can resume its search for the optimum. NOTE: you are not “dragging” the goal with the mouse - you merely click once in the location where you want the goal to move to. This can be a slow process, because the values of the patches need to be recomputed.
Particles are initialized in a random cluster in the upper right portion of the graphics pane, and self-organize into a triangular lattice. The default position of the optimum is at the lower left portion of the graphics pane.
Try running the simulation with the default settings. You will see that the formation quickly locates the optimum. The center of mass of the formation is drawn at each time step so you can see the trajectory that the formation follows. The function looks like a simple bowl, with the bottom of the bowl located at the lower left. Brighter yellow patches are “higher,” whereas darker patches are “lower.”
If you move the NOISE slider while the simulation is running, the simulation will pause as the fitness values of all the patches are recomputed.
Start with the default zoom level of one. Run the simulation several times. Gradually increase the amount of noise. Look at the depiction of the function in the graphics pane (the yellow patches). Is the global structure of the function still discernable? Can the robots find the optimum? How would you describe the trajectory? You might want to move the speed slider to the right when you increase the noise.
Zoom in once, so that the amount of zoom is 10. When there is no noise you will start to see some structure in the function. Run the simulation several times, gradually increasing the amount of noise. Is the optimization task more difficult? Why? An alternative way of viewing the “zoom” factor is that the function has not really changed, but that the spacing between particles has been reduced by a factor of 10. Does that viewpoint help explain the results?
Zoom in once more, so that the amount of zoom is 100. Now you should see definite structure (when there is no noise). Can the robots find the global minimum? If not, what happens and why?
Does increasing the number of robots and/or increasing the DESIRED_SEPARATION help the search process? If it does, explain why.
Regardless of the amount of zoom, you can always use the mouse to move the optimum. When you click in the graphic pane, the optimum moves to that location. Will the formation move towards the new location? This task is referred to as “tracking.”
The noise used in this simulation is drawn from a uniform distribution. Try alternative forms of noise, such as Gaussian.
This simulation uses the “Rastrigan” function to optimize. Try other functions - Chapters 18 and 19 provide a suite of commonly used functions.
The natural extension of this model is to higher dimensions. This is hard to do with NetLogo. See Chapter 19 for details about how to do this in C, where a comparison is made with particle swarm optimization (PSO) on nonstationary environments.
Note, in order to change any NetLogo simulation, you must have the source code (i.e., “apo.nlogo”) downloaded to your computer, as well as NetLogo itself. You can not change the code when you are running the simulation with your browser.
This simulation makes use of a NetLogo feature that allows the toroidal nature of the world to be turned off, creating an environment with strict boundaries.
This simulation also creates the variable “fitness” for every patch. This variable is used to store the function value f(x,y) at every patch.
Also, mouse events are used to move the global optimum of the function, allowing you to “translate” the function in the Cartesian coordinate system.
Chapter 18 of the book provides an alternative approach towards using artificial physics for function optimization. You might be interested in modifying this version of APO to run as described in Chapter 18.
If you mention this model in an academic publication, we ask that you include these citations for the model itself and for the NetLogo software:
- Spears, W. M. and Spears, D. F. (eds.) Physics-based Swarm Intelligence: From Theory to Practice, Springer-Verlag, (2011).
- Wilensky, U. (1999). NetLogo. http://ccl.northwestern.edu/netlogo/. Center for Connected Learning and Computer-Based Modeling, Northwestern University, Evanston, IL.
Copyright 2011 William M. Spears. All rights reserved.
Permission to use, modify or redistribute this model is hereby granted, provided that both of the following requirements are followed:
a) this copyright notice is included.
b) this model will not be redistributed for profit without permission from William M. Spears. Contact William M. Spears for appropriate licenses for redistribution for profit.
; William M. Spears, September 2011
; Artificial Physics Optimization Tutorial Code
; For research and educational use only
breed [apo-robots apo-robot] ; Introduce the "apo-robot" breed
globals [center_of_mass_x center_of_mass_y Noise_Level
G FR D FMAX DeltaT apoF APOMIN OPTX OPTY zoom numbots]
turtles-own [hood r F Fx Fy v vx vy dvx dvy deltax deltay]
patches-own [fitness] ; Store fitness landscape for visualization
;; This code assumes APO is minimizing
to setup
clear-all ; Clear everything
set numbots Number_of_Robots
set Noise_Level Noise ; Used to detect if the Noise slider changes
update-info
set OPTX -125 ; Location of optimum
set OPTY -125
set apoF 0.10 ; See theory in Chapter 6
set zoom 1
; Create and initialize particles
create-apo-robots Number_of_Robots [setup-apo-robots]
; Computes center of mass and displays location
set center_of_mass_x (sum [xcor] of apo-robots) / numbots
set center_of_mass_y (sum [ycor] of apo-robots) / numbots
ask patch (round center_of_mass_x) (round center_of_mass_y) [set pcolor magenta]
; Distance of center of mass to optimum
set APOMIN sqrt((center_of_mass_x - OPTX) * (center_of_mass_x - OPTX) +
(center_of_mass_y - OPTY) * (center_of_mass_y - OPTY))
calculate-patches
reset-ticks
end
to run-and-monitor
if (count turtles < 1) [user-message "Please click HALT and then SETUP first" stop]
update-info
set G (0.9 * FMAX * (D ^ 2) / (2 * (sqrt 3))) ; Gravitational constant set by theory
ask apo-robots [update-apo-robots]
ask turtles [move]
if (mouse-down? and mouse-inside?) [ ; Use mouse click to move optimum
ask patch OPTX OPTY [ask patches in-radius 6 [set pcolor black]]
set OPTX mouse-xcor ; Reset x-coordinate of optimum
set OPTY mouse-ycor ; Reset y-coordinate of optimum
calculate-patches
]
; Draw the location of the optimum
ask patch OPTX OPTY [ask patches in-radius 6 [set pcolor blue]]
; Computes center of mass and displays location
set center_of_mass_x (sum [xcor] of apo-robots) / numbots
set center_of_mass_y (sum [ycor] of apo-robots) / numbots
ask patch (round center_of_mass_x) (round center_of_mass_y) [set pcolor magenta]
; Distance of center of mass to optimum
set APOMIN sqrt((center_of_mass_x - OPTX) * (center_of_mass_x - OPTX) +
(center_of_mass_y - OPTY) * (center_of_mass_y - OPTY))
do-plot ; Update graph
tick
end
to setup-apo-robots ; Set up the robots
setxy ((world-width / 3) + (random-normal 0 1))
((world-height / 3) + (random-normal 0 1)) ; Place at upper right
set vx 0 set vy 0 ; Start with no motion (and no linear momentum)
set shape "circle" set color magenta set size 6
end
to update-apo-robots ; Run artificial physics
set Fx 0 set Fy 0 ; Initialize force components to zero
set vx (1 - FR) * vx ; Slow down according to friction
set vy (1 - FR) * vy
set hood [who] of other apo-robots ; Get the IDs of your neighbors
foreach hood [
set deltax (([xcor] of apo-robot ?) - xcor)
set deltay (([ycor] of apo-robot ?) - ycor)
set r sqrt (deltax * deltax + deltay * deltay)
if (r < 1.5 * D) [ ; The generalized split Newtonian law
set F (G / (r ^ 2))
if (F > FMAX) [set F FMAX] ; Bounds check on force magnitude
ifelse (r > D)
[set Fx (Fx + F * (deltax / r)) ; Attractive force, x-component
set Fy (Fy + F * (deltay / r))] ; Attractive force, y-component
[set Fx (Fx - F * (deltax / r)) ; Repulsive force, x-component
set Fy (Fy - F * (deltay / r))] ; Repulsive force, y-component
; The modification to AP that performs optimization!
ifelse ((rastrigan xcor ycor Noise_Level) < ; Move towards fitness minimum
(rastrigan ([xcor] of apo-robot ?)
([ycor] of apo-robot ?) Noise_Level))
[set Fx (Fx - (apoF * (deltax / r))) ; Repulsive force, x-component
set Fy (Fy - (apoF * (deltay / r)))] ; Repulsive force, y-component
[set Fx (Fx + (apoF * (deltax / r))) ; Attractive force, x-component
set Fy (Fy + (apoF * (deltay / r)))] ; Attractive force, y-component
]
]
set dvx DeltaT * Fx
set dvy DeltaT * Fy
set vx (vx + dvx) ; The x-component of velocity
set vy (vy + dvy) ; The y-component of velocity
set deltax DeltaT * vx
set deltay DeltaT * vy
if ((deltax != 0) or (deltay != 0))
[set heading (atan deltax deltay)]
end
to move
fd (sqrt (deltax * deltax + deltay * deltay)) ; Move the robot
end
to update-info ; Update information from the sliders
set FMAX Force_Maximum
set FR Friction
set DeltaT Time_Step
set D Desired_Separation
if (Noise != Noise_Level) [ ; If the amount of noise changes,
set Noise_Level Noise ; note this fact, and
calculate-patches ; recalculate the fitness values of the patches
]
end
to-report rastrigan [x y n] ; Computes the Rastrigan function (n = noise level)
set x (x / zoom)
set y (y / zoom)
let tempx ((x - OPTX / zoom) ^ 2 - (10 * cos (360 * (x - OPTX / zoom))) + 10)
let tempy ((y - OPTY / zoom) ^ 2 - (10 * cos (360 * (y - OPTY / zoom))) + 10)
report (tempx + tempy + (n * ((random-float 1.0) - 0.5)))
end
to zoom-in ; Increase the zoom by a factor of 10
if (zoom = 0) [user-message "Please click HALT and then SETUP first" stop]
set zoom 10 * zoom
calculate-patches ; Recalculate the fitness values of the patches
end
to reset-zoom ; Reset the zoom back to 1
set zoom 1
calculate-patches ; Recalculate the fitness values of the patches
end
to calculate-patches ; Calculate the fitness values of the patches
ask patches [set fitness (rastrigan pxcor pycor Noise_Level)]
let lower min [fitness] of patches ; Find the minimum value
let upper max [fitness] of patches ; Find the maximum value
; Color the patches shades of yellow, scaled by the min and max
ask patches [set pcolor (scale-color yellow fitness lower upper)]
; Draw the location of the optimum
ask patch OPTX OPTY [ask patches in-radius 6 [set pcolor blue]]
end
to do-plot
set-current-plot "Distance to Optimum" ; Select the "Distance to Optimum" plot
set-current-plot-pen "apopen" ; Select the apo pen
plot APOMIN ; Plot the distance of the center of mass from the optimum
end