Why do I get this assignment-size error only in parfor, but not in for loops?

Hi

I am trying to do something fairly simple. I want to use ode45 for solving an ode for a large number of initial conditions.

The function that I pass to ode45 itself calls fairly simple functions, but I think that is not relevant. At least I hope so.

Running the loop over the ininital conditions in a normal for loop runs works perfectly, but when I replace for with parfor, I get the error

"Unable to perform assignment because the left and right sides have a different number of elements."

What I've tried (as you can see in the code) is to replace the initial condition with a constant [1,1,1,1,1]. The error stays the same.

p_initialCond = [IC.x0(:,:)',IC.phi0(:),IC.k0(:,:)']; %<--- not used
cells1 = cell(npar,1);
cells2 = cell(npar,1);
parfor i = 1:npar
    [cells1{i},cells2{i}]= ode45(@f_dXdt_annular_current,...
                                [0,1],...
                                [1,1,1,1,1],...
                                options...
                                );
end

The error occurs in line 26 of the following code

function [dXdt] = f_dXdt_annular_current(t,X)
% rename for reuse
dXdt = zeros(5,1);
x   = X(1); 
y   = X(2);
kx  = X(4);
ky  = X(5);

% velocity field 
global r1; global r2; global umax;
[U,V] = f_annular_current(x,y,umax,r1,r2);
% velocity field derivatives 
[Ux,Uy,Vx,Vy] = f_annular_current_derivatives(x,y,umax,r1,r2);

[dOmdkx,dOmdky] = f_Omega_k(kx,ky,U,V);
[dOmdx,dOmdy] = f_Omega_x(kx,ky,Ux,Uy,Vx,Vy);

% calculate derivative vector
dXdt(1) = dOmdkx; %   ERROR LINE<-------------------- line 26
dXdt(2) = dOmdky;
dXdt(3) = kx*dOmdkx+ky*dOmdky;
dXdt(4) = -dOmdx;
dXdt(5) = -dOmdy;
end

The two relevant functions called are

function [u,v,u_theta] = f_annular_current(x,y,umax,r0,R0)
r = sqrt(x.^2+y.^2);
u_theta = -umax*4*(r-r0).*(r-R0)/(R0-r0)^2;
outside = r>R0|r<r0;
u_theta(outside) = 0;

u = u_theta .* -y./r;
v = u_theta .*  x./r;

u(r==0) = 0;
v(r==0) = 0;
end

and

function [Om_kx,Om_ky] = f_Omega_k(kx,ky,U,V)
g = 9.81;
k = sqrt(kx^2+ky^2);

Om_kx = 0.5   * sqrt(g/k)*  kx/k + U;
Om_ky = 0.5   * sqrt(g/k)*  ky/k + V;
end

both of which are fairly basic functions, I assume, so I don't understand how calling everything in a parfor loop causes this.

The full error message I have received is:

Error using f_dXdt_annular_current (line 26)

Unable to perform assignment because the left and right sides have a different number of elements.

Error in odearguments (line 90)

f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.

Error in ode45 (line 115)

odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);

Error in First_Try_phase_field (line 86)

parfor i = 1:npar

It would be awesome if someone knows what's going on. Thanks!

[Alvaro]

The line that defines U

u = u_theta .* -y./r;

and the line that defines the future value of dOmdkx

Om_kx = 0.5   * sqrt(g/k)*  kx/k + U;

seem to the suggest that you are expecting dOmdkx to be a multidimensional matrix which then you are trying to store inside of a vector with

dXdt(1) = dOmdkx;

This could be the source of the error but it is hard to tell without knowing what the other variables are. Two immediate suggestions would be to first try getting rid of those global variables which tend to cause issues with parfor.

https://www.mathworks.com/matlabcentral/answers/?term=global+variables+parfor

You can instead pass those as broadcast variables. If that doesn't work, the second suggestion would be to remove the problematic assignment and write a dummy function to return the value dOmdkx to the base workspace and see whether it is what you expect.