%   Esse programa resolve o sistema de equações não-lineares descrito
%   no m-file mod_fiscal_04_ss.

% Make a starting guess at the solution
x0 = [1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1; 1];

% Option to display output
options=optimset('Display','iter','MaxFunEvals',200000,'TolFun',1e-8);  
%options=optimset('Display','iter');  

% Call solver
[x,fval] = fsolve(@mod_fiscal_04_ss,x0,options);  

%   Dicionário:
%   x(1)=wc, x(2)=lc, x(3)=cc; 
%   x(4)=c, x(5)=cl;
%   x(6)=rk; 
%   x(7)=ynt, x(8)=y, x(9)=ptilnt;
%   x(10)=yt, x(11)=ptilt;
%   x(12)=ktilt, x(13)=yt_star, x(14)=w, x(15)=ptx;
%   x(16)=ktilnt;
%   x(17)=lt, x(18)=lnt;
%   x(19)=qt, x(20)=qnt;
%   x(21)=il;
%   x(22)=er;
%   x(23)=wl, x(24)=ll, x(25)=l;
%   x(26)=q;
%   x(27)=b_star, x(28)=sp, x(29)=b;
%   x(30) = pib
%   x(31) = def
%   x(32) = llt
%   x(33) = lg
%   x(34) = gt
%   x(35) = wg
%   x(36) = CCT  % ADICIONADO
%   x(37) = WMT  % ADICIONADO

%   Número de variáveis do sistema:
N = 37;

%   Valores assumidos por alguns parâmetros do modelo:
par_pop = nr_i;
par_pop2 = nr2_i;
beta = beta_i;
h = h_i;
sigma_c = sigma_c_i;
sigma_l = sigma_l_i;
eta_t_k = eta_t_k_i;
eta_t_l = eta_t_l_i;
eta_nt_k = eta_nt_k_i;
eta_nt_l = eta_nt_l_i;
delta = delta_i;
delta_til = delta_til_i;
gamma = gamma_i;
eta = eta_i;
eta_t_g = eta_t_g_i;
eta_nt_g = eta_nt_g_i;
mu_t = mu_t_i;
mu_nt = mu_nt_i;
tao_star = tao_star_i;
mu_w = mu_w_i;
phi = phi_i;
theta = theta_i;
a = a_i;

big = big_i;
bag = bag_i;

%   Valores atribuídos às variáveis exógenas no steady-state:

%   Setor externo
y_star = y_star_i;
pt_star = pt_star_i;
taoc_star = taoc_star_i;
r_star = r_star_i;

%   Variáveis fiscais
g = g_ss_i;
ig = ig_ss_i;
tl = t_ss_i;
wm = wm_ss_i;
taoc = taoc_ss_i;
taow = taow_ss_i;
taok = taok_ss_i;
taox = taox_ss_i;

uwg_ss = uwg_ss_i;
lg_llt_ss = lg_llt_ss_i;

uwt_ss = uwt_ss_i;

for jjj = 1:N
    if (imag(x(jjj)) > 0.0001)
        disp(['A raiz ...', num2str(jjj) ,'... do sistema é imaginária!'])
    else
        x(jjj) = real(x(jjj));
    end;
end;

disp(['Parâmetros necessários para a calibragem:'])
disp(' ')

kl_ss_i = x(21)/delta;
kg_ss_i = ig/delta_til;
d_y_ss_i = x(29)/x(8);
b_star_y_ss_i = x(27)/x(8);
sp_y = x(28)/x(8);

wc_ss_i = x(1);
lc_ss_i = x(2);
cc_ss_i = x(3);
c_ss_i = x(4);
cl_ss_i = x(5);
y_ss_i=x(8);
p_t_til_ss_i = x(11);
p_nt_til_ss_i = x(9);
wl_ss_i = x(23);
w_ss_i = x(14);
rer_ss_i = x(22);
l_ss_i = x(25);
pib_ss_i = x(31);
def_ss_i = x(30);
llt_ss_i = x(32);
lg_ss_i = x(33);
gt_ss_i = x(34);
wg_ss_i = x(35);

cct_ss_i = x(36);
wmt_ss_i = x(37);

aig_i = 1+big_i*exp(-x(33));
uig_lg_ss_i = big_i*exp(-x(33))*x(33)/(1+big_i*exp(-x(33)));
aag_i = 1+bag_i*exp(-x(33));
uag_lg_ss_i = bag_i*exp(-x(33))*x(33)/(1+bag_i*exp(-x(33)));

ll_llt_ss_i = x(24)/x(32);
g_gt_ss_i = g/x(34);
wglg_gt_ss_i = x(33)*x(35)/x(34);

y_ywglg_ss_i = x(8) / (x(8)+x(33)*x(35));
wglg_ywglg_ss_i = x(33)*x(35)/(x(8)+x(33)*x(35));
wglg_y_ss_i = x(33)*x(35)/x(8);
wlt_y_ss_i = (x(14)*x(25)+x(33)*x(35))/x(8);
wglg_ss_i = x(33)*x(35);

disp(' ')
disp(['ll_llt = ' num2str(ll_llt_ss_i)])
disp(['lg_llt = ' num2str(lg_llt_ss_i)])
disp(['g_gt = ' num2str(g_gt_ss_i)])
disp(['wglg_gt = ' num2str(wglg_gt_ss_i)])
disp(['y_ywglg = ' num2str(y_ywglg_ss_i)])
disp(['wglg_ywglg = ' num2str(wglg_ywglg_ss_i)])
disp(['wglg_y = ' num2str(wglg_y_ss_i)])
disp(['wlt_y = ' num2str(wlt_y_ss_i)])
disp(['wglg = ' num2str(wglg_ss_i)])
disp(['uig_lg = ' num2str(uig_lg_ss_i)])
disp(['aig = ' num2str(aig_i)])
disp(['uag_lg = ' num2str(uag_lg_ss_i)])
disp(' ')


disp(['Razão entre as exportações de bens tradables e a produção total de bens tradables...'])
y_y_ss_i = x(13)/(x(10)+x(13));
disp(['y_y = ' num2str(y_y_ss_i)])
disp(' ')
disp(['Razão entre o consumo dos indivíduos do tipo ricardiano e o consumo agregado...'])
cl_c_ss_i = x(5)/x(4);
disp(['cl_c = ' num2str(cl_c_ss_i)])
disp(' ')
disp(['Razão entre o consumo dos indivíduos do tipo não-ricardiano que trabalha e o consumo agregado...'])
cc_c_ss_i = x(3)/x(4);
disp(['cc_c = ' num2str(cc_c_ss_i)])
disp(' ')
disp(['Razão entre o consumo dos indivíduos do tipo não-ricardiano que NAO trabalha e o consumo agregado...'])
cct_c_ss_i = x(36)/x(4);
disp(['cc_c = ' num2str(cct_c_ss_i)])
disp(' ')
disp(['Razão entre o consumo e o produto...'])
c_y_ss_i = x(4)/x(8);
disp(['c_y = ' num2str(c_y_ss_i)])
disp(' ')
disp(['Razão entre o consumo dos indivíduos do tipo ricardiano e o produto...'])
cl_y_ss_i = cl_c_ss_i*c_y_ss_i;
disp(['cl_c = ' num2str(cl_y_ss_i)])
disp(' ')
disp(['Razão entre o consumo dos indivíduos do tipo não-ricardiano que trabalha e o produto...'])
cc_y_ss_i = cc_c_ss_i*c_y_ss_i;
disp(['cc_y = ' num2str(cc_y_ss_i)])
disp(' ')
disp(['Razão entre o consumo dos indivíduos do tipo não-ricardiano que NAO trabalha e o produto...'])
cct_y_ss_i = cct_c_ss_i*c_y_ss_i;
disp(['cc_y = ' num2str(cct_y_ss_i)])
disp(' ')
disp(['Razão entre a tributação lump-sum e o produto...'])
t_y = tl/x(3);
disp(['t_y = ' num2str(t_y)])
disp(' ')
disp(['Razão entre as transferências aos não ricardianos que trabalham e o produto...'])
wm_y_ss_i = (wm_ss_i*par_pop)/x(8);
disp(['wmtot_y = ' num2str(wm_y_ss_i)])
disp(' ')
disp(['Razão entre as transferências aos não ricardianos que nao trabalham e o produto...'])
wmt_y_ss_i = (wmt_ss_i*par_pop2)/x(8);
disp(['wmtot_y = ' num2str(wmt_y_ss_i)])
disp(' ')
disp(['Razão entre a renda do trabalho no setor privado e o produto...'])
wl_y_ss_i = w_ss_i*l_ss_i/y_ss_i;
disp(['wl_y = ' num2str(wl_y_ss_i)])
disp(' ')
disp(['Razão entre a quantidade de insumos importados e o produto...'])
q_y_ss_i = x(26)/x(8);
disp(['q_y = ' num2str(q_y_ss_i)])
disp(' ')
disp(['Razão entre a balanca comercial e o produto...'])
bc_y_ss_i = (1/(1+taoc)*x(11)*x(13)-(1+taox)/(1+taoc_star)*x(22)*pt_star*x(26))/x(8);
disp(['bc_y = ' num2str(bc_y_ss_i)])
disp(' ')

bc_ss_i = (1/(1+taoc)*x(11)*x(13)-(1+taox)/(1+taoc_star)*x(22)*pt_star*x(26));
rendak_ss_i = x(8)/(1+taoc)-x(14)*x(25)+bc_ss_i;
gg_ss_i = (g_ss_i+ig_ss_i)/(1+taoc_ss_i)+nr_i*wm_ss_i+nr2_i*wmt_ss_i+wg_ss_i*lg_ss_i;

ctb_y_ss_i = ( ( taow*(x(14)*x(25)+x(33)*x(35)) + taok*rendak_ss_i) + ( (taoc/(1+taoc))*(x(4)+g_ss_i)+taox*x(15)*x(26)) + (t_ss_i) ) / (x(30)*x(31));
ctl_y_ss_i = ctb_y_ss_i - (wm_ss_i*par_pop+wmt_ss_i*par_pop2)/(x(30)*x(31));

disp(['Parâmetros da aproximação de 1a ordem da condição BP = 0...'])
yt_star_y_ss_i = x(13)/x(8);
p_t_til_star_ss_i = pt_star/(1+taoc_star);
disp(['ytstar/y = ' num2str(yt_star_y_ss_i)])
disp(' ')

disp(['Parâmetros da aproximação de 1a ordem da condição usos = recursos...'])
disp(' ')
g_y_ss_i = g_ss_i/x(8);
disp(['Razão consumo do governo em B+S/produto = ' num2str(g_y_ss_i)])
disp(' ')
kl_y_ss_i = kl_ss_i/x(8);
disp(['Razão capital mantido pelos agentes ricardianos/produto = ' num2str(kl_y_ss_i)])
disp(' ')
kg_y_ss_i = (1/delta_til)*(ig_ss_i/x(8));
disp(['Razão capital mantido pelo governo/produto = ' num2str(kg_y_ss_i)])
disp(' ')
ig_y_ss_i = ig_ss_i/x(8);
disp(['Razão investimento público/produto = ' num2str(ig_y_ss_i)])
disp(' ')
il_y_ss_i = delta*kl_y_ss_i;
disp(['Razão investimento privado/produto = ' num2str(il_y_ss_i)])
disp(' ')

disp(['Parâmetros relativos à utilização de insumos em ss...'])
disp(' ')
kt_k_ss_i = x(12)/(x(12)+x(16));
knt_k_ss_i = x(16)/(x(12)+x(16));
disp(['Razão capital tradables/capital privado = ' num2str(kt_k_ss_i)])
disp(['Razão capital non-tradables/capital privado = ' num2str(knt_k_ss_i)])
disp(' ')

lt_l_ss_i = x(17)/x(25);
lnt_l_ss_i = x(18)/x(25);
disp(['Razão trabalho tradables/trabalho = ' num2str(lt_l_ss_i)])
disp(['Razão trabalho non-tradables/trabalho = ' num2str(lnt_l_ss_i)])
disp(' ')

qt_q_ss_i = x(19)/x(26);
qnt_q_ss_i = x(20)/x(26);
disp(['Razão importados tradables/importados = ' num2str(qt_q_ss_i)])
disp(['Razão importados non-tradables/importados = ' num2str(qnt_q_ss_i)])
disp(' ')

disp(['Parâmetros relativos ao PIB/deflator do PIB...'])
disp(' ')
y_pib_i = (x(8)+x(33)*x(35))/(x(30)*x(31));
disp(['Absorcao domestica/PIB = ' num2str(y_pib_i)])
exp_pib_i = (1/(1+taoc))*x(11)*x(13)/(x(30)*x(31));
disp(['Export/PIB = ' num2str(exp_pib_i)])
imp_pib_i = x(15)*x(26)/(x(30)*x(31));
disp(['Import/PIB = ' num2str(imp_pib_i)])
disp(' ')

%   Informa alguns resultados.
disp(' ')
disp(['Alguns resultados do steady-state:'])
disp(' ')
disp(['Produto/PIB = ' num2str(x(8)/((x(30)*x(31))))])
disp(['Consumo total/PIB = ' num2str((x(4)+x(34))/(x(30)*x(31)))])
disp(['  Consumo privado total/PIB = ' num2str(x(4)/(x(30)*x(31)))])
disp(['    Consumo dos indivíduos ricardianos/PIB = ' num2str(x(5)/(x(30)*x(31)))])
disp(['    Consumo dos indivíduos nao ricardianos que trabalham/PIB = ' num2str(x(3)/(x(30)*x(31)))])
disp(['    Consumo dos indivíduos nao ricardianos que nao trabalham/PIB = ' num2str(x(36)/(x(30)*x(31)))])
disp(['  Consumo total do governo /PIB = ' num2str(x(34)/(x(30)*x(31)))])
disp(['    Consumo do governo em B+S/PIB = ' num2str(g_ss_i/(x(30)*x(31)))])
disp(['    Consumo do governo em salarios/PIB = ' num2str(x(33)*x(35)/(x(30)*x(31)))])
disp(['Investimento total/PIB = ' num2str((x(21)+ig_ss_i)/(x(30)*x(31)))])
disp(['  Investimento privado/PIB = ' num2str(x(21)/(x(30)*x(31)))])
disp(['  Investimento do governo/PIB = ' num2str(ig_ss_i/(x(30)*x(31)))])
disp(['Absorcao domestica/PIB = ' num2str(y_pib_i)])
disp(['Export/PIB = ' num2str(exp_pib_i)])
disp(['Import/PIB = ' num2str(imp_pib_i)])
disp(['Renda do trabalho/PIB = ' num2str((x(14)*x(25)+x(33)*x(35))/(x(30)*x(31)))])
disp(['  Renda do trabalho no s.privado/PIB = ' num2str((x(14)*x(25))/(x(30)*x(31)))])
disp(['  Renda do trabalho no s.publico/PIB = ' num2str((x(33)*x(35))/(x(30)*x(31)))])
disp(['  Renda do trabalho dos ricardianos    /PIB = ' num2str((x(23)*x(24)+x(33)*x(35))/(x(30)*x(31)))])
disp(['  Renda do trabalho dos nao ricardianos/PIB = ' num2str((x(1)*x(2))/(x(30)*x(31)))])
disp(['Renda do capital/PIB = ' num2str( (rendak_ss_i) /((x(30)*x(31))))])
disp(['Gastos totais do governo/PIB = ' num2str( (gg_ss_i) /((x(30)*x(31))))])
disp(['    Consumo do governo em B+S/PIB = ' num2str(g_ss_i/(x(30)*x(31)))])
disp(['    Salarios no s.publico/PIB = ' num2str((x(33)*x(35))/(x(30)*x(31)))])
disp(['    Transferencias/PIB = ' num2str((wm_ss_i*par_pop+wmt_ss_i*par_pop2)/(x(30)*x(31)))])
disp(['Carga Tributaria Bruta/PIB = ' num2str(ctb_y_ss_i)])
disp(['Carga Tributaria Liquida/PIB = ' num2str(ctl_y_ss_i)])
disp(['   Tributação sobre renda do trabalho/PIB = ' num2str(taow*(x(14)*x(25)+x(33)*x(35))/(x(30)*x(31)))])
disp(['   Tributação sobre renda do capital/PIB = ' num2str( (taok*rendak_ss_i) /((x(30)*x(31))))])
disp(['   Impostos sobre produtos/PIB = ' num2str( ( (taoc/(1+taoc))*(x(4)+g_ss_i)+taox*x(22)*pt_star*x(26)) /((x(30)*x(31))))])
disp(['   Tributos lump-sum/PIB = ' num2str( (t_ss_i) /((x(30)*x(31))))])
disp(['Trabalho s.privado/trabalho total = ' num2str(x(25)/(x(32)+x(2)))])
disp(['Trabalho s.publico/trabalho total = ' num2str(x(33)/(x(32)+x(2)))])
disp(['Trabalho ricardianos s.privado/trabalho total ricardianos = ' num2str(x(24)/x(32))])
disp(['Trabalho ricardianos s.publico/trabalho total ricardianos = ' num2str(x(33)/x(32))])
disp(['Trabalho ricardianos s.privado    /trabalho total s.privado = ' num2str(x(24)/x(25))])
disp(['Trabalho nao ricardianos s.privado/trabalho total s.privado = ' num2str(x(2)/x(25))])
disp(['Produção de bens tradables/nao-tradables = ' num2str(x(11)*x(10)/(x(9)*x(7)))])
disp(['Salário ricardianos/salarios nao ricardianos = ' num2str(x(23)/x(1))])
disp(['Estoque de capital privado/PIB = ' num2str(x(21)/delta/(4*x(30)*x(31)))])
disp(['Estoque de capital publico/PIB = ' num2str(ig/delta_til/(4*x(30)*x(31)))])
disp(['Dívida interna/PIB em ss = ' num2str(x(29)/(4*x(30)*x(31)))])
disp(['Dívida externa/PIB em ss = ' num2str(x(27)/(4*x(30)*x(31)))])
disp(['Superávit primário/PIB em ss = ' num2str(x(28)/(x(30)*x(31)))])
disp(' ')

