using lsqnonlin () for solution of an nonlinear equation
Show older comments
Please consider the following code:
function diff = myfun(D,Y)
t = 0:30:600;
l = length(t);
a = 0.5*0.75*10^(-4);
for i = 1:l
Yp(i)=0;
for n = 1:10
y(i,n) = exp(-D*t(i)*(pi^2)*(n^2)/(a^2))/n^2;
Yp(i) = Yp(i) + y(i,n);
end
end
diff = 1- (6*Yp/pi^2)-Y;
%function
Y = [0.0566 0.4432 0.5539 0.6783 0.7303 0.3569 0.4001 0.4278 0.4499 0.4720 ...
0.4500 0.5157 0.5237 0.5492 0.5590 0.5799 0.5890 0.6000 0.6150 0.6300...
0.6450];
D0 = 2.0*10^(-12);
options=optimset('LargeScale','on','Display','iter','MaxFunEvals',1e20,'TolFun',2e-50,'TolX',2e-50,'LevenbergMarquardt','on');
[D,resnorm,residual,exitflag,output,lambda]=lsqnonlin(@myfun,D0,[],[],options);
I am trying to get D value using the lsqnonlin() but I am getting the message:
Error using ==> feval
Undefined function or method 'diffusion' for input arguments of type 'double'.
Error in ==> lsqnonlin at 200
initVals.F = feval(funfcn{3},xCurrent,varargin{:});
Caused by:
Failure in initial user-supplied objective function evaluation. LSQNONLIN cannot continue.
Please help me with the code why I am getting this error message. What could be the better possible way to obtain the solution for D value.
Accepted Answer
I suspect that the code you are running is not the code you have shown. The posted myfun(D,Y) takes 2 arguments, but you are not passing a value for Y to myfun in any way. MATLAB should be complaining that Y is missing, unless your call to lsqnonlin looks something like this,
[D,resnorm,residual,exitflag,output,lambda]=lsqnonlin(@(D) myfun(D,Y), D0,[],[],options);
Also, if D is a scalar, it would be probably be much simpler to use fminsearch instead of lsqnonlin to minimize norm(diff).
17 Comments
KB
on 26 Jul 2014
Matt J
on 26 Jul 2014
It is running now
Now that what? What is your actual code?
Matt J
on 26 Jul 2014
No, this new version of the code still does not run and gives the error below. LevenbergMarquardt is indeed not an optimset option listed in the lsqnonlin documentation,
Error using optimset (line 215)
Unrecognized parameter name 'LevenbergMarquardt'. Please see the optimset reference page in the
documentation for a list of acceptable option parameters. Link to reference page.
KB
on 27 Jul 2014
Your objective function consists of exponential terms with awkwardly scaled time constants,
exp(-D*t(i)*(pi^2)*(n^2)/(a^2))
The quantities t(i)*(pi^2)*(n^2)/(a^2) are on the order of 1e+11, so the gradient will be numerically 0 except near D that are scaled very small. When I initialize at D0=0 , I reach a final D=1.4010e-13 after about 90 iterations, and the Resnorm moves from 1.9655 to 2.8453e-04. It looks like 'a' is meant to be much larger, though.
This was on R2013b, running with
options=optimset('Display','iter','MaxFunEvals',1e20,'TolFun',2e-50,'TolX',2e-50);
KB
on 27 Jul 2014
Matt J
on 27 Jul 2014
The finite difference calculations done by lsqnonlin aren't smart enough to know how small D needs to be. You could play with DiffMinChange, I suppose, but scaling the function like below also works. Note also that you can verify by plotting that you've found a valid minimum.
fun=@(D)diffusion(D,Y);
[Dsc,Resnorm]=lsqnonlin(@(Dsc) fun((1e-10)*Dsc),D0,0,[],options);
D=Dsc*1e-10
D_interval=linspace(0,2*D,1000);
plot(D_interval, arrayfun(@(D)norm(fun(D))^2,D_interval));
xlabel 'D'
ylabel 'norm(diff)^2'
KB
on 30 Jul 2014
Matt J
on 30 Jul 2014
Looks like your diffusion function is being confused with one of MATLAB's
You might want to rename yours, or at least move it higher on your path.
KB
on 30 Jul 2014
LSQNONLIN requires all values returned by user functions to be of data type double.
Make sure your function is returning type double.
lsqnonlin gives D around 8.9*10^-9 whereas scaling gives 2.8*106(-13).
Not for me. lsqnonlin gave me the 2.7*1e(-13) value.
What if it is totally unknown i.e. at what range the D value should fall.
As you can see from my code, I was able to find it by giving nothing nower than a positivity bound D>=0 to lsqnonlin. That's only because the function is unimodal, however. When you have local minima, you would need to be able to provide an accurate initial guess D0.
KB
on 30 Jul 2014
Your code, and its output, are difficult to read. You should use the
button to format it distinctly from your text (like my posted code appears). In any case, the code you've shown is not what I ran. I'm still applying lsqnonlin to the scaled version of the problem:
fun=@(D)diffusion(D,Y);
[Dsc,Resnorm]=lsqnonlin(@(Dsc) fun((1e-10)*Dsc),D0,0,[],options);
D=Dsc*1e-10;
Thanks Matt. You have been really awesome. Could you help me out with this data set:
and I am using similar code:
function diff = kani(D,Y)
t = 0:120:7200;
l=length(t);
a = 0.5*0.75*10^(-4);
for i = 1:l
Yp(i)=0;
for n = 1:10
y(i,n) = exp(-D*t(i)*(pi^2)*(n^2)/(a^2))/n^2;
Yp(i) = Yp(i) + y(i,n);
end
end
y;
Yp;
diff = 1- (6*Yp/pi^2)-Y;
Y = [ 0 0.0035 0.0066 0.0091 0.0121 0.0143 0.0150 0.0161 0.0167 0.0174 0.0188 0.0197 0.0199 0.0210 0.0218 0.0226 0.0229 0.0226 0.0237 0.0247 0.0244 0.0262 0.0249 0.0252 0.0248 0.0267 0.0285 0.0272 0.0279 0.0277 0.0292 0.0294 0.0289 0.0301 0.0282 0.0285 0.0299 0.0285 0.0304 0.0306 0.0309 0.0315 0.0317 0.0310 0.0310 0.0304 0.0339 0.0326 0.0320 0.0323 0.0335 0.0312 0.0332 0.0333 0.0316 0.0311 0.0315 0.0307 0.0313 0.0312 0.0328];
fun=@(D)kani(D,Y);
D0 = 0;
options = optimset('Display','iter','MaxFunEvals',1e20,'TolFun',2e-50,'TolX',2e-50);
[Dsc,Resnorm]=lsqnonlin(@(Dsc) fun((1e-10)*Dsc),D0,0,[],options);
D=Dsc*1e-10
Resnorm
D_interval=linspace(0,2*D,1000);
plot(D_interval, arrayfun(@(D)norm(fun(D))^2,D_interval));
xlabel 'D'
ylabel 'norm(diff)^2' }
and this is what I got:
I have tried to solve this by scaling the D value very low e^(-20), but still not getting it.
More Answers (0)
Categories
Find more on Get Started with Optimization Toolbox in Help Center and File Exchange
