concatenation of solutions of ODE

I want to plot the solutions of simple SIR model under two time intervals (tspan1 and tspan2). After geting the solution of the model within this time intervals, I want to concatenate and plote the solution curve. But, it is not runing due to the error ''imensions of arrays being concatenated are not consistent''. I am wondering any one who can help me by fixing this error. This is my matlab code!

Thank you

function SS

tspan1=[0 10];

t2f=20;

tspan2=[10 20];

S0=1000;

I0=1;

R0=0;

Initial0=[S0, I0, R0];

%%% ODE

function dx = tem(t, x)

dx = zeros(3,1);

beta = 0.003;

delta = 1;

dx(1) = -beta * x(1) * x(2);

dx(2) = beta * x(1) * x(2) - delta * x(2);

dx(3) = delta * x(2);

end

[t1,xx1] = ode45(@(t,x)tem(t,x),tspan1, Initial0);

S1=xx1(:,1);

I1=xx1(:,2);

R1=xx1(:,3);

Initial1=[S1(t2f),I1(t2f),R1(t2f)];

[t2,x2] = ode45(@(t,x)tem(t,x),tspan2, Initial1);

S2=x2(:,1);

I2=x2(:,2);

R2=x2(:,3);

T=[t1,t2];

S = [S1; S2];

I = [I1; I2];

R = [R1; R2];

%plot

plot(T,S,'k', 'Linewidth', 1.75);

% plot(t,S,'k',t,I,'r', t,R,'c', 'Linewidth', 1.75);

% xlabel('Time')

% ylabel('Population')

% title('Dynamics of SIR model')

% legend('S', 'I', 'R');

end

Answer 1

[Star Strider]

I would create one ‘tspan1’ vector with a specific number of elemensts, and the same for ‘tspan2’ although they are only required to have the same numbers of elements if you want to horizontally concatenate them. Since they are column vectors, they can have different numbers of elements if you want to vertically concatenate them. The problem is that unless you specify the elements of the ‘tspan’ argument, the MATLAB ODE integrators will produce time vectors of whatever lengths they believe is appropriate, and in a situation such as this, they might not always be the same.

One approach —

tspan1 = linspace(0, 10, 25);                           % Use 'linspace' To Generate 'tspan'
t2f=20;
tspan2 = linspace(10, 20, 25);

    S0=1000;
    I0=1;
    R0=0;
Initial0=[S0, I0, R0];

        [t1,xx1] = ode45(@(t,x)tem(t,x),tspan1, Initial0);
        S1=xx1(:,1);
        I1=xx1(:,2);
        R1=xx1(:,3);
      Initial1=[S1(t2f),I1(t2f),R1(t2f)];
       [t2,x2] = ode45(@(t,x)tem(t,x),tspan2, Initial1);
        S2=x2(:,1);
        I2=x2(:,2);
        R2=x2(:,3);
    T=[t1; t2];
     S = [S1; S2];
     I = [I1; I2];
     R = [R1; R2];
        %plot
        figure
 plot(T,S,'k', 'Linewidth', 1.75);
 hold on
 plot(T,I,'c', 'Linewidth', 1.75);
 plot(T,R,'m', 'Linewidth', 1.75);
 hold off
    % plot(t,S,'k',t,I,'r', t,R,'c', 'Linewidth', 1.75);
    % xlabel('Time')
    % ylabel('Population')
    % title('Dynamics of SIR model')
    legend('S', 'I', 'R', 'Location','best');

    % end

    
function SS
    tspan1=[0 10];
    t2f=20;
    tspan2=[10 20];
    S0=1000;
    I0=1;
    R0=0;
    Initial0=[S0, I0, R0];
end
%%% ODE
function dx = tem(t, x)
  dx = zeros(3,1);
        beta = 0.003; 
        delta = 1;
  dx(1) = -beta * x(1) * x(2);
  dx(2) = beta * x(1) * x(2) - delta * x(2);
  dx(3) = delta * x(2);
end 

Make appropriate changes to get the result you want.

EDIT — (7 Aug 2023 at 00:24)

Corrected typographical errors, clarified discussion.

.

[Temesgen]

Dear Star Strider , thank you very much for your help. I want to sketch the solution for one of the state variables for example for infected (I) in the time interval tspan1. And then vary the parameter values and sketch the output in the same plot that will continue from the first plot (that means the plot after some changes will be plotted in tspan2). My model is so huge and complicated so that i am using sample model. If I can manage this i will customise into my own.

So could you update by considering this?

Thank you!

[Star Strider]

I am not certain what you are actually doing. If your function has parameters that you want to change (your SIR model does not), your current approach is correct in stopping the integration, saving the last integrated values to use as initial conditions for the next integration. The only change you need to do is to change the parameters for your ODE function aftter the end of the first integration (or any subsequent integration if you want to do this more than once).

Sketching that out —

tspan = [0 10];

p = [0.003 1;                                                   % Parameter Matrix
     0.006 2;
     0.009 3];

S0=1000;
I0=1;
R0=0;
Initial0=[S0, I0, R0];

for k = 1:size(p,1)
    [t,x]  = ode45(@(t,x)tem(t,x,p(k,:)),tspan, Initial0);
    tc{k,:}= t;                                                 % Store Results As Cell Arrays (For Convenioence)
    xc{k,:} = x;
    Initial0 = x(end,:);
    tspan = tspan+t(end);
end

tv = cell2mat(tc);                                              % Conovert To Numeric Arrays For Plotting
xm = cell2mat(xc);

figure
plot(tv, xm)
grid
xlabel('Time')
ylabel('Population')
title('Dynamics of SIR model')
legend('S', 'I', 'R', 'Location','best');

function dx = tem(t, x, p)
dx = zeros(3,1);
beta = p(1); 
delta = p(2);
dx(1) = -beta * x(1) * x(2);
dx(2) = beta * x(1) * x(2) - delta * x(2);
dx(3) = delta * x(2);
end 

I am certain you can follow this, however if you have any problems, I will do my best to resolve them.

.

[Temesgen]

Thank you so much!!!

[Star Strider]

My pleasure!

If my Answer helped you solve your problem, please Accept it!

.