function output = commitment(Omega,A,B,SigmaE,C,Ci,W,beta,n1,n2,X10,Theta1,Theta2);
%-----------------------------------------------------------------------------------------
% PURPOSE: performs optimal commitment under certainty for Svensson and 
%          Woodford (2002), "Indicator Variables for Optimal Policy"           
%-----------------------------------------------------------------------------------------
% USAGE: output  = commitment(Omega,A,B,SigmaE,C,Ci,W,beta,n1,n2,X10,Theta1,Theta2)
% where:  Omega  = n2xn2 (possibly singular) matrix structural coefficients
%         A      = nxn matrix of structural coefficients for [X(t)' x(t)']
%         B      = nxk matrix of structural coefficients for i(t)
%         SigmaE = n1xn1 covariance matrix of structural shocks e(t)
%         C      = nyxn matrix of target variable coefficients for [X(t)' x(t)']
%         Ci     = nyxk matrix of target variable coefficients for i(t)
%         W      = nyxny matrix of policy weights
%         beta   = discount rate in the central bank's loss function
%         n1     = number of predetermined state variables
%         n2     = number of forward looking state variables
%         X10    = initial value of the state vector X(t) for t = 0
%         Theta1 = ntxn matrix of extra-model coefficients for [X(t)' x(t)']
%         Theta2 = ntxk matrix of extra-model coefficients for i(t)
%----------------------------------------------------------------------------------------- 
%  model setup:
%  [X(t+1)' Omega*x(t+1|t)'] = A*[X(t)' x(t)'] + B*i(t) + [e(t+1)' 0']
%  Y(t) = C*[X(t)' x(t)'] + Ci*i(t)
%  L(t) = Y(t)'*W*Y(t)
%  T(t) = Theta1*[X(t)' x(t)'] + Theta2*i(t)
%
%  model solution:
%  i(t)    = F*X(t) + Phi*lam2(t-1)
%  x(t)    = G*X(t) + Gam*lam2(t-1)
%  lam1(t) = T*X(t) + Theta*lam2(t-1)
%  X(t+1)  = H*X(t) + Psi*lam2(t-1) + e(t+1)
%  lam2(t) = S*X(t) + Sig*lam2(t-1)
%-----------------------------------------------------------------------------------------
% RETURNS: a structure
% output.F       = F     output.Phi   = Phi
% output.G       = G     output.Gam   = Gam
% output.T       = T     output.Theta = Theta
% output.H       = H     output.Psi   = Psi
% output.S       = S     output.Sig   = Sig
% output.numb    = number of unstable eigenvalues in the linear system
% output.Loss    = Loss 
% output.weights = central bank policy weights
% output.CovT    = SigmaT (Covariance matrix of extra-model variables T(t))
% output.stdT    = stdT (Standard deviations of extra-model varialbes T(t))
%-----------------------------------------------------------------------------------------
% written by: Gregory Givens
%             5 May 2003
%-----------------------------------------------------------------------------------------

n = n1 + n2; nN = size(A,1); nn = size(SigmaE,1);
k = size(B,2); h = size(Omega,1); 

%check matrix dimensions%
%----------------------------------------------------------------
if nN ~= n;
   error('wrong number of predetermined and forward variables');
end;
if nn ~= n1;
   error('covariance matrix dimensions do not agree');
end;
if h ~= n2;
   error('dimensions of Omega do not agree');
end;
%----------------------------------------------------------------

%construct structural matrices and reorder%
%----------------------------------------------------------------
A_11 = A(1:n1,1:n1); A_12 = A(1:n1,n1+1:n);
A_21 = A(n1+1:n,1:n1); A_22 = A(n1+1:n,n1+1:n);
B_1 = B(1:n1,:); B_2 = B(n1+1:n,:);
Ihat = [eye(n1) zeros(n1,n2);zeros(n2,n1) Omega];
Q = [C';Ci']*W*[C Ci];
Lyy = Q(1:n,1:n);     Lyi = Q(1:n,n+1:n+k);
Liy = Q(n+1:n+k,1:n); Lii = Q(n+1:n+k,n+1:n+k);
R = A - B*inv(Lii)*Liy; U = B*inv(Lii)*B'; V = -Lyy + Lyi*inv(Lii)*Liy;
Dd = [zeros(n,n) R';Ihat U]; Gg = [V (1/beta)*Ihat';R zeros(n,n)];
Dd1 = Dd(:,1:n1); Dd2 = Dd(:,n1+1:n); Dd3 = Dd(:,n+1:n+n1); Dd4 = Dd(:,n+n1+1:2*n);
Gg1 = Gg(:,1:n1); Gg2 = Gg(:,n1+1:n); Gg3 = Gg(:,n+1:n+n1); Gg4 = Gg(:,n+n1+1:2*n);
Ddd = [Dd1 Dd4 Dd2 Dd3]; Ggg = [Gg1 Gg4 Gg2 Gg3];
%-----------------------------------------------------------------

%solve the linear system using Klein (2000)%
%-----------------------------------------------------------------
[mat1,mat2] = klein(Ddd,Ggg,n,n);
G = mat2(1:n2,1:n1); Gam = mat2(1:n2,n1+1:n);
T = mat2(n2+1:n,1:n1); Thet = mat2(n2+1:n,n1+1:n);
H = mat1(1:n1,1:n1); Psi = mat1(1:n1,n1+1:n); 
S = mat1(n1+1:n,1:n1); Sig = mat1(n1+1:n,n1+1:n);
FF = -(inv(Lii)*Liy*[eye(n1) zeros(n1,n2);G Gam] + ...
   inv(Lii)*B'*[T Thet;zeros(n2,n1) eye(n2)]*mat1);
F = FF(:,1:n1); Phi = FF(:,n1+1:n);
z = abs(eig(Ggg,Ddd)); numb = sum((z>=1));
%-----------------------------------------------------------------

%calculate loss%
%-----------------------------------------------------------------
Lambda = C*[eye(n1) zeros(n1,n2);G Gam] + [Ci*F Ci*Phi];
vec_Vtil = inv(eye(n^2)-beta*kron(mat1',mat1'))*vec(Lambda'*W*Lambda);
Vtil = reshape(vec_Vtil,n,n);
delta = (beta/(1-beta))*trace(Vtil*[SigmaE zeros(n1,n2);zeros(n2,n1) zeros(n2,n2)]);
Loss = [X10' zeros(n2,1)']*Vtil*[X10;zeros(n2,1)] + delta;
%-----------------------------------------------------------------

%calculate 2nd moments for extra-model variables%
%-----------------------------------------------------------------
Lambda_D = Theta1*[eye(n1) zeros(n1,n2);G Gam] + [Theta2*F Theta2*Phi];
vec_SigmaX = inv(eye(n^2)-kron(mat1,mat1))*vec([SigmaE zeros(n1,n2);zeros(n2,n1) zeros(n2,n2)]);
SigmaX = reshape(vec_SigmaX,n,n);
SigmaT = Lambda_D*SigmaX*Lambda_D';
stdT = sqrt(diag(SigmaT)');
%-----------------------------------------------------------------

%input results into 'output' structure%
%-----------------------------------------------------------------
output.meth = 'optimal commitment: complete information';
output.F = F; output.Phi = Phi;
output.G = G; output.Gam = Gam;
output.T = T; output.Theta = Thet;
output.H = H; output.Psi = Psi;
output.S = S; output.Sig = Sig;
output.numb = numb;
output.Loss = Loss; output.weights = diag(W)';
output.CovT = SigmaT; output.stdT = stdT;
%-----------------------------------------------------------------