% Programa Final

clear   %limpia memoria
close all %cierra figuras
clc %limpia la pantalla de comandos

% leer datos
Eco_n1=load("Eco_n1.txt");

% calibrar el modelo
[ns, y, a,  beta, alpha, theta, LT, KT ] = Calibrar(Eco_n1);


% prueba de la calibracion
% w=1
% semilla r=p(1,ns)=y(1,ns)=1

xx=zeros(1,ns*2+1);

xx(1)=log(1); %semilla de r
for i =1:ns
    xx(1+i)=log(1); %semilla de p
    xx(1+ns+i)=log(y(i)); %semilla de y(dato)
end

ff0=ff_eq(xx,ns,beta,alpha,a,LT,KT,theta);

if sum(abs(ff0)) > 10^(-6)
    disp('error calibracion')
end

%%% control del newton
niter=1000;
EPSILON=10^(-11);
eps=10^(-5);

xx_sol=Resuelve_con_Newton(niter,EPSILON,eps,...
    xx,ns,beta,alpha,a,LT,KT,theta);

%% Evalua el equilibrio
[r, p, y, k, l, pinv, Ycon, c, x, inv]=...
    evalua_eq(xx_sol,ns,beta,alpha,a,LT,KT,theta);



LT2=LT+1;
%KT2=KT-1;
KT2=KT
xx_sol2=Resuelve_con_Newton(niter,EPSILON,eps,...
    xx,ns,beta,alpha,a,LT2,KT2,theta);

[r2, p2, y2, k2, l2, pinv2, Ycon2, c2, x2, inv2]=...
    evalua_eq(xx_sol2,ns,beta,alpha,a,LT2,KT,theta);

ipc=1
ipc2=(p2*c')./(p*c')

w2=1/ipc2
r2=r2/ipc2
p2=p2/ipc2

p(1)./p(2)   % unidades del bien 2 por unidad del bien 1 pais base

p2(1)./p2(2) % unidades del bien 2 por unidad del bien 1 pais 2