Event function and ODE.

5 次查看(过去 30 天)
Arnab
Arnab 2024-11-20
评论: Arnab 2024-11-22
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 个评论
Arnab
Arnab 2024-11-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 2024-11-21
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.

请先登录,再进行评论。

采纳的回答

Torsten
Torsten 2024-11-20
移动:Torsten 2024-11-20
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 个评论
Torsten
Torsten 2024-11-21
编辑:Torsten 2024-11-22
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 2024-11-22
yes, i wanted something like this.
Thank you.

请先登录,再进行评论。

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Ordinary Differential Equations 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by