Event function and ODE.

46 views (last 30 days)
Arnab
Arnab on 20 Nov 2024 at 13:17
Commented: Arnab on 22 Nov 2024 at 13:34
I was solving an one dimensional ODE using ODE45.
I have two events in event function, where the execution of second event is dependent on the time instant of the first event.
How could I do that? Please help!
Elaborately...
1) Second event will comes into existence only after the execution of first event and is dependent on the time instant of first event.
for example..
Let first event instant be t=te. Then second instant will be at t=te+1.
  3 Comments
Arnab
Arnab on 21 Nov 2024 at 5:21
Thank you sir.
I am attaching the code below...which i tried
In the "events1" function I want second event occurs after some time at which first event occurs. That is , if "te" is the time aat which first event occurs then second event will occurs at time t=te+0.01.
Also I need y(te+0.01) as per defined in the else part i.e., y0 =z-0.5*ye(end);
Can I do such stuff or i have to do something else, kindly help. Thank you
function mtry2
y0 =2;
yplot = []; % store the solution
T=1; % constant
tp=[]; %collection of step size of time
tstar = 0;
z=[];
tau = 0.01; %delay
te=0.2652;
options1 = odeset('Events',@events1 );
while 1
tspan = linspace(tstar, tstar + 0.05, 50);
[t, y, te, ye, ie] = ode45(@event1, tspan, y0, options1); % solving the ode
%ie = intersection count of occuring of the events
%te = store the collection of intersecting oints
event1 =find(ie == 1);
event2 =find(ie == 2);
if isempty(event1) && isempty(event2)
% if isempty(ie) % Extend the interval
tstar = t(end);
y0 = y(end);
yplot = [yplot; y(:, 1)];
tp = [tp; t];
elseif ~isempty(event1) && isempty(event2)
tstar = te(end);
y0 = ye(end);
yplot = [yplot; y(:, 1)];
tp = [tp; t];
z=ye(end);
% event2 =find(ie == 2, 1);
% elseif isempty(event2)
% tstar = t(end); % Introduce the delay by adding tau
%
% y0 = y(end, :);
% % [t, y] = ode45(@event1, [tstar te(end)+tau], y0);
% yplot = [yplot; y(:, 1)];
% tp = [tp; t];
else
tstar=te(end);
y0 =z-0.5*ye(end);
yplot = [yplot; y(:, 1)];
tp = [tp; t];
end
% end
if (tstar >= 1)
break;
end
end
plot(tp,yplot,'b','LineWidth',2);
% hold on
% plot(tp,2*exp(-log(40)*tp/T));
grid on
ylabel('y(t)','FontSize',20,...
'FontWeight','bold','Color','k');
xlabel('t','FontSize',20,...
'FontWeight','bold','Color','k');
set(gca,'FontSize',25,'LineWidth',2.5);
grid on;
end
%==================================================
function dy = event1(~,y) % intermediate dunction to solve the ode
dy1=((y(1))^1.5)*sign(y(1)); % given ode
dy=dy1;
end
function [value,isterminal,direction] = events1(t,y) % when to stop the integration
% y is symbolic variable
% t is a vector
M=2;
tau = 0.01; %delay
te=0.2652;
ep=0.1;
T=1;
V=abs(y(1));
u=V/(1+V); %black curve
B=-log(M*M/ep)/T;
value =[u-M*exp(B*t);t-te-tau]; %blue curve
isterminal = [1;1]; %flag to stop the integration
direction = [0;0]; %dummy variable
end
Steven Lord
Steven Lord on 21 Nov 2024 at 22:54
Can you show the mathematical form of the ODE you're trying to solve, not the code? Show it and describe what it's computing in text, please.

Sign in to comment.

Accepted Answer

Torsten
Torsten on 20 Nov 2024 at 13:24
Moved: Torsten on 20 Nov 2024 at 13:29
Don't use both events simultaneously in your event function.
Start the solver with the first event function and integrate until the first event happens.
Return control to the calling program.
Restart the solver with the value of the solution variable obtained so far and the second event function (or up to te+1).
  7 Comments
Torsten
Torsten on 21 Nov 2024 at 22:20
Edited: Torsten on 22 Nov 2024 at 1:32
I'm not sure if this is what you want.
tfinish = 1.0;
tstart = 0.0;
treached = tstart;
tend = tfinish;
tspan = [tstart, tend];
y0 = 2;
options = odeset('Events',@event);
T = [];
Y = [];
while 1
[t, y] = ode45(@fun, tspan, y0, options); % solving the ode
T = [T;t];
Y = [Y;y];
treached = t(end);
if treached >= tfinish
break
end
tstart = treached;
tend = min(tfinish,treached + 0.01);
tspan = [tstart, tend];
y0 = y(end,1);
[t, y] = ode45(@fun, tspan, y0); % solving the ode
T = [T;t];
Y = [Y;y];
treached = t(end);
if treached >= tfinish
break
end
tstart = treached;
tend = tfinish;
tspan = [tstart,tend];
y0 = y0-0.5*y(end,1);
end
plot(T,Y)
function [value,isterminal,direction] = event(t,y) % when to stop the integration
% y is symbolic variable
% t is a vector
M=2;
ep=0.1;
T=1;
V=abs(y(1));
u=V/(1+V); %black curve
B=-log(M*M/ep)/T;
value =u-M*exp(B*t); %blue curve
isterminal = 1; %flag to stop the integration
direction = 0; %dummy variable
end
function dy = fun(~,y) % intermediate dunction to solve the ode
dy=((y(1))^1.5)*sign(y(1)); % given ode
end
Arnab
Arnab on 22 Nov 2024 at 13:34
yes, i wanted something like this.
Thank you.

Sign in to comment.

More Answers (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!