% Perform a matrix vector multiplication b = A*x with a matrix A being a
% Kronecker product given by A = kron( kron(...,As{2}), As{1} ).
function b = kronmvm(As,x,transp)
if nargin>2 && ~isempty(transp) && transp   % transposition by transposing parts
  for i=1:numel(As)
    if isnumeric(As{i})
      As{i} = As{i}';
    else
      As{i}.mvm = As{i}.mvmt;
      As{i}.size = [As{i}.size(2),As{i}.size(1)];
    end
  end
end
m = zeros(numel(As),1); n = zeros(numel(As),1);                  % extract sizes
for i=1:numel(n)
  if isnumeric(As{i})
    [m(i),n(i)] = size(As{i});
  else
    m(i) = As{i}.size(1); n(i) = As{i}.size(2);
  end
end
d = size(x,2);
b = x;
for i=1:numel(n)
  a = reshape(b,[prod(m(1:i-1)), n(i), prod(n(i+1:end))*d]);    % prepare  input
  if isnumeric(As{i})                                % apply along one dimension
    tmp = reshape(permute(a,[1,3,2]),[],n(i))*As{i}';
    b = permute(reshape(tmp,[size(a,1),size(a,3),m(i)]),[1,3,2]);
  else
    b = zeros([prod(m(1:i-1)), m(i), prod(n(i+1:end))*d]);      % prepare output
    for j=1:size(b,3), b(:,:,j) = As{i}.mvm(a(:,:,j)')'; end
  end
end
b = reshape(b,prod(m),d);                        % bring result in correct shape
