clear
cd D:\4.00\userexample
pkg load symbolic
syms Omegabar Qbar
Ybar=1
pbeta=.99
epsilon=-1 
sigma=.03 
theta=.3 
k=13
delta=.025
rho=1.1
%  
Ybar=1;
Lambdabar=pbeta
Rlbar=1+1/Qbar
Rbar=1/pbeta
Rebar=1/pbeta
munbar=Lambdabar*Omegabar*Rbar
musbar=Lambdabar*Omegabar*(Rlbar-Rbar)
muebar=Lambdabar*Omegabar*(1-pbeta)
epersbar=musbar/muebar
xbar=-epersbar+(epersbar^(2)+(2/k)*(1-epsilon*epersbar))^(1/2)
Thetabar=theta*(1+epsilon*xbar+(k/2)*xbar^(2))
phibar=Omegabar/(Thetabar-musbar-xbar*muebar)
Ebar=xbar*Qbar
Nbar=Qbar/phibar
Dbar=Qbar-Nbar-Ebar
Nnbar=delta*(Qbar+Ybar)
Pibar=(sigma)*(Qbar+Ybar-Rebar*Ebar-Rbar*Dbar)+Nnbar
Cbar=Pibar+(Rbar-1)*Dbar+(Rebar-1)*Ebar

eq1=Omegabar==(sigma+(1-sigma)*(munbar+phibar*musbar+phibar*xbar*muebar))
eq2=Nbar==(1-sigma)*(Qbar+Ybar-Rebar*Ebar-Rbar*Dbar)+Nnbar
[c1 c2]=solve(eq1,eq2,Omegabar,Qbar)