Search code examples
matlabfor-loopnested-loopsparfor

Using a for loop inside a parallel for loop in Matlab


I am trying to use a for loop inside of a parfor loop in Matlab.
The for loop is equivalent to the ballode example in here.
Inside the for loop a function ballBouncing is called which is a system of 6 differential equations.

So, what I am trying to do is to use 500 different sets of parameter values for the ODE system and run it, but for each parameter set, a sudden impulse is added, which is handled through the code in 'for' loop.

However, I don't understand how to implement this using a parfor and a for loop as below.
I could run this code by using two for loops but when the outer loop is made to be a parfor it gives the errors,
the PARFOR loop cannot run due to the way variable results is used,
the PARFOR loop cannot run due to the way variable y0 is used and
Valid indices for results are restricted in PARFOR loops

results=NaN(500,100);
x=rand(500,10);

parfor j=1:500

    bouncingTimes=[10,50];%at time 10 a sudden impulse is added
    refine=2;
    tout=0;
    yout=y0;%initial conditions of ODE system
    paras=x(j,:);%parameter values for the ODE 
    for i=1:2
        tfinal=bouncingTimes(i);
        [t,y]=ode45(@(t,y)ballBouncing(t,y,paras),tstart:1:tfinal,y0,options);
        nt=length(t);
        tout=[tout;t(2:nt)];
        yout=[yout;y(2:nt,:)];

        y0(1:5)=y(nt,1:5);%updating initial conditions with the impulse
        y0(6)=y(nt,6)+paras(j,10);

        options = odeset(options,'InitialStep',t(nt)-t(nt-refine),...
                                 'MaxStep',t(nt)-t(1));
        tstart =t(nt);
    end

    numRows=length(yout(:,1));
    results(1:numRows,j)=yout(:,1);

end
results;

Can someone help me to implement this using a parfor outer loop.


Solution

  • Fixing the assignment into results is relatively straightforward - what you need to do is ensure you always assign a whole column. Here's how I would do that:

    % We will always need the full size of results in dimension 1
    numRows = size(results, 1);
    parfor j = ...
        yout = ...; % variable size
        yout(end:numRows, :) = NaN; % Expand if necessary
        results(:, j) = yout(1:numRows, 1); % Shrink 'yout' if necessary
    end
    

    However, y0 is harder to deal with - the iterations of your parfor loop are not order-independent because of the way you're passing information from one iteration to the next. parfor can only handle loops where the iterations are order-independent.