l=0.65;
ksi=7.5;
chi=1;
psi =8.5;

eta =1;
phiw =0.75;
phip =0.75;
epsp =5.5;
epsw =5.5;

psip=1.5;


alpha =0.35;
beta =0.99;
delta =0.025;
zetap =0.7;
phi =2;
mu =0.21;
gamma =0.12;
theta =0.6;
Fss=0.09;
sigmass=0.25;
gss=0.2;
Rss=1.02;
psiy=0.5;
rhoa =0.9;
rhoj =0.9;
rhos =0.9;
rhog=0.9;
rhor=0.9;
rhoh=0.9;



A=1;
q=1;
J=1;
h=0.001;
vp=1;
pi=1;
pis=1;
mc=(epsw-1)/epsw;
g=0.2;


sigma=sigmass;
omegabar=exp(norminv(Fss)*sigmass-sigmass^2/2);

F = normcdf((log(omegabar)+sigma^2/2)/sigma);
Gl = normcdf((log(omegabar)+sigma^2/2)/sigma-sigma);
GAMMA = omegabar*(1-F)+Gl;
Fpl = normcdf((log(omegabar)+sigma^2/2)/sigma);
Gpl = normcdf((log(omegabar)+sigma^2/2)/sigma-sigma);
GAMMApl = omegabar*(1-Fpl)+Gpl;
dFpl = 1/(omegabar*sigma)*normpdf((log(omegabar)+sigma^2/2)/sigma);


Rd=1/(beta*(1-gamma));
Re = solve_Re(1.5, Rd, Fss, sigmass, mu);
Leverage = 1/(1-(GAMMA-mu*Gpl)*Re/Rd);

SP=Re/Rd;
Rk=Re-1+delta;
WS = (1-alpha)*(mc*alpha^alpha/Rk^alpha)^(1/(1-alpha));
W=WS;

N = solve_labor(2, WS, epsw, psi, eta, alpha, Rk, delta, mc);
C = WS*(epsw-1)/(epsw*psi*N^eta);
K_N = (Rk/(alpha*mc))^((1/(alpha-1)));
K = K_N*N;
Y = K^alpha*N^(1-alpha);
I=delta*K;
g = 0.2;
G=g*Y;

x1 = 1/(1-phip*beta)*(mc*Y/C);
x2 = 1/(1-phip*beta)*(Y/C);

H1 = psi*N^(1+eta)/(1-phiw*beta);
H2 = N/C/(1-phiw*beta);
V = q*K/Leverage;


We = V-theta*(1-GAMMA)*Re*q*K
S=V/q;
D=Re*q*K*h/S;
B=J*beta*D/(1-beta);
b=B/q;



SS = [C Rd pi B J D WS W H1 H2 N pis A K vp Rk mc x1 x2 Y Re q I Leverage V...
omegabar b SP S G g h sigma F Gl GAMMA Fpl Gpl GAMMApl dFpl];










