Error using oddiyxy1 (line 23) Error: The variable a in a parfor cannot be classified. See Parallel for Loops in MATLAB, "Overview".

n=91;         h=1/(n-1);        nt=1080;        tau=0.001; 
for i=1:n 
    x(i)=(i-1)*h;       y(i)=(i-1)*h;  
end  

for j=1:n
    for i=1:n 
        p0(i,j)=1;    
    end   
end

for jt=1:nt 
    t(jt)=jt*tau;    
    for j=1:n  
        for i=1:n 
            pt(i,j)=exp(-2*x(i)*y(j)*t(jt));
        end      
    end   
end

tic
parfor jt=1:nt        
    for j=2:n-1     
        for i=2:n-1  
            a(i)=1;       c(i)=1;       b(i)=2+2*h*h/tau;
            ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));   
            d(i)=2*h*h/tau*p0(i,j)+(p0(i,j-1)-2*p0(i,j)+p0(i,j+1))+ff;  
        end
        aa1(1)=0; bb1(1)=1;    
        for i=2:n-1
        aa1(i)=c(i)/(b(i)-aa1(i-1)*a(i));      bb1(i)=(d(i)+bb1(i-1)*a(i))/(b(i)-aa1(i-1)*a(i));
        end
        p(n,j)=exp(-2*y(j)*t(jt));  
        for i=n-1:-1:1 
            p(i,j)=aa1(i)*p(i+1,j)+bb1(i);
        end                   
    end     
    
    for i=1:n 
        p(i,n)=exp(-2*x(i)*t(jt));          p(i,1)=1;
    end      
    
    for j=1:n 
        for i=1:n      
            p0(i,j)=p(i,j);  
        end  
    end         
    
    for i=2:n-1  
        for j=2:n-1  
            a(j)=1;       c(j)=1;       b(j)=2+2*h*h/tau;
            ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));   
            d(j)=2*h*h/tau*p0(i,j)+(p0(i-1,j)-2*p0(i,j)+p0(i+1,j))+ff; 
        end
        aa1(1)=0;       bb1(1)=1;  
        for j=2:n-1 
        aa1(j)=c(j)/(b(j)-aa1(j-1)*a(j));      bb1(j)=(d(j)+bb1(j-1)*a(j))/(b(j)-aa1(j-1)*a(j));
        end
        p(i,n)=exp(-2*x(i)*t(jt));  
        for j=n-1:-1:1 
            p(i,j)=aa1(j)*p(i,j+1)+bb1(j);
        end    
    end          
    
    for j=1:n  
        p(n,j)=exp(-2*y(j)*t(jt));        p(1,j)=1;
    end  
    
    for i=1:n   
        for j=1:n      
            p0(i,j)=p(i,j);  
        end  
    end       
end 
toc

for j=1:n
    px(j)=p0(21,j);        pxt(j)=pt(21,j);
end             

x=0:h:1;
y=0:h:1;
figure ('Position',[600 160 700 600]); 
subplot(3,2,1);
meshc(x,y,pt) 
title('АНИК ЕЧИМ 3D ГРАФИГИ')

subplot(3,2,2);
meshc(x,y,p0) 
title('CОНЛИ ЕЧИМ 3D ГРАФИГИ')

subplot(3,2,3);
contour(x,y,pt,'ShowText','on','LineWidth',2); 
title('АНИК ЕЧИМ КОНТУР ГРАФИГИ')

subplot(3,2,4);
contour(x,y,p0,'ShowText','on','LineWidth',2); 
title('CОНЛИ ЕЧИМ КОНТУР ГРАФИГИ')

subplot(3,2,5);
plot(x,px,x,pxt); 
title('КЕСИМДА АНИК ВА CОНЛИ ЕЧИМ ЎЗГАРИШИ X БЎЙИЧА')

Answer 1

[Raymond Norris]

I'm not sure how long it took you, but your code as it's written only took 1-2s, which wouldn't make a good problem to parallelize. I'm assuming n is much larger and that the code takes much longer to run.

With that said, I'm not sure this can be parallelized. To begin with, I'll show you how I would resolve the classification of a (and there are others, such as b and c). In doing so, you'll see how the parfor has a loop dependence on p0.

I've also added a few suggestions to cleanup/vectorize your code (listed as "% MW") -- there's more that could be done, but this gives you the idea.

function ismailov
n=91;         h=1/(n-1);        nt=1080;        tau=0.001;

% MW: Preallocate x & y
x = nan(n,1);
y = nan(n,1);
for i=1:n
    x(i)=(i-1)*h;       y(i)=(i-1)*h;
end

% MW: No need for double for-loop
% for j=1:n
%     for i=1:n
%         p0(i,j)=1;
%     end
% end
p0 = ones(n);

t = nan(nt,1);
pt = nan(n);
for jt=1:nt
    t(jt)=jt*tau;
    for j=1:n   % Anlitik yechimni hisoblash
        for i=1:n
            pt(i,j)=exp(-2*x(i)*y(j)*t(jt));
        end
    end
end

tic
parfor jt=1:nt    % Vaqt bo`yicha tsikl}
    p0 = unit_of_work(n, h, tau, jt, t, x, y, p0);
end % jt
toc

% MW: No need for double for-loop
% for j=1:n
%     px(j)=p0(21,j);
%     pxt(j)=pt(21,j);
% end
px = p0(21,:);
pxt = pt(21,:);

x=0:h:1;
y=0:h:1;
figure ('Position',[600 160 700 600]);
subplot(3,2,1);
meshc(x,y,pt) %Uch o'lchovli 3D grafika
title('АНИК ЕЧИМ 3D ГРАФИГИ')

subplot(3,2,2);
meshc(x,y,p0) %Uch o'lchovli 3D grafika
title('CОНЛИ ЕЧИМ 3D ГРАФИГИ')

subplot(3,2,3);
contour(x,y,pt,'ShowText','on','LineWidth',2); %kontur chizish
title('АНИК ЕЧИМ КОНТУР ГРАФИГИ')

subplot(3,2,4);
contour(x,y,p0,'ShowText','on','LineWidth',2); %kontur chizish
title('CОНЛИ ЕЧИМ КОНТУР ГРАФИГИ')

subplot(3,2,5);
plot(x,px,x,pxt); %Ikki o'lchovli grafika
title('КЕСИМДА АНИК ВА CОНЛИ ЕЧИМ ЎЗГАРИШИ X БЎЙИЧА')
end

function p0 = unit_of_work(n, h, tau, jt, t, x, y, p0)
% MW: Preallocate a, b, c, d, aa1, bb1, p
a = ones(1, n-1);
c = ones(1, n-1);
b = rand(1, n-1);
d = rand(1, n-1);
aa1 = rand(1, n-1);
bb1 = rand(1, n-1);
p = rand(n);
for j=2:n-1     % Tenglamalar koeffitsientini hisoblash
    for i=2:n-1
        % a(i)=1;       c(i)=1;       b(i)=2+2*h*h/tau;
        b(i)=2+2*h*h/tau;
        ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
        d(i)=2*h*h/tau*p0(i,j)+(p0(i,j-1)-2*p0(i,j)+p0(i,j+1))+ff;
    end
    aa1(1)=0; bb1(1)=1;
    for i=2:n-1
        aa1(i)=c(i)/(b(i)-aa1(i-1)*a(i));
        bb1(i)=(d(i)+bb1(i-1)*a(i))/(b(i)-aa1(i-1)*a(i));
    end
    p(n,j)=exp(-2*y(j)*t(jt));
    for i=n-1:-1:1
        p(i,j)=aa1(i)*p(i+1,j)+bb1(i);
    end
end %j

for i=1:n % Chegaraviy qiymatlarni hisoblash
    %   p(i,n)=(4.0*p(i,n-1)-p(i,n-2)+h)/3;     %   p(i,1)=(4.0*p(i,2)-p(i,3)+2*h)/3;
    p(i,n)=exp(-2*x(i)*t(jt));          p(i,1)=1;
end

% MW: No need for double for-loop
% for j=1:n  % k+1 qatlam uchun boshlangich qiymatni hisoblash
%     for i=1:n
%         p0(i,j)=p(i,j);
%     end
% end
p0 = p;

for i=2:n-1    % Tenglamalar koeffitsientini hisoblash
    for j=2:n-1
        a(j)=1;       c(j)=1;       b(j)=2+2*h*h/tau;
        ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
        d(j)=2*h*h/tau*p0(i,j)+(p0(i-1,j)-2*p0(i,j)+p0(i+1,j))+ff;
    end
    aa1(1)=0; bb1(1)=1;
    for j=2:n-1
        aa1(j)=c(j)/(b(j)-aa1(j-1)*a(j));
        bb1(j)=(d(j)+bb1(j-1)*a(j))/(b(j)-aa1(j-1)*a(j));
    end
    p(i,n)=exp(-2*x(i)*t(jt));
    for j=n-1:-1:1
        p(i,j)=aa1(j)*p(i,j+1)+bb1(j);
    end
end  %i

for j=1:n  % Chegaraviy qiymatlarni hisoblash
    %p(n,j)=(4.0*p(n-1,j)-p(n-2,j)+2*h)/3;      %p(1,j)=(4.0*p(2,j)-p(3,j)+2*h)/3;
    p(n,j)=exp(-2*y(j)*t(jt));        p(1,j)=1;
end

% MW: No need for double for-loop
% for i=1:n   % keyingi qatlam uchun boshlangich qiymatni hisoblash
%     for j=1:n
%         p0(i,j)=p(i,j);
%     end
% end
p0 = p;
end

Let's look closer at the parfor loop

parfor jt=1:nt    % Vaqt bo`yicha tsikl}
    p0 = unit_of_work(n, h, tau, jt, t, x, y, p0);
end % jt

I've refactored the for-loop to call unit_of_work. The only variable modified in your for-loop is p0. But notice that p0 is not dependent on the loop variable jt. For instance

parfor jt=1:nt
    p0(jt,:) = unit_of_work(n, h, tau, jt, t, x, y);
end % jt

In this case, we would be building up each row independent of any other row, which then allows us to assign each row to a worker. Instead, the entire matrix p0 is dependent on the previous state of the matrix -- what I'd call a cascading algorithm -- and therefore can't be parallelized.

Answer 2

[Walter Roberson]

Minimal example:

parfor jt=1:1
    a(1) = 1;
end 

Error: Unable to classify the variable 'a' in the body of the parfor-loop. For more information, see Parallel for Loops in MATLAB, "Solve Variable Classification Issues in parfor-Loops".

For whatever reason, parfor is unable to determine that a is acting as a local variable. If you assign to all of a within the parfor loop, such as if you had

a = zeros(1,n-1);

then MATLAB would be able to figure it out. As it is, MATLAB worries that there might be an existing a or that a might be used after the loop.