ERROR IN PARALLEL COMPUTING HEAT EQUATION

Hello guys, so i've been working to this project and i've made the serial code for heat equation but when i try to do parallel computing it doesnt work and keep saying

'Unable to classify the variable 'C' in the body of the parfor-loop.'

I tried parallel computing on two codes

1.

tic;
lx = 200;
ly = 200;
h = 50;
x0 = 75;
xn = 50;
y0 = 0;
yn = 100;
%memanggil fungsi array batas
    nx = 4;
    ny = 4;
    S =ones(1,nx+1);
    R =ones(1,ny+1);
    
    %mendefinisikan baris dan kolom batas distribusi suhu
    bx0 = S*x0; bx0(1) = (x0+y0)/2; bx0(nx+1) = (x0+yn)/2;
    bxn = S*xn; bxn(1) = (xn+y0)/2; bxn(nx+1) = (xn+yn)/2;
    by0 = R*y0; by0(1) = (x0+y0)/2; by0(ny+1) = (xn+y0)/2;
    byn = R*yn; byn(1) = (x0+yn)/2; byn(ny+1) = (xn+yn)/2;
 
%memanggil elemen matriks perhitungan persamaan aljabar yang diperlukan
    nn = (nx-1)*(ny-1);
    ni = zeros(ny+1, nx+1);
    b = zeros(nn,1);
    C = zeros(nn,nn);
    parfor i = 2:ny
        for j = 2:nx
            ij = (i-2)*(nx-1)+j-1;
            ni(i,j) = ij;
        end
    end
    
    parfor i = 2:ny
        for j = 2:nx
            ij = (i-2)*(nx-1)+j-1;
            C(ij,ij) = -4;
            %Sisi Atas
            if ni(i-1,j) ~= 0
                ic = ni(i-1,j);
                C(ij,ic) = 1;
            else
                b(ij) = b(ij) - by0(j);
            end
            
            %Sisi Bawah
            if ni(i+1,j) ~= 0
                ic = ni(i+1,j);
                C(ij,ic) = 1;
            else
                b(ij) = b(ij) - byn(j);
            end
            
            %Sisi Kiri
            if ni(i,j-1) ~= 0
                ic = ni(i, j-1);
                C(ij,ic) = 1;
            else
                b(ij) = b(ij) - bx0(i)
            end
            
            %Sisi Kanan
            if ni(i,j+1) ~= 0
                ic = ni(i,j+1);
                C(ij,ic) = 1;
            else
                b(ij) = b(ij)-bxn(i);
            end
        end
    end


%memanggil fungsi penyelesaian persamaan linier aljabar gauss-seidel
    n = length(b)
    m = length(C)
    
    parfor i = 1:n
        b(i) = b(i)/C(i,i);
        C(i,:) = C(i,:)./C(i,i)
        C(i,i) = 0 %membuat diagonal matriks menjadi nol
    end
    
    G = zeros(m,1); %membuat matriks awal nol
    iter = 1;
    while(1)
        Golb = G;
        parfor j=1:n
            G(j) = b(j) - (C(j,:)*G);
            e = (G(j) - Golb(j))/G(j);
            ea(j) = abs(e);
        end
        iter = iter+1;
        if(min(ea)<eps), break, end
    end

%memanggil matriks kontur/solusi
H = zeros(ny+1, nx+1);
H(:,1) = bx0;
H(:, nx+1) = bxn;
H(1,:) = by0;
H(ny+1,:) = byn;

nn = (nx-1)*(ny-1);
k = [1:(nx-1):nn];
parfor i = 2:nx
    j = k(i-1);
    for p=2:ny
        H(i,p) = G(j+p-2);
    end
end
H
t1 = toc;
disp(['proses: ' num2str(t1)])

2. The simpler code

'Unable to classify the variable 'U' in the body of the parfor-loop.'

tic
clc;
clear all;
a =1;
b =1;
c =1;
xgrid =16;
ygrid =16;
tgrid =16;
deltax =a/(xgrid-1);
deltay =b/(ygrid-1);
deltat =c/(tgrid-1);
U = zeros(ygrid,xgrid);
U(:,:) =0; %syarat awal
U(:,1) =75; %left
U(:,xgrid) =50; %right
U(1,:) =100; %bottom
U(ygrid,:) =0; %top
alpha =0.1;
Uold =U;
x = linspace(0,a,xgrid);
y = linspace(0,b,ygrid);
t = linspace(0,c,tgrid);
parfor k=1:tgrid
    t1 = deltat*(k);
    for j=2:ygrid-1
        for i=2:xgrid-1
            U(j,i)=(alpha^2*deltat/(deltax^2*deltay^2))*...
                (deltay^2*(Uold(j,i+1)-2*Uold(j,i)+Uold(j,i-1))...
                +deltax^2*(Uold(j+1,i)-2*Uold(j,i)+Uold(j-1,i)))...
                +Uold(j,i);
        end
    end
    figure(1)
    image(x,y,U);
    colorbar
    title(['t= ', num2str(t1), 'detik']);
    Uold = U;
end
U
t1 = toc;
disp(['proses: ' num2str(t1)])

[Iqamatulhaq Rakayama]

pls help!

Answer 1

[Altaïr]

Hey @Iqamatulhaq Rakayama,

The use of a parallel pool is suitable when each iteration in a loop is independent and can be executed in any order. More information on this can be accessed using the command:

web(fullfile(docroot, 'parallel-computing/decide-when-to-use-parfor.html'))

The outer for-loop with "k=1:tgrid" as the control variable cannot be parallelized due to dependencies. However, in the second implementation the inner for-loops can be parallelized as demonstrated below:

U_sliced = U(2:ygrid-1, 2:xgrid-1);
    parfor j=2:ygrid-1
        temp = zeros([1 xgrid-2]);
        for i=2:xgrid-1
            temp(1,i-1)=(alpha^2*deltat/(deltax^2*deltay^2))*...
                (deltay^2*(Uold(j,i+1)-2*Uold(j,i)+Uold(j,i-1))...
                +deltax^2*(Uold(j+1,i)-2*Uold(j,i)+Uold(j-1,i)))...
                +Uold(j,i);
        end
        U_sliced(j-1,:) = temp;
    end
U(2:ygrid-1, 2:xgrid-1) = U_sliced;

Slicing the matrix U and using a temporary array, temp, is necessary. This is because, if a nested for-loop is used to index into a sliced array, that array cannot be used elsewhere in the parfor loop. More details can be found in the documentation:

web(fullfile(docroot, 'parallel-computing/troubleshoot-variables-in-parfor-loops.html'))

For further details on the classification of variables in a parfor loop, the following command can be used:

web(fullfile(docroot, 'coder/ug/classification-of-variables-in-parfor-loops.html'))

The full code is attached for reference.