Having Trouble Using For-Loop Over Distributed Range with ODE45
2 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
Thomas Veith
il 16 Lug 2019
Commentato: Thomas Veith
il 17 Lug 2019
Hi all,
So I'm trying to solve an ODE where the values for two of my inputs change over time.
For the first half of my tspan, delta1 should equal 1, and delta0 should equal 0. Then, in the second half, this gets reversed (i.e., delta1=0 and delta0=1). However, the way I've written my code, it only takes the second set of values. That is, it's solving the ODE as is delta1=0 and delta0=0 for the entire time - and I'm not sure why?
Here is my code:
% assign number of parameter sets and create parameters vector
n = 1000;
parameters = zeros(n,6);
% initialize parameters
for k=1:n
alpha0=rand;psa0=29.00;gamma=rand;psi=0.00;beta=rand;mse=0;
tspan=[0 2272];
parameters(k,:) = [alpha0,psa0,gamma,psi*gamma,beta,mse];
end
% solve ODE
for k = 1:size(parameters,1)
options=odeset('RelTol', 1e-8, 'AbsTol', 1e-8);
tspan=[0 2272];
for t=drange(0:1200)
delta1=1;delta0=0;
end
for t=drange(1201:2272)
delta1=0;delta0=1;
end
ODE = @(t,x,delta1,delta0,gamma,psi,beta) [delta1*-gamma*x(1) + delta0*psi*x(1); beta*x(2) - delta1* x(1)*x(2)];
sol = ode45(@(t,x) ODE(t,x,delta1,delta0,parameters(k,3),parameters(k,4),parameters(k,5)), tspan, [parameters(k,1) parameters(k,2)], options);
% calculate MSE and update parameters matrix
a = deval(sol,patientdata1(:,1));
a = a';
cost = (a(:,2)-patientdata1(:,2))./patientdata1(:,2);
mse = cost'*cost;
if min(a(:,2))<= 0.3
mse = mse * 100;
end
parameters(k,6) = mse;
end
Any help would be appreciated. Thanks in advance!
0 Commenti
Risposta accettata
Torsten
il 17 Lug 2019
Modificato: Torsten
il 17 Lug 2019
delta1 = @(t) (t<=1200)
delta0 = @(t) (t>1200)
ODE = @(t,x,gamma,psi,beta) [delta1(t)*(-gamma)*x(1) + delta0(t)*psi*x(1); beta*x(2) - delta1(t)* x(1)*x(2)];
But a better solution would be to integrate from t=0 to t=1200, store the solution, and call ODE45 for a second time from t=1200 to t=2272 with delta0 and delta1 interchanged.
5 Commenti
Più risposte (0)
Vedere anche
Categorie
Scopri di più su Ordinary Differential Equations in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!