Change input at each time step of the ODE solver 'ode45'

17 vues (au cours des 30 derniers jours)
KC
KC le 5 Déc 2015
Commenté : Jan le 5 Mai 2017
I am not sure how to change an input parameter 'β' at each time step. My code is below - which gives me an error. Can anybody help please!
t = [7 14 21 28 35 42 49 56 63 70 77 84];
for i=1:12;
beta(i) = 0.43e-08 + (4.28e-08 - 0.43e-08)*exp(-0.20*t(i));
end
f = @(t,x) [3494-0.054*x(1)-beta*x(1)*x(3); beta*x(1)*x(3) - 0.41*x(2); ...
50000*x(2) - 23*x(3)];
[t,xa1] = ode45(f,t,[64700 0 0.0033],beta);
  1 commentaire
Jan
Jan le 5 Déc 2015
And the error message is:
Error using vertcat
Dimensions of matrices being concatenated are not consistent.
Error in @(t,x)[3494-0.054*x(1)-beta*x(1)*x(3);beta*x(1)*x(3)-0.41*x(2);50000*x(2)-23*x(3)]

Connectez-vous pour commenter.

Réponse acceptée

Jan
Jan le 6 Déc 2015
Please consider, that Matlab's ODE integrators cannot handle non-smooth functions sufficiently. See http://www.mathworks.com/matlabcentral/answers/59582#answer_72047 .
The only reliable method to run the integration is a loop over the time intervals:
function yourIntegration
tResult = [];
xResult = [];
tStep = [7 14 21 28 35 42 49 56 63 70 77 84];
x0 = [64700 0 0.0033];
for index = 2:numel(tStep)
% Integrate:
beta = 0.43e-08 + (4.28e-08 - 0.43e-08) * exp(-0.20*t(index - 1))
af = @(t,x) f(t, x, beta);
t = tStep(index-1:index);
[t, x] = ode45(af, t, x0);
% Collect the results:
tResult = cat(1, tResult, t);
xResult = cat(1, xResult, x);
% Final value of x is initial value for next step:
x0 = x(end, :);
end
function dx = f(t,x, beta)
dx = [3494-0.054*x(1)-beta*x(1)*x(3); ...
beta*x(1)*x(3) - 0.41*x(2); ...
50000*x(2) - 23*x(3)];
  7 commentaires
Saiprasad Gore
Saiprasad Gore le 5 Mai 2017
Thanks a lot, I had a similar problem. I wanted to switch the eqn depending on condition after every step. I hope this will work in my case too. Can you tell me how to give ode45 just 1 step without intermediate adaptive steps?
Jan
Jan le 5 Mai 2017
@Saiprasad Gore: This is not possible with ode45.

Connectez-vous pour commenter.

Plus de réponses (1)

Walter Roberson
Walter Roberson le 6 Déc 2015
f = @(T,x) [3494-0.054*x(1)-interp1(t,beta,T,'linear','extrap')*x(1)*x(3); interp1(t,beta,T,'linear','extrap')*x(1)*x(3) - 0.41*x(2); ...
50000*x(2) - 23*x(3)];
  2 commentaires
KC
KC le 12 Déc 2015
Thanks Walter!
sam
sam le 15 Juin 2016
Modifié(e) : sam le 16 Juin 2016
@Walter Roberson
Hi Walter,
Why do we have to do interpolation if we already know the exact expression of the variables? Couldnt we just input the exact expression of the variables into the Matlab ode45 solver? If we could, could you kindly tell me how to do this? Thanks.

Connectez-vous pour commenter.

Tags

Community Treasure Hunt

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

Start Hunting!

Translated by