How to use if statement with ode23

I have a formula in the form of:

dy/dt = f - a*y
f(t) = x-k*t

I can use the ode23 function to solve this for both:

x = 2;
a = 3;
x1 = 4;
a1 = 5;
Tb = 0;
Te = 0.6;
Ti = 0.15;
F = @(t,y) (f*t)-(a*y);
[t,y] = ode23(F,[0 0.6],0);
F1 = @(t,y) (f*t2)-(k*y2);
[t2,y2] = ode23(F1,[0 0.6],0);

When f(t) is between min en max use formula 1, if f(t) is between max and min use formula 2.

Hope someone can help me out here.

[Walter Roberson]

Your F1 uses t2 and y2 before they are defined.

Answer 1

[Jan]

ode23 can handle smooth functions only. See Answers: Discontinuities in ODE45 . The only reliable solution is an event function, which stops the integration and restarts it with the current values and the new function. See https://www.mathworks.com/help/matlab/math/ode-event-location.html

F1  = @(t,y) (f*t)-(a*y);
F2  = @(t,y) (f*t2)-(k*y2);
opt = odeset('Events', @(t,y) DetectEvent(t, y, minValue, maxValue));
t   = Tb;
y0  = 0;
while t < Te
  if "y is inside the interval"   % (This is not clear in your description)
    F = F1;
  else
    F = F2;
  end
  [tt, yy] = ode23(F, [t, Tb], y0, opt);
  t  = tt(end);
  y0 = yy(end);
  ... some code to collect the yy values
end

[m4 Chrennikov]

Hi! When solving ode like this:

sol = ode45(___)

Is there any way to concat parts of solution(sol variable returned after each stop event) together?

[Walter Roberson]

[T1, Y1] = ode45(...);
[T2, Y2] = ode45(...);
T = [T1; T2];
Y = [Y1; Y2];

Note: if you restarted from exactly where the first one left off this will give you a duplicate point; you can eliminate that with

T = [T1; T2(2:end,:)];
Y = [Y1; Y2(2:end,:)];

[m4 Chrennikov]

My question is about sol struct.

There are two ways to get results:

[t,y,te,ye,ie] = ode45(odefun,tspan,y0,options)

and

sol = ode45(___)

It is clear how to concat parts of solution for first case but not for second one. As I understand by now there is no automatic way to concat sol struct, but maybe I haven't searched well enough.

[Walter Roberson]

sol1 = ode45(...);
sol2 = ode45(...);
sol.x = [sol1.x, sol2.x(2:end)];
sol.y = [sol1.y; sol2.y(2:end,:)];
sol.solver = sol1.solver;
xe = []; ye = []; ie = [];
if isfield(sol1, 'xe')
  xe = sol1.xe;
  ye = sol1.ye;
  ie = sol1.ie;
end
if isfield(sol2, 'xe')
  xe = [xe, sol2.xe];
  ye = [ye; sol2.ye];
  ie = [ie, sol2.ie];
end
if ~isempty(xe)   %sol only contains event info if events were detected
  sol.ie = ie;
  sol.ye = ye;
  sol.ie = ie;
end

Note: this assumes that you used the same event function list for both calls.

[m4 Chrennikov]

Hmmm, my fault. The only reason I want to use sol struct is that I use deval function to process ode solution. If you inspect its source code you can see that it utilizes not only data mentioned above but some 3d matrix sol.idata.f3d. My idea for next day was to try to concat this this matrix if there's no out-of-box solution. Hope I was clear enough)

[m4 Chrennikov]

works for me:

function sol1 = concatSol( sol1, sol2 )
%CONCATSOL Summary of this function goes here
%   Detailed explanation goes here
if ~strcmp(sol1.solver,'ode45')
    error('ode45 solver required')
end
if ~strcmp(sol1.solver,sol2.solver)
    error('sol1 & sol2 sovers are not the same')
end

sol1.x = [sol1.x, sol2.x(2:end)];
sol1.y = [sol1.y, sol2.y(:,2:end)];

sol1.xe = [sol1.xe, sol2.xe];
sol1.ye = [sol1.ye, sol2.ye];
sol1.ie = [sol1.ie, sol2.ie];

sol1.stats.nsteps = sol1.stats.nsteps + sol2.stats.nsteps;
sol1.stats.nfailed = sol1.stats.nfailed + sol2.stats.nfailed;
sol1.stats.nfevals = sol1.stats.nfevals + sol2.stats.nfevals;

sol1.idata.f3d = cat(3,sol1.idata.f3d,sol2.idata.f3d(:,:,2:end));
end