1function k = commutation(n, m, sparseflag)
2% k = commutation(n, m, sparseflag)
3% -------------------------------------------------------------------------
4% Returns Magnus and Neudecker's commutation matrix of dimensions n by m,
5% that solves k*vec(X)=vec(X')
6% =========================================================================
7% INPUTS
8%   n:          [integer] row number of original matrix
9%   m:          [integer] column number of original matrix
10%   sparseflag: [integer] whether to use sparse matrices (=1) or not (else)
11% -------------------------------------------------------------------------
12% OUTPUTS
13%   k:          [n by m] commutation matrix
14% -------------------------------------------------------------------------
15% This function is called by
16%   * get_perturbation_params_derivs.m (previously getH.m)
17%   * get_identification_jacobians.m (previously getJJ.m)
18%   * pruned_state_space_system.m
19% -------------------------------------------------------------------------
20% This function calls
21%   * vec (embedded)
22% =========================================================================
23% Copyright (C) 1997 Tom Minka <minka@microsoft.com>
24% Copyright (C) 2019-2020 Dynare Team
25%
26% This file is part of Dynare.
27%
28% Dynare is free software: you can redistribute it and/or modify
29% it under the terms of the GNU General Public License as published by
30% the Free Software Foundation, either version 3 of the License, or
31% (at your option) any later version.
32%
33% Dynare is distributed in the hope that it will be useful,
34% but WITHOUT ANY WARRANTY; without even the implied warranty of
35% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
36% GNU General Public License for more details.
37%
38% You should have received a copy of the GNU General Public License
39% along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
40% =========================================================================
41% Original author: Thomas P Minka (tpminka@media.mit.edu), April 22, 2013
42
43if nargin < 2
44    m = n(2);
45    n = n(1);
46end
47if nargin < 3
48    sparseflag = 0;
49end
50
51if sparseflag
52    k = reshape(kron(vec(speye(n)), speye(m)), n*m, n*m);
53else
54    k = reshape(kron(vec(eye(n)), eye(m)), n*m, n*m);
55end
56
57function V = vec(A)
58    V = A(:);
59end
60
61end
62