function A = coeffmat(n)

B = 4*diag(ones(1, n)) - diag(ones(1, n-1), 1) - diag(ones(1, n-1), ...
    -1);
O = zeros(n);
nI = -diag(ones(1,n));

A = [B nI zeros(n, (n-2)*n)];

for i=2:(n-1)
  A = [A; zeros(n, (i-2)*n) nI B nI zeros(n, (n-i-1)*n)];
end

A = [A; zeros(n, (n-2)*n) nI B];
