Trouble using integral in a nested f solve
1 view (last 30 days)
Show older comments
I'm running an fsolve of a function which involves evaluating the integral of another function. Everything is a scalar in this problem. Every line prior to the integral line runs fine and spits out a scalar answer. However when I try to evaluate the integral it says "Subscripted assignment dimension mismatch". Any help would be appreciated! Here's the code: I'm running it as a Monte Carlo but we can just set all parameter and variable values: gamma=-1; beta=1; c=7;
x1= 1;
x2= 1;
xi1=0.2;
xi2=0.21;
mc1=7;
mc2=7;
xxi=zeros(1,2);
xxi0 = [0,0];
pd0 = [0,0];
pm1=8;
pm2=8;
pm0=0;
The main function:
function F = func(xxi, gamma, beta, mc1, mc2, pm1, pm2, x1, x2, pd1, pd2, xi1, xi2)
xxi1 = xxi(1); xxi2 = xxi(2);
fun3 = @(a, b)funpd21(gamma, beta, mc1, mc2, x1, x2, a, b); fun4 = @(c, d)funpd22(gamma, beta, mc1, mc2, x1, x2, c, d);
fun1 = @(xxi1,e)(fun3(xxi1,e)-mc1).*exp(x1*beta+fun3(xxi1,e).*gamma+xxi1)./(1+exp(x1.*beta+fun3(xxi1,e).*gamma+xxi1)+exp(x2.*beta+fun4(xxi1,e).*gamma+e)).*normpdf(e,0.2,1); fun2 = @(f,xxi2)(fun4(f,xxi2)-mc2).*exp(x2*beta+fun4(f,xxi2).*gamma+xxi2)./(1+exp(x1.*beta+fun3(f,xxi2).*gamma+f)+exp(x2.*beta+fun4(f,xxi2).*gamma+xxi2)).*normpdf(f,0.2,1);
F(1) = normcdf(xxi(2),0.2,1).*(pm1-mc1).*exp(x1*beta+pm1.*gamma+xxi1)./(1+exp(x1.*beta+pm1.*gamma+xxi1))+(1-normcdf(xxi(2),0.2,1))*integral(@(g)fun1(xxi1,g),xxi2,5)/(1-normcdf(xxi2,0.2,1)); F(2) = normcdf(xxi(1),0.2,1).*(pm2-mc2).*exp(x2*beta+pm2.*gamma+xxi2)./(1+exp(x2.*beta+pm2.*gamma+xxi2))+(1-normcdf(xxi(1),0.2,1))*integral(@(h)fun2(h,xxi2),xxi1,5)/(1-normcdf(xxi1,0.2,1)); end
And the functions earlier in the nest: function pd91 = funpd21(gamma, beta, mc1, mc2, x1,x2,xi1,xi2)
options = optimset('Display','off');
f21 = @(pd)funpd(pd, gamma, beta, mc1, mc2, x1,x2,xi1,xi2); pd9 = fsolve(f21, [0,0], options);
pd91 = pd9(1);
end
function pd92 = funpd22(gamma, beta, mc1, mc2, x1,x2,xi1,xi2)
options = optimset('Display','off');
f20 = @(pd)funpd(pd, gamma, beta, mc1, mc2, x1,x2,xi1,xi2); pd9 = fsolve(f20, [0,0], options);
pd92 = pd9(2);
end
function F = funpd(pd, gamma, beta, mc1, mc2, x1, x2, xi1, xi2) pd1 = pd(1); pd2 = pd(2); F(1) = (pd1-mc1).*gamma+1./(1-exp(x1.*beta+pd1.*gamma+xi1)./(1+exp(x1.*beta+pd1.*gamma+xi1)+exp(x2.*beta+pd2.*gamma+xi2))); F(2) = (pd2-mc2).*gamma+1./(1-exp(x2.*beta+pd2.*gamma+xi2)./(1+exp(x1.*beta+pd1.*gamma+xi1)+exp(x2.*beta+pd2.*gamma+xi2))); end
0 Comments
Answers (0)
See Also
Categories
Find more on Special Functions in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!