config const n = 6, epsilon = 1.0e-5, p = 2.0, t = 10.0;

const BigD = {0..n+1, 0..n+1},
D = BigD[1..n, 1..n];

var A, Temp : [BigD] real;
var red_odd: sparse subdomain(D);
var red_even: sparse subdomain(D);
var red: sparse subdomain(D);

red_odd+= {1..n, 1..n} by 2 align (1,1);
red_even+= {1..n, 1..n} by 2 align (2,2);
red = red_odd + red_even;

var black_odd: sparse subdomain(D);
var black_even: sparse subdomain(D);
var black: sparse subdomain(D);

black_odd+= {1..n, 1..n} by 2 align (1,2);
black_even+= {1..n, 1..n} by 2 align (2,1);
black= black_odd + black_even;

do {
    forall (i, j) in red do
        A[i,j] = (A[i-1,j] + A[i+1,j] + A[i,j-1] + A[i,j+1] + p/t ) / 4;
    forall (i, j) in black do
        A[i,j] = (A[i-1,j] + A[i+1,j] + A[i,j-1] + A[i,j+1] + p/t ) / 4;

    const delta = max reduce abs(A[D] - Temp[D]);
    Temp[D] = A[D];

} while (delta > epsilon);

writeln(A);