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,:)          ....    endBut, 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.
How to implement the timestepping in TR-BDF2 using matlab ode solver as end-points (upper bound)?
matlab