I want to solve the Shortest Path problem using PSO in MATLAB.I encoded the path using priority encoding [1], and I'm using constriction and velocity clamping [2].
The problem I'm facing is that the code is extremely slow compared to Dijkstra. I, first, test using Dijkstra to get a benchmark time, and then run PSO to find the lowest cost it can achieve in that time. The result of PSO is always much higher.
If I check how fast each iteration finishes, I find that it takes few seconds for a path with 1000+ nodes on an Intel Core i3-2120 processor.
In the code below, you need to run data.m first to initialise the cost matrix, and then run Dijkstra to get a time benchmark. After that, modify the allowedTime variable in pso.m in seconds.
Parameters:
data.m
pso.m
dimensions range)dimensions range)dijkstra.m
costMatrix, <start_node_id>, <end_node_id>)I'm sorry about the messy code and not using functions, but I needed to make everything inline and see all values after the code is done, or when I break.
data.m
% <<<<<<<<<<<<<<<<<< data definition >>>>>>>>>>>>>>>> %
clear;
clc;
fprintf('Generating data ...\n\n');
dimensions = 5000;
costMatrix = randi(dimensions, dimensions);
fprintf('DONE!\n\n');
pso.m
%% initialization
clc;
fprintf('Initialising swarm ...\n\n');
% parameters
% >>>>>>>>>>>>>>>>>>> edit <<<<<<<<<<<<<<<<<<< %
allowedTime = 15 * 60;
swarm_size = 50;
% SP algorithm specific.
startNode = 1;
endNode = 2;
% velocity equation params.
correction_factor_p = 2.05;
correction_factor_g = 2.05;
contrictionFactor = 0.72984;
% ^^^^^^^^^^^^^^^^^^^ end ^^^^^^^^^^^^^^^^^^^ %
gbest = 1;
oldFitness = 1000000001;
iterations = 0;
% pre-allocate arrays.
swarmPos = zeros(swarm_size, dimensions);
swarmVel = zeros(swarm_size, dimensions);
swarmBestPos = zeros(swarm_size, dimensions);
swarmBestPath = cell(swarm_size, 1);
swarmBestFit = zeros(1, swarm_size);
result = zeros(1, swarm_size);
upperBound = zeros(1, dimensions);
lowerBound = zeros(1, dimensions);
% set bounds.
for i = 1 : dimensions
upperBound(i) = 100;
lowerBound(i) = -100;
end
% ---- initiate swarm -----
for i = 1 : swarm_size
for j = 2 : dimensions
swarmPos(i,j) = lowerBound(j) + rand * (upperBound(j) - lowerBound(j));
swarmVel(i,j) = rand * (upperBound(j) - lowerBound(j)) / 2;
swarmBestPos(i,j) = swarmPos(i,j);
swarmBestPath{i}(j) = -1;
swarmBestFit(i) = 1000000000; % best fitness so far
end
end
% set starting node to avoid on accidental access.
for i = 1 : swarm_size
swarmPos(i,1) = -99999999;
swarmVel(i,1) = -99999999;
swarmBestPos(i,1) = -99999999;
swarmBestPath{i}(1) = startNode;
end
% ^^^^^^^^^^^^^^^^ END: initialisation ^^^^^^^^^^^^^^^^ %
% >>>>>>>>>>>>>>>>>>> START: swarming <<<<<<<<<<<<<<<<<<< %
clc;
fprintf('Swarming ...\n\n');
tic;
%% iterations
while true
% reset results to allow summing.
parfor i = 1 : swarm_size
result(i) = 0;
end
% <<<<<<<<<<<<<<<<< START: movement and fitness >>>>>>>>>>>>>>>>> %
for i = 1 : swarm_size
for j = 2 : dimensions
swarmPos(i,j) = swarmPos(i,j) + swarmVel(i,j); % update x position
if (swarmPos(i,j) > upperBound(j))
swarmPos(i,j) = swarmPos(i,j) - (swarmPos(i,j) - lowerBound(j)) / 2;
elseif (swarmPos(i,j) < lowerBound(j))
swarmPos(i,j) = swarmPos(i,j) + (lowerBound(j) - swarmPos(i,j)) / 2;
end
end
%tic;
% <<< inline fitness function >>> %
tempPath = [];
tempPath(1) = startNode;
invalidBuild = false;
tempPos = swarmPos(i,:);
for j = 2 : dimensions
for k = 2 : dimensions
[discard, maxPos] = max(tempPos);
cost = costMatrix(tempPath(j - 1), maxPos);
tempPos(maxPos) = -9999999 - k;
if (cost < 100000000)
tempPath(j) = maxPos;
result(i) = result(i) + cost;
break;
elseif (k == dimensions)
invalidBuild = true;
end
end
if (invalidBuild)
result(i) = 1000000000;
break;
elseif (tempPath(j) == endNode)
break;
end
end
% ^^^ END: fitness function ^^^ %
% if new position is better
if result(i) < swarmBestFit(i)
for d = 1 : dimensions
swarmBestPos(i,d) = swarmPos(i,d); % update best x,
end
swarmBestPath{i} = tempPath;
swarmBestFit(i) = result(i); % and best value
end
end
% ^^^^^^^^^ END: movement and fitness ^^^^^^^^^ %
% <<<<<<<<<<<<<<<<< update global best >>>>>>>>>>>>>>>>> %
for i = 1 : swarm_size
if swarmBestFit(i) < swarmBestFit(gbest)
gbest = i; % global best i.
took = toc; % time taken to reach this best.
end
end
% <<<<<<<<<<<<<<<<< update velocity >>>>>>>>>>>>>>>>> %
for i = 1 : swarm_size
for j = 2 : dimensions
swarmVel(i,j) = contrictionFactor * (swarmVel(i,j) ...
+ correction_factor_p * rand * (swarmBestPos(i,j) - swarmPos(i,j)) ...
+ correction_factor_g * rand * (swarmBestPos(gbest,j) - swarmPos(i,j)));
if (swarmVel(i,j) > (upperBound(j) - lowerBound(j)) / 2)
swarmVel(i,j) = (upperBound(j) - lowerBound(j)) / 2;
end
end
end
% <<<<<<<<<<<<<<<<< print global bests if changed >>>>>>>>>>>>>>>>> %
if ( oldFitness ~= swarmBestFit(gbest) )
oldFitness = swarmBestFit(gbest);
% update display
clc
fprintf('Best particle position:\n');
sizeTemp = size(swarmBestPath{gbest}, 2);
for i = 1 : sizeTemp
if (swarmBestPath{gbest}(i) ~= 0)
fprintf('%d\n', swarmBestPath{gbest}(i));
end
end
fprintf('\nBest fitness: %d\n\n', swarmBestFit(gbest));
end
iterations = iterations + 1;
% end on timeout
elapsedTime = toc;
if (elapsedTime > allowedTime)
break;
end
end
clc;
fprintf('>>>>>>>>>>>>>>> FINISHED <<<<<<<<<<<<<<<\n\n\n');
fprintf('Best path:\n');
sizeTemp = size(swarmBestPath{gbest}, 1);
for i = 1 : sizeTemp
if (swarmBestPath{gbest}(i) ~= 0)
fprintf('%d\n', swarmBestPath{gbest}(i));
end
end
fprintf('\nBest cost: %d\n\n', swarmBestFit(gbest));
fprintf('\nTook: %d iterations, and %.2f seconds.\n\n', iterations, took);
dijkstra.m
function dijkstra(matriz_costo, s, d)
% This is an implementation of the dijkstra´s algorithm, wich finds the
% minimal cost path between two nodes. It´s supoussed to solve the problem on
% possitive weighted instances.
% the inputs of the algorithm are:
%farthestNode: the farthest node to reach for each node after performing
% the routing;
% n: the number of nodes in the network;
% s: source node index;
% d: destination node index;
%For information about this algorithm visit:
%http://en.wikipedia.org/wiki/Dijkstra%27s_algorithm
%This implementatios is inspired by the Xiaodong Wang's implememtation of
%the dijkstra's algorithm, available at
%http://www.mathworks.com/matlabcentral/fileexchange
%file ID 5550
%Author: Jorge Ignacio Barrera Alviar. April/2007
n=size(matriz_costo,1);
S(1:n) = 0; %s, vector, set of visited vectors
dist(1:n) = inf; % it stores the shortest distance between the source node and any other node;
prev(1:n) = n+1; % Previous node, informs about the best previous node known to reach each network node
dist(s) = 0;
iterations = 0;
tic;
while sum(S)~=n
candidate=[];
for i=1:n
if S(i)==0
candidate=[candidate dist(i)];
else
candidate=[candidate inf];
end
end
[u_index u]=min(candidate);
S(u)=1;
for i=1:n
if(dist(u)+matriz_costo(u,i))<dist(i)
dist(i)=dist(u)+matriz_costo(u,i);
prev(i)=u;
end
end
iterations = iterations + 1;
end
sp = [d];
while sp(1) ~= s
if prev(sp(1))<=n
sp=[prev(sp(1)) sp];
else
error;
end
end;
spcost = dist(d);
took = toc;
fprintf('Best path:\n');
fprintf('%d\n', sp);
fprintf('\nBest cost: %d\n\n', spcost);
fprintf('\nTook: %d iterations, and %.2f seconds.\n\n', iterations, took);
(1) A Nondominated Sorting Genetic Algorithm for SP Routing Problem
(2) Constriction factors and Parameters
I'm a newbie to the world of Artificial Intelligence. I have tried my level best to help you.
PSO is a metaheuristic approach. Problems with incomplete or imperfect information or limited computation capacity can be solved using metaheuristics. Such problems cannot be solved using Dijkstra's as it needs full topology details. Hence use of the algorithm depends on the problem domain.
Since PSO is a stochastic process, initial result will not be the optimal. As the iterations are performed, the cost function value reduces. Where as Dijkstra's finds the shortest path by one iteration.
So time consumption for PSO will be more compared to Dijkstra. Use of these algorithms are problem specific.
Hope this will help you!