Implicit function in a for loop

12 次查看（过去 30 天）
Peter Galbacs
Peter Galbacs 2023-1-10
评论： Torsten 2024-12-17
Hi all,
I think I need your help. I have a model of three difference equations for three variables. Two initial values are needed, then the system is supposed to work through a for loop. Unfortunately, one of the main equations is implicit, so it must be solved numerically in every round. Here is what I have now:
clear
syms alpha beta delta theta A % parameters
alpha=1 % parameter values
beta=0.95
delta=0.1
theta=0.5
A=2
k1=10 % inital value for capital
h1=0.5 % initial value for labor supply
k(1)=k1 % 1st element of capital vector
h(1)=h1 % 1st element of labor supply vector
for t=1:100
c(t)=((1-theta)/alpha)*(A^theta)*(k(t)^theta)/(h(t)^theta)*(1-h(t)) % endogenous consumption
k(t+1)=(1-delta)*k(t)+(A^(1-theta))*((h(t))^(1-theta))*(k(t)^theta)-c(t)
eqn=beta*(1-delta)+beta*theta*A^(1-theta)*(h(t+1)^(1-theta))/(k(t+1)^(1-theta))==(k(t+1)^theta)/(h(t+1)^theta)*(1-h(t+1))/(1-h(t))*(h(t)^theta)/(k(t)^theta)
h(t+1)=solve(eqn,h(t+1))
end
Running ends up in this error message: Index exceeds the number of array elements. Index must not exceed 1.
How to proceed, guys? Any idea is warmly welcome. And thanks a lot.

请先登录，再回答此问题。

采纳的回答

Walter Roberson
Walter Roberson 2023-1-10
Ran in:
syms alpha beta delta theta A % parameters
syms h_tplus1
Q = @(v) sym(v);
alpha = Q(1); % parameter values
beta = Q(0.95);
delta = Q(0.1);
theta = Q(0.5);
A = Q(2);
k1 = Q(10); % inital value for capital
h1 = Q(0.5); % initial value for labor supply
k(1) = double(k1); % 1st element of capital vector
h(1) = double(h1); % 1st element of labor supply vector
for t=1:100
c(t) = double(((1-theta)/alpha)*(A^theta)*(k(t)^theta)/(h(t)^theta)*(1-h(t))); % endogenous consumption
k(t+1) = double((1-delta)*k(t)+(A^(1-theta))*((h(t))^(1-theta))*(k(t)^theta)-c(t));
eqn = beta*(1-delta)+beta*theta*A^(1-theta)*(h_tplus1^(1-theta))/(k(t+1)^(1-theta))==(k(t+1)^theta)/(h_tplus1^theta)*(1-h_tplus1)/(1-h(t))*(h(t)^theta)/(k(t)^theta);
thissol = double(vpasolve(eqn, h_tplus1, [-inf inf])); %constrain to real
h(t+1) = thissol;
end
plot(h)
  5 个评论
显示 3更早的评论
Peter Galbacs
Peter Galbacs 2023-1-11
Thanks a lot, Walter. An elegant solution. I will study your code to fully absorb the technical details. You helped a lot!
Mendes Pereira
Mendes Pereira 2024-12-17
The answer above does not represent the dynamics of h(t) in the original model. See further details below.

请先登录，再进行评论。

更多回答（1 个）

Mendes Pereira
Mendes Pereira 2024-12-17
The numerical simulation above is the simulation of something dynamic, but it is not the simulation of any element of the model described in the initial problem. If I am not making some mistake, the three equations in the original model are as follows:
Using the same parameter values as above, this system has a steady state given by: . However, the figure above shows that h(t) converges to 1 after a large number of iterations, while its actual steady state value is 0.426471.
I have a similar problem, but I do not know how to run a for loop with an implicit equation as part of the model. Some help would be very much appreciated.
  1 个评论
Torsten
Torsten 2024-12-17
If h(t) converges to 1 like in the code from above, the behaviour of the term (1-h(t+1))/(1-h(t)) is unpredictable.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Programming 的更多信息

标签

产品


版本

R2022b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by