function A = laplacian2d_sparse2(k)
% LAPLACIAN2D_SPARSE2 : discrete Laplace operator in 2 dimensions
%
% A = laplacian2d_sparse2(k) returns a k^2-by-k^2 sparse matrix

m = k^2;

IJV = zeros(5*m,3);  % enough space for "nonzero triples"
nz = 0;

for x = 1:k
    for y = 1:k
        i = (y-1)*k + x;
        
        % A(i,i) = 4;
        nz = nz+1; 
        IJV(nz,:) = [i i 4];
        
        if y > 1
            % A(i,i-k) = -1;
            nz = nz+1;
            IJV(nz,:) = [i i-k -1];
        end;
        if y < k
            % A(i,i+k) = -1;
            nz = nz+1;
            IJV(nz,:) = [i i+k -1];
        end;
        if x > 1
            % A(i,i-1) = -1;
            nz = nz+1;
            IJV(nz,:) = [i i-1 -1];
        end;
        if x < k
            % A(i,i-k) = -1;
            nz = nz+1;
            IJV(nz,:) = [i i+1 -1];
        end;     
    end;
end;

I = IJV(1:nz,1);
J = IJV(1:nz,2);
V = IJV(1:nz,3);
A = sparse(I,J,V,m,m);
