Matrix input in fzero solver
3 次查看(过去 30 天)
显示 更早的评论
Hello,
I am stuck with my code. I am trying to vectorize loops for performance issues. However, I have a problem when the vectorization involves fzero. Here I post an example of what I want to get. First code line shows the script with the loops, and the second one the function file.
MAIN SCRIPT
clear all
clc
grida=10;
gridz=60;
w=1.296;
mu1=0.85;
theta=0.1;
s = rng;
bbar_AQ = zeros(grida,gridz);
pistar_n_f = zeros(1,gridz);
k_star_f = zeros(1,gridz);
auxiliar = zeros(grida,gridz);
debt = zeros(grida,gridz);
A=randn(10,60);
aa=0.001:(w-0.002)/(grida-1):w-0.001;
z = linspace(1, 12, gridz);
opt=optimset('display', 'off');
for j=1:grida
for t=1:gridz
k_star_f(t) =(z(t))^(1/(1-mu1-theta));
pistar_n_f(t) = (1-mu1)*(mu1/w)^(mu1/(1-mu1))*(z(t))^(1/(1-mu1))*(k_star_f(t)^(theta/(1-mu1)));
bbar_AQ(j,t) = (pistar_n_f(t) + aa(j))/0.75;
auxiliar(j,t) =k_star_f(t);
if bbar_AQ(j,t)>=auxiliar(j,t)
debt(j,t)=A(j,t);
else
if j==1
bbar_0 = 0.1;
else
bbar_0 = 0.11;
end
debt(j,t)=fzero(@(x)b_max(x,z(t),aa(j),theta,mu1,w),bbar_0,opt);
end
end
end
FUNCTION FILE
function f=b_max(x,z,a,theta,mu,w)
f=(((1-mu)*(mu/w)^(mu/(1-mu))*(z*(x+a)^theta)^(1/(1-mu))))-x;
The above code runs fine. Howeve, I am trying to vectorize it but when I arrive to the fzero line of code, I am stuck about how to proceed. This is what I did:
clear all
clc
grida = 10;
gridz = 60;
w = 1.296;
mu1 = 0.85;
theta = 0.1;
aa = 0.001:(w-0.002)/(grida-1):w-0.001;
z = linspace(1, 12, gridz);
rng(s)
A=randn(10,60);
% Compute k_star_f and pistar_n_f using vectorized operations
k_star_f = z.^(1./(1-mu1-theta));
pistar_n_f = (1-mu1)*(mu1/w).^(mu1./(1-mu1)).*(z.^(1./(1-mu1))).*(k_star_f.^(theta./(1-mu1)));
% Compute bbar_AQ using vectorized operations
expanded_aa = repmat(aa(:), 1, gridz);
bbar_AQ = (pistar_n_f + expanded_aa) / 0.75;
% Compute auxiliar using vectorized operations
auxiliar = k_star_f;
auxiliar=repmat(auxiliar,grida,1);
% Compute debt using vectorized operations based on the condition
debt2 = A2 .* (bbar_AQ2 >= auxiliar2) + fzero(@(x)b_max(x,z(t),aa(j),theta,mu1,w),bbar_0,opt).* (bbar_AQ2 < auxiliar2);
Any help would be really appreciated.
Thank in advanced!
8 个评论
Torsten
2023-9-8
编辑:Torsten
2023-9-8
So, it will not possible to use parfor loop?
No. The loop as written has to be processed sequentially because you made the initial values of the i-th call to "fzero" depend on the result of the (i-1)th call.
You might try to set bbar_0 outside and only leave the line
bbar_constrained(j,t)=fzero(@(x)b_max_2_OnlyPi(x,z(t),eta,aa(j),theta,mu1,w),bbar_0,opt);
within the double loop.
This loop should then be possible to be processed in parallel.
回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Loops and Conditional Statements 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!