valid indices in a parfor loop

Can someone please explain why the following code does not work in the parfor loop? The error is about invalid indices on the matrix 'mat'.

inds = 1:ptsPerCore*nCores;

inds = reshape(inds,ptsPerCore,nCores);

mat = zeros(ptsPerCore*nCores,n);

parfor ii=1:nCores

par_inds = inds(:,ii);

for jj=1:ptsPerCore

cur_ind = par_inds(jj);

vec = randn(1,n);

mat(cur_ind,:) = vec;

end

end

Answer 1

[Matt J]

You are trying to use the parfor loop to fill a pre-existing matrix mat. The only type of parfor variable for which this is allowed is a Sliced Variable. Indexing into sliced variables must obey certain rules:

https://www.mathworks.com/help/parallel-computing/sliced-variable.html

[Neil]

Thanks @Matt J, but I still do not understand. Which rule is being violated? Is there a workaround?

[Matt J]

Which rule is being violated?

"Form of Indexing. Within the first-level of indexing for a sliced variable, exactly one indexing expression is of the form i, i+k, i-k, or k+i. The index i is the loop variable and k is a scalar integer constant or a simple (non-indexed) broadcast variable. Every other indexing expression is a positive integer constant, a simple (non-indexed) broadcast variable, a nested for-loop index variable, colon, or end"

Is there a workaround?

The whole thing looks equivalent to

mat=randn(ptsPerCore*nCores,n);

Not sure why you don't just do that.

[Neil]

Ok. I think I am getting it now, and your workaround is indeed very simple. But I should have specified that the above code was just an illustrative example.

The main goal is to have a parfor loop that accesses a subset of indices 'par_inds' out of an entire unordered set 'inds' and then to have a nested for loop that runs through all the indices of the subset.

Is that possible?

[Matt J]

It might be possible, but it would be important that the par_inds are unique and don't overlap from worker to worker. Otherwise, the tasks are not truly parallel.

[Neil]

In my case, I know that the par_inds are definitely unique and don't overlap. But how do I convince Matlab that this is true?

[Matt J]

If so, use inds to re-order the mat matrix before the loop:

matReordered=reshape( mat(inds,:)  ,[ptsPerCore,nCores,n]);

parfor i=1:nCores
    
    matReordered(:,i,:)=rand(ptsPerCore,1,n); %assign something
        
end

mat(inds,:)=matReordered(:,:);

[Neil]

Thanks @Matt J. Reordering before the loop along with the 3d output proposed by @Edric Ellis solved my problem as I needed to keep the nested loop.

Answer 2

[Edric Ellis]

There are a couple of ways to work around this. You can either unroll the inner loop, or else make the output matrix 3d so that it is sliceable (and then use reshape to get back).

%% Original code
ptsPerCore = 7;
nCores = 5;
n = 3;
inds = 1:ptsPerCore*nCores;
inds = reshape(inds,ptsPerCore,nCores);
mat = zeros(ptsPerCore*nCores,n);
for ii=1:nCores
    par_inds = inds(:,ii);
    
    for jj=1:ptsPerCore
        cur_ind = par_inds(jj);
        vec = cur_ind.*(1:n);
        mat(cur_ind,:) = vec;
    end
end

%% Alternative 1 - run as a single loop
mat2 = zeros(ptsPerCore*nCores,n);
parfor iijj = 1:(nCores * ptsPerCore)
    vec = iijj.*(1:n);
    mat2(iijj,:) = vec;
end

Starting parallel pool (parpool) using the 'local' profile ...
Connected to the parallel pool (number of workers: 2).

assert(isequal(mat2, mat));

%% Alternative 2 - use a 3d output
mat3 = zeros(ptsPerCore, nCores, n);
parfor ii=1:nCores
    par_inds = inds(:,ii);
    
    for jj=1:ptsPerCore
        cur_ind = par_inds(jj);
        vec = cur_ind.*(1:n);
        mat3(jj, ii, :) = vec;
    end
end
% Need to reshape back to original size
mat3 = reshape(mat3, nCores*ptsPerCore, n);
assert(isequal(mat3, mat));

[Neil]

Thanks @Edric Ellis. The 3d output worked when I reordered before the loop as proposed by @Matt J.