Azzera filtri
Azzera filtri

ode45 change function index when something happends

7 visualizzazioni (ultimi 30 giorni)
I have a proble with one of my codes, I want to change index of a function when y value reach a given value of y
I will give an easy example:
A = 1;
B = 2;
tspan = [0 5];
y0 = [0 0.01];
[t,y] = ode45(@(t,y) odefcn(t,y,A,B), tspan, y0);
function dydt = Fun(t,y,A,B)
if y(1)>0.6
%Then A=3(example) change forever
end
dydt = zeros(2,1);
dydt(1) = y(2);
dydt(2) = (A/B)*t.*y(1);
end
%If you do this
function dydt = Fun(t,y,A,B)
if y(1)>0.6
A=3
end
dydt = zeros(2,1);
dydt(1) = y(2);
dydt(2) = (A/B)*t.*y(1);
end
%Then @ Command Window
A =
1
A =
3
A =
1
A =
3
%I want it to be all 3s.
  1 Commento
Ameer Hamza
Ameer Hamza il 24 Apr 2020
It appears that you didn't paste the odefcn you are trying to solve, but the thing printed in the command window is the correct behavior of ODE. It does not increase the variables monotonically forward. It needs to take forward and backward steps according to its internal algorithm. From the code you wrote, I will say that your current logic is correct.

Accedi per commentare.

Risposta accettata

darova
darova il 24 Apr 2020
Try persistent variable
function main
clc % clear command window
clear functions % clear persistent variable
A = 1;
B = 2;
tspan = [0 5];
y0 = [0 0.01];
[t,y] = ode45(@(t,y) Fun(t,y,A,B), tspan, y0);
function dydt = Fun(t,y,A,B)
persistent flag % declare 'flag' as persisten variable
if isempty(flag) % if 'flag' is empty define 'false'
flag = false;
end
if y(1)>0.6
flag = true; % if y(1) > 0.6 define flag 'true'
end
if flag == true
A = 3;
end
sprintf('y(1) = %f; A = %d',y(1),A)
dydt = zeros(2,1);
dydt(1) = y(2);
dydt(2) = (A/B)*t.*y(1);
See more about persistent
  3 Commenti
Ameer Hamza
Ameer Hamza il 24 Apr 2020
Using persistent in ODE is not correct. The ODE solver needs to take forward and backward steps while solving the ODE. Using persistent in this way will invalidate the ode45 algorithm, and the final answer will be wrong. The OP is confused because the ode45 solver moves forward and backward at the boundary of y(1)=0.6, and therefore in the command window, it prints A=1 and A=3 alternatively. But that is the correct behavior of ode45.
darova
darova il 24 Apr 2020
  • I want to change index of a function when y value reach a given value of y
He wanted A=3 after first y(1)>0.3
So A=3 even if y(1)<0.3 after condition y(1)>0.3 first time achieved

Accedi per commentare.

Più risposte (1)

Steven Lord
Steven Lord il 24 Apr 2020
Rather than changing your ODE function "on the fly" like this, you probably want to use a restart-based approach.
  1. Solve the system of ODEs until your condition to switch is satisfied. If the switch was time-based, use an appropriate tspan input to ode45. If it's based on the solution values use an events function.
  2. Use the final results from that call to ode45 to generate a new initial condition vector.
  3. Make whatever modifications you need to your ODE function.
  4. Call ode45 again with the new function, a new time span (starting where the previous call ended), and the new initial condition vector.
See the ballode example which does this using an events function to detect when the bouncing ball reaches the ground.

Tag

Community Treasure Hunt

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

Start Hunting!

Translated by