Sunday, April 4, 2010

Particle Swarm Optimization (PSO) Sample Code using Java

Yes, I'm still coding and I'm proud of it :)

This post assume that the reader has already known about Particle Swarm Optimization (PSO) method, and hence I wouldn't spare a space to explain about it. But if you'd like to ask the method's component on the code I provide below then I'll be gladly explain it in greater length. [This post is written by Gandhi Manalu - gandhim.wordpress.com]

In this post, I'll describe and provide a sample code of PSO method for solving a very simple function optimization. Let say that the function to be minimized is as follow:

$latex f(x,y) = (2.8125-x+xy^4)^2+(2.25-x+xy^2)^2+(1.5-x+xy)^2$

With the following constraints: $latex 1<=x<=4; -1<=y<=1$

In order to solve this problem using PSO, we'll need these classes:

  • Position: to represent Position part of the particle
  • Velocity: to represent Velocity part of the particle
  • Particle: the particle itself [This post is written by Gandhi Manalu - gandhim.wordpress.com]
  • SimplePSO: the main control of the program
  • PSOConstants: an interface to define parameters used in the PSO
Since we're going to solve two-variable function optimization, we'll need to provide two-dimensional position and velocity. For the Position we have:

package org.gandhim.pso;

/**
*
* @author gandhim
*/
public class Position {
private double x;
private double y;

public Position(double x, double y) {
this.x = x;
this.y = y;
}

public double getX() {
return x;
}

public void setX(double x) {
this.x = x;
}

public double getY() {
return y;
}

public void setY(double y) {
this.y = y;
}
}

For the Velocity we have:

package org.gandhim.pso;

/**
*
* @author gandhim
*/
public class Velocity {
private double x;
private double y;

public Velocity(double x, double y) {
this.x = x;
this.y = y;
}

public double getX() {
return x;
}

public void setX(double x) {
this.x = x;
}

public double getY() {
return y;
}

public void setY(double y) {
this.y = y;
}
}
And the Particle is as follow:[This post is written by Gandhi Manalu - gandhim.wordpress.com]

package org.gandhim.pso;

/**
*
* @author gandhim
*/
public class Particle {
private Position location;
private Velocity velocity;
private double fitness;

public double getFitness() {
calculateFitness();
return fitness;
}

public void calculateFitness() {
double x = this.location.getX();
double y = this.location.getY();

fitness = Math.pow((2.8125 - x + x * Math.pow(y, 4)), 2) +
Math.pow((2.25 - x + x * Math.pow(y, 2)), 2) +
Math.pow((1.5 - x + x * y), 2);
}

public Position getLocation() {
return location;
}

public void setLocation(Position location) {
this.location = location;
}

public Velocity getVelocity() {
return velocity;
}

public void setVelocity(Velocity velocity) {
this.velocity = velocity;
}
}

Pay attention to the calculateFitness() method, it is where we put the function evaluation. Now, we're ready for the main process of the PSO. In this class, we'll need several methods:
  • initializeSwarm() - to initialize the swarm used in the method
  • execute() - the main part of the process[This post is written by Gandhi Manalu - gandhim.wordpress.com]
Here's the code (partial code):

private void initializeSwarm() {
Particle p;
Random generator = new Random();

for (int i = 0; i < SWARM_SIZE; i++) {
p = new Particle();
double posX = generator.nextDouble() * 3.0 + 1.0;
double posY = generator.nextDouble() * 2.0 - 1.0;
p.setLocation(new Position(posX, posY));

double velX = generator.nextDouble() * 2.0 - 1.0;
double velY = generator.nextDouble() * 2.0 - 1.0;
p.setVelocity(new Velocity(velX, velY));

swarm.add(p);
}
}

public void execute() {
Random generator = new Random();
initializeSwarm();

evolutionaryStateEstimation();

int t = 0;
double w;

while (t < MAX_ITERATION) {
// calculate corresponding f(i,t) corresponding to location x(i,t)
calculateAllFitness();

// update pBest
if (t == 0) {
for (int i = 0; i < SWARM_SIZE; i++) {
pBest[i] = fitnessList[i];
pBestLoc.add(swarm.get(i).getLocation());
}
} else {
for (int i = 0; i < SWARM_SIZE; i++) {
if (fitnessList[i] < pBest[i]) {
pBest[i] = fitnessList[i];
pBestLoc.set(i, swarm.get(i).getLocation());
}
}
}

int bestIndex = getBestParticle();
// update gBest
if (t == 0 || fitnessList[bestIndex] < gBest) {
gBest = fitnessList[bestIndex];
gBestLoc = swarm.get(bestIndex).getLocation();
}

w = W_UP - (((double) t) / MAX_ITERATION) * (W_UP - W_LO);

for (int i = 0; i < SWARM_SIZE; i++) {
// update particle Velocity
double r1 = generator.nextDouble();
double r2 = generator.nextDouble();
double lx = swarm.get(i).getLocation().getX();
double ly = swarm.get(i).getLocation().getY();
double vx = swarm.get(i).getVelocity().getX();
double vy = swarm.get(i).getVelocity().getY();
double pBestX = pBestLoc.get(i).getX();
double pBestY = pBestLoc.get(i).getY();
double gBestX = gBestLoc.getX();
double gBestY = gBestLoc.getY();

double newVelX = (w * vx) + (r1 * C1) * (pBestX - lx) + (r2 * C2) * (gBestX - lx);
double newVelY = (w * vy) + (r1 * C1) * (pBestY - ly) + (r2 * C2) * (gBestY - ly);
swarm.get(i).setVelocity(new Velocity(newVelX, newVelY));

// update particle Location
double newPosX = lx + newVelX;
double newPosY = ly + newVelY;
swarm.get(i).setLocation(new Position(newPosX, newPosY));
}

t++;
}
}

And the last is the interface for storing the constants:

package org.gandhim.pso;

/**
*
* @author gandhim
*/
public interface PSOConstants {
int SWARM_SIZE = 30;
int DIMENSION = 2;
int MAX_ITERATION = 300;
double C1 = 2.0;
double C2 = 2.0;
double W_UP = 1.0;
double W_LO = 0.0;
}

I made the interface just for the sake of flexibility. A sample of running the program is as follow: PSO Result

As we can see from the result, the program found the solution of the problem for (x=3.0 and y=0.5). You might have noticed that the program is not optimized yet, for example it could have been stopped when it already found the solution.[This post is written by Gandhi Manalu - gandhim.wordpress.com]

Actually, there are other PSO's components that have not been implemented in this program, for example constraints handling. But since the problem we solved here is a very simple one it doesn't really need the constraint handling.

I hope that this post is useful for you.

1 comment: