parameter estimation with lsqcurvefit and boundary value ode (BVP4c)

7 views (last 30 days)
I am trying to find out the value of parameter theta for a boundary value problem. Attaching the code and expermental data with which I want to fit my model. But it is showing error as function value and YDATA are not equal. Is it possible to uselsqcurvefit for bvp?
function sol1=burial(theta,x)
options = bvpset('stats','on') ;
solinit =bvpinit (x, [1 0])
sol = bvp4c (@main, @bvpbcs, solinit,options);
function dydx= main(x,y)
phi= 0.708+((0.744-0.708)*exp(-0.7*x))
D_phi= -0.0252*exp(-0.7*x)
D0= (((0.438*4)+9.60)*10^(-6))/3.17098e-8;
Ds= D0./(1-(2*log(phi)));
D_Ds=-(396770375294112771.*exp(-(7.*x)/10))./(21990232555520000.*(2.*log((9.*exp(-(7.*x)./10))./250 + 177/250) - 1).^2*((9.*exp(-(7.*x)/10))./250 + 177/250));
wf=0.03;
w=((wf*0.704)-(theta*0.744))./phi;
D_w=-(63*exp(-(7*x)/theta)*((93*theta)/125 - 66/3125))/(2500*((9*exp(-(7*x)/theta))/250 + 177/250)^2);
dydx = [ y(2);
((w*y(2))+(y(1)*D_w)-(Ds*y(2)*D_phi)-(phi*D_Ds*y(2)))/(Ds*phi)];
end
function res = bvpbcs (ya , yb)
% provides the boundary conditions at the end points a and b
res = [ ya(1)-551
yb(1)-381];
end
sol1=sol.y(1,:);
end
----------------------
x=[0
0.974025974
3.051948052
5.064935065
7.012987013
8.896103896
10.97402597
12.85714286
14.87012987];
Cl_val=[550.9202454
550.9202454
528.8343558
494.4785276
462.5766871
438.0368098
415.9509202
386.5030675
381.595092];
theta0=1;
[theta,Rsdnrm,Rsd,ExFlg,OptmInfo,Lmda,Jmat]=lsqcurvefit(@burial,theta0,x,Cl_val);
parameter= theta;
Error using lsqcurvefit (line 263)
Function value and YDATA sizes are not equal.

Accepted Answer

Star Strider
Star Strider on 13 Aug 2020
The lsqcurvefit function needs to have the output of the objective function the same size (with respect to row and column dimensions) as the data it is to be compared to.
If you replace the current ‘sol1’ assignment in ‘burial’ with these three assignments:
solx = sol.x;
soly=sol.y(1,:);
sol1 = interp1(solx(:), soly(:), x);
the output will be the required size and the parameter estimation will work, producing:
parameter =
2.9515
Note that lsqcurvefit also issues a ‘Local minimum possible.’ notice. It could be worthwhile to try different initial values for ‘theta0’ to see if better estimates are possible, or perhaps using ga or another Global Optimization Toolbos function to see if a better estimate is possible.
I was almost finished writing an earlier version of this when Windows 10 crashed (again), after yet another downgrade. It would really be great if Micro$oft put out an OS as stable and as reliable as Windows 7, and that actually allowed me to configure my ccomputer the way I want and not the way Micro$oft wants me to. I am looking at a good learning reference for Linux.

More Answers (0)

Community Treasure Hunt

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

Start Hunting!