ODE solver-using initial conditions in events function

I want to solve a second-order ODE (time derivatives on x) using numerical solver like ode45:

mx''+cx'+kx+F1=0

where F1 is a nonlinear function of x and it's initial value x(0). But when x becomes x(0)-5, I want F1 to stop functioning and switch the ODE to

mx''+cx'+kx+F2=0

where F2 is just a constant and I will simulate it until x' is zero.

I know how to set event functions to make the integration stop at x'=0, but how can I switch between F1 and F2? Event functions seem not to support using initial conditions in the event describing.

One background is that F1 and F2 can be merged. When x becomes x(0)-5 and further, F1 will be a complex number while it's real part is just F2, so I can solve one equation

mx''+cx'+kx+real(F1)=0

until x' is zero. But then the simulation will not give a fine resolution at the switching point of F1 and F2. So if anyone can suggest a way to make fine resolution near the switching point of F1 and F2, that also solves my issue. Thank you!

Answer 1

[Torsten]

If you want to use the event option, you can pass x(0) as an additional argument to the event function function and return control to the calling program when x(0)-5 is reached.

I suggest you just use mx''+cx'+kx+real(F1)=0 and choose tspan fine enough to resolve the switching point.

[Haocheng Yang]

Yes I understand the second way, thought it will lower the overall simulation speed as I don't need very fine steps outside the switching point.

I also tried to include x(0) in the event function but not successful. I converted that second order equations to two first order equations, so I define two initial conditions

x_initial=[x0 v0] and options=odeset('Events',@onedof_events)

The part running the solver is

[tout,xout]=ode45(EOM,[t0 T],x_initial,options)

and my event function is

function [value isterminal direction]=onedof_events(t,x,x_initial)

value=x_initial(1)-5-x(1);

isterminal=1;

direction=0;

end

But this gives me error like "Not enough input arguments." and the full messages are here

Would you have any suggestions on this issue?

[Torsten]

options=odeset('Events',@(t,x)onedof_events(t,x,x_initial)

instead of

options=odeset('Events',@onedof_events)

But you don't need to care how long the tspan vector is - the steps the solver takes are independent from how many elements tspan comprises.

[Haocheng Yang]

Great! Your code works for the force-switching event. And I figured out that for both methods, better resolution can be achieved by setting lower relative and absolute tolerance. Thank you for the help!