Matlab - skrypt z metodą jacobiego

Witam, mam problem z tym zadaniem. Mam do wyznaczenia wartości w węzłach środkowych metodą prostej iteracji. Model ma kształt nieregularny i zadane warunki brzegowe jak poniżej. Z góry dziękuję za pomoc!

%Model parameters

Lx=15;
Ly=15;

%Number of nodes

nx=Lx+1;
ny=Ly+1;

h=zeros(ny,nx);

%Boundary conditions

%Right Nodes

for i = 1:5
for j= 12:16
h(i,j) = NaN;
end
end

%Upper pressure1

h(1,1:11)=linspace(7,7,11);

%Down pressure

h(ny,1:nx)=linspace(6,6,nx);

%Left side

h(2:ny-1,1)=linspace(4,4,ny-2);

%Right side1

h(2:5,11)=linspace(13,13,4);

%Upper pressure2

h(6,11:16)=linspace(5,5,6);

%Right side2

h(7:15,nx)=linspace(9,9,9);

tol=1d-6;
err=1;
m=0;
hkp1=h;

while err > tol

m = m+1;

%Left region

for i = 2:ny-1
    for j = 2:10
        
        hkp1(i,j) = 0.25*(h(i+1,j)+h(i,j+1)+h(i-1,j)+h(i,j-1));
        
    end
end

%Right region

for i = 7:15
    for j = 11:15
        
        hkp1(i,j) = 0.25*(h(i+1,j)+h(i,j+1)+h(i-1,j)+h(i,j-1));
        
    end
end

err = sqrt(sum(sum((hkp1-h).^2)));

h = hkp1;

end

[X,Y] = meshgrid(0:1:Lx,0:-1:-Ly);
ch = contour(X,Y,h,10, 'color', 'r');