Efficient way to perform elementary integration step with NDSolve internal method

Question

I'm trying to tweak the NDSolve function to perform one elementary integration step (using some explicitly selected stepping algorithm via Method option). Such a possibility is crucial for me, since after each performed step I'd like to check for some condition and modify state vector, then continue the integration. To the best of my knowledge such an operation is not possible with the 'pure' NDSolve, therefore I'm looking for an alternative solution.

So far, I managed to come up with the following solution (please note that the order of a method is selected to mimic the case with a lot of internal steps, this is equivalent to the higher order and longer integration period - tmax):

eq = {{x''[t] == -x[t], x[0] == 0, x'[0] == 1}, x, t};
tmax = 1;
method = {"ExplicitRungeKutta", "DifferenceOrder" -> 2};

Block[{state, options, timer},
  options = {MaxSteps -> Infinity, Method -> method};
  state = First@NDSolve`ProcessEquations[Sequence @@ eq, options];
  timer = First@Timing[NDSolve`Iterate[state, tmax]];
  {timer, state@"CurrentTime"["Forward"], state@"TimeStepsUsed"["Forward"], state@"SolutionVector"["Forward"]}
]

(* ==> {0.048003, 1., 7032, {0.841471, 0.540302}} *)


Block[{state, options, ctime, ic, dt, count, timer},
  options = {MaxSteps -> 1, Method -> method};
  state = First@NDSolve`ProcessEquations[Sequence @@ eq, options];
  count = 0;
  ctime = state@"CurrentTime"["Forward"];
  timer = First@Timing@While[ctime < tmax,
      Quiet[NDSolve`Iterate[state, tmax*1.000001], NDSolve::mxst];
      ctime = state@"CurrentTime"["Forward"];
      ic = state@"SolutionVector"["Forward"];
      dt = state@"TimeStep"["Forward"];
      state = First@NDSolve`Reinitialize[state, {x[ctime] == ic[[1]], x'[ctime] == ic[[2]]}, StartingStepSize -> dt];
      count++;
     ];
  {timer, state@"CurrentTime"["Forward"], count, state@"SolutionVector"["Forward"]}
]

(* ==> {3.1722, 1., 2948, {0.841471, 0.540302}} *)`

Timing results are not optimistic. The situation is worse when one want use step size from the previous step (with StartingStepSize -> dt commented in Reinitialize).

Is there a better (more efficient) way to to this?

Answer 1

Not sure if this does what you want, but you can use StepMonitor or EvaluationMonitor and the trick to make it stop is to do a Throw[] first thing. You can obtain the state of the solution in there also

Catch[NDSolve[{y'[t] == 0, y[0] == 0}, y[t], {t, 0, 1}, 
    StepMonitor :> {Print[Row[{"t=", t, " y[t]=", y[t]}]]; Throw["one step done"]}]]

this prints the time of the one step and the value of y is then

t=0.000102648 y[t]=0.