# How to implement the timestepping in TR-BDF2 using matlab ode solver as end-points (upper bound)?

Even though i know that i cannot advance the solver ode23tb, which uses TR-BDF2 itself, in just one timestep, i want to set an upper bound on the integration period after which I do something to output and feed the modified output into ode solver again. The upper bounds serve for limiting the timesteps taken in the integration interval tspan = [t_prev ti[j]], (see below), to be as low as possible, hopefully just one. To implement it, i try to do:

    y_prev = rho0;
t_prev = 0
ti = [list of timestep in  TR-BDF2]
j = 1;
while t_prev < tf
#do something to the output
...
# integrate up to adaptive ti in  TR-BDF2 even though inside the interval ode45 has its own timestep
tspan = [t_prev ti[j]]
j = j + 1;
[T, Y] = ode23tb(dydt, tspan, y_prev);
t_prev = T(end,:)
y_prev = Y(end,:)
....
end


But, how can I know what ti is beforehand? I look at https://math.la.asu.edu/~gardner/TRBDF.pdf , but do not have an idea on how to do this.

• It is unclear to me what you are trying to achieve. Why can't you do whatever operations you need to perform in your dydt function? Or alternatively, if you want to run ode23tb in a loop, why can't you simply break tf into some predetermined number of steps? If neither of those options is acceptable, why not do your own implementation of trapezoidal rule (apparently you don't want to take advantage of the automatic time step selection in ode23tb)? Dec 3, 2016 at 13:13
• @Bill Greene Yes neither of the options work for me. Ideally I should implement my own TR-BDF2 with adaptive timestep, which uses the same timestepping as ode23tb. In this way I do not have to specify the tspan. After each timestep, I do something to the output and use this to feed into TR-BDF2 and advance in another timestep.... Dec 3, 2016 at 23:42
• @Bill Greene However, the adaptive timestep version of TR-BDF2 is a lot difficult to implement. I wonder if you know where I can find this? Dec 3, 2016 at 23:43
• @Bill Greene In dydt function, the dydt function is evaluated at every increment. In a single timestep, there are several increments, for example, in Runge–Kutta methods. So the architecture here is: increment $\subset$ timestep $\subset$ timespan Dec 4, 2016 at 3:51
• You might be able to use the event mechanism of the ode solvers to do what you want: mathworks.com/help/matlab/math/ode-event-location.html Dec 4, 2016 at 15:44