Overview

Particle swarm optimization (PSO) is a type of swarm intelligence inspired by the behavior of swarms of animals. In this article, I will introduce a simple example of particle swarm optimization.

Parabolic formula

The paraboloidal equation is given in the following form.

\begin{aligned}
z = x^2+y^2
\end{aligned}

Obviously, the minimum value is $ z = 0 $ when $ (x, y) = (0, 0) $. This is calculated using the particle swarm optimization method.

Explanation of particle swarm optimization method

See the article below. (Is the subscript in §2 partly wrong?) Particle Swarm Optimization and Nonlinear System

Let's implement §2 of this article.

Source code

`main.py`


# -*- coding: utf-8 -*-

import numpy as np
import random

#Evaluation function: z = x^2 + y^2
def criterion(x, y):
    z = x * x + y * y
    return z

#Function to update the position of particles
def update_position(x, y, vx, vy):
    new_x = x + vx
    new_y = y + vy
    return new_x, new_y

#A function that updates the velocity of particles
def update_velocity(x, y, vx, vy, p, g, w=0.5, ro_max=0.14):
    #Parameter ro is given randomly
    ro1 = random.uniform(0, ro_max)
    ro2 = random.uniform(0, ro_max)
    #Update the particle velocity
    new_vx = w * vx + ro1 * (p["x"] - x) + ro2 * (g["x"] - x)
    new_vy = w * vy + ro1 * (p["y"] - y) + ro2 * (g["y"] - y)
    return new_vx, new_vy


def main():
    N = 100  #Number of particles
    x_min, x_max = -5, 5
    y_min, y_max = -5, 5
    #Particle position,speed,Personal vest,Initialize the global best
    ps = [{"x": random.uniform(x_min, x_max), 
        "y": random.uniform(y_min, y_max)} for i in range(N)]
    vs = [{"x": 0.0, "y": 0.0} for i in range(N)]
    personal_best_positions = list(ps)
    personal_best_scores = [criterion(p["x"], p["y"]) for p in ps]
    best_particle = np.argmin(personal_best_scores)
    global_best_position = personal_best_positions[best_particle]
    
    T = 30  #Time limit(Number of loops)
    for t in range(T):
        for n in range(N):
            x, y = ps[n]["x"], ps[n]["y"]
            vx, vy = vs[n]["x"], vs[n]["y"]
            p = personal_best_positions[n]
            #Update the position of the particles
            new_x, new_y = update_position(x, y, vx, vy)
            ps[n] = {"x": new_x, "y": new_y}
            #Update the velocity of particles
            new_vx, new_vy = update_velocity(
                new_x, new_y, vx, vy, p, global_best_position)
            vs[n] = {"x": new_vx, "y": new_vy}
            #Find the evaluation value,Update your personal vest
            score = criterion(new_x, new_y)
            if score < personal_best_scores[n]:
                personal_best_scores[n] = score
                personal_best_positions[n] = {"x": new_x, "y": new_y}
        #Update the global best
        best_particle = np.argmin(personal_best_scores)
        global_best_position = personal_best_positions[best_particle]
    #Optimal solution
    print(global_best_position)
    print(min(personal_best_scores))

if __name__ == '__main__':
    main()

result

`result`


{'y': 0.00390598718159734, 'x': -0.0018420875049243782}
1.86500222386e-05

Visualization

It looks like $ N (= 100) $ particles are concentrated on $ (x, y) = (0, 0) $.

Other

How do you determine some parameters such as the number of particles? (Trial and Error?)

The order of position update, speed update, personal best, global best is out of order depending on the literature ... In this article, we will update the position ⇒ update the speed ⇒ update the personal best ⇒ update the global best.