code real RlRes=46, RlOut=48;
def  S = 10;

proc SetBoundary(MV, MF);
real MV; int MF;
[MF(1,1):= 1;  MV(1,1):= 1.0;
 MF(6,7):= -1; MV(6,7):= -1.0;
];

func real CalcDiff(MV, MF, DV, W, H);
real MV; int MF; real DV; int W, H;
int  I, J, N; real V, Total;
[Total:= 0.0;
for I:= 0 to H-1 do
    for J:= 0 to W-1 do
        [V:= 0.0;  N:= 0;
        if I then [V:= V + MV(I-1,J);  N:= N+1];
        if J then [V:= V + MV(I,J-1);  N:= N+1];
        if I+1 < H then [V:= V + MV(I+1,J);  N:= N+1];
        if J+1 < W then [V:= V + MV(I,J+1);  N:= N+1];
        V:= MV(I,J) - V/float(N);  DV(I,J):= V;
        if MF(I,J) = 0 then Total:= Total + V*V;
        ];
return Total;
];

func real Iter(MV, MF, W, H);
real MV; int MF, W, H;
real DV, Diff, Cur; int I, J;
[DV:= RlRes(W);  for I:= 0 to W-1 do DV(I):= RlRes(H);
Diff:= 1E10;
Cur:= [0.0, 0.0, 0.0];
while Diff > 1E-24 do
        [SetBoundary(MV, MF);
        Diff:= CalcDiff(MV, MF, DV, W, H);
        for I:= 0 to H-1 do
            for J:= 0 to W-1 do
                MV(I,J):= MV(I,J) - DV(I,J);
        ];
for I:= 0 to H-1 do
    for J:= 0 to W-1 do
        Cur(MF(I,J)+1):= Cur(MF(I,J)+1) +
                DV(I,J) * float(-(I>0) - (J>0) - (I<H-1) - (J<W-1));
                                \middle=4; side=3; corner=2
return (Cur(2)-Cur(0))/2.0;
];

real MeshV(S,S); int MeshF(S,S);
RlOut(0, 2.0 / Iter(MeshV, MeshF, S, S))