ode45 error Unable to perform assignment because the indices on the left side are not compatible with the size of the right side

4 次查看（过去 30 天）
Christian Gärtner
Christian Gärtner 2020-3-19
评论： Christian Gärtner 2020-3-26
Hi Guys,
I am trying to solve a set of 3 differential equations using ode45. But i can't get it to work because I am getting the following error:
%Unable to perform assignment because the indices on the left side are not compatible with the size of the right side.
%Error in freylaw_stat (line 29)
%dy(1,1)=-Ct+(R+r)*y(3); %Mr the constant r is introduced here according to Looft et al.(2018) and is an intended variation of the
%original equation
%Error in Gesamtskript_stat>@(t,y)freylaw_stat(t,y,t_min,Fext_M_P,MVC_Verlauf_P) (line 119)
%[t,y]=ode45(@(t,y) freylaw_stat(t,y,t_min,Fext_M_P,MVC_Verlauf_P),timerange,IC); %PROZENT MUs
%Error in odearguments (line 90)
%f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
%Error in ode45 (line 115)
%odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
%Error in Gesamtskript_stat (line 119)
%[t,y]=ode45(@(t,y) freylaw_stat(t,y,t_min,Fext_M_P,MVC_Verlauf_P),timerange,IC); %PROZENT MUs
following is the code I am using for ode45:
timerange=[T0_min tEnd_min];
IC=[100; 0; 0];
MVC_Verlauf_P=100;
[t,y]=ode45(@(t,y) freylaw_stat(t,y,t_min,Fext_M_P,MVC_Verlauf_P),timerange,IC);
with Fext_M_P being of the type 741979x1 double and using the following function:
function dy=freylaw_stat(t,y,t_min,Fext_M_P,MVC_Verlauf)
L=10; %Modellspezifischer Faktor
r=15; %rest-recovery-multiplier nach 2018#Looft et al.
F=0.00912; %Ermüdungsfaktor nach Frey Law (2012)
R=0.00094; %Erholungsfaktor nach Frey Law (2012)
Fext_I=Fext_M_P;
MVC=MVC_Verlauf-y(3);
TL0=(Fext_I/MVC)*100;
if (y(2)<TL0)&(y(1)>(TL0-y(2)))
Ct=L*(TL0-y(2)); %if Ma<TL and Mr>(TL-Ma)
elseif (y(2)<TL0)&(y(1)<(TL0-y(2)))
Ct=L*y(1); %if Ma<TL and Mr<(TL-Ma)
else
Ct=L*(TL0-y(2)); %if Ma>=TL
end
dy(1,1)=-Ct+(R+r)*y(3); %Mr
dy(2,1)=Ct-F*y(2); %Ma
dy(3,1)=F*y(2)-R*y(3); %Mf
end
Can someone tell me what I am doing wrong here, and how to change it? Your help is greatly appreciated!
Thank you and best regards,
Christian
  8 个评论
显示 6更早的评论
Ameer Hamza
Ameer Hamza 2020-3-20
To solve them as a differential equation with ode45, you need dy to be 3*1, not 741979*3. The equation you pasted also show that C(t) is a function of time. So you need to specify, how the time t is related to the elements of vector Ct. For example, if t=0, C(0) is the first element of Ct, i.e., Ct(1), then at the next time step, say t=0.1, which elemnt of Ct is equal to the value of C(0.1).
Christian Gärtner
Christian Gärtner 2020-3-23
thanks for your answer Ameer!
okay, I did not know that. So I would need to use the following code in the function?
dy(1,:)=-Ct+(R+r)*y(3); %Mr the constant r is introduced here according to Looft et al.(2018) and is an intended variation of the original equation
dy(2,:)=Ct-F*y(2); %Ma
dy(3,:)=F*y(2)-R*y(3); %Mf
And for C(t): it has the same length as the time vector so the corresponding value to t=0 is the first element of C. In the next time step the second element of C is the corresponding value, and so on. C(t) is calculated with the following equation:
There LD/LR=10=const, TL is a vector with 741979 elements, where the first cell corresponds to t=0, the second cell corresponds to the next time step and so on. Ma and Mr are the values which are calculated with this differential equation at the corresponding time steps.
I can also send you the whole file with the vectors, if that makes things easier for you.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

Ameer Hamza
Ameer Hamza 2020-3-24
Check the following code. I created a random vector Fext_M_P to test the output of my code. Also the time range is from 0 to 10.
Fext_M_P = rand(741979, 1);
timerange=linspace(0, 10, numel(Fext_M_P));
IC=[100; 0; 0];
MVC_Verlauf_P=100;
[t,y]=ode45(@(t,y) freylaw_stat(t,y,timerange,Fext_M_P,MVC_Verlauf_P),timerange,IC);
function dy=freylaw_stat(t,y,timerange,Fext_M_P,MVC_Verlauf)
L=10; %Modellspezifischer Faktor
r=15; %rest-recovery-multiplier nach 2018#Looft et al.
F=0.00912; %Ermüdungsfaktor nach Frey Law (2012)
R=0.00094; %Erholungsfaktor nach Frey Law (2012)
index = floor(interp1(timerange([1 end]), [1 numel(Fext_M_P)], t));
Fext_I=Fext_M_P(index);
MVC=MVC_Verlauf-y(3);
TL0=(Fext_I/MVC)*100;
if (y(2)<TL0)&(y(1)>(TL0-y(2)))
Ct=L*(TL0-y(2)); %if Ma<TL and Mr>(TL-Ma)
elseif (y(2)<TL0)&(y(1)<(TL0-y(2)))
Ct=L*y(1); %if Ma<TL and Mr<(TL-Ma)
else
Ct=L*(TL0-y(2)); %if Ma>=TL
end
dy(1,1)=-Ct+(R+r)*y(3); %Mr
dy(2,1)=Ct-F*y(2); %Ma
dy(3,1)=F*y(2)-R*y(3); %Mf
end
  1 个评论
Christian Gärtner
Christian Gärtner 2020-3-26
Hello there, sorry for the late answer, I am quite busy lately.
Thank you a lot for your answer, it worked great to get a result for the ODE-system!
Still, I was expecting another result, for dy(3,1). It is not behaving as expected, although dy(1,1) and dy(2,1) are perfectly reasonable results.
I have attached a picture of the plotted results where MR is dy(1,1); MA is dy(2,1) and MF is dy(3,1). Also I have drawn a yellow line, as dy(3,1) is expected to behave.
Maybe you have got any idea, what this different behaviour is related to? I will check units and stuff, but I think that they are right.
What I think may be the issue is that MVC is not decreasing over the timerange? Actually it should decrease in the form of:
MVC=MVC_Verlauf-y(3)
But I am unsure how to check MVC, because I do not know how to take a look at it.
Thank you!

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Calculus 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by