How do I stop ode solver when output is complex number?

I am currently working on a coupled ode problem, and using ode45 solver to solve it.

My code is something like this:

bubble.m

function rdot = f(t, r)
%Some mathematical expressions
....
...
...
rdot(1) = r(2); % first derivative of r1
rdot(2) = (X11 - X21 - (X31*(X41 - X51)))/X81; % second derivative of r1
rdot(3) = r(4); % first derivative of r2
rdot(4) = (X12 - X22 - (X32*(X42 - X52)))/X82; % second derivative of r2
rdot = rdot';

And,

bubble_plotter.m

clc;
clear all;
%close all;
time_range = [0 1d-4];
r1_eq = 10d-6;
r2_eq = 1d-6;
initial_conditions = [r1_eq 0.d0 r2_eq 0.d0];
[t, r] = ode45('bubble', time_range, initial_conditions);
r1_us=1000000*r(:,1);
r1dot=r(:,2);
r2_us=1000000*r(:,3);
r2dot=r(:,4);

%% Recovery of double derivatives
for i=1:numel(t)
  t_actual = t(i);
  r_actual = r(i,:);
  rdot(i,:) = bubble_mettin(t_actual,r_actual);
end

r1ddot = rdot(:,2);
r2ddot = rdot(:,4);
figure(101);
plot(normalized_time, r1_us, 'b-', normalized_time, r2_us, 'k-.');

For certain input parameter values, ode solver successfully runs only for first few t-values and outputs corresponding real positive values for r, which in the code is r(1) and r(3). But after that the solver starts outputing negative or complex values (that is, imaginary numbers) for r(1) and r(3). I wish the program stops right before it starts yielding complex or negative values. Also the so far calculated values of r(1), r(2), r(3), r(4), rdot(2) and rdot(4) must be receoverable.

Another condition to stop the solver could be that r(1) + r(3) should not excced a given value.

[Reha Lalsab]

Did you find a solution to this?

[Vikash Pandey]

Yes, i solved it myself.

[Juan]

Would you please share how did you solve this issue?

[Umar Abdurrehman]

please share how you solved it

Answer 1

[madhan ravi]

see

https://www.mathworks.com/matlabcentral/answers/282449-how-to-stop-ode45-when-one-of-states-reach-certain-value#answer_220745

http://www.mathworks.com/help/matlab/ref/odeset.html#input_argument_d0e709348

use conditions as

isreal(r(1))==0 | r(1)<0 %an example isreal tells whether the number is real returns 0 if it's not

[Vikash Pandey]

Madhan, that looks quite complicated. Any other easier way?

[madhan ravi]

Vikash it's the only way , otherwise what you can do is compute everything and filter the datas after the calculation?

[Vikash Pandey]

Actually I did not understood what si there in that link? How do I implement that in my code?

[madhan ravi]

Opt    = odeset('Events', @myEvent); %add this before ode45
[t, r] = ode45('bubble', time_range, initial_conditions);

save the below function separately:
function [value, isterminal, direction] = myEvent(t, r)
value      = isreal(r(1))==0 | r(1)<0 | isreal(r(3))==0 | r(3)<0;
isterminal = 1;   % Stop the integration
direction  = 0;
end

[Vikash Pandey]

Let me check