special case of using parallel computing toolbox in order to solve decomposed problems

Hello,

How can I have a for loop with workers of parallel computing doing "1-some functions with different input for each worker, 2-gather the result of iteration from workers 3-doing computation with the result of all workers, 4-update the global input to each worker until stopping criteria and so on?

pseudonym code can be like this:

a=cell(number of workers,1); %first input in cell

b=cell(number of workers,1); %first input in cell

for i=1:maxiter % a loop in the client

tell the workers to find "X" based on some identical function and computations by different inputs (the cells of a and b)

do a global computation to find z with the resut of X gathered form each worker and update the input parameters (a and b)

if fcn(a,b,X,z)<eps

stop

end

end

Answer 1

[Walter Roberson]

perhaps use spmd. As each iteration calculates results, labSend to the other workers. You might need to barrier() to synchronize. Depending on how the results are combined it might make sense to select one lab (such as lab 1) to send all the results to and have it combine the results and send out the updated parameters; if you use labSendRecieve on all the labs except the coordinator then they willeeffectively block until they get a response without any barrier(). The tricky part might be ensuring that the coordinator has results from all of the other workers.

[Bill Masters]

Hi Walter

The code using SPMD is right below but in the first step when spmd starts to work I have the error below:

"Brace indexing is not supported for variables of this type."

I did a simple test with this structure and there was no error, but in the main problem, confronted this error.

You can check my code as below:

%% problem definition (a simple problem to check if the code works)
delete(gcp('nocreate'));
clear;
close all;
clc;

N=4;
c=cell(1,N);
for i=1:4
    c{i}=rand(2,1)+0.5;
end
cc=[c{1}; c{2}; c{3}; c{4}];

x0=cell(1,N);
for i=1:4
    x0{i}=randn(2,1);
end
x00=[x0{1}; x0{2}; x0{3}; x0{4}];

N=4; %number of workers
A=cell(1,N);
A{1}=[2 3;1 5];
A{2}=[4 3;2 6];
A{3}=[5 2;3 7];
A{4}=[6 1;7 2];
AA=[A{1} A{2} A{3} A{4}];
b=cell(1,N);
x0=cell(1,N);
bb=AA*x00;
for i=1:N
    b{i}=bb/N;
end
% Defining a sample problem finishes here
% ----------------------------------------------------------------
%% admm call (the code written to solve generated problem)
rho=1.5;
alpha=1.5;
%% PADMM
QUIET    = 0;
MAX_ITER = 1000;
ABSTOL   = 1e-4;
RELTOL   = 1e-2;
%[m n] = size(A);
x=cell(4,1);
for i=1:4
    x{i}=zeros(2,1);
end
z=zeros(2,1);
u=cell(N,1);
for i=1:4
    u{i}=zeros(2,1);
end
if ~QUIET
    fprintf('%3s\t%10s\t%10s\t%10s\t%10s\t%10s\n', 'iter', ...
      'r norm', 'eps pri', 's norm', 'eps dual', 'objective');
end
parpool(4);
spmd   %affromentioned error happens in this line
    cw=c{labindex}; %define input parameters of workers
    Aw=A{labindex}; %define input parameters of workers
    bw=b{labindex}; %define input parameters of workers
    xw=x{labindex}; %define input parameters of workers
    uw=u{labindex}; %define input parameters of workers
    for k = 1:MAX_ITER
    % x-update
    z_old=z; % the initial z vector that would be used in ahead lines of code
   
    tmp = [ rho*eye(2), transpose(Aw{labindex}); Aw{labindex}, zeros(2) ] \ [ rho*(z - uw{labindex}) - cw{labindex}; bw{labindex} ];
    xw{labindex} = tmp{labindex}(1:n);
    labBarrier;
    x_bar=gop(@plus,xw{labindex})/4; %x_bar the average of the xw that is calculated for each worker
    u_bar=gop(@plus,u{labindex})/4; %u_bar the average of the uw that is calculated for each worker
    x_hat{labindex} = alpha*x_bar + (1 - alpha)*z_old; %x_hat variable calculated for each worker
    z = subplus(x_bar + u_bar);    % variable z must be calculated in each iteration by the results of workers
    uw{labindex} = uw{labindex} + (x_hat{labindex} - z);  % the new amount of uw, calculated in this iteration usable in the rest of the iterations
    labBarrier;
    history.objval(k)  = gop(@plus,objective(cw{labindex}, xw{labindex}));

    history.r_norm(k)  = norm(x_bar - z);
    history.s_norm(k)  = norm(-rho*(z - zold));

    history.eps_pri(k) = sqrt(2)*ABSTOL + RELTOL*max(norm(x_bar), norm(-z));
    history.eps_dual(k)= sqrt(2)*ABSTOL + RELTOL*norm(rho*u_bar);

    if ~QUIET
        fprintf('%3d\t%10.4f\t%10.4f\t%10.4f\t%10.4f\t%10.2f\n', k, ...
            history.r_norm(k), history.eps_pri(k), ...
            history.s_norm(k), history.eps_dual(k), history.objval(k));
    end

    if (history.r_norm(k) < history.eps_pri(k) && ...
       history.s_norm(k) < history.eps_dual(k))
         break;
    end
    end
        
end

[Bill Masters]

I've found that the for loop inside the SPMD block causes the error and the question is how to implement for loops in SPMD the result of iterations being saved? searched for any answer to this question nothing found.

[Walter Roberson]

for i=1:4
    c{i}=rand(2,1)+0.5;
end

So c is a cell containing a 2 x 1 numeric vector, and there is one cell for each lab.

spmd   %affromentioned error happens in this line
    cw=c{labindex}; %define input parameters of workers

so cw is the current lab's 2 x 1 numeric vector.

        tmp = [ rho*eye(2), transpose(Aw{labindex}); Aw{labindex}, zeros(2) ] \ [ rho*(z - uw{labindex}) - cw{labindex}; bw{labindex} ];

Notice the cw{labindex} . But cw is a numeric vector already extracted from c, so you are trying to use {} indexing on a 2 x 1 numeric vector.

uw=u{labindex}; %define input parameters of workers

uw is this lab's u. You use uw{labindex} in that same assignment to tmp and then in the next line

    xw{labindex} = tmp{labindex}(1:n);

so even though xw was already extracted, you are assigning to it based on lab index.

You do not assign to any of the variables A, b, c, x, u that you extract from so it is not obvious why you are using {labindex} at any point after you copy out of A, b, c, x, u.

Exactly which variables are you expecting to be outputs? The only thing I can see is that it looks like maybe you are expecting your output to be history and that you intend to access it as a "composite" after the spmd.

Reminder: it is more efficient to write into local variables in for loops, and to assign the completed results into the composite after the loop.

[Bill Masters]

My purpose of do this is to assign the workers their own initial vector and then start the for loop seperately. in other words, I want the for loop to run in each worker with their own inputs. If I don't write the code like this, how can each worker call their own inputs?

[Walter Roberson]

%% problem definition (a simple problem to check if the code works)
delete(gcp('nocreate'));
clear;
close all;
clc;

N=4;
c=cell(1,N);
for i=1:4
    c{i}=rand(2,1)+0.5;
end
cc=[c{1}; c{2}; c{3}; c{4}];

x0=cell(1,N);
for i=1:4
    x0{i}=randn(2,1);
end
x00=[x0{1}; x0{2}; x0{3}; x0{4}];

N=4; %number of workers
A=cell(1,N);
A{1}=[2 3;1 5];
A{2}=[4 3;2 6];
A{3}=[5 2;3 7];
A{4}=[6 1;7 2];
AA=[A{1} A{2} A{3} A{4}];
b=cell(1,N);
x0=cell(1,N);
bb=AA*x00;
for i=1:N
    b{i}=bb/N;
end
% Defining a sample problem finishes here
% ----------------------------------------------------------------
%% admm call (the code written to solve generated problem)
rho=1.5;
alpha=1.5;
%% PADMM
QUIET    = 0;
MAX_ITER = 1000;
ABSTOL   = 1e-4;
RELTOL   = 1e-2;
%[m n] = size(A);
x=cell(4,1);
for i=1:4
    x{i}=zeros(2,1);
end
z=zeros(2,1);
u=cell(N,1);
for i=1:4
    u{i}=zeros(2,1);
end
if ~QUIET
    fprintf('%3s\t%10s\t%10s\t%10s\t%10s\t%10s\n', 'iter', ...
      'r norm', 'eps pri', 's norm', 'eps dual', 'objective');
end
parpool(4);
spmd
    cw = c{labindex}; %define input parameters of workers
    Aw = A{labindex}; %define input parameters of workers
    bw = b{labindex}; %define input parameters of workers
    xw = x{labindex}; %define input parameters of workers
    uw = u{labindex}; %define input parameters of workers
    for k = 1:MAX_ITER
    % x-update
    z_old = z; % the initial z vector that would be used in ahead lines of code
   
    tmp = [ rho*eye(2), transpose(Aw); Aw, zeros(2) ] \ [ rho*(z - uw) - cw; bw ];
    xw = tmp(1:n);
    labBarrier;
    x_bar = gop(@plus, xw)/4; %x_bar the average of the xw that is calculated for each worker
    u_bar = gop(@plus, u)/4; %u_bar the average of the uw that is calculated for each worker
    x_hat = alpha*x_bar + (1 - alpha)*z_old; %x_hat variable calculated for each worker
    z = subplus(x_bar + u_bar);    % variable z must be calculated in each iteration by the results of workers
    uw = uw + (x_hat - z);  % the new amount of uw, calculated in this iteration usable in the rest of the iterations
    labBarrier;
    history.objval(k)  = gop(@plus, objective(cw, xw));

    history.r_norm(k)  = norm(x_bar - z);
    history.s_norm(k)  = norm(-rho*(z - zold));

    history.eps_pri(k) = sqrt(2)*ABSTOL + RELTOL*max(norm(x_bar), norm(-z));
    history.eps_dual(k)= sqrt(2)*ABSTOL + RELTOL*norm(rho*u_bar);

    if ~QUIET
        fprintf('%3d\t%10.4f\t%10.4f\t%10.4f\t%10.4f\t%10.2f\n', k, ...
            history.r_norm(k), history.eps_pri(k), ...
            history.s_norm(k), history.eps_dual(k), history.objval(k));
    end

    if (history.r_norm(k) < history.eps_pri(k) && ...
       history.s_norm(k) < history.eps_dual(k))
         break;
    end
    end
        
end

%bring composite values into local workspace
%reverse loop order to get pre-allocation effects
%without needing to define structure ahead of time

for i = N : -1 : 1
    xhats(i) = x_hat{i};     
    histories(i) = history{i};
end

Answer 2

[Raymond Norris]

Hi Bill,

There might be a couple of options

1. parfeval*
2. createJob/createTask

parfor won't do the trick because it doesn't allow for early exit (i.e. if fcn<eps then stop). parfeval provides the ability to run a collection of 'futures'. Each future could be the same function with different input arguments or different functions entirely. Once the futures are submitted, loop through and query their results. If you've found your threshhold, delete the remaining futures.

I'm putting an * next to parfeval because I don't fully under how you'd update the input parameters to something that's already been submitted, which would could cause issues for parfeval. I think you're suggesting something like:

a=cell(number of workers,1);  %first input in cell
b=cell(number of workers,1);  %first input in cell
for i=1:maxiter % a loop in the client
    X(i) = identical_fcn(a{i},b{i})
    z = X ...
    % Would updating a and b effect how you're calling identical_fcn at the next iteration?
    a{i} = z ...   update some element of a, based on z before calling the next iteration of the identical_fcn
    b{i} = z ...   ditto
    if fcn(a,b,X,z)<eps
        break
    end
end

A bit more code would help refine this.

Raymond

[Bill Masters]

Hi Raymond,

thanks for your reply,

here the question is how I can make workers to call their own input parameter in each iteration?

I'm a little confused in using indices for main loop in client and the indices of workers that refers to a cell parameter that each cell belongs to a specific worker.

[Bill Masters]

%global input parameters such as z,rho and alpha defined before. they are all the same for all of the workers
x0=cell(1,4) % the amount of x0 for each worker (I have 4 workers)
u=cell(1,4) % the amount of u for each worker (I have 4 workers)
for k = 1:MAX_ITER
% x-update and calculation of u
    f(k)=parfeval(p,@tmp,1,MAX_ITER);
    [i,v]=fetchNext(f);
    x{k}=v(1:2);
    % z-update
    z_old = z; %a variable that is the same for all the workers
    x_hat{k} = alpha*x{k} + (1 - alpha)*z_old; %this variable differs from one worker to other
    z = subplus(x_hat{k} + u{k});
    u{k} = u{k-1} + (x_hat{k} - z); % the value of u in the last iter is needed
    % termination checks
    if norm(x{k}-z)<eps  % there are other termination conditions that I don't think any necessity here to increase loc in this post
    break;
    end

function f=tmp(i)
tmp = [ rho*eye(2), A'; A, zeros(2) ] \ [ rho*(z - u) - c; b];  %A,u,c,b are input parameters that differ from one worker to another and I don't know how to call the right one for workers
end  

[Raymond Norris]

Hi Bill,

So I don't make any assumptions, a couple of suggestions/thoughts

1. I think you're missing an end for the for loop (before the definition of the tmp function)
2. the tmp function should be returning (I assume f) not tmp
3. the tmp function doesn't use the input argument i
4. the call to tmp in the parfeval only passes in MAX_ITER (which maps to i). Nowhere are any of the variables rho, etc. defined. Even if they are global they won't be passed to each of the workers (for starters, you’d need to declare them as global in tmp as well, but this still won’t work).
5. parfeval and fetchNext are (typically) called from within different for loops, because fetchNext is blocking until MATLAB fetches the results. This, along with updating z_old, x_hat, and u, all indicate that you'll have issues parallelizing this. This is a cascading waterfall algorithm, u(k) is dependent on u(k-1).

Based on #5, perhaps the best way to think of parallezing this is if you need to run this entire block of code multiple times. Then each worker can be given the entire code you have listed above and you can run 4 of these (for example). They wouldn't be tied together, but you could run multiple sets, all at the same time.

But perhaps I've misconstrued some things here.

Raymond

[Bill Masters]

Hi Raymond,

thanks for your suggestions

I think using SPMD makes sense for my problem. I have written the code using spmd but yet I have errors. I would appreciate If you could check it in the below answer's comment.