Piecewise Plot for ODE

Hi,

I am trying to do a piecewise plot but I have no idea (or maybe a little but not enough) to do this. For a period over 6 seconds, the first four seconds plot one function and the next 2 seconds plot another function.

Below if my attempt.

 if true
    % code
  function rhs=pacemaker(t,E)
  R= 1;
  C= 1;
  rhs = -(1./(R*C))*E;

 function rhs=pacemakerRecharge(t,E)
 rhs = 0;

hold
% set axis 
axis([0 28 -1 15])
% define a set of initial conditions
Eo1=0;
Eo2 = 12;
% solve equation for different initial conditions
ts= [1:1:6]
if  (0<=ts<4)
[t,E] = ode45(@pacemakerRecharge,ts,Eo1);
  plot (t, E)
else if (4<=ts<6)
  [t,E] = ode45(@pacemaker,ts,Eo2);
  plot (t, E)
end
% label plots
title('Voltage across a heart in a pacemaker')
xlabel('t')
ylabel('Voltage')
end

Thanks. and apologies for the awful code.

Answer 1

[Star Strider]

See if this does what you want:

R= 1;
C= 1;
pacemaker = @(t,E) -(1./(R*C))*E;               % Anonymous Function Version
pacemakerRecharge = @(t,E) 0;

Eo1=0;
Eo2 = 12;
ts1= [1:0.1:4];
[t1,E1] = ode45(pacemakerRecharge,ts1,Eo1);
ts2 = [4.1:0.1:6];
[t2,E2] = ode45(pacemaker,ts2,Eo2);

T = [t1; t2];
E = [E1; E2];

figure(1)
plot(T, E)
grid
title('Voltage across a heart in a pacemaker')
xlabel('t')
ylabel('Voltage')

I re-wrote your two functions as Anonymous Functions (easier for me), then divided your time vector into two separate vectors, ‘ts1’ and ‘ts2’, integrated each function, then concatenated the resulting time and voltage vectors into the ‘T’ and ‘E’ vectors.

your code isn‘t ‘awful’, it’s just a bit inefficient. Your code will go better as you gain experience. You also did not need the if block, since all you needed to do was to divide the time vector into two parts, as I did.

[Hollie Bell]

Thanks heaps.

that works perfect. If I wanted t add another time period ts3, how would I change this:

if true
  % code
T = [t1; t2];
E = [E1; E2];
end

I tried

if true
  % code
T = [t1; t3];
E = [E1; E3];
end

but that didn't work for me, Thanks very much for your time.

[Star Strider]

My pleasure!

What about it didn’t work? That’s not enough information for me to help you. If it threw an error, please copy and paste the entire red text from the Command Window to a comment here, including the line that threw the error.

I don’t know what ‘t3’ and ‘E3’ are, or their sizes, so I can’t comment on them.

[Star Strider]

Hollie Bell’s Answer moved here...

Your code worked perfectly thanks.

Can you show me how I could add another time interval the same as the first one .. for instance like this:

if true
  % code
ts3= [6:0.1:9.9];
[t3,E3] = ode45(pacemakerRecharge,ts3,Eo1);
end

and when I've added that in, how do I change this bit of your code to reflect the additional time period:

 if true
  % code
 T = [t1; t2];
 E = [E1; E2];
 end

[Star Strider]

My pleasure.

If you want to go from 1 to 4 and then from 4 to 9.9, your ‘ts2’ vector becomes:

ts2 = [4.1:0.1:9.9];

and your ‘ts1’ vector remains the same. Since you have one variable, your initial value for ‘E2’ is a scalar, and can be whatever you want it to be.

You’re correct in describing your integration as ‘piecewise’, noting that the ‘pieces’ are contiguous and not overlapping. I suspect that was the problem with your ‘t3’, if you defined it as going from 1 to 9.9. (I don’t know those details, so that is simply a guess.)