%------------------------------%
% compute the momentum matrix  %
%------------------------------%

vecsigns= sign(uxold(1:J)+uxold(2:J+1));
vecsignsleft =sign(uxleft);
vecsignsright=sign(uxright);


B(1:J  ,1) =	-0.5*lambda*uxold(1:J).*(ones(J,1)+ ...
				vecsigns(1:J));
B(:    ,2) =  0.5*lambda*uxold.*([vecsigns;vecsignsright]+ ...
					    [vecsignsleft;vecsigns]);
B(2:J+1,3) =	 0.5*lambda*uxold(2:J+1).*(ones(J,1)- ...
				vecsigns(1:J));

mmat=spdiags(B,[-1 0 1],J+1,J+1);

