Question about ODE45

function Math462hw2Q5partcmodel
kG=0.002;
kAV=0.004;
k0=0.3;
k12=1;
k21=0.7;
T=1500;
tspan=0:1:168; %We have to do this for the rest intervals of same length from 0 to 168.
x1=0;
x2=0;
T=0;
miut=15;
y=[x1,x2,T];

[t,soln]=ode45(@(t,y)dotfunctions(t,y,k0,k12,k21,kG,kAV),tspan,y);

%Plot solution
plot(t,soln);
xlabel('time');
ylabel('value');
title('Progress of treatment');

    function ddt=dotfunctions(~,xt,k0,k12,k21,kG,kAV)
        x1=xt(1);
        x2=xt(2);
        T=xt(3);
        
        x1dot=miut-(k0+k12)*x1+k21*x2;
        x2dot=k12*x1-k21*x2;
        Tdot=kG*T-kAV*T*x1;
        ddt=[x1dot;x2dot;Tdot];
    end
end

For this code, how am I supposed to change the code if I only want miut to be 15 at t=0,24,48,72,96,120,144? Thanks!

[Walter Roberson]

What do you want it to be at other times?

Are you asking for an impulse equal to 15 at those particular times ?

[Zhukun Wang]

I want it to be zero at all other times Yeah, I want an impulse function.

Answer 1

[Walter Roberson]

The ode*() functions cannot handle impulses. More completely, they are only designed to handle objectives in which the second derivative of the function you supply is continuous. That rules out impulses, and also rules out most cases of abs() or min() or max() or if statements.

What you need to do is stop integrating and re-start from where you left off. That would look like

breakpoints = 0:24:24:144;
yinit = y0;
for K = 1 : length(breakpoints)
    [t{K}, y{K}] = ode45(fun, breakpoints(K):breakpoints(K+1), yinit);
    yinit = y{K}(end,:);
    yinit(1) = yinit(1) + miut;   %apply the impulse
end

Then afterwards, stitch together all of the t{:} and y{:}