Conditional with time in ode45

Hi, how would you change a variable in a function when 't' in ODE45 reaches a certain value? I have already tried formatting it like the code below, but that doesn't work I have also tried using the if statement in the function but then it simply won't read the 't' value. So if I want to change upsilon after a certain time has passed in ODE45 how would I approach this?

clear all, close all, clc

mu = 6.25e10^-3;
eta = 6.25e10^-3;
beta = 0.03813350836;
gamma = 0.01165323939;
delta = 1/7;
alpha = 0.00513392178;


tspan = [1:1:500];
y0 = [0.99 0 0.01 0 0];
[t,y] = ode45(@(t,y)odefcn(y,mu,eta,upsilon,beta,gamma,delta,alpha),tspan,y0);

if t>300
    upsilon = 0.001;
else
    upsilon = 0;
end

plot(t,y)
legend('S','E','I','R','D')

function dydt = odefcn(y,mu,eta,upsilon,beta,gamma,delta,alpha)

dydt = zeros(5,1);
dydt(1) = mu - y(1)*(beta*y(3)+mu+upsilon);
dydt(2) = beta*y(1)*y(3)-y(2)*(mu+delta);
dydt(3) = delta*y(2)-y(3)*(mu+gamma+alpha);
dydt(4) = gamma*y(3)+upsilon*y(1)-eta*y(4);
dydt(5) = alpha*y(3);
end

Answer 1

[Alan Stevens]

Include upsilon within the odefcn

mu = 6.25e10^-3;
eta = 6.25e10^-3;
beta = 0.03813350836;
gamma = 0.01165323939;
delta = 1/7;
alpha = 0.00513392178;


tspan = [1:1:500];
y0 = [0.99 0 0.01 0 0];
[t,y] = ode45(@(t,y)odefcn(t, y,mu,eta,beta,gamma,delta,alpha),tspan,y0);


plot(t,y)
legend('S','E','I','R','D')

function dydt = odefcn(t,y,mu,eta,beta,gamma,delta,alpha)

if t>300
    upsilon = 0.001;
else
    upsilon = 0;
end


dydt = zeros(5,1);
dydt(1) = mu - y(1)*(beta*y(3)+mu+upsilon);
dydt(2) = beta*y(1)*y(3)-y(2)*(mu+delta);
dydt(3) = delta*y(2)-y(3)*(mu+gamma+alpha);
dydt(4) = gamma*y(3)+upsilon*y(1)-eta*y(4);
dydt(5) = alpha*y(3);
end

[Brilliant Purnawan]

I see, it works. Cheers!

Answer 2

[Steven Lord]

Solve the system from time t = 0 to time t = 300 with the first value of upsilon. Use the value of the solution at t = 300 as the initial condition for a second call to the ODE solver with the updated value of upsilon.

f = @(t, y, k) k*y;
[t1, y1] = ode45(@(t, y) f(t, y, 2), [0 2], 1);
[t2, y2] = ode45(@(t, y) f(t, y, 1), [2 4], y1(end));
plot(t1, y1, 'r-', t2, y2, 'k--')

If you don't have a fixed time when the value of the parameter should change (you instead have a condition you can check to decide when the parameter should change) use an Events function to detect the change. See the ballode example for an illustration of this technique.