 function y = Mprod(x)
% y = Mx where x is a column n-vector, M = eye(n,n) + S_n,
% and S_n is the n-by-n shift matrix
  n = length(x);
  y = x + [x(2:n);x(1)];