Parallelization of Floyd-Warshall algorithm in matlab

Question

I know this question has been answered several times before but I checked all the previous answers to remedy my situation but it didn't help.

What i need to do is to parallelize my loop so that each city (or the inner loop) is processed parallel. But while using parfor i get the error "The variable A in a parfor cannot be classified". The 2d matrix has a fixed size of n X n. i don't know see the problem. kindly help me out...

the c implementation that i was provided was done using mpi.h . using mpicc. what i need to achieve is that there should be n processes, each responsible to find the shortest paths from its local city to all other cities.

Every case if different. In my case:

my_first_city=2;
my_last_city=n;
parpool(n-1);

parfor (int_city=2:n,n-1)
% broadcast all --  create threads to handle each city
  for local_city=my_first_city:n
    for city2=2:n
          A(local_city,city2)=min(A(local_city,city2),A(local_city,int_city)+A(int_city,city2));
    end
  end
end

Here's my function to compute shortest paths:

function [ A,init ] = floydWarshall(input_matrix, n )
%% Floyd_Warshall algorithm is an analytical algorithm for finding shortest paths in weighted graph ,
%  for example an adjacency matrix or a map graph.
% Floyd_Warshall algorithm compares all possible paths through a graph between each pair of vertices,
% The complexity of this algorithm is O(n^3) where n is the number of vertices, or nodes.
%% Floyd_Warshall
% inputs : 
%       n          = number of vertices to initialize an adjacency matrix.
%    input_matrix  = the input matrix of initial weights or path costs. a nXn matrix
% outputs: 
%       A  = the matrix after floydWarshall algo is applied and the matrix contains the shortest
%            paths from each node to each other node
%     init = The original matrix with all the costs from each node to each other node.
if(nargin<2)
  n=size(input_matrix);
elseif (nargin<1)
   n=size(input_matrix);
   A=magic(n);
end


for i=1:n    % marking the border rows and columns with cities
    A(i,1)=i-1;
    A(1,i)=i-1;
end

for i=1:n    % making sure that the distance from a city i to city i is 0
    A(i,i)=0;
end

for i=2:n   
    for j=2:n
        A(i,j)=input_matrix(i,j);  % input matrix, values 
    end
end


init=A;   % temp variable to store the old matrix
for int_city=2:n
    for city1=2:n
        for city2=2:n
            A(city1,city2)=min(A(city1,city2),A(city1,int_city)+A(int_city,city2));
        end  % floyd-warshall
    end
end

Answer 1

Here is a code in MATLAB which implements Floyd-Warshall algorithm

%%%%% Step 0: Initialization and Parameters %%%%%

N = size(D,1); 
INF =  1000*max(max(D))*N;  %% effectively infinite distance
Y.coords = cell(length(dims),1); 
R = zeros(1,length(dims)); 

%%%%% Step 1: Construct neighborhood graph %%%%%
disp('Constructing neighborhood graph...'); 

K = n_size;
if n_fcn == 'k'
     [~, ind] = sort(D); % For matrices, sort(X) sorts each column of X in ascending order.
     for i=1:N % N is number of points
          D(i,ind((2+K):end,i)) = INF; %setting to infinity distances that more far number of neighborhoods
     end                               %two because we need to have at least one neighbor???
elseif n_fcn == 'epsilon'
     warning off    %% Next line causes an unnecessary warning, so turn it off
     D =  D./(D<=epsilon); 
     D = min(D,INF); 
     warning on
end

D = min(D,D');    %% Make sure distance matrix is symmetric


% Finite entries in D now correspond to distances between neighboring points. 
% Infinite entries (really, equal to INF) in D now correspond to 
%   non-neighoring points. 

%%%%% Step 2: Compute shortest paths %%%%%
disp('Computing shortest paths...'); 

tic; 
for k=1:N
     D = min(D,  repmat(D(:,k),[1 N])  +  repmat(D(k,:),[N 1])  ); % compares each matrix element
     if rem(k,100) == 0 % rem    Remainder after division.
          disp([' Iteration: ', num2str(k), '     Estimated time to completion: ', num2str((N-k)*toc/k/60), ' minutes']); 
     end
end

%%%%% Remove outliers from graph %%%%%
disp('Checking for outliers...'); 
n_connect = sum(~(D==INF));        %% number of points each point connects to 1xN
[~, firsts] = min(D==INF);       %% first point each point connects to
[comps, ~, ~] = unique(firsts);    %% represent each connected component once
size_comps = n_connect(comps);     %% size of each connected component
[~, comp_order] = sort(size_comps);  %% sort connected components by size
comps = comps(comp_order(end:-1:1)); %%move upside down   
size_comps = size_comps(comp_order(end:-1:1)); 
n_comps = length(comps);               %% number of connected components
comp=1;                              %% default: use largest component

disp(['  Number of connected components in graph: ', num2str(n_comps)]); 
disp(['  Embedding component ', num2str(comp), ' with ', num2str(size_comps(comp)), ' points.']); 
Y.index = find(firsts==comps(comp)); 

D = D(Y.index, Y.index); 
N = length(Y.index); 

I know you'll probably get some questions You can ask. By the way the fastest algorithm for finding distances is Dijkstra's Algorithm(but you should implement it in C).

Hope it will help you! I am newbie as well, but I had similar task before and code 100% correct.

Answer 2

I believe that your problem is the fact that you are not "slicing" A matrix.

parfor construct in Matlab creates processes. That means that all the processes will compete to update variable A, simultaneously. I don't know whether Matlab implements a shared memory between the processes and proper synchronization. It looks like it doesn't.

If Matlab was creating threads, then it would have been easier to synchronize accesses to A because all threads whould have had access to it, being in a single process. With processes is more complicated.

So your problem is that A is shared between processes. To avoid this problem you can split the A matrix into n variables (equal to the number of processes), give each process a slice and then reconstruct A with the output given by the n processes.

Best would be to do the slicing manually, by assigning n sub matrices (or an array of n sub-matrices). Actually the compiler gives you the error because it cannot slice it on its own. It needs you to do it.

See the post here, describing a similar problem.

Good luck.