function A=delsq3d(m,n,p);
% A=delsq3d(m,n,p)
% makes 3 dimensional Laplacian over m*n*p grid
%
Tm=spdiags(ones(m,1)*[-1 2 -1],-1:1,m,m);
Tn=spdiags(ones(n,1)*[-1 2 -1],-1:1,n,n);
Tp=spdiags(ones(p,1)*[-1 2 -1],-1:1,p,p);

A=kron(Tm,speye(n*p))+kron(speye(m),kron(Tn,speye(p)))+kron(speye(m*n),Tp);