function foc=rintang_resid(X,param);
%RIntang_resid  Residuals for the 3-sector model with intangible capital.
%

%               Ellen R. McGrattan, 3-18-05
%               Revised, ERM, 8-10-05

%______________________________________________________________________________
%
% PARAMETERS

beta       = param(1);
deltI      = param(2);
deltT      = param(3);
eta        = param(4);
gamma      = param(5);
phi        = param(6);
psi        = param(7);
sigma      = param(8);
theta      = param(9);
kT0        = param(10);
kI0        = param(11);
T          = param(12);
inputs     = zeros(T,15);
inputs(:)  = param(13:12+15*T);
A          = inputs(:,1);
lwedge     = inputs(:,2);
xwedge     = inputs(:,3);
g          = inputs(:,4);
tauc       = inputs(:,5);
gtauh      = inputs(:,6);
gtaux      = inputs(:,7);
tauk       = inputs(:,8);
taup       = inputs(:,9);
taud       = inputs(:,10);
taug       = inputs(:,11);
chi        = inputs(:,12);
yn         = inputs(:,13);
xn         = inputs(:,14);
hn         = inputs(:,15);
bderiv     = param(13+15*T);
uderiv     = param(14+15*T);

bhat       = beta*(1+gamma)^(-sigma);
lshare     = 1-theta-phi;
grate      = (1+eta)*(1+gamma);
tauh       = 1-(1-gtauh).*lwedge;
taux       = (1+gtaux)./xwedge-1;

%______________________________________________________________________________
%
% VARIABLES APPEARING IN RESIDUALS

x          = zeros(T,3);
x(:)       = X;
kTn        = x(:,1);
kIn        = x(:,2);
hb         = x(:,3);
kT         = [kT0;kTn(1:T-1)];
kI         = [kI0;kIn(1:T-1)];
xT         = grate*kTn-(1-deltT)*kT;
xI         = grate*kIn-(1-deltI)*kI;
h          = hb+hn;
tyb        = A.*kT.^theta.*kI.^phi.*hb.^lshare;
yb         = tyb-xI;
c          = yb+yn-xT-xn-g;
l          = 1-h;
uc         = c.^(-sigma).*l.^(psi*(1-sigma));
rT         = theta*tyb./kT;
rI         = phi*tyb./kI;
xtau       = (1+taux).*(1-taud)./(1+tauc);
dtau       = (1-taud).*taup./(1+tauc);
ptau       = (1-taup).*(1-taud)./(1+tauc);
gtau       = (taug+(1-chi).*(1-gtauh)+chi.*(1-taup).*(1-taud))./(1+tauc);
t          = 1:T-1;
t1         = 2:T;

foc        = [[xtau(t).*uc(t)- ...
                   bhat*uc(t1).*(ptau(t1).*(rT(t1)-tauk(t1))+ ...
                            (1-deltT)*xtau(t1)+deltT*dtau(t1));
                kTn(T)-bderiv*kTn(T-1)];
              [gtau(t).*uc(t)- ...
                   bhat*uc(t1).*((1-deltI)*gtau(t1)+ptau(t1).*rI(t1));
                kIn(T)-uderiv*kIn(T-1)];
              psi*(1+tauc).*c.*hb./tyb-lshare*(1-tauh).*l]; 
