Error in ode45 (Index exceeds matrix dimensions)
2 次查看(过去 30 天)
显示 更早的评论
Can someone please explain what I have to fix for the AF_function_in(t) of my code?
My function is as followed:
function dx=myeqn(t,x,p)
global tu;
function output=AF_function_in(t)
output=interp1(tu(:,1),tu(:,2),t);
end
dx(1,1)= (p(1)/p(2))*(AF_function_in(t)-x(1))- (p(4)/p(2))*(x(1)-x(2));
dx(2,1)= (p(4)/p(3)) *(x(1)-x(2));
end
When I run:
tv=1x225 matrix;
[t,x]=ode45(@(t,x) myeqn (t,x,[0.1 0.2 0.3 0.4],tv,[0 0]));
I get an error of:
Error in myeqn/AF_function_in (line 4)
output=interp1(tu(:,1),tu(:,2),t);
Error in myeqn (line 6)
dx(1,1)= (p(1)/p(2))*(AF_function_in(t)-x(1))- (p(4)/p(2))*(x(1)-x(2));
Error in @(t,x)myeqn(t,x,[0.1,0.2,0.3,0.4])
采纳的回答
Star Strider
2017-6-1
First, ‘tu’ has to be an (Nx2) matrix.
Second, do not use global variables! They make your code much more difficult to troubleshoot. Pass ‘p’ and ‘tu’ as extra parameters instead.
Try this:
function dx=myeqn(t,x,p,tu)
output=interp1(tu(:,1),tu(:,2),t);
dx(1,1)= (p(1)/p(2))*(AF_function_in(t)-x(1))- (p(4)/p(2))*(x(1)-x(2));
dx(2,1)= (p(4)/p(3)) *(x(1)-x(2));
end
and then call it in your ODE solver as:
[t,x] = ode45(@(t,x) myeqn(t,x,p,tu), tspan, x0);
with ‘x0’ being your (2x1) initial conditions vector.
18 个评论
Arbol
2017-6-1
编辑:Walter Roberson
2017-6-3
Thank you, I have been working too long on this and i get stuck ahha.
and my final code is as followed:
function dx=myeqn(t,x,p,tu)
output=interp1(tu(:,1),tu(:,2),t);
dx(1,1)= (p(1)/p(2))*(output-x(1))- (p(4)/p(2))*(x(1)-x(2));
dx(2,1)= (p(4)/p(3)) *(x(1)-x(2));
end
Arbol
2017-6-2
@Star Strider Also, I see that you understand how to use lsqnonlin or lsqcurvefit by following this: https://www.mathworks.com/matlabcentral/answers/254566-how-do-i-fit-coupled-differential-equations-to-experimental-data.
However, I couldn't able to get my code work. Hopefully you can enlighten me. I have tried all possible way. Using this code I made, the function stopped prematurely, and the First-order optimality converge to one value. I use lsqnonlin function to fit my model as followed:
Fp= 0.1; fp=0.1; fis=0.1; PS=0.1;
param0 = [Fp,fp,fis,PS];
fun=@(p) odefnc_1(p,tu(:,1),tu(:,3),tu);
%Then call the lsqnonlin using:
[newparam,resnorm,RESIDUAL,EXITFLAG,OUTPUT,LAMBDA,Jacobian]=...
lsqnonlin(fun,param0,lb,ub,options);
function diff= odefnc_1(p,texp,x,tu)
F=p(1);
fp=p(2);
fis=p(3);
PS=p(4);
init = [0 0];
[tv,y]=ode45(@myeqn, texp, init);
P=y(:,1);
I=y(:,2);
diff = (I+P)-x;
function dx=myeqn(t,x)
output=interp1(tu(:,1),tu(:,2),t);
dx(1,1)= (p(1)/p(2))*(output-x(1))- (p(4)/p(2))*(x(1)-x(2));
dx(2,1)= (p(4)/p(3)) *(x(1)-x(2));
end
end
Star Strider
2017-6-2
If lsqcurvefit converged to a local minimum (non-global solution), your code works. It simply did not converge on the optimal (global) minimum. If you know about what the optimal parameter set should be, use those as your starting values.
The Global Optimization Toolbox has functions such as patternsearch that can search a wider area of your parameter space, and may be able to help you find a global optimum. Another option is a genetic algorithm (the ga function). These are not fast, but they will often provide a global optimum if one exists.
Otherwise, keep experimenting with different initial parameter estimates until you get one that are close enough so that lsqcurvefit finds your global minimum.
Arbol
2017-6-3
it basically returns the initial parameter that put in. For example, if I have x0=[1 1 1 1], it will return [1 1 1 1]. I don't think it is going anywhere.
Star Strider
2017-6-3
I cannot determine what the problem is. Your code appears to be correct.
One possibility is to vectorize these lines, as:
dx(1,1) = (p(1)./p(2)) .* (output-x(1)) - (p(4)./p(2)).*(x(1)-x(2));
dx(2,1) = (p(4)./p(3)) .* (x(1)-x(2));
Arbol
2017-6-3
I have tried that, but no luck... I even tried lsqnonlin function, or fit function. Or played around with the options of Tolfun and TolX But no luck on my end. When I run examples of other people's work, their seem to work though @.@.
Walter Roberson
2017-6-3
I had to invent some numbers to run. I do not see the behaviour you describe:
tu = sort( rand(225,3) );
lb = zeros(1,4); ub = 10000 * ones(1,4);
options = [];
Fp= 0.1; fp=0.1; fis=0.1; PS=0.1;
param0 = [Fp,fp,fis,PS];
fun=@(p) odefnc_1(p,tu(:,1),tu(:,3),tu);
%Then call the lsqnonlin using:
[newparam,resnorm,RESIDUAL,EXITFLAG,OUTPUT,LAMBDA,Jacobian]=...
lsqnonlin(fun,param0,lb,ub,options);
which I used together with
function diff= odefnc_1(p,texp,x,tu)
F=p(1);
fp=p(2);
fis=p(3);
PS=p(4);
init = [0 0];
[tv,y]=ode45(@myeqn, texp, init);
P=y(:,1);
I=y(:,2);
diff = (I+P)-x;
function dx=myeqn(t,x)
output=interp1(tu(:,1),tu(:,2),t);
dx(1,1)= (p(1)/p(2))*(output-x(1))- (p(4)/p(2))*(x(1)-x(2));
dx(2,1)= (p(4)/p(3)) *(x(1)-x(2));
end
end
I notice, by the way, that pass four variables into odefnc_1, but you never use the last of them, tu
Arbol
2017-6-3
I have the same behavior when I run your code. I reedited my code as your request. As followed:
tu = sort( rand(225,3) );
lb = zeros(1,4); ub = 10000 * ones(1,4);
options = [];
Fp= 0.1; fp=0.1; fis=0.1; PS=0.1;
param0 = [Fp,fp,fis,PS];
fun=@(p) odefnc_test(p,tu(:,1),tu(:,3),tu);
%Then call the lsqnonlin using:
[newparam,resnorm,RESIDUAL,EXITFLAG,OUTPUT,LAMBDA,Jacobian]=...
lsqnonlin(fun,param0,lb,ub,options);
[tv,y]=ode45(@(t,x) myeqn(t,x,newparam,tu), tu(:,1), [0 0]);
P=y(:,1);
I=y(:,2);
plot(tu(:,1),tu(:,3),'o', tv,I+P);
I used it with:
function diff= odefnc_test(p,texp,x,tu)
F=p(1);
fp=p(2);
fis=p(3);
PS=p(4);
init = [0 0];
[tv,y]=ode45(@(t,x) myeqn(t,x,p,tu), texp, init);
P=y(:,1);
I=y(:,2);
diff = (I+P)-x;
function dx=myeqn(t,x,p,tu)
output=interp1(tu(:,1),tu(:,2),t);
dx(1,1)= (p(1)/p(2))*(output-x(1))- (p(4)/p(2))*(x(1)-x(2));
dx(2,1)= (p(4)/p(3)) *(x(1)-x(2));
end
end
Star Strider
2017-6-3
I just now noticed that you are calling the lsqnonlin function, not lsqcurvefit. The lsqnonlin function is appropriate here. The problem is that I do not see where you are comparing the results of your objective function to your dependent variable. (If you are, and I am not seeing it because I do not know what your dependent variable is, please post that function, and your call to it, specifically.)
If ‘y’ is your dependent variable, I would anticipate that your cost function (that you pass to lsqnonlin to optimise) would be something like this (for lsqnonlin):
fun=@(p) y - odefnc_test(p,tu(:,1),tu(:,3),tu);
Try something like that to see if you get the result you want.
Otherwise, it might be easier to use lsqcurvefit instead of lsqnonlin, since it incorporates the cost function in its code. As I mentioned previously, if your regression hypersurface has many local minima, lsqcurvefit may have a very difficult task in getting the ‘best’ fit global minimum. If that emerges as a problem, I would consider using a genetic algorithm (the ga function, although simple but effective ones are not difficult to write yourself) or the patternsearch function.
Arbol
2017-6-3
Sorry, I kind of messed up with the coding there. I rewrote my function for lsqnonlin as:
function diff= odefnc_1(p,texp,tu)
init = [0 0];
[tv,y]=ode45(@myeqn, texp, init);
P=y(:,1);
I=y(:,2);
diff = (I+P);
function dx=myeqn(t,x)
output=interp1(tu(:,1),tu(:,2),t);
dx(1,1)= (p(1)/p(2))*(output-x(1))- (p(4)/p(2))*(x(1)-x(2));
dx(2,1)= (p(4)/p(3)) *(x(1)-x(2));
end
end
And ran:
fun=@(p) odefnc_1(p,tu(:,1),tu)-tissData;
% %Then call the lsqnonlin using:
[newparam,resnorm,RESIDUAL,EXITFLAG,OUTPUT,LAMBDA,Jacobian]=...
lsqnonlin(fun,param0,lb,ub,options);
And got the same result, but 'newparam' is different than 'param0'. I also tried this for lsqcurvefit, and the same approach. LOL. Bad luck for me. Just one comment, my odefunction will return this attached graph
Where tissue is the result of 'param0', and Data is from the data. Do you think the big gap between the Data and Tissue is what cause this?
Star Strider
2017-6-3
I am not certain what I am seeing in the plot. If the parameters are now changing, that means the code is working as you want it to. The remaining issue seems to be a scaling problem. and that is a property of your differential equations. The gradient-descent algorithms may be getting ‘stuck’ in a very flat area of your regression space. Different initial parameter estimates could help. You will likely get a better fit with ga or patternsearch in this situation.
Arbol
2017-6-5
Hey @Star Strider, I got it to fit! I just had to pick out a good starting point (param0) for the fit to work! So that is something to note for. Have to pick out a good initial parameter!
更多回答(0 个)
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!发生错误
由于页面发生更改,无法完成操作。请重新加载页面以查看其更新后的状态。
您也可以从以下列表中选择网站:
如何获得最佳网站性能
选择中国网站(中文或英文)以获得最佳网站性能。其他 MathWorks 国家/地区网站并未针对您所在位置的访问进行优化。
美洲
- América Latina (Español)
- Canada (English)
- United States (English)
欧洲
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom(English)
亚太
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)