parfor with nontrivial indexing

Suppose I want to use parfor to do something in parallel with 2D array. The problem is that I cannot use complex expressions as index. So if anyone can help me figure out this example, I will be grateful!

M=zeros(5);
parfor k=1:5^2
     i=ceil(k/5);
     j=mod(k-1,5)+1;
     M(i,j)=k^2;
end

In my case there are complex expressions dependent on other arrays (PDE solution) instead of k^2.

The error is "The variable M in parfor cannot be classified."

Another problem with the same error:

M=cell(2);
M{1}=zeros(5);
parfor k=1:5
     M{1}(k)=k;
end

For some reason parpool does not allow to address specific element of array in a cell. Please, help!

Answer 1

[Walter Roberson]

In the case of your first example, you can transform it to linear indexing with

   parfor k = 1 : 5^2
     M(k) = k^2;
   end
   M = reshape(M,5,5).';

the transpose being there as your indexing with mod is equivalent to filling sequentially by rows instead of by columns.

I know that you indicated that your computations are more complex than that; on the other hand, this illustrates that sometimes you can reframe the problem.

If you have a computation that is guaranteed to fill every matrix location exactly once, but not sequentially in your controlling variable, then you can store the results in a vector indexed by the controlling variable (which parfor will be able to see clearly is a unique address), and you can also store a companion vector which is the linear index of where the element should go in the final array. Then after the parfor, you can index the value array at the companion vector to get the elements in linear order, and then reshape() it to the desired dimension. You might find it more efficient to use a companion array, the rows of which are the indices you would like to store at, as afterwards you can do a mass conversion of indices to linear index:

companionlin = sub2idx(DesiredSize, companionidx(:,1), companionidx(:,2));
finalarray = reshape(valuevector(companionlin),DesiredSize);

If there might be duplicates in the companion index list then the implication is that you are attempting to store multiple times to the same location, and the result of that would depend upon the order of operations.

Unless, that is, that you are doing a summation or multiplication reduction in those locations, in which case the result would be the same up to round-off error, provided you operated sequentially; if you try to do a reduction in parallel then unless MATLAB can clearly deduce that you are doing a reduction (in which case it would be able to classify the variable), MATLAB would not be able to promise that different parallel workers are not trying to write to the same location at the same time.

If you are trying to do a parallel summation or multiplication reduction into individual locations that are not sequential, then the way to handle it is to make the value list and companion index array I referred to before (possibly recording the individual indices instead of converting to linear, if that is what is most convenient), and then use accumarray() to do the reduction.