-
Notifications
You must be signed in to change notification settings - Fork 6
/
SolveStationaryProblem.m
53 lines (46 loc) · 1.42 KB
/
SolveStationaryProblem.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
function [res] = SolveStationaryProblem(xState,Para,xInit)
s=1;
u2bdiff=xState(1);
RR=xState(2);
options=optimset('Display','off');
[x,fvec,exitflag]=fsolve(@(x) SteadyStateCharacterization(x,u2bdiff,RR,Para,1) ,xInit,options);
Para.theta=[Para.theta_1 Para.theta_2];
Para.alpha=[Para.alpha_1 Para.alpha_2];
Par=Para;
u2btild=u2bdiff;
R=RR;
s_=s;
n1=Para.n1;
n2=Para.n2;
ctol=Para.ctol;
%% GET THE Policy Rules
psi= Par.psi;
beta = Par.beta;
P = Par.P;
theta_1 = Par.theta(1);
theta_2 = Par.theta(2);
g = Par.g;
alpha = Par.alpha;
sigma = 1;
c1_1=x(1);
c1_2=x(2);
c2_1=x(3);
%compute components from unconstrained guess
[c2_2 grad_c2_2] = computeC2_2(c1_1,c1_2,c2_1,R,s_,P,sigma);
[l1 l1grad l2 l2grad] = computeL(c1_1,c1_2,c2_1,c2_2,grad_c2_2,...
theta_1,theta_2,g,n1,n2);
[btildprime grad_btildprime] = computeBtildeprime(c1_1,c1_2,c2_1,c2_2,grad_c2_2,l1,l2,l1grad,l2grad,...
u2btild,s_,psi,beta,P);
% x' - definition
u2btildprime=psi*[c2_1^(-1) c2_2^(-1)].*btildprime;
% State next period
X(1,:) = [psi*c2_1^(-1)*btildprime(1),c2_1^(-1)/c1_1^(-1)];%state next period
X(2,:) = [psi*c2_2^(-1)*btildprime(2),c2_2^(-1)/c1_2^(-1)];%state next period
Vobj = P(s_,1)*(alpha(1)*uBGP(c1_1,l1(1),psi)+alpha(2)*uBGP(c2_1,l2(1),psi));
Vobj = Vobj + P(s_,2)*(alpha(1)*uBGP(c1_2,l1(2),psi)+alpha(2)*uBGP(c2_2,l2(2),psi));
%Vobj=-Vobj;
PolicyRules=[c1_1 c1_2 c2_1];
res.PolicyRules=PolicyRules;
res.Value=Vobj/(1-Para.beta);
res.exitflag=exitflag;
end